整理并分析了里面一些错误和不适之处，已经在代码中修改，原因是Python版本的更新，这次代码比较全，有学习相关章节的读者可以参考下面的代码，博主学习时，已经测试过，如有问题，请留言交流哈。
# -*-coding:utf-8 -*-
import random
from numpy import *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):j=i #we want to select any J not equal to iwhile (j==i):j = int(random.uniform(0,m))return jdef clipAlpha(aj,H,L):if aj > H:aj = Hif L > aj:aj = Lreturn ajdef smoSimple(dataMatIn, classLabels, C, toler, maxIter):dataMatrix = mat(dataMatIn); labelMat = mat(classLabels).transpose()b = 0; m,n = shape(dataMatrix)alphas = mat(zeros((m,1)))iter = 0while (iter < maxIter):alphaPairsChanged = 0for i in range(m):fXi = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + bEi = fXi - float(labelMat[i])#if checks if an example violates KKT conditionsif ((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)) + bEj = fXj - float(labelMat[j])alphaIold = alphas[i].copy()alphaJold = alphas[j].copy()if (labelMat[i] != labelMat[j]):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 = 2.0 * dataMatrix[i,:]*dataMatrix[j,:].T -\dataMatrix[i,:]*dataMatrix[i,:].T -\dataMatrix[j,:]*dataMatrix[j,:].Tif eta>= 0: print("eta>=0"); continuealphas[j] -= labelMat[j]*(Ei - Ej)/etaalphas[j] = clipAlpha(alphas[j],H,L)if (abs(alphas[j] - alphaJold) < 0.00001): print("j not moving enough"); continuealphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])#update i by the same amount as j#the update is in the oppostie directionb1 = b - Ei- labelMat[i]*(alphas[i]-alphaIold)*\dataMatrix[i,:]*dataMatrix[i,:].T - \labelMat[j]*(alphas[j]-alphaJold)*\dataMatrix[i,:]*dataMatrix[j,:].Tb2 = b - Ej- labelMat[i]*(alphas[i]-alphaIold)*\dataMatrix[i,:]*dataMatrix[j,:].T - \labelMat[j]*(alphas[j]-alphaJold)*\dataMatrix[j,:]*dataMatrix[j,:].Tif (0 < alphas[i]) and (C > alphas[i]): b = b1elif (0 < alphas[j]) and (C > alphas[j]): b = b2else: b = (b1 + b2)/2.0alphaPairsChanged += 1print("iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged))if (alphaPairsChanged == 0):iter += 1else:iter = 0print("iteration number: %d" % iter)return b,alphas#以上测试，没有问题。class optStruct:def __init__(self,dataMatIn,classLabels,C,toler,kTup):self.X=dataMatInself.labelMat=classLabelsself.C=Cself.tol=tolerself.m=shape(dataMatIn)[0]self.alphas=mat(zeros((self.m,1)))self.b=0self.eCache=mat(zeros((self.m,2)))#changed for kernelself.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):# changed for kernelfXk=float(multiply(oS.alphas,oS.labelMat).T*(oS.K[:,k]+oS.b))#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=-1maxDeltaE=0Ej=0oS.eCache[i]=[1,Ei]validEcaheList=nonzero(oS.eCache[:,0].A)[0]if (len(validEcaheList))>1:for k in validEcaheList:if k==i:continueEk=calcEk(oS,k)deltaE=abs(Ei-Ek)if(deltaE>maxDeltaE):maxK=kmaxDeltaE=deltaEEj=Ekreturn maxK,Ejelse: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) #this has been changed from selectJrandalphaIold = 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 0eta = 2.0 * oS.K[i,j] - oS.K[i,i] - oS.K[j,j] #changed for kernel#eta=2.0 * oS.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].Tif eta >= 0: print("eta>=0"); return 0oS.alphas[j] -= oS.labelMat[j]*(Ei - Ej)/etaoS.alphas[j] = clipAlpha(oS.alphas[j],H,L)updateEk(oS, j) #added this for the Ecacheif (abs(oS.alphas[j] - alphaJold) < 0.00001): print("j not moving enough"); return 0oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])#update i by the same amount as jupdateEk(oS, i) #added this for the Ecache #the update is in the oppostie direction#changed for kernelb1 = 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]#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]if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2else: oS.b = (b1 + b2)/2.0return 1else: return 0def smoP(dataMatIn,classLabels,C,toler,maxIter,kTup=('lin',0)):oS=optStruct(mat(dataMatIn),mat(classLabels).transpose(),C,toler,kTup)iter=0entireSet=TruealphaParisChanged=0while (iter< maxIter) and ((alphaParisChanged>0) or (entireSet)):alphaParisChanged=0if entireSet:for i in range(oS.m):alphaParisChanged+=innerL(i,oS)print("fullset,iter: %d i: %d, pairs changed %d" %(iter,i,alphaParisChanged))iter+=1else:nonBoundIs=nonzero((oS.alphas.A>0)* (oS.alphas.A <C))[0]for i in nonBoundIs:alphaParisChanged+=innerL(i,oS)print("non-bound,iter: %d i:%d,pairs changed %d" %(iter,i,alphaParisChanged))iter+=1if entireSet:entireSet=Falseelif (alphaParisChanged==0):entireSet=Trueprint("iteration number: %d" %iter)return oS.b,oS.alphasdef 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 wdef kernelTrans(X,A,kTup):'''The first argument in the tuple is a string describing what type of kernel should be used.the other arguments are optional arguments that may be needed for a kernel.In the case of the linear kernel, a dot product is taken between the two inputs,which are the full dataset and a row of the dataset.In the case of the radial bias function,the Gaussian function is evaluated for every elementin the matrix in the for loop. After the for loop is finished, you apply the calculations over the entire vector.'''#kTup tuple contain the kernel function information:sgma.m,n=shape(X)K=mat(zeros((m,1)))#linear kernelif kTup[0]=='lin':K=X*A.T#radial bias functionelif kTup[0]=='rbf':for j in range(m):deltaRow=X[j,:]-AK[j]=deltaRow*deltaRow.TK=exp(K/(-1*kTup[1]**2))else:raise NameError('Houstion We Have a Problem -- \That Kernel is not recognized')return Kdef testRbf(k1=1.3):dataArr,labelArr=loadDataSet('data/testSetRBF.txt')b,alphas=smoP(dataArr,labelArr,200,0.001,10000,('rbf',k1))dataMat=mat(dataArr)labelMat=mat(labelArr).transpose()#martix.A means transform martix into arraysvInd=nonzero(alphas.A>0)[0]sVs=dataMat[svInd]labelSV=labelMat[svInd]print("There are %d Support Vectors" %shape(sVs)[0])m,n=shape(dataMat)errorCount=0for i in range(m):kernelEval=kernelTrans(sVs,dataMat[i,:],('rbf',k1))predict=kernelEval.T *multiply(labelSV,alphas[svInd])+bif sign(predict)!=sign(labelArr[i]):errorCount+=1print("the training error rate is: %f" %(float(errorCount)/m))dataArr,labelArr=loadDataSet('data/testSetRBF2.txt')errorCount=0dataMat=mat(dataArr)labelMat=mat(labelArr).transpose()m,n=shape(dataMat)for i in range(m):kernelEval=kernelTrans(sVs,dataMat[i,:],('rbf',k1))predict=kernelEval.T *multiply(labelSV,alphas[svInd])+bif sign(predict)!=sign(labelArr[i]):errorCount+=1print("The test error rate is: %f" %(float(errorCount)/m))def img2vector(filename):returnVect=zeros((1,1024))fr=open(filename)for i in range(32):lineStr=fr.readline()for j in range(32):returnVect[0,32*i+j]=int(lineStr[j])return returnVectdef loadImages(dirName):from os import listdirhwlabels=[]trainingFileList=listdir(dirName)m=len(trainingFileList)trainingMat=zeros((m,1024))for i in range(m):fileNameStr=trainingFileList[i]fileStr=fileNameStr.split('.')[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,hwlabelsdef testDigits(kTup=('rbf',10)):#same as testRbf.dataArr,labelArr=loadImages('digits/trainingDigits')b,alphas=smoP(dataArr,labelArr,200,0.0001,10000,kTup)dataMat=mat(dataArr)labelMat=mat(labelArr).transpose()#martix.A means transform martix into arraysvInd=nonzero(alphas.A>0)[0]sVs=dataMat[svInd]labelSV=labelMat[svInd]print("There are %d Support Vectors" %shape(sVs)[0])m,n=shape(dataMat)errorCount=0for i in range(m):kernelEval=kernelTrans(sVs,dataMat[i,:],kTup)predict = kernelEval.T * multiply(labelSV, alphas[svInd]) + bif sign(predict) != sign(labelArr[i]):errorCount += 1print("the training error rate is: %f" % (float(errorCount) / m))dataArr, labelArr = loadImages('digits/testDigits')errorCount = 0dataMat = mat(dataArr)labelMat = mat(labelArr).transpose()m, n = shape(dataMat)for i in range(m):kernelEval = kernelTrans(sVs, dataMat[i,:], kTup)predict = kernelEval.T * multiply(labelSV, alphas[svInd]) + bif sign(predict) != sign(labelArr[i]):errorCount += 1print("The test error rate is: %f" % (float(errorCount) / m))
'''
if __name__=='__main__':dataArr,labelArr=loadDataSet('data/testSet.txt')b,alphas=smoP(dataArr,labelArr,0.6,0.001,40)ws=calcWs(alphas,dataArr,labelArr)testRbf()testDigits()'''