前言
前面我們進行了很多的理論性研究，下面我們開始用代碼進行實現(xiàn)。
編程求解線性SVM
- 可視化數(shù)據(jù)集
import matplotlib.pyplot as pltimport numpy as np#讀取數(shù)據(jù)def loadDataSet(fileName): dataMat = []; labelMat = [] fr = open(fileName) for line in fr.readlines(): #逐行讀取，濾除空格等 lineArr = line.strip().split(”) dataMat.append([float(lineArr[0]),float(lineArr[1])]) #添加數(shù)據(jù) labelMat.append(float(lineArr[2])) #添加標簽 return dataMat,labelMat#數(shù)據(jù)可視化def showDataSet(dataMat,labelMat): data_plus = [] data_minus = [] for i in range(len(dataMat)): if labelMat[i] > 0: data_plus.append(dataMat[i]) else: data_minus.append(dataMat[i]) data_plus_np = np.array(data_plus) #轉(zhuǎn)換為numpy矩陣 data_minus_np = np.array(data_minus)#轉(zhuǎn)換為numpy矩陣 plt.scatter(np.transpose(data_plus_np)[0],np.transpose(data_plus_np)[1]) #正樣本 plt.scatter(np.transpose(data_minus_np)[0],np.transpose(data_minus_np)[1]) #負樣本 plt.show()dataMat,labelMat = loadDataSet(‘testSet.txt’)showDataSet(dataMat,labelMat)
這個數(shù)據(jù)集顯然線性可分。
應(yīng)用簡化版SMO算法處理小規(guī)模數(shù)據(jù)集
import randomdef selectJrand(i,m): j=i while(j==i): #選擇一個不等于i的j j = int(random.uniform(0,m)) return jdef clipAlpha(aj,H,L): if aj > H: aj = H if L > aj: aj = L return ajdataArr,labelArr =loadDataSet(‘testSet.txt’)labelArr
[-1.0, -1.0, 1.0, -1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, 1.0, -1.0, -1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0, -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0]
data = np.mat(dataArr)data[2,:]matrix([[ 7.55151, -1.58003]])
可以看出來，這里使用的類別標簽是-1和1
SMO算法的偽代碼：
創(chuàng)建一個alpha向量并將其初始化為 0 向量 當(dāng)?shù)?span id="29a4o8x" class="wpcom_tag_link">次數(shù)小于最大迭代次數(shù)時 (外循環(huán)) 對數(shù)據(jù)集中的每個數(shù)據(jù)向量 (內(nèi)循環(huán)): 如果該數(shù)據(jù)向量可以被優(yōu)化: 隨機選擇另外一個數(shù)據(jù)向量 同時優(yōu)化這兩個向量 如果兩個向量都不能被優(yōu)化, 退出內(nèi)循環(huán)如果所有向量都沒被優(yōu)化, 增加迭代數(shù)目, 繼續(xù)下一次循環(huán)#簡化版SMO算法def smoSimple(dataMatIn,classLabels,C,toler,maxIter): dataMatrix = np.mat(dataMatIn); labelMat = np.mat(classLabels).transpose() b = 0; m,n = np.shape(dataMatrix) alphas = np.mat(np.zeros((m,1)))#初始化alpha參數(shù)，設(shè)置為0 iterSmo = 0 #初始化迭代次數(shù) while(iterSmo < maxIter): alphaPairsChanged = 0#用于記錄alpha是否已經(jīng)進行優(yōu)化 #步驟1. 計算誤差Ei for i in range(m): fXi = float(np.multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + b Ei = fXi – float(labelMat[i]) if ((labelMat[i]*Ei 0)): j = selectJrand(i,m) #步驟1. 計算誤差Ej fXj = float(np.multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[j,:].T)) + b Ej = fXj – float(labelMat[j]) #保存更新前的alpha值，使用淺拷貝 alphaIold = alphas[i].copy() alphaJold = alphas[j].copy() #步驟2：計算上界H和下界L 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 #步驟3：計算eta 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 #步驟4：更新alpha_j alphas[j] -= labelMat[j]*(Ei-Ej)/eta #步驟5：修剪alpha_j alphas[j] = clipAlpha(alphas[j],H,L) if (abs(alphas[j] – alphaJold) < 0.00001) : print("j not moving enough") ; continue #步驟6：更新alpha_i alphas[i] += labelMat[j]*labelMat[i]*(alphaJold – alphas[j]) #步驟7：更新b_1和b_2 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 #步驟8：根據(jù)b_1和b_2更新b if (0 alphas[i]) : b = b1 elif (0 alphas[j]): b = b2 else: b = (b1 + b2)/2.0 #統(tǒng)計優(yōu)化次數(shù) #如果程序執(zhí)行到for循環(huán)的最后一行都不執(zhí)行continue語句，那么就已經(jīng)成功地改變了一對alpha，同時可以增加alphaPairsChanged的值 alphaPairsChanged += 1 #打印統(tǒng)計信息 print("第%d次迭代 樣本:%d, alpha優(yōu)化次數(shù):%d" % (iterSmo,i,alphaPairsChanged)) #更新迭代次數(shù) #在for循環(huán)之外，需要檢查alpha值是否做了更新，如果有更新則將iterSmo設(shè)為0后繼續(xù)運行程序。只有在所有數(shù)據(jù)集上遍歷maxIter次，且不再發(fā)生任何alpha修改之后，程序才會停止并退出while循環(huán) if (alphaPairsChanged == 0): iterSmo += 1 else: iterSmo = 0 print("迭代次數(shù)：%d" % iterSmo) return b,alphasb,alphas = smoSimple(dataArr,labelArr,0.6,0.001,500)
L==H第0次迭代 樣本:1, alpha優(yōu)化次數(shù):1第0次迭代 樣本:3, alpha優(yōu)化次數(shù):2第0次迭代 樣本:5, alpha優(yōu)化次數(shù):3L==H第0次迭代 樣本:8, alpha優(yōu)化次數(shù):4L==Hj not moving enoughj not moving enoughL==HL==Hj not moving enoughL==H第0次迭代 樣本:30, alpha優(yōu)化次數(shù):5第0次迭代 樣本:31, alpha優(yōu)化次數(shù):6L==HL==H第0次迭代 樣本:54, alpha優(yōu)化次數(shù):7L==HL==H第0次迭代 樣本:71, alpha優(yōu)化次數(shù):8L==HL==HL==H第0次迭代 樣本:79, alpha優(yōu)化次數(shù):9L==H第0次迭代 樣本:92, alpha優(yōu)化次數(shù):10j not moving enoughL==H迭代次數(shù)：0第0次迭代 樣本:1, alpha優(yōu)化次數(shù):1j not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enoughL==HL==Hj not moving enoughj not moving enoughj not moving enoughj not moving enoughL==Hj not moving enoughL==HL==Hj not moving enoughj not moving enough第0次迭代 樣本:37, alpha優(yōu)化次數(shù):2第0次迭代 樣本:39, alpha優(yōu)化次數(shù):3第0次迭代 樣本:52, alpha優(yōu)化次數(shù):4j not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enough第0次迭代 樣本:71, alpha優(yōu)化次數(shù):5j not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enough迭代次數(shù)：0j not moving enoughj not moving enoughj not moving enough第0次迭代 樣本:8, alpha優(yōu)化次數(shù):1L==Hj not moving enough第0次迭代 樣本:23, alpha優(yōu)化次數(shù):2L==Hj not moving enoughj not moving enoughL==Hj not moving enoughj not moving enoughj not moving enough第0次迭代 樣本:39, alpha優(yōu)化次數(shù):3L==Hj not moving enough第0次迭代 樣本:52, alpha優(yōu)化次數(shù):4j not moving enough第0次迭代 樣本:55, alpha優(yōu)化次數(shù):5L==HL==HL==HL==HL==Hj not moving enough第0次迭代 樣本:79, alpha優(yōu)化次數(shù):6第0次迭代 樣本:92, alpha優(yōu)化次數(shù):7迭代次數(shù)：0j not moving enoughL==Hj not moving enoughj not moving enoughL==Hj not moving enough第0次迭代 樣本:23, alpha優(yōu)化次數(shù):1j not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enoughL==HL==H第0次迭代 樣本:51, alpha優(yōu)化次數(shù):2j not moving enoughj not moving enoughj not moving enoughj not moving enoughL==H第0次迭代 樣本:69, alpha優(yōu)化次數(shù):3L==Hj not moving enough第0次迭代 樣本:94, alpha優(yōu)化次數(shù):4j not moving enoughj not moving enough迭代次數(shù)：0j not moving enoughj not moving enoughj not moving enoughj not moving enoughj not moving enough
…
迭代次數(shù)：497j not moving enoughj not moving enoughj not moving enough迭代次數(shù)：498j not moving enoughj not moving enoughj not moving enough迭代次數(shù)：499j not moving enoughj not moving enoughj not moving enough迭代次數(shù)：500
bmatrix([[-3.83785102]])alphas[alphas>0]matrix([[0.1273855 , 0.24131542, 0.36872064]])np.shape(alphas[alphas>0])(1, 3)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#分類結(jié)果可視化def showClassifer(dataMat, w, b): #繪制樣本點 data_plus = [] #正樣本 data_minus = [] #負樣本 for i in range(len(dataMat)): if labelMat[i] > 0: data_plus.append(dataMat[i]) else: data_minus.append(dataMat[i]) data_plus_np = np.array(data_plus) #轉(zhuǎn)換為numpy矩陣 data_minus_np = np.array(data_minus) #轉(zhuǎn)換為numpy矩陣 plt.scatter(np.transpose(data_plus_np)[0], np.transpose(data_plus_np)[1], s=30, alpha=0.7) #正樣本散點圖 plt.scatter(np.transpose(data_minus_np)[0], np.transpose(data_minus_np)[1], s=30, alpha=0.7) #負樣本散點圖 #繪制直線 x1 = max(dataMat)[0] x2 = min(dataMat)[0] a1, a2 = w b = float(b) a1 = float(a1[0]) a2 = float(a2[0]) y1, y2 = (-b- a1*x1)/a2, (-b – a1*x2)/a2 plt.plot([x1, x2], [y1, y2]) #找出支持向量點 for i, alpha in enumerate(alphas): if abs(alpha) > 0: x, y = dataMat[i] plt.scatter([x], [y], s=150, c=’none’, alpha=0.7, linewidth=1.5, edgecolor=’red’) plt.show()#計算wdef get_w(dataMat, labelMat, alphas): alphas, dataMat, labelMat = np.array(alphas), np.array(dataMat), np.array(labelMat) w = np.dot((np.tile(labelMat.reshape(1, -1).T, (1, 2)) * dataMat).T, alphas) return w.tolist()w = get_w(dataMat,labelMat,alphas)showClassifer(dataMat,w,b)
