本文采用smo算法计算svm
程序有点问题,开始才用的libsvm的代码,准备将其java代码写成python的,后面发现用libsvm的数据格式老是出问题。就参考了
机器学习实战的代码。
程序有很多要优化的地方
1)核函数要完善,这里只写了线性核函数。但是整个程序中没有用核函数进行计算。
2)一些异常状况的处理。
整个迭代公式可以参考
http://blog.csdn.net/macyang/article/details/38782399/
个人觉得非常棒,就是后面的smo要各种计算,推导。
其实最后迭代也是比较简单的:
1)找出误差ei
2)找出2k(i,j)-k(i,i)-k(j,j)
3)在ai的约束条件下,更新ai
4)更新aj
5)更新b
6)返回结果
#coding=utf-8 #!/usr/bin/python import pprint from numpy import * import matplotlib.pyplot as plt def loadDataSet(fileName): dataMat = []; labelMat = [] fr = open(fileName) for line in fr.readlines(): lineArr = line.strip().split('\t') dataMat.append([float(lineArr[0]), float(lineArr[1])]) labelMat.append(float(lineArr[2])) #只有一列 return dataMat, labelMat def load_data(path): ''' :param path:这里将libsvm中的数据格式进行读取 :return: ''' data_set=[] label_set=[] file_object=open(path) for line in file_object.readlines(): lineArr = line.strip().split(" ") label_set.append(int(lineArr[0])) #第一列作为lable tempdataline=[] for i in range(1,len(lineArr)): tempdataline.append(float(lineArr[i].split(':')[1])) data_set.append(tempdataline) print shape(mat(data_set)) return data_set,label_set class mySVM(object): def __init__(self): ''' :这里先初始化数据和权重 ''' self.data_set=[] self.label_set=[] self.b=[] #最后结果中的位移偏量 self.a=[] #最后结果中的拉格朗日乘子 self.kernel=[] #核函数 self.c=0.6 #设置惩罚因子 self.tol=0.01 #设置容忍极限 def k(self,x1,x2): ''' :param x1: 传入两个向量 :param x2: :return: 返回核函数 ''' linek=mat(x1).transpose()*mat(x2) #这里先采用线性核函数,后面可以改成高斯核函数 return linek[0][0] #因为矩阵乘法最后返回的是矩阵,所以要取第一个元素 def train(self,data_set,label_set): ''' :param data_set:传入数据集 :param label_set: 传入标签集 :return: 返回训练好的a和b ''' data_set=mat(data_set) label_set=mat(label_set).transpose() m,n = shape(data_set) pprint.pprint(data_set) print(m,n) a=mat(zeros((m,1))) max_iter=100 #迭代次数 b=0 toler=0.1 c=1 print a print lable_set while(max_iter>=0): max_iter=max_iter-1 for i in range(m):#遍历a #计算 ei=ui-yi tempu=multiply(a,label_set).T*data_set*data_set[i,:].T ui=tempu[0][0] + b ei=ui-label_set[i][0] ''' * 把违背KKT条件的ai作为第一个 * 满足KKT条件的情况是: * yi*f(i) >= 1 and alpha == 0 (正确分类) * yi*f(i) == 1 and 0<alpha < C (在边界上的支持向量) * yi*f(i) <= 1 and alpha == C (在边界之间) * * 这里借用了另一个网友的话 * * ri = y[i] * Ei = y[i] * f(i) - y[i]^2 >= 0 * 如果ri < 0并且alpha < C 则违反了KKT条件 * 因为原本ri < 0 应该对应的是alpha = C * 同理,ri > 0并且alpha > 0则违反了KKT条件 * 因为原本ri > 0对应的应该是alpha =0 ''' if((label_set[i]*ei<-toler) and (a[i]<c)) or ((label_set[i]*ei>toler) and (a[i]>0)): j=0 #这里应该是寻找MAX|E1 - E2| ,为了简单点,就随机选择了一个.只是速度要慢些 while(j==i): j=int(random.uniform(0, m)) uj=float(multiply(a, label_set).T*data_set*data_set[j, :].T)+b ej=uj-float(label_set[j]) aj_old=a[j] ai_old=a[i] if(label_set[i]<>label_set[j]): #开始计算L和H L=max(0, label_set[j]-a[i]) H=min(c, c+a[j]-a[i]) else: L=max(0, a[j]+a[i]-c) H=min(c, a[j]+a[i]) #这里计算2k(i,j)-k(i,i)-k(j,j) ,这里应该用核函数,先这样将就用了.不先核化 eta=2.0*data_set[i, :]*data_set[j, :].T-data_set[i, :]*data_set[i, :].T- data_set[j, :]*data_set[j, :].T if (eta >= 0):#如果eta等于0或者大于0 则表明a最优值应该在L或者U上 continue ''' 这里更新a[j] ''' a[j]=a[j]-label_set[j]*(ei-ej)/eta if(a[j]<L): a[j]=L elif(a[j]>H): a[j]=H #更新ai a[i]+=label_set[j]*label_set[i]*(aj_old-a[j]) #更新b ''' 根据公式 bnew1 =bold −E1 −y1 (α1new −α1old)K(x1 ,x1 )−y2 (α2new −α2old)K(x2 ,x2 ) bnew2 =bold −E2 −y1 (α1new −α1old)K(x1 ,x1 )−y2 (α2new −α2old)K(x2 ,x2 ) ''' b1=b-ei-label_set[i]*(a[i]-ai_old)*data_set[i, :]*data_set[i, :].T- \ label_set[j]*(a[j]-aj_old)* data_set[i, :]*data_set[j, :].T b2=b-ej-label_set[i]*(a[i]-ai_old)*data_set[i, :]*data_set[i, :].T- \ label_set[j]*(a[j]-aj_old)* data_set[i, :]*data_set[j, :].T #这里选择合适的b if(0<a[i]) and (c>a[i]): b=b1 elif(0<a[j]) and (c>a[j]): b=b2 else: b=(b1+b2)/2.0 return a,b if __name__ == '__main__': print("-------start load data-----") c=0.6 #定义系数c path="./testSet.txt" data_set,lable_set=loadDataSet(fileName=path) mysvm= mySVM() a,b=mysvm.train(data_set=data_set,label_set=lable_set) print a,b