摘要：本篇內容為機器學習實戰第章支持向量機部分程序清單。支持向量機優點泛化錯誤率低，計算開銷不大，結果易解釋。注以上給出的僅是簡化版算法的實現，關于完整的算法加速優化并應用核函數，請參照機器學習實戰第頁。

本篇內容為《機器學習實戰》第 6 章 支持向量機部分程序清單。所用代碼為 python3。

支持向量機
優點：泛化錯誤率低，計算開銷不大，結果易解釋。 
缺點：對參數調節和核函數的選擇敏感，原始分類器不加修改僅適用于處理二分類問題。
適用數據類型：數值型和標稱型數據。

1996 年，John Platt 發布了一個稱為SMO的強大算法，用于訓練 SVM。SMO表示序列最小優化 (Sequential Minimal Optimization)。

SMO算法的工作原理是：每次循環中選擇兩個alpha進行優化處理。一旦找到一對合適的alpha，那么就增大其中一個同時減小另一個。這里的“合適”是指兩個alpha必須要符合一定的條件，第一個條件是這兩個alpha必須要在間隔邊界之外，第二個條件是這兩個alpha還沒有進行過區間化處理或者不在邊界上。

下面給出簡化版的SMO算法程序清單。

該SMO函數的偽代碼如下：
創建一個alpha向量并將其初始化為 0 向量
當迭代次數小于最大迭代次數時（外循環）
···對數據集中的每個數據向量（內循環）：
······如果該數據向量可以被優化：
·········隨機選擇另外一個數據向量
·········同時優化這兩個向量
·········如果兩個向量都不能被優化，退出內循環
···如果所有向量都沒被優化，增加迭代數目，繼續下一次循環

# coding=utf-8
# import sys
from numpy import *

def loadDataSet():
    dataMat = []
    labelMat = []
    fr = open("testSet.txt")

    for line in fr.readlines():
        lineArr = line.strip().split("	")
        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
    while (j == i):
        j = int(random.uniform(0, m))
    return j

# 用于調整大于 H 或小于 L 的 alpha 值
def clipAlpha(aj, H, L):
    if aj > H:
        aj = H
    if L > aj:
        aj = L
    return aj

在 python 提示符下，執行代碼并得到結果：

>>> import svmMLiA
>>> dataArr, labelArr = svmMLiA.loadDataSet()
>>> 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]

可以看出，這里采用的類別標簽是 -1 和 1。

# 參數：數據集，類別標簽，常數C，容錯率，退出前最大的循環次數
def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
    dataMatrix = mat(dataMatIn)
    # 由于轉置了類別標簽，我們得到的是一個列向量而不是列表
    labelMat = mat(classLabels).transpose()
    b = 0
    m,n = shape(dataMatrix)
    # 構建一個 alpha 列矩陣，矩陣中元素都初始化為0
    alphas = mat(zeros((m, 1)))
    # iter 變量存儲的是在沒有任何 alpha 改變的情況下便利數據集的次數
    # 當這個變量達到輸入值 maxIter 時，函數結束運行并退出
    iter = 0

    while(iter < maxIter):

        # 每次循環當中，將 alphaPairsChanged 先設為0，在對整個集合順序遍歷
        # 變量 alphaPairsChanged 用于記錄 alpha 是否已經進行優化
        alphaPairsChanged = 0
        for i in range(m):
            # 計算 fXi，即我們預測的類別
            fXi = float(multiply(alphas, labelMat).T * (dataMatrix*dataMatrix[i,:].T) + b)
            # 與真實值比對，計算誤差 Ei
            Ei = fXi - float(labelMat[i])
            # 如果誤差很大，可以對該數據實例所對應的 alpha 值進行優化

            # 不論正間隔還是負間隔都會被測試
            # 檢查 alpha 值，保證其不能等于 0 或 C
            if((labelMat[i]*Ei < -toler) and (alphas[i] < C)
               or (labelMat[i]*Ei > toler) and (alphas[i] > 0)):

                # 用輔助函數 selectJrand 隨機選擇第二個 alpha 值，即 alpha[j]
                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()

                # 計算 L 和 H，調整 alpha 到 0 與 C 之間
                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 是 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

                # 計算出一個新的 alpha[j]，并進行調整
                alphas[j] -= labelMat[j] * (Ei - Ej) / eta
                alphas[j] = clipAlpha(alphas[j], H, L)

                # 檢查 alpha[j] 是否有輕微改變，是的話則退出 for 循環
                if(abs(alphas[j] - alphaJold) < 0.00001):
                    print("j not moving enough")
                    continue

                # 對 alpha[i] 進行和 alpha[j] 同樣的改變
                # 改變的大小一樣，方向正好相反
                alphas[i] += labelMat[j] * labelMat[i] * (alphaJold - alphas[j])

                # 對 alpha[i] 和 alpha[j] 進行優化之后，給它們設置一個常數項 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

        if(alphaPairsChanged == 0):
            iter += 1
        else:
            iter = 0
        print("iteration number: %d" % iter)
    return b, alphas

在 python 提示符下，執行代碼并得到結果：

b, alphas = svmMLiA.smoSimple(dataArr, labelArr, 0.6, 0.001, 40)

再執行：

>>> for i in range(100):
...     if alphas[i] > 0.0:
...             print(dataArr[i], labelArr[i])
... 
[3.542485, 1.977398] -1.0
[7.108772, -0.986906] 1.0
[4.658191, 3.507396] -1.0
[7.40786, -0.121961] 1.0
[3.457096, -0.082216] -1.0
[5.286862, -2.358286] 1.0
[6.080573, 0.418886] 1.0
[6.543888, 0.433164] 1.0
[1.966279, -1.840439] -1.0

所輸出的數據點即為支持向量。

注：以上給出的僅是簡化版SMO算法的實現，關于完整的SMO算法加速優化并應用核函數，請參照《機器學習實戰》第 99 頁。

不足之處，歡迎指正。