svm 核函数

#coding=utf-8
import re
from numpy import *
def load_data():
data=[];label=[]
o=open('test.txt')
for line in o.readlines():
line_arr=re.split(r'(-\d*|\d*)',line.strip())
data.append([1.0,float(int(line_arr[1])),float(int(line_arr[3]))])
label.append(int(line_arr[5]))
return data,label
def kernel_trans(x,a,k_tup):
m,n=shape(x)
k=mat(zeros((m,1)))
if k_tup[0]=='lin':k=x*a.T
elif k_tup[0]=='rbf':
for j in range(m):
delta_row=x[j,:]-a
k[j]=delta_row*delta_row.T
k=exp(k/(-1*k_tup[1]**2))
else:
raise NameError('Houston we have a problem that kernel is not recognized')
return k

class self_var:
"""函数间传递变量用，省得定义全局变量"""
def __init__(self,data_matrix,label_matrix,c,toler,k_tup):
self.x=data_matrix
self.y=label_matrix
self.c=c
self.toler=toler
self.m,self.n=shape(data_matrix)
self.alphas=mat(zeros((self.m,1)))
self.b=0
self.ecache=mat(zeros((self.m,2)))
self.k=mat(zeros((self.m,self.m)))
for i in range(self.m):
self.k[:,i]=kernel_trans(self.x,self.x[i,:],k_tup)

def select_jrand(i,m):
j=i
while j==i:
j=int(random.uniform(0,m))
return j
def clip_alpha(aj,h,l):
if aj>h:
aj=h
if l>aj:
aj=l
return aj
def clac_ek(sv,k):
fxk=float(multiply(sv.alphas,sv.y).T*sv.k[:,k]+sv.b)
ek=fxk-float(sv.y[k])
return ek
def select_j(i,sv,ei):
max_k=-1;max_delta_e=0;ej=0
sv.ecache[i]=[1,ei]
valid_ecache_list=nonzero(sv.ecache[:,0].A)[0]
if len(valid_ecache_list)>1:
for k in valid_ecache_list:
if k==i:continue
ek=clac_ek(sv,k)
delta_e=abs(ei-ek)
if delta_e>max_delta_e:
max_k=k;max_delta_e=delta_e;ej=ek
return max_k,ej
else:
j=select_jrand(i,sv.m)
ej=clac_ek(sv,j)
return j,ej
def update_ek(sv,k):
ek=clac_ek(sv,k)
sv.ecache[k]=[1,ek]
def inner_loop(i,sv):
ei=clac_ek(sv,i)
if ((sv.y[i]*ei<-sv.toler)and(sv.alphas[i]<sv.c)) or \
((sv.y[i]*ei>sv.toler)and(sv.alphas[i]>0)):
j,ej=select_j(i,sv,ei)
alpha_iold=sv.alphas[i].copy()
alpha_jold=sv.alphas[j].copy()
if sv.y[i] != sv.y[j]:
l=max(0,sv.alphas[j]-sv.alphas[i])
h=min(sv.c,sv.c+sv.alphas[j]-sv.alphas[i])
else:
l=max(0,sv.alphas[j]+sv.alphas[i]-sv.c)
h=min(sv.c,sv.alphas[j]+sv.alphas[i])
if l==h:print 'l==h';return 0
eta=2.0*sv.k[i,j]-sv.k[i,i]-sv.k[j,j]
if eta>=0:print 'eta>=0';return 0
sv.alphas[j]-=sv.y[j]*(ei-ej)/eta
sv.alphas[j]=clip_alpha(sv.alphas[j],h,l)
update_ek(sv,j)
if abs(sv.alphas[j]-alpha_jold)<0.0001:
print 'j not moving enough';return 0
sv.alphas[i]+=sv.y[j]*sv.y[i]*(alpha_jold-sv.alphas[j])
update_ek(sv,i)
b1=sv.b-ei-sv.y[i]*(sv.alphas[i]-alpha_iold)*sv.k[i,i]-sv.y[j]*\
(sv.alphas[j]-alpha_jold)*sv.k[i,j]
b2=sv.b-ej-sv.y[i]*(sv.alphas[i]-alpha_iold)*sv.k[i,j]-sv.y[j]*\
(sv.alphas[j]-alpha_jold)*sv.k[j,j]
if 0<sv.alphas[i] and sv.c>sv.alphas[i]:sv.b=b1
elif 0<sv.alphas[j] and sv.c>sv.alphas[j]:sv.b=b2
else:sv.b=(b1+b2)/2.0
return 1
else:
return 0
def smop(data_matrix,label_matrix,c,toler,max_iter,k_tup=('lin',0)):
sv=self_var(data_matrix,label_matrix.T,c,toler,k_tup)
ite=0
entire_set=True;alpha_pairs_changed=0
while ite<max_iter and (alpha_pairs_changed>0 or entire_set):
alpha_pairs_changed=0
if entire_set:
for i in range(sv.m):
alpha_pairs_changed+=inner_loop(i,sv)
print 'fullset,iter:%d i:%d,pairs changed %d'%(ite,i,alpha_pairs_changed)
ite+=1
else:
non_boundis=nonzero((sv.alphas.A>0)*(sv.alphas.A<c))[0]
for i in non_boundis:
alpha_pairs_changed+=inner_loop(i,sv)
print 'non-bound,iter:%d i:%d,paris changed %d'%(ite,i,alpha_pairs_changed)
ite+=1
if entire_set:entire_set=False
elif alpha_pairs_changed==0:entire_set=True
print 'iteration number: %d'%ite
return sv.b,sv.alphas
def test(k1=1.3):
data,label=load_data()
data_matrix=mat(data);label_matrix=mat(label)
m,n=shape(data_matrix)
b,alphas=smop(data_matrix,label_matrix,200,0.0001,10000,('rbf',k1))
sv_ind=nonzero(alphas.A>0)[0]
svs=data_matrix[sv_ind]
label_sv=label_matrix.T[sv_ind]
print 'there are %d support vectors'%shape(svs)[0]
error_count=0
for i in range(m):
kernel_eval=kernel_trans(svs,data_matrix[i,:],('rbf',k1))
predict=kernel_eval.T*multiply(label_sv,alphas[sv_ind])+b
if sign(predict)!=sign(label[i]):error_count+=1
print 'the training error rate is : %f'%(float(error_count)/m)

test.txt:

[2, 4, '1']
[5, 8, '1']
[8, 5, '1']
[9, 9, '1']
[4, 1, '1']
[-1, -1, '1']
[5, 8, '1']
[9, 3, '1']
[1, 8, '1']
[7, 3, '1']
[9, 3, '1']
[3, 8, '1']
[4, 6, '1']
[9, 7, '1']
[7, 1, '1']
[5, 2, '1']
[9, 6, '1']
[6, 9, '1']
[9, 8, '1']
[7, -1, '1']
[4, 5, '1']
[9, 8, '1']
[-1, 4, '1']
[4, 3, '1']
[6, -1, '1']
[9, 9, '1']
[-1, 3, '1']
[9, 8, '1']
[1, 7, '1']
[5, 8, '1']
[7, 8, '1']
[1, 5, '1']
[-1, 7, '1']
[1, 9, '1']
[7, 8, '1']
[2, 5, '1']
[7, 4, '1']
[2, 1, '1']
[6, 1, '1']
[-1, 1, '1']
[2, 4, '1']
[6, -1, '1']
[8, -1, '1']
[4, 9, '1']
[8, 3, '1']
[9, 8, '1']
[8, 9, '1']
[5, 9, '1']
[9, 6, '1']
[4, 2, '1']
[8, 7, '1']
[1, 9, '1']
[3, 8, '1']
[-1, 1, '1']
[1, 1, '1']
[-1, 9, '1']
[-1, 6, '1']
[1, 5, '1']
[2, 6, '1']
[9, 5, '1']
[5, -1, '1']
[2, 4, '1']
[5, 9, '1']
[9, 5, '1']
[6, 3, '1']
[9, 3, '1']
[3, 6, '1']
[8, 6, '1']
[7, 7, '1']
[-1, -1, '1']
[5, 4, '1']
[2, 9, '1']
[5, 7, '1']
[3, 9, '1']
[6, 9, '1']
[8, 2, '1']
[8, 3, '1']
[8, -1, '1']
[2, 4, '1']
[9, 2, '1']
[-1, 3, '1']
[6, 8, '1']
[5, 4, '1']
[5, -1, '1']
[5, 3, '1']
[7, 6, '1']
[-1, 4, '1']
[3, 9, '1']
[7, 5, '1']
[8, 3, '1']
[9, 7, '1']
[8, 3, '1']
[3, 5, '1']
[2, 6, '1']
[1, 9, '1']
[6, 2, '1']
[3, 5, '1']
[9, 7, '1']
[5, 6, '1']
[7, 2, '1']
[11, -1, '-1']
[13, 7, '-1']
[16, 5, '-1']
[11, -1, '-1']
[17, 6, '-1']
[9, 5, '-1']
[15, 1, '-1']
[13, 7, '-1']
[12, 6, '-1']
[9, 5, '-1']
[17, 4, '-1']
[10, 8, '-1']
[10, 8, '-1']
[9, 5, '-1']
[13, 2, '-1']
[13, 6, '-1']
[9, -1, '-1']
[11, 9, '-1']
[17, 2, '-1']
[9, 7, '-1']
[16, 4, '-1']
[12, 1, '-1']
[11, 8, '-1']
[10, 1, '-1']
[17, 7, '-1']
[12, -1, '-1']
[16, 5, '-1']
[18, 2, '-1']
[15, 6, '-1']
[9, 5, '-1']
[13, 4, '-1']
[13, 2, '-1']
[11, 2, '-1']
[17, 7, '-1']
[16, 1, '-1']
[15, -1, '-1']
[9, 4, '-1']
[16, 7, '-1']
[13, 1, '-1']
[17, -1, '-1']
[18, 4, '-1']
[12, 3, '-1']
[11, 7, '-1']
[14, 6, '-1']
[9, 5, '-1']
[11, 9, '-1']
[12, 4, '-1']
[17, 8, '-1']
[11, 2, '-1']
[12, 5, '-1']
[13, -1, '-1']
[12, 2, '-1']
[11, 1, '-1']
[14, 1, '-1']
[17, -1, '-1']
[18, 3, '-1']
[11, 5, '-1']
[18, 2, '-1']
[12, 4, '-1']
[15, 8, '-1']
[17, 9, '-1']
[18, 5, '-1']
[14, 9, '-1']
[16, 9, '-1']
[18, 5, '-1']
[9, 1, '-1']
[14, 4, '-1']
[13, 2, '-1']
[12, 9, '-1']
[16, 8, '-1']
[15, 4, '-1']
[12, -1, '-1']
[16, 9, '-1']
[14, 3, '-1']
[12, 9, '-1']
[17, 5, '-1']
[11, 4, '-1']
[13, 6, '-1']
[16, 3, '-1']
[16, 2, '-1']
[11, 3, '-1']
[11, 1, '-1']
[17, 9, '-1']
[18, 2, '-1']
[11, 8, '-1']
[14, 3, '-1']
[11, -1, '-1']
[18, 6, '-1']
[12, 6, '-1']
[11, -1, '-1']
[14, -1, '-1']
[16, 5, '-1']
[12, 7, '-1']
[15, -1, '-1']
[15, 1, '-1']
[18, 9, '-1']
[9, -1, '-1']
[18, -1, '-1']
[18, 6, '-1']
[9, 3, '-1']