機器學習——支持向量機(SVM)之拉格朗日乘子法，KKT條件以及簡化版SMO算法分析

SVM有很多實現，現在只關注其中最流行的一種實現，即序列最小優化（Sequential Minimal Optimization，SMO）算法，然后介紹如何使用一種核函數（kernel）的方式將SVM擴展到更多的數據集上。

 

1.基於最大間隔分隔數據

幾個概念：

1.線性可分（linearly separable）：對於圖6-1中的圓形點和方形點，如果很容易就可以在圖中畫出一條直線將兩組數據點分開，就稱這組數據為線性可分數據

2.分隔超平面（separating hyperplane）：將數據集分隔開來的直線稱為分隔超平面

3.如果數據集是1024維的，那么就需要一個1023維的超平面來對數據進行分隔

4.間隔（margin）：數據點到分隔面的距離稱為間隔

5.支持向量（support vector）:離分隔超平面最近的那些點

 

支持向量機的優點：泛化錯誤率低，計算開銷不大，結果易解釋

支持向量機的缺點：對參數調節和核函數的選擇敏感，原始分類器不加修改僅適用於處理二類問題

適用數據類型：數值型和標稱型數據

 

2.尋找最大間隔

如何求解數據集的最佳分隔直線？

分隔超平面的形式可以寫成  

其中 w = (w1,w2,w3...wd)為法向量，決定了超平面的方向，其中d等於數據的維度，

這很好理解，假設二維的(x1,x2)點可以被 ax+b=0 分隔，這里面直線 ax+b=0 是一維的，但是這里面a和x都是二維的

b為位移項，決定了超平面與原點之間的距離

對於圖6-3中A點到分隔直線的距離為 

表示向量的模，，w與w共軛的內積再開方

假設超平面（w,b）能將訓練樣本正確分類，即對於 ,

有 

則兩個異類支持向量到超平面的距離之和為

欲找到具有“最大間隔（maximum margin）”的划分超平面，也就是要找到能滿足 中約束的參數w和b，使得最大，即

  　　  ，

其中約束條件為　s.t.  ，其實這個約束條件就是把兩個不等式合並成了一個

顯然，為了最大化間隔，僅需最大化 ，這等價於最小化 ，於是上式可重寫為

　　，

其中約束條件為　s.t. 

這就是支持向量機（Support Vector Machine，簡稱SVM）的基本型

 

對於這類帶有不等式約束的最優化問題，可以使用拉格朗日乘子法（Lagrange Multiplier）對其進行求解。

首先先了解一下最優化問題的分類，最優化問題可以分為一下三類：

<1>無約束的優化問題，可以寫成：

      

對於第<1>類的優化問題，常常使用的方法就是Fermat定理，即使用求取f(x)的導數，然后令其為零，可以求得候選最優值，再在這些候選值中驗證；如果是凸函數，可以保證是最優解。

<2>有等式約束的優化問題，可以寫成：

　　

　　約束條件 

對於第<2>類的優化問題，常常使用的方法就是拉格朗日乘子法（Lagrange Multiplier) ，即把等式約束 用一個系數與目標函數f(x)寫為一個式子，稱為拉格朗日函數，而系數稱為拉格朗日乘子。拉格朗日函數的形式如下，u即拉格朗日乘子：

　　

通過拉格朗日函數對各個變量求導，令其為零，可以求得候選值集合，然后驗證求得最優值。

<3>有不等式約束的優化問題，可以寫成：

　　

　　約束條件  

　　　　　　　　　

對於第<3>類的優化問題，常常使用的方法就是KKT條件。同樣地，我們把所有的等式、不等式約束與f(x)寫為一個式子，也叫拉格朗日函數，系數也稱拉格朗日乘子，通過一些條件，可以求出最優值的必要條件，這個條件稱為KKT條件。

 

拉格朗日乘子法（Lagrange Multiplier）：

 拉格朗日乘子法可以應用於有等式約束的優化問題和有不等式約束的優化問題，

對於有等式約束的優化問題， ，約束條件  ，其中 f(x) 被稱為目標函數

<1>當沒有約束條件的時候，我們通過求導的方法，來尋找最優點

設 是這個最優點，即此時有

　　 

此外，如果 f(x) 是一個實值函數， 是一個n維向量的話，那么 f(x) 對向量 的導數被定義為

　　

 

<2>當有一個等式約束條件 的時候，舉例設

　　

　　

從幾何的角度看，可以看成是在一個曲面 上尋找函數 的最小值

設目標函數  ，當 z 取不同的值的時候，相當於可以投影在曲面 上，即形成等高線，如下圖

當 取d1或者d2的時候，形成了虛線的等高線，此時我們要求的是最小的z

現在我們約束 的時候，在曲面 上形成一條曲線，即圖中紅色的曲線

假設形成的紅色曲線與等高線相交，交點就是同時滿足等式約束條件和目標函數的可行域的值，但肯定不是最優值，因為相交意味着肯定還存在其它的等高線在該條等高線的內部或者外部，使得新的等高線與目標函數的交點的值更大或者更小。

只有到等高線與目標函數的曲線相切的時候，可能取得最優值，如上圖中的情況，即等高線和目標函數的曲線在該點的法向量必須有相同方向，所以最優值必須滿足：f(x)的梯度 = * h(x)的梯度， 是常數，表示左右兩邊同向。這個等式就是拉格朗日函數 對參數x求導后，令其等於0的結果，即，最后求得最優解是 。

 

以上的約束條件只是 ，這時約束條件是一個等式。

<3>那么當約束條件是多個等式的時候，拉格朗日函數應該改寫成

，是需要求其最小值的目標函數，是拉格朗日乘子，是約束條件

 

<4>如果約束條件是等式和不等式，且有多個的時候，就需要若干個KKT(Karush-Kuhn-Tucker)條件，即對於

　　

　　約束條件  

　　　　　　　　　

想取得最優解，要滿足以下條件（KKT條件）：

　　1. L(x, a, b)對x求導為零；

　　2. h(x) =0;

　　3. a*g(x) = 0;

 

一步步來說明，此時我們要求f(x)的最小值，可以構建拉格朗日函數

　　

為什么這么構建呢？因為當 g(x)<=0，h(x)=0，a>=0 的時候，（注意此處的h(x)=0是KKT條件的第2個）

a*g(x)<=0，所以 在取得最大值的時候，即a*g(x)=0的時候，就等於f(x)，（注意此處的a*g(x)=0是KKT條件的第3個）

所以需要求解的目標函數的最小值寫成表達式是

　　

上式的對偶表達式是

　　

由於我們要求的最優化是滿足強對偶的，即對偶式子的最優值等於原式子的最優解

設當 的時候有最優解，此時

　　

　　　　　

　　　　　

即函數 在 處取得最小值

用fermat定理，對於函數 求取導數后令其等於0后的結果，（注意此處的導數=0是KKT條件的第1個）

即f(x)的梯度 + a*g(x)的梯度 + b*h(x)的梯度 = 0

 

所以對於多個不等式和等式混合的約束條件的時候，拉格朗日函數可以寫成

　　

　　約束條件  

　　　　　　　　　

從而上面提到的支持向量機（Support Vector Machine，簡稱SVM）的基本型

　　，其中約束條件為　s.t. 

就可以重新寫成

　　

其中 ，

接下來由第1個KKT條件（導數等於0），令 對 和 b 的偏導等於0，可得

　　

把上面兩個等式帶入 中，得

　　  的對偶表達式等於

　　

　　

約束條件

在最后推出的公式中，y表示的是標簽向量，u表示的是拉格朗日乘子，x表示的是待分類的點的坐標，也是一個向量

現在約束條件變成了等式，所以我們的任務變成了

　　

約束條件

 

以上的推導有一個假設條件，就是數據是線性可分的，而實際上，大多數情況是線性不可分的，對於這種情況：

我們把約束條件從 改成 ，允許出現例外

所以支持向量機（Support Vector Machine，簡稱SVM）的基本型就需要改寫成，加上了松弛變量（這里面有泛函的數學原理，參考《神經網絡與機器學習 第3版》）

　　，其中約束條件為

通過和上面一樣的對 w,b和求偏導的過程，進行簡化之后，也是得到

　　

但是約束條件變成了

 

現在可以應用上面的等式約束條件的方法，求出 u 后，在求出w與b即可得到模型

　　 

 

2.SMO高效優化算法

接下來分析SMO算法，因為有約束條件 ，所以一次不能只更新一個 ，因為如果一次只更新一個 的話，那么累加之后的結果就不能等於0，這就不滿足約束條件了，所以SMO算法采取的是挑選一對 進行更新，這樣就能滿足約束條件。

舉例，當更新 的時候，一次同時更新 和 ，根據約束條件 ，

那么 ，即1和2相加，等於3到m相加的負數。

圖來自網上

當y1和y2異號的時候，即y1=1,y2=-1或者y1=-1,y2=1的時候，我們先固定一個 的值，然后求 的范圍

<1>設y1=1,y2=-1

則如上圖所示，橫坐標軸是 ，縱坐標軸是 ，因為 ，且  是大於等於0的，

所以當 和 有可能滿足紅線的情況，即y1=1,y2=-1時，，那么 的范圍是 ，即

所以當 和 也有可能滿足藍線的情況，即y1=1,y2=-1時，，那么 的范圍是 ，即 

由於不知道是紅線還是藍線，所以 的取值范圍為：

　　

 <2>設y1=-1,y2=1

同樣可能通過畫圖來說明，計算出來的 的取值范圍仍然是：

　　

綜上所述，當y1和y2異號的時候， 的取值范圍為：

　　

 

同理，當y1和y2同號的時候， 的取值范圍為：

　　

 

以上內容的參考文獻：

支持向量機入門系列-2：等式約束極小的最優性條件

深入理解拉格朗日乘子法（Lagrange Multiplier) 和KKT條件

《機器學習》——周志華

《神經網絡與機器學習（第3版）》

《機器學習實戰》

 

以下內容參考的網址：

支持向量機（SVM），SMO算法原理及源代碼剖析

《機器學習實戰》

 

我們需要求解的問題是

約束條件為 

現在把求最大值轉換成求最小值，其實是一樣的，現在問題變成了

約束條件為 

 

因為

所以

接下來將 和 代入

得

 

以下推導內容是CSDN上一篇博文的內容，直接貼圖片，符號表示有些許區別

 

對最小化的目標函數和約束條件進行優化，以下是Python代碼對應的偽代碼流程

 

 

 

源代碼及代碼注釋：

 

from numpy import *
from time import sleep
import time

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

#輔助函數，用於在某一區間范圍內隨機選擇一個整數
#i是第一個alpha的下標，m是所有alpha的數目，只要函數值不等於輸入值i，函數就會進行隨機選擇
def selectJrand(i,m):	
    j=i 		#we want to select any J not equal to i
    while (j==i):
        j = int(random.uniform(0,m))
    return j

def clipAlpha(aj,H,L):	#輔助函數，數值太大或者太小時對其進行調整
    if aj > H: 
        aj = H
    if L > aj:
        aj = L
    return aj

#Platt SMO算法中的外循環確定要優化的最佳aplha集合，而簡化版卻會跳過這一部分
#首先在數據集上遍歷每一個alpha，然后在剩下的alpha集合中隨機選擇另一個alpha，構建alpha對
#簡化版SMO算法，數據集，類別標簽，常數C，容錯率，取消前最大的循環次數
def smoSimple(dataMatIn, classLabels, C, toler, maxIter):	
    dataMatrix = mat(dataMatIn); labelMat = mat(classLabels).transpose()	#把列表dataMatIn和classLabels轉換成矩陣dataMatrix和labelMat
    b = 0; m,n = shape(dataMatrix)	#初始的b=0,m=100,n=2
    alphas = mat(zeros((m,1)))		#生成一個100×1的零矩陣
    iter = 0
    while (iter < maxIter):	#在iter小於maxIter的時候循環
        alphaPairsChanged = 0	#用於記錄alpha是否已經進行優化
        for i in range(m):	#循環所有數據集，總共有100個
            fXi = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + b	#fXi是預測的類別，i是alpha的下標，把第i個x代入公式計算
            Ei = fXi - float(labelMat[i])						#Ei是第i數據向量的計算誤差
            #判斷每一個alpha是否被優化過，如果誤差很大，就對該alpha值進行優化，toler是容錯率
            if ((labelMat[i]*Ei < -toler) and (alphas[i] < C)) or ((labelMat[i]*Ei > toler) and (alphas[i] > 0)):	
                j = selectJrand(i,m)								#隨機選擇另外一個數據向量j
                print "隨機選擇的j下標是",j
                fXj = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[j,:].T)) + b	#同時優化這兩個向量，如果都不能被優化，退出內循環
                Ej = fXj - float(labelMat[j])							#Ej是第j數據向量的計算誤差
                alphaIold = alphas[i].copy();alphaJold = alphas[j].copy();
                #保證alpha在0與C之間
                if (labelMat[i] != labelMat[j]):		#當y1和y2異號，計算alpha的取值范圍 	
                    L = max(0, alphas[j] - alphas[i])
                    H = min(C, C + alphas[j] - alphas[i])
                else:						#當y1和y2同號，計算alpha的取值范圍
                    L = max(0, alphas[j] + alphas[i] - C)
                    H = min(C, alphas[j] + alphas[i])
                #如果L和H相等，就不做任何改變，本次循環結束直接運行下一次for循環
                if L==H: print "L==H"; continue			
                #eta是alpha[j]的最優修改量，eta=K11+K22-2*K12,也是f(x)的二階導數，K表示內積
                eta = 2.0 * dataMatrix[i,:]*dataMatrix[j,:].T - dataMatrix[i,:]*dataMatrix[i,:].T - dataMatrix[j,:]*dataMatrix[j,:].T
                #如果二階導數-eta <= 0，說明一階導數沒有最小值，就不做任何改變，本次循環結束直接運行下一次for循環
                if eta >= 0: print "eta>=0"; continue
                alphas[j] -= labelMat[j]*(Ei - Ej)/eta		#利用公式更新alpha[j]，alpha2new=alpha2-yj(Ei-Ej)/eta
                alphas[j] = clipAlpha(alphas[j],H,L)		#再判斷一次alpha[j]的范圍
                #如果alphas[j]沒有調整，就忽略下面語句，本次循環結束直接運行下一次for循環
                if (abs(alphas[j] - alphaJold) < 0.00001): print "j改變的不夠"; continue
                alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])	#調整alphas[i]，修改量與j相同，但是方向相反
                #已經計算出了alpha，接下來根據模型的公式計算b
                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
                #根據公式確定偏移量b，理論上可選取任意支持向量來求解，但是現實任務中通常使用所有支持向量求解的平均值，這樣更加魯棒
                if (0 < alphas[i]) and (C > alphas[i]): b = b1
                elif (0 < alphas[j]) and (C > alphas[j]): b = b2
                else: b = (b1 + b2)/2.0
                alphaPairsChanged += 1
                print "總共100次迭代的次數: %d 當前i的下標:%d, 成對改變的alpha的數量 %d" % (iter,i,alphaPairsChanged)
                #print "輸出現在的alphas",alphas.T
                
	#只有滿足alphaPairsChanged為0的條件，即100個alpha的誤差都在范圍之內，才能進入下一個遞歸，不然只能在100個alpha再次循環優化
        if (alphaPairsChanged == 0): iter += 1
        else: iter = 0
        print "iteration number: %d" % iter
        plotBestFit(alphas,dataArr,labelArr,b)
    return b,alphas

 

 

 

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

def plotBestFit(alphas,dataArr,labelArr,b):	#畫出數據集和超平面
    import matplotlib.pyplot as plt
    dataMat,labelMat=loadDataSet('testSet.txt')	#數據矩陣和標簽向量
    dataArr = array(dataMat)			#轉換成數組
    n = shape(dataArr)[0] 
    xcord1 = []; ycord1 = []			#聲明兩個不同顏色的點的坐標
    xcord2 = []; ycord2 = []
    for i in range(n):
        if int(labelMat[i])== 1:
            xcord1.append(dataArr[i,0]); ycord1.append(dataArr[i,1])
        else:
            xcord2.append(dataArr[i,0]); ycord2.append(dataArr[i,1])
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(xcord1, ycord1, s=30, c='red')
    ax.scatter(xcord2, ycord2, s=30, c='green', marker='s')
    x = arange(2.0, 6.0, 0.1)
    ws = calcWs(alphas,dataArr,labelArr)
    #最佳擬合曲線，因為0是兩個分類（0和1）的分界處（Sigmoid函數）
    #圖中y表示x2,x表示x1
    y = [(i*ws[0]+array(b)[0])/-ws[1] for i in x]
    ax.plot(x, y)
    plt.xlabel('X1'); plt.ylabel('X2');
    plt.show()

 

if __name__ == '__main__':
    dataArr,labelArr = loadDataSet('testSet.txt')
    b,alphas = smoSimple(dataArr,labelArr,0.6,0.001,40)
    plotBestFit(alphas,dataArr,labelArr,b)

 

用超平面對數據集進行划分