今天看了Python語言寫的使用SVM中的SMO進行優化，使用RBF函數進行手寫體識別，下面簡單整理一下整個過程及思路，然后詳細介紹各個部分。
整個過程:
(1)獲取訓練數據集trainingMat和labelMat;
(2)利用SMO進行優化獲得優化參數alphas和b，這一步即是進行訓練獲得最優參數;
(3)使用alphas和b帶入RBF高斯核函數計算訓練集輸出並計算訓練錯誤率;
(4)獲取測試數據集testMat和labelMat1;
(5)使用(2)的參數alphas和b帶入RBF高斯核函數計算輸出，從而計算分類錯誤率;
Python代碼:
def testDigits(kTup=('rbf',10)):
dataArr,labelArr=loadImages('D://softwareTool/Python/python_exerciseCode/Chap6_SVM//trainingDigits');
b,alphas=smoP(dataArr,labelArr,200,0.0001,1000,kTup);
dataMat=mat(dataArr);
labelMat=mat(labelArr).transpose();
# 取得支持向量的索引
svInd=nonzero(alphas.A>0)[0];
sVs=dataMat[svInd];
labelSV=labelMat[svInd];
print("there are ",shape(sVs)[0],' Support Vectors');
m,n=shape(dataMat);
errorCount=0.0;
for i in range(m):
kernelEvl=kernelTrans(sVs,dataMat[i,:],kTup);
# 計算輸出公式
predict=kernelEvl.T*multiply(labelSV,alphas[svInd])+b;
if sign(predict)!=sign(labelMat[i]):
errorCount+=1.0;
print("the training error rate is:",errorCount/(len));
dataArr,labelArr=loadImages('D://softwareTool/Python/python_exerciseCode/Chap6_SVM//testDigits');
dataMat=mat(dataArr);
labelMat=mat(labelArr).transpose();
m,n=shape(dataMat);
errorCount=0.0;
for i in range(m):
kernelEval=kernelTrans(sVs,dataMat[i,:],kTup);
predict=kernelEval*multiply(labelSV,alphas[svInd])+b;
if sign(predict)!=sign(labelMat[i]):
errorCount+=1.0;
print("the test error rate is: ",errorCount/float(ms));
上面是整個主框架和主程序。
下面分模塊介紹各個部分：
(1)獲取訓練數據集和訓練標簽:
如給的訓練數據存放在文件trainingDigits中,其里面有多個.txt子文件,每個.txt文件存放的是一幅32*32的的圖像，每一幅圖像表示0-9的一個數字；
如這幅圖顯示的為數字3，最后將每幅圖像轉化為一個32*32=1024的列向量，如果訓練樣本數為m，則dataMat為m*1024的矩陣，labelMat為1*m的列向量。
下面先說說在Python中怎樣將一幅圖像（例32*32）轉化為一個列向量（1*1024）:函數輸入為這幅圖像的文件名‘3_177.txt’
偽代碼: 初始化列向量returnVec為zeors((1,1024));
遍歷每行：
讀取每行（1*32的列向量）內容;
遍歷每列:
將每列內容強轉化為int型依次添加到returnVec中;
返回returnVec;
Python代碼:
# 將一幅圖像轉化為一個向量 def img2vector(filename): # 將一幅32*32的圖像轉化為一個1*1024的圖像 returnVect=zeros((1,1024)); fr=open(filename);
# 遍歷每行 for i in range(32): # 讀取文件每行內容 lineStr=fr.readline(); for j in range(32): returnVect[0,i*32+j]=int(lineStr[j]); return returnVect;
注意:fr.readline()和fr.readlines()區別;
下面寫一下調用上面的函數img2vector(filename)怎樣讀取所有文件內容得到訓練數據矩陣trainingMat(m*1024)和類標簽l向量hwLabels(m*1):
偽代碼: 取得所有子文件名(如'0_0.txt','0_1.txt',...),將其存放在一個列表里;
計算子文件個數(m),即統計樣本個數;
初始化trainingMat和hwLabels;
遍歷每個子文件:
取得該文件名(如'3_177.txt');
使用split('.')將文件名字符串切分為2部分,取前面部分(如'3_177');
使用split('_')將字符串切分為2部分,取前面部分('3'),同時強制轉化為int類型;
由前面得到的整數獲得對應類標簽添加到hwLabels里(由於本例使用2分類,將數字0-8分為1類,數字9分為-1類,所以數字3的類標簽為1);
調用函數img2vector(filename)將該文件的一幅圖像(32*32)轉化為一個列向量(1*1024),存放到trainingMat里;
返回trainingMat,hwLabels;
Python代碼:
def loadImages(dirName):
from os import listdir; hwLabels=[]; # 將所有子文件名存放在一個列表中 trainingFileList=listdir(dirName); # 統計子文件個數 m=len(trainingFileList); # 初始化訓練數據矩陣 trainingMat=zeros((m,1024)); # 遍歷每個子文件 for i in range(m): # 獲得子文件名如'0_2.txt' fileNameStr=trainingFileList[i]; # 用'.'將上面獲得的字符串切分為2部分,取前面部分'0_2' fileStr=fileNameStr.split('.')[0]; # 用'_'將上面獲得的字符串切分為2部分,取前面部分'0'且強轉為int型得到類標簽 classNumStr=int(fileStr.split('_')[0]); if classNumStr==9: hwLabels.append(-1); else: hwLabels.append(1); # 將每幅圖像轉化為一個1*1024的列向量存放在trainingMat中 trainingMat[i,:]=img2vector(dirName+'/'+fileNameStr);
return trainingMat,hwLabels;
(2)使用SMO(序列最小最優化)算法進行參數優化訓練取得最優參數b和alphas
支持向量機的學習問題可以形式化為求解凸二次規划問題,這樣的凸二次規划問題具有全局最優解,並且有許多最優化算法可以用於這一問題的求解.但當訓練樣本容量很大時,這些算法往往變得非常低效,以至無法使用.而SMO是一種快速實現算法.
SMO算法要解決如下凸二次規划的對偶問題:
在這個問題中,變量是拉格朗日乘子,一個變量alphai對應於一個樣本點(xi,yi);變量的總數等於訓練樣本數m.
SMO算法是一種啟發式算法,其基本思路是:如果所有變量的解都滿足此最優化問題的KKT條件,那么這個最優化問題的解就得到了.因為KKT條件是該最優化問題的充分必要條件.否則,選擇兩個變量,固定其他變量,針對這兩個變量構建一個二次規划問題.這個二次規划問題關於這兩個變量的解應該更接近原始二次規划問題的解,因為這使得原始二次規划問題的目標函數值變得更小,重要的是,這時子問題可以通過解析方法求解,這樣就可以大大提高整個算法的計算速度.子問題有兩個變量,一個是違反KKT條件最嚴重的那一個,另一個由約束條件自動確定.如此,SMO算法將原問題不斷分解為子問題並對子問題求解,進而達到求解原問題的目的.
注意,子問題的兩個變量中只有一個是自由變量.假設alpha1,alpha2為兩個變量,alpha3,alpha4,...,alphaN固定,那么由等式約束可知:
如果alpha2確定,那么alpha1也隨之確定.所以子問題中同時更新兩個變量.
整個SMO算法包括兩個部分:求解兩個變量二次規划的解析方法和選擇變量的啟發式方法.
簡化版的SMO算法處理小規模數據集:Platt SMO算法的完整實現需要大量代碼.Platt SMO算法中的外循環確定要優化的最佳alpha對.而簡化版會跳過這一部分,首先在數據集上遍歷每一個alpha,然后在剩下的alpha集合中隨機選擇另一個alpha,從而構成alpha對.就是要同時改變兩個alpha.為此,將構建兩個輔助函數.一個輔助函數用於在某個區間范圍內隨機選擇一個整數;一個輔助函數用於在數值太大時進行調整.下面是兩個輔助函數的Python代碼:
# 隨機選擇一個j,函數輸入i,m;i為alpha的下標,m為alpha的個數 def selectJrand(i,m): j=i; while(j==i): j=int(random.uniform(0,m)); return j; # 用於調整大於H或者小於L的alpha,使其介於[L,H]內
# 輸入參數aj,H,L def clipAlpha(aj,H,L): if aj>H: aj=H; if aj<L: aj=L; return aj;
簡化版的SMO偽代碼大致如下:
偽代碼: 創建一個alpha向量並將其初始化為0向量,初始化b為0;
當迭代次數小於最大迭代次數時(外循環):
遍歷每個樣本,即對數據集中的每個數據向量(內循環):
計算alpha1的函數值和誤差;
誤差和alpha1滿足優化條件(如果該數據向量可以被優化):
隨機選擇另外一個數據向量alpha2;
計算alpha2的函數值和誤差;
計算L和H;
計算eta;
更新alpha1和alpha2;
更新b1,b2;
計算b;
同時優化這兩個向量;
如果兩個向量都不能被優化,退出內循環;
如果所有向量都沒被優化,增加迭代數目,繼續下一次循環.
Python源代碼:
def smoSimple(dataMatIn,classLabels,C,toler,maxIter): # 輸入5個參數. dataMatrix=mat(dataMatIn); # m*n labelMat=mat(classLabels).transpose(); m,n=shape(dataMatrix); b=0; alphas=mat(zeros((m,1)));# 創建一個alpha向量並將其初始化為0向量 iterNum=0; while(iterNum<maxIter): alphasPairsChanged=0; for i in range(m): # 計算猜測的類別,P143式(7.104)g(x)=sigma(ai*yi*K(x,xi))+b fXi=float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T))+b; # 計算誤差 Ei=fXi-float(labelMat[i]); # 如果該隨機變量可以被優化(滿足約束條件),隨機選擇另外一個數據向量;隨機選擇j # 如果誤差很大,對alpha值進行優化 # 在if語句中,不過是正間隔還是負間隔都會被測試,同時檢查alpha值不能為0或C,等於這兩個值時,表示在邊界上 if ((Ei*labelMat[i]<-toler) and (alphas[i]<C))\ or ((Ei*labelMat[i]>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,alphaJold;alphaIold,alphaJold為初始可行解 alphaIold=alphas[i].copy(); alphaJold=alphas[j].copy(); # 計算L和H,用於將alphas[j]調整到0-C之間 if(labelMat[i]!=labelMat[j]): L=max(0,alphaJold-alphaIold); H=min(C,C+alphaJold-alphaIold); else: L=max(0,alphaJold+alphaIold-C); H=min(C,alphaJold+alphaIold); if L==H: print('L==H'); continue; # 結束本次循環 # eta是alpha[j]的最優修改量eta=K11+K22-2*K12 eta=2*dataMatrix[i,:]*dataMatrix[j,:].T-\ dataMatrix[i,:]*dataMatrix[i,:].T-\ dataMatrix[j,:]*dataMatrix[j,:].T; if eta>=0: print("eta>=0"); continue; # 結束本次循環 # 同時改變alphas[i]和alphas[j],一個增大,一個減小 # 最優化問題沿着約束方向未經剪輯時的解 alphas[j]-=labelMat[j]*(Ei-Ej)/eta; # 調整大於H或者小於L的alpha alphas[j]=clipAlpha(alphas[j],H,L); # 檢查alphas[j]是否有輕微改變;如果是的話,退出for循環 if(abs(alphas[j]-alphaJold)<0.0001): print("j not moving enough"); continue; # 由alpha2new求得alpha1new alphas[i]+=labelMat[i]*labelMat[j]*(alphaJold-alphas[j]); # 每次完成兩個變量的優化后,都要重新計算閾值b b1=-Ei-labelMat[i]*dataMatrix[i,:]*dataMatrix[i,:].T*\ (alphas[i]-alphaIold)-labelMat[j]*dataMatrix[i,:]*\ dataMatrix[j,:].T*(alphas[j]-alphaJold)+b; b2=-Ej-labelMat[i]*dataMatrix[i,:]*dataMatrix[j,:].T*\ (alphas[i]-alphaIold)-labelMat[j]*dataMatrix[j,:]*\ dataMatrix[j,:].T*(alphas[j]-alphaJold)+b; if alphas[i]>0 and alphas[i]<C: b=b1; elif alphas[j]>0 and alphas[j]<C: b=b2; else: b=(b1+b2)/2; alphasPairsChanged+=1; print("iter: ",iterNum," i: ",i," pairs changed: ",alphasPairsChanged); # 如果所有向量(遍歷完所有樣本后)都沒被優化,增加迭代數目,繼續下一次循環 if(alphasPairsChanged==0): iterNum+=1; else: iterNum=0; print("iteration number: ",iterNum); return b,alphas;
Platt SMO算法是通過一個外循環來選擇第一個alpha值的,並且其選擇過程會在兩種方式之間進行交替:一種方式是在所有數據集上進行單邊掃描,另一種方式是在非邊界alpha中實現單邊掃描.而所謂非邊界alpha指的是那些不等於邊界0或C的alpha值.
本文：http://www.cnblogs.com/yuzhuwei/p/4141713.html