机器学习实战-支持向量机
来源:互联网 发布:mirrorlink已认证软件 编辑:程序博客网 时间:2024/04/18 05:22
支持向量机
优点:泛化错误率低,计算开销不大,结果易于理解
缺点:对参数调节和核函数的选择敏感,原始分类器不加修改仅适用于处理二类问题
适用数据类型:数值型和标称型数据
如果数据集是1024维,需使用一个1023维的某对象来对数据进行分隔,改对象称作超平面。
#6-1 SMO算法辅助函数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,labelMatdef selectJrand(i,m):#防止alpha下标和alpha的数目相同 j=i while(j==i): j = int(random.uniform(0,m)) return jdef clipAlpha(aj,H,L): if aj > H: aj = H if L > aj: aj = L return aj
实验这部分代码
import svmMLiAdataArr,labelArr = svmMLiA.loadDataSet("testSet.txt")labelArrOut[7]: [-1.0, -1.0,#略过若干结果... -1.0, -1.0]
下面开始SMO算法的第一个版本,伪代码大致如下:
创建一个alpha向量并将其初始化为0向量当迭代次数小于最大迭代次数时(外循环) 对数据集中的每个数据向量(内循环): 如果该数据向量可以被优化: 随机选择另外一个数据向量 同事优化这两个向量 如果两个向量都不能被优化,退出内循环如果所有向量都没被优化,增加迭代数目,继续下一次循环
#6-2 简化版SMO算法def smoSimple(dataMatIn,classLabels,C,toler,maxIter):#数据集,类别标签,常数C,容错率,退出前最大的循环次数 dataMatrix = mat(dataMatIn); labelMat = mat(classLabels).transpose() b = 0; m,n = shape(dataMatrix) alphas =mat(zeros((m,1))) iter = 0 while(iter < maxIter): alphaPairsChanged = 0 for i in range(m): fXi = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T))+b Ei = fXi - float(labelMat[i]) if ((labelMat[i]*Ei < -toler) and (alphas[i] < C)) or ((labelMat[i]*Ei > toler) and (alphas[i] > 0)): j = selectJrand(i,m) fXj = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[j,:].T)) +b Ej = fXj - float(labelMat[j]) alphaIold = alphas[i].copy(); alphaJold = alphas[j].copy();#深度拷贝 if (labelMat[i]!=labelMat[j]):#保证alpha在0和C之间 L = max(0,alphas[j]-alphas[i]) H = min(C,C+alphas[j]-alphas[i]) else: L = max(0,alphas[j]+alphas[i]-C) H = min(C,alphas[j]+alphas[i]) if L==H:print "L==H"; continue #eta是alpha[j]的最优修改量 eta = 2.0*dataMatrix[i,:]*dataMatrix[j,:].T - dataMatrix[i,:]*dataMatrix[i,:].T - dataMatrix[j,:]*dataMatrix[j,:].T if eta>=0:print "eta>=0";continue alphas[j]-=labelMat[j]*(Ei-Ej)/eta alphas[j]=clipAlpha(alphas[j],H,L) if (abs(alphas[j]-alphaJold)<0.00001):print"j not moving enough"; continue alphas[i]+=labelMat[j]*labelMat[i]*(alphaJold-alphas[j])#修改方向相反 b1 = b - Ei- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[i,:].T - labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[i,:]*dataMatrix[j,:].T b2 = b - Ej- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[j,:].T - labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[j,:]*dataMatrix[j,:].T if (0 < alphas[i]) and (C > alphas[i]): b = b1#设置常数项B elif (0 < alphas[j]) and (C > alphas[j]): b = b2 else: b = (b1 + b2)/2.0 alphaPairsChanged += 1 print "iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged) if (alphaPairsChanged == 0): iter += 1 else: iter = 0 print "iteration number: %d" % iter return b,alphas
#测试实际效果import svmMLiAdataArr,labelArr = svmMLiA.loadDataSet("testSet.txt")b,alphas = svmMLiA.smoSimple(dataArr,labelArr,0.6,0.001,40)
#略过部分iteration number: 28iter: 28 i:29, pairs changed 1iteration number: 0j not moving enoughiteration number: 1j not moving enoughj not moving enoughiteration number: 40
对结果进行观察
In [27]: bOut[27]: matrix([[-3.83810926]])
我们可以直接观察alpha本身,但是其中的0元素过多,为了观察大于0的元素的数量,可以
alphas[alphas>0]#适用于NumPy类型Out[28]: matrix([[ 0.12749752, 0.24132585, 0.36882337]])
由于SMO算法的随机性,读者运行后的结果可能不同。
#获得支持向量的格式In [30]: shape(alphas[alphas>0])Out[30]: (1L, 3L)#了解哪些数据点是支持向量In [32]: for i in range(100): ...: if alphas[i]>0.0:print dataArr[i],labelArr[i] ...: [4.658191, 3.507396] -1.0[3.457096, -0.082216] -1.0[6.080573, 0.418886] 1.0
下面开始讨论完整版Platt SMO算法。他通过一个外循环来选择第一个alpha值,并且选择过程会在两种方式之间进行交替:一种是在所有数据集是进行单遍扫描,另一种方式则是非边界alpha中实现单遍扫描。
class optStruct: def __init__(self,dataMatIn,classLabels,C,toler): self.X = dataMatIn self.labelMat = classLabels self.C = C self.tol = toler self.m = shape(dataMatIn)[0] self.alphas = mat(zeros((self.m,1)))#是否有效的标志位和实际的E值 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] = kernelTrans(self.X,self.X[i,:],kTup)def calcEk(oS,k): fXk = float(multiply(oS.alphas,oS.labelMat).T*(oS.X*oS.X[k,:].T)) + oS.b Ek = fXk - float(oS.labelMat[k]) return Ekdef selectJ(i,oS,Ei): maxK = -1; maxDeltaE=0; Ej=0 oS.eCache[i] = [1,Ei] validEcacheList=nonzero(oS.eCache[:,0].A)[0] if(len(validEcacheList))>1: for k in validEcacheList: if k==i: continue Ek = calcEk(oS,k) deltaE=abs(Ei-Ek) if(deltaE>maxDeltaE):#选择最大步长 maxK=k;maxDeltaE=deltaE;Ej=Ek return maxK,Ej else: j=selectJrand(i,oS.m) Ej=calcEk(oS,j) return j,Ejdef updateEk(oS, k): Ek = calcEk(oS, k) oS.eCache[k] = [1,Ek]#寻找决策边界的优化例程def innerL(i,oS): Ei=calcEk(oS,i) if((oS.labelMat[i]*Ei<-oS.tol)and(oS.alphas[i]< oS.C)) or ((oS.labelMat[i]*Ei > oS.tol) and (oS.alphas[i] > 0)): j,Ej = selectJ(i, oS, Ei) alphaIold = oS.alphas[i].copy();alphaJold = oS.alphas[j].copy(); if (oS.labelMat[i] != oS.labelMat[j]): L = max(0,oS.alphas[j] - oS.alphas[i]) H = min(oS.C, oS.C + oS.alphas[j] - oS.alphas[i]) else: L = max(0,oS.alphas[j] + oS.alphas[i] - oS.C) H = min(oS.C, oS.alphas[j] + oS.alphas[i]) if L==H: print "L==H"; return 0 eta = 2.0*oS.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].T if eta >= 0: print "eta>=0"; return 0 oS.alphas[j] -= oS.labelMat[j]*(Ei - Ej)/eta oS.alphas[j] = clipAlpha(oS.alphas[j],H,L) updateEk(oS, j)#更新缓存误差 if (abs(oS.alphas[j] - alphaJold) < 0.00001):print "j not moving enough"; return 0 oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j]) updateEk(oS, i) b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.X[i,:]*oS.X[i,:].T - oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.X[i,:]*oS.X[j,:].T b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.X[i,:]*oS.X[j,:].T- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.X[j,:]*oS.X[j,:].T if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1 elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2 else: oS.b = (b1 + b2)/2.0 return 1 else: return 0def smoP(dataMatIn, classLabels, C, toler, maxIter): oS = optStruct(mat(dataMatIn),mat(classLabels).transpose(),C,toler)#构建数据结构容纳所有数据 iter = 0 entireSet = True; alphaPairsChanged = 0 while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)): alphaPairsChanged = 0 if entireSet: for i in range(oS.m): alphaPairsChanged += innerL(i,oS) print "fullSet, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged) iter += 1 else: nonBoundIs = nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0] for i in nonBoundIs: alphaPairsChanged += innerL(i,oS) print "non-bound, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged) iter += 1 if entireSet: entireSet = False elif (alphaPairsChanged == 0): entireSet = True print "iteration number: %d" % iter return oS.b,oS.alphas
下面观察上述执行结果:
import svmMLiAlidataArr,labelArr = svmMLiAli.loadDataSet("testSet.txt")b,alphas = svmMLiAli.smoP(dataArr,labelArr,0.6,0.001,40)
fullSet, iter: 0 i:0, pairs changed 1fullSet, iter: 0 i:1, pairs changed 1fullSet, iter: 0 i:2, pairs changed 2fullSet, iter: 0 i:3, pairs changed 3L==HfullSet, iter: 0 i:4, pairs changed 3L==HfullSet, iter: 0 i:5, pairs changed 3L==HfullSet, iter: 0 i:6, pairs changed 3fullSet, iter: 0 i:7, pairs changed 3fullSet, iter: 0 i:8, pairs changed 4fullSet, iter: 0 i:9, pairs changed 4j not moving enoughfullSet, iter: 0 i:10, pairs changed 4fullSet, iter: 0 i:11, pairs changed 4fullSet, iter: 0 i:12, pairs changed 4fullSet, iter: 0 i:13, pairs changed 4fullSet, iter: 0 i:14, pairs changed 4fullSet, iter: 0 i:15, pairs changed 4fullSet, iter: 0 i:16, pairs changed 4j not moving enough#省略部分结果
下面开始计算w
def calcWs(alphas,dataArr,classLabels): X = mat(dataArr); labelMat = mat(classLabels).transpose() m,n = shape(X) w = zeros((n,1)) for i in range(m): w += multiply(alphas[i]*labelMat[i],X[i,:].T) return w
In [2]: ws = calcWs(alphas,dataArr,labelArr)In [3]: wsOut[3]: array([[ 0.65307162], [-0.17196128]])
现在对数据进行分类处理,比如对第一个数据分类,可以:
In [4]: datMat = mat(dataArr)In [5]: datMat[0]*mat(ws) + bOut[5]: matrix([[-0.92555695]])
如果该数字大于0,则属于1类,小于则属于-1类:
In [6]: labelArr[0]Out[6]: -1.0In [7]: datMat[2]*mat(ws) + bOut[7]: matrix([[ 2.30436336]])In [8]: labelArr[2]Out[8]: 1.0In [9]: datMat[1]*mat(ws) + bOut[9]: matrix([[-1.36706674]])In [10]: labelArr[1]Out[10]: -1.0
下面研究不能线性处理的情况,需要使用核函数(kernel)。我们需要将数据从一个特征空间转换到另外一个特征空间。
#6-6 核转换函数def kernelTrans(X, A, kTup):#kT是核函数信息,第一个是类型,另外两个是可选参数 m,n = shape(X) K = mat(zeros((m,1))) if kTup[0]=="lin":K = X * A.T elif kTup[0] == "rbf": for j in range(m): deltaRow = X[j,:] - A K[j] = deltaRow*deltaRow.T K = exp(K/(-1*kTup[1]**2))#元素除,NumPy中指矩阵元素展开计算 else: raise NameError("Houston We Have a Problem -- That Kernel is not recognized") return Kclass optStruct: def __init__(self,dataMatIn,classLabels,C,toler, kTup):#增加kTup self.X = dataMatIn self.labelMat = classLabels self.C = C self.tol = toler self.m = shape(dataMatIn)[0] self.alphas = mat(zeros((self.m,1)))#是否有效的标志位和实际的E值 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] = kernelTrans(self.X,self.X[i,:],kTup)
另外还需要修改:
def innerL(i,oS): #eta = 2.0*oS.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].T eta = 2.0*oS.K[i,j] - oS.K[i,i] - oS.K[j,j] #b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.X[i,:]*oS.X[i,:].T - oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.X[i,:]*oS.X[j,:].T #b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.X[i,:]*oS.X[j,:].T- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.X[j,:]*oS.X[j,:].T b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,i] - oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[i,j] b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.K[i,j]- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.K[j,j]def calcEk(oS,k): #fXk = float(multiply(oS.alphas,oS.labelMat).T*(oS.X*oS.X[k,:].T)) + oS.b fXk = float(multiply(oS.alphas,oS.labelMat).T*oS.K[:,k] + oS.b)
#6-8 利用核函数进行分类的径向基测试函数def testRbf(k1=1.3): dataArr,labelArr = loadDataSet("testSetRBF.txt") b,alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, ("rbf", k1)) datMat = mat(dataArr); labelMat = mat(labelArr).transpose() svInd = nonzero(alphas.A>0)[0] #返回数组中值不为零的元素的下标 sVs = datMat[svInd] labelSV = labelMat[svInd] print "there are %d Support Vectors" % shape(sVs)[0] m,n = shape(datMat) errorCount = 0 for i in range(m): kernelEval = kernelTrans(sVs,datMat[i,:],("rbf",k1)) predict=kernelEval.T * multiply(labelSV,alphas[svInd]) + b if sign(predict)!=sign(labelArr[i]): errorCount += 1 print "the training error rate is: %f" % (float(errorCount)/m) dataArr,labelArr = loadDataSet('testSetRBF2.txt') errorCount = 0 datMat=mat(dataArr); labelMat = mat(labelArr).transpose() m,n = shape(datMat) for i in range(m):#测试数据集 kernelEval = kernelTrans(sVs,datMat[i,:],('rbf', k1)) predict=kernelEval.T * multiply(labelSV,alphas[svInd]) + b if sign(predict)!=sign(labelArr[i]): errorCount += 1 print "the test error rate is: %f" % (float(errorCount)/m)
测试
In [19]: svmMLiA.testRbf()#省略部分数据L==HfullSet, iter: 4 i:91, pairs changed 0L==HfullSet, iter: 4 i:92, pairs changed 0fullSet, iter: 4 i:93, pairs changed 0fullSet, iter: 4 i:94, pairs changed 0fullSet, iter: 4 i:95, pairs changed 0L==HfullSet, iter: 4 i:96, pairs changed 0fullSet, iter: 4 i:97, pairs changed 0fullSet, iter: 4 i:98, pairs changed 0fullSet, iter: 4 i:99, pairs changed 0iteration number: 5there are 29 Support Vectorsthe training error rate is: 0.130000the test error rate is: 0.150000
支持向量的数目存在一个最优值。SVM的优点在于它能对数据进行高效分类。如果支持向量太少,就可能会得到一个很差的决策边界;如果支持向量太多,也就相当于每次都利用整个数据集进行分类,这种分类情况称作k近邻。
先加入第二章knn算法中的img2vector()函数,然后加入如下代码:
#6-9 基于SVM的手写数字识别def loadImages(dirName): from os import listdir hwLabels = [] trainingFileList = listdir(dirName) m = len(trainingFileList)#总文件个数 trainingMat = zeros((m,1024)) for i in range(m): fileNameStr = trainingFileList[i] fileStr = fileNameStr.split('.')[0]#按“.”分开,取第0行 classNumStr = int(fileStr.split('_')[0]) if classNumStr == 9: hwLabels.append(-1) else: hwLabels.append(1) trainingMat[i,:] = img2vector('%s/%s' % (dirName, fileNameStr)) return trainingMat, hwLabels def testDigits(kTup=('rbf', 10)): dataArr,labelArr = loadImages('trainingDigits') b,alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, kTup) datMat=mat(dataArr); labelMat = mat(labelArr).transpose() svInd=nonzero(alphas.A>0)[0] sVs=datMat[svInd] labelSV = labelMat[svInd]; print "there are %d Support Vectors" % shape(sVs)[0] m,n = shape(datMat) errorCount = 0 for i in range(m): kernelEval = kernelTrans(sVs,datMat[i,:],kTup) predict=kernelEval.T * multiply(labelSV,alphas[svInd]) + b if sign(predict)!=sign(labelArr[i]): errorCount += 1 print "the training error rate is: %f" % (float(errorCount)/m) dataArr,labelArr = loadImages('testDigits') errorCount = 0 datMat=mat(dataArr); labelMat = mat(labelArr).transpose() m,n = shape(datMat) for i in range(m): kernelEval = kernelTrans(sVs,datMat[i,:],kTup) predict=kernelEval.T * multiply(labelSV,alphas[svInd]) + b if sign(predict)!=sign(labelArr[i]): errorCount += 1 print "the test error rate is: %f" % (float(errorCount)/m)
In [35]: svmMLiA.testDigits(("rbf",20))#省略部分结果L==HfullSet, iter: 3 i:397, pairs changed 0L==HfullSet, iter: 3 i:398, pairs changed 0L==HfullSet, iter: 3 i:399, pairs changed 0fullSet, iter: 3 i:400, pairs changed 0j not moving enoughfullSet, iter: 3 i:401, pairs changed 0iteration number: 4there are 51 Support Vectorsthe training error rate is: 0.000000the test error rate is: 0.016129
根据课本,σ取10的时候可以得到最小的错误率。可以观察到一个有趣的现象,即最小的训练错误率并不对应于最小的支持向量数目。另外,线性核函数并不是特别的糟糕,可以以牺牲线性核函数的错误率来换取分类速度的提高。
- 机器学习实战之支持向量机
- 机器学习实战-支持向量机
- 机器学习实战 支持向量机
- 机器学习实战_06支持向量机
- 机器学习实战笔记-支持向量机
- 机器学习实战-6支持向量机
- 《机器学习实战》学习笔记:支持向量机
- Python学习-机器学习实战-ch06 支持向量机
- 《机器学习实战》学习笔记-[6]-支持向量机SMO
- 机器学习实战——第六章:支持向量机
- 【机器学习实战】第6章 支持向量机
- python机器学习实战5:支持向量机
- 机器学习实战-第六章(支持向量机)
- PYTHON机器学习实战——SVM支持向量机
- 机器学习实战【5】(SVM-支持向量机)
- 机器学习-支持向量机
- 机器学习--支持向量机
- 机器学习-支持向量机
- 机器学习算法原理与编程学习(1)
- Eclipse Java注释模板设置详解以及版权声明
- c c++ 混用
- android:process验证记录
- 在运行时请求权限 android
- 机器学习实战-支持向量机
- Unity Transform
- for...in...循环搭配对象解构赋值的易错点
- Android7.0中文文档(API)--- BaseExpandableListAdapter
- IOS进程间通信
- Swift3之细致理解Optional(可选类型)
- 解决android studio3.0 preview一直卡在Gradle:Resolve dependencies':classpath'的问题
- DDoS攻击与CC攻击的区别
- web网页 页面布局的几种方式