新的代碼修改了下代碼風格鲫骗，盡量參照了Google的Python風格規(guī)范洪囤。之前寫的沒去改，可能有點不一致消玄。
SVM看了好長時間，卡了很久，數學知識嚴重欠缺。以下是個人的一些拙見采转，很長很丑，希望對你們也有些幫助瞬痘，寫的不對的地方敬請指正故慈。

支持向量機（support vector machines，SVM）是一種二分類模型图云，它的基本模型是定義在特征空間上的間隔最大化的線性分類器惯悠，間隔使它有別于感知機，核技巧使它成為實質上的非線性分類器竣况。

SVM包括線性可分支持向量機、線性支持向量機筒严、非線性支持向量機丹泉。當訓練數據線性可分時，通過硬間隔最大化（hard margin maximization）學習線性可分支持向量機鸭蛙；當訓練數據近似線性可分時摹恨，通過軟間隔最大化（soft margin maximization）學習線性支持向量機；當訓練數據線性不可分時娶视，通過核技巧及軟間隔最大化學習非線性支持向量機晒哄。

與感知機區(qū)別

一、SVM推導

二肪获、SMO算法

貼上第二節(jié)全部代碼

from numpy import *
import numpy as np
import matplotlib.pyplot as plt


def loadDataSet(filename):
    """
        Function:
            讀取數據
        Parameters:
            fileName - 文件名
        Returns:
            dataMat - 數據矩陣
            labelMat - 數據標簽
        Modify:
            2018-09-21
    """
    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[-1]))
    return dataMat, labelMat


def selectJrand(i, m):
    """
        Function:
            隨機選擇alphaJ
        Parameters:
            i - alpha下標
            m - alpha參數個數
        Returns:
            j - 隨機選擇alpha另一個下標
        Modify:
            2018-09-21
    """
    j = i
    while (j == i):
        j = int(random.uniform(0, m))
    return j


def clip_alpha(aj, H, L):
    """
        Function:
            修剪alpha
        Parameters:
            aj - alpha值
            H - alpha上限
            L - alpha下限
        Returns:
            aj - alpah值
        Modify:
            2018-09-21
    """
    if aj > H:
        aj = H
    if L > aj:
        aj = L
    return aj


def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
    """
        Function:
            簡化版SMO算法
        Parameters:
            dataMatIn - 數據矩陣
            classLabels - 數據標簽
            C - 松弛變量
            toler - 容錯率
            maxIter - 最大迭代次數
        Returns:
            無
        Modify:
            2018-09-21
    """
    # 轉換為numpy的mat存儲
    dataMatrix = np.mat(dataMatIn)
    # 類別標簽列向量
    labelMat = np.mat(classLabels).transpose()
    # 初始化b參數寝凌，統(tǒng)計dataMatrix的維度
    b = 0
    m, n = np.shape(dataMatrix)
    # 初始化alpha參數，設為0
    alphas = np.mat(np.zeros((m, 1)))
    # 初始化迭代次數
    iter_num = 0
    # 最多迭代matIter次
    while (iter_num < maxIter):
        alphaPairsChanged = 0
        for i in range(m):
            # 步驟1：計算誤差Ei
            fXi = float(np.multiply(alphas, labelMat).T * (dataMatrix * dataMatrix[i, :].T)) + b
            Ei = fXi - float(labelMat[i])
            # 優(yōu)化alpha孝赫，更設定一定的容錯率较木。
            if ((labelMat[i] * Ei < -toler) and (alphas[i] < C)) or ((labelMat[i] * Ei > toler) and (alphas[i] > 0)):
                # 隨機選擇另一個與alpha_i成對優(yōu)化的alpha_j
                j = selectJrand(i, m)
                # 步驟1：計算誤差Ej
                fXj = float(np.multiply(alphas, labelMat).T * (dataMatrix * dataMatrix[j, :].T)) + b
                Ej = fXj - float(labelMat[j])
                # 保存更新前的aplpha值，使用深拷貝
                alphaIold = alphas[i].copy()
                alphaJold = alphas[j].copy()
                # 步驟2：計算上下界L和H
                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] = clip_alpha(alphas[j], H, L)
                if (abs(alphas[j] - alphaJold) < 0.00001):
                    print('alpha_j變化太小')
                    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：根據b_1和b_2更新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
                # 統(tǒng)計優(yōu)化次數
                alphaPairsChanged += 1
                # 打印統(tǒng)計信息
                print("第%d次迭代樣本:%d, alpha優(yōu)化次數:%d" % (iter_num, i, alphaPairsChanged))
        # 更新迭代次數
        if (alphaPairsChanged == 0):
            iter_num += 1
        else:
            iter_num = 0
        print("迭代次數: %d" % iter_num)
    return b, alphas


def show_classifier(dataMat, class_labels,w, b, alphas):
    """
        Function:
            可視化分類結果
        Parameters:
            dataMat - 數據矩陣
            w - 直線法向量
            b - 直線截距
            alphas - alphas值
        Returns:
            無
        Modify:
            2018-09-22
    """
    # 繪制樣本點
    data_plus = []
    data_minus = []
    for i in range(len(dataMat)):
        if class_labels[i] > 0:
            data_plus.append(dataMat[i])
        else:
            data_minus.append(dataMat[i])
    data_plus_np = np.array(data_plus)
    data_minus_np = np.array(data_minus)
    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 alpha > 0:
            x, y = dataMat[i]
            plt.scatter([x], [y], s=150, c='none', alpha=0.7, linewidth=1.5, edgecolor='red')
    plt.show()


def calc_ws(alphas, data_arr, class_labels):
    """
        Function:
            計算w
        Parameters:
            data_arr - 數據矩陣
            class_labels - 數據標簽
            alphas - alphas值
        Returns:
            w - 計算得到的w
        Modify:
            2018-09-22
    """
    X = np.mat(data_arr)
    label_mat = np.mat(class_labels).transpose()
    m, n = np.shape(X)
    w = np.zeros((n, 1))
    for i in range(m):
        w += np.multiply(alphas[i] * label_mat[i], X[i, :].T)
    return w


class OptStruct(object):
    def __init__(self, dataMatIn, classLabels, C, toler):
        """
            Parameters:
                dataMatIn - 數據矩陣
                classLabels - 數據標簽
                C - 松弛變量
                toler - 容錯率
        """
        self.X = dataMatIn
        self.label_mat = classLabels
        self.C = C
        self.tol = toler
        # 數據矩陣行數
        self.m = np.shape(dataMatIn)[0]
        # 根據矩陣行數初始化alpha參數為0
        self.alphas = np.mat(np.zeros((self.m, 1)))
        self.b = 0
        # 根據矩陣行數初始化誤差緩存青柄，第一列為是否有效的標志位伐债，第二列為實際的誤差E的值
        self.eCache = np.mat(np.zeros((self.m, 2)))


def calc_Ek(oS, k):
    """
        Function:
            計算誤差
        Parameters:
            oS - 數據結構
            k - 下標為k的數據
        Returns:
             Ek - 下標為k的數據誤差
        Modify:
            2018-09-22
    """
    fXk = float(np.multiply(oS.alphas, oS.label_mat).T * (oS.X * oS.X[k, :].T) + oS.b)
    Ek = fXk - float(oS.label_mat[k])
    return Ek


def select_J(i, oS, Ei):
    """
        Function:
            內循環(huán)啟發(fā)方式
        Parameters:
            i - 下標為i的數據的索引值
            oS - 數據結構
            Ei - 下標為i的數據誤差
        Returns:
            j, maxK - 下標為j或maxK的數據的索引值
            Ej - 下標為j的數據誤差
        Modify:
            2018-09-22
    """
    maxK = -1
    max_deltaE = 0
    Ej = 0
    oS.eCache[i] = [1, Ei]
    # 返回Ei不為0的索引值數組
    valid_eCache_list = nonzero(oS.eCache[:, 0].A)[0]
    # 有不為0的Ei
    if (len(valid_eCache_list)) > 1:
        for k in valid_eCache_list:
            if k == i: continue
            Ek = calc_Ek(oS, k)
            deltaE = abs(Ei - Ek)
            if (deltaE > max_deltaE):
                maxK = k
                max_deltaE = deltaE
                Ej = Ek
        return maxK, Ej
    # 沒有不為0的Ei
    else:
        j = selectJrand(i, oS.m)
        Ej = calc_Ek(oS, j)
    return j, Ej


def update_Ek(oS, k):
    """
        Function:
            計算Ek,并更新誤差緩存
        Parameters:
            oS - 數據結構
            k - 下標為k的數據的索引值
        Returns:
            無
        Modify:
            2018-09-22
    """
    Ek = calc_Ek(oS, k)
    # 更新Ei緩存
    oS.eCache[k] = [1, Ek]


def inner_L(i, oS):
    """
        Function:
            優(yōu)化的SMO算法
        Parameters:
            i - 下標為i的數據的索引值
            oS - 數據結構
        Returns:
            1 - 有任意一對alpha值發(fā)生變化
            0 - 沒有任意一對alpha值發(fā)生變化或變化太小
        Modify:
            2018-09-22
    """
    # 步驟1：計算誤差Ei
    Ei = calc_Ek(oS, i)
    if ((oS.label_mat[i] * Ei < -oS.tol) and (oS.alphas[i] < oS.C)) or ((oS.label_mat[i] * Ei > oS.tol) and (oS.alphas[i] > 0)):
        # 使用內循環(huán)啟發(fā)方式2選擇alpha_j,并計算Ej
        j, Ej = select_J(i, oS, Ei)
        # 保存更新前的aplpha值预侯，使用深拷貝
        alpha_i_old = oS.alphas[i].copy()
        alpha_j_old = oS.alphas[j].copy()
        # 步驟2：計算上下界L和H
        if (oS.label_mat[i] != oS.label_mat[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 0
            # 步驟3：計算eta
        eta = 2.0 * oS.X[i, :] * oS.X[j, :].T - oS.X[i, :] * oS.X[i, :].T - oS.X[j, :] * oS.X[j, :].T
        if eta >= 0:
            print("eta>=0")
            return 0
        # 步驟4：更新alpha_j
        oS.alphas[j] -= oS.label_mat[j] * (Ei - Ej) / eta
        # 步驟5：修剪alpha_j
        oS.alphas[j] = clip_alpha(oS.alphas[j], H, L)
        # 更新Ei至誤差緩存
        update_Ek(oS, j)
        if (abs(oS.alphas[j] - alpha_j_old) < 0.00001):
            print('alpha_j變化太小')
            return 0
        # 步驟6：更新alpha_i
        oS.alphas[i] += oS.label_mat[j] * oS.label_mat[i] * (alpha_j_old - oS.alphas[j])
        # 步驟7：更新b1和b2
        b1 = oS.b - Ei - oS.label_mat[i] * (oS.alphas[i] - alpha_i_old) * oS.X[i, :] * oS.X[i, :].T - \
             oS.label_mat[j] * (oS.alphas[j] - alpha_j_old) * oS.X[i, :] * oS.X[j, :].T
        b2 = oS.b - Ej - oS.label_mat[i] * (oS.alphas[i] - alpha_i_old) * oS.X[i, :] * oS.X[j, :].T - \
             oS.label_mat[j] * (oS.alphas[j] - alpha_j_old) * oS.X[j, :] * oS.X[j, :].T
        # 步驟8：根據b1和b2更新b
        if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]):
            oS.b = b1
        elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]):
            oS.b = b2
        else:
            oS.b = (b1 + b2) / 2.0
        return 1
    else:
        return 0


def smoP(data_arr, class_labels, C, toler, max_iter):
    """
    Function:
        完整的線性SMO算法
    Parameters：
        dataMatIn - 數據矩陣
        classLabels - 數據標簽
        C - 松弛變量
        toler - 容錯率
        maxIter - 最大迭代次數
    Returns:
        oS.b - SMO算法計算的b
        oS.alphas - SMO算法計算的alphas
    Modify:
        2018-09-23
    """
    oS = OptStruct(np.mat(data_arr), np.mat(class_labels).transpose(), C, toler)
    iter = 0
    entire_set = True
    alpha_pairs_changed = 0
    # 遍歷整個數據集alpha也都沒有更新或者超過最大迭代次數則退出循環(huán)
    while (iter < max_iter) and ((alpha_pairs_changed > 0) or (entire_set)):
        alpha_pairs_changed = 0
        # 遍歷整個數據集
        if entire_set:
            for i in range(oS.m):
                # 使用優(yōu)化的SMO算法
                alpha_pairs_changed += inner_L(i, oS)
                print("全樣本遍歷:第%d次迭代 樣本:%d, alpha優(yōu)化次數:%d" % (iter, i, alpha_pairs_changed))
            iter += 1
        # 遍歷非邊界值
        else:
            # 遍歷不在邊界0和C的alpha
            non_bound_is = np.nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
            for i in non_bound_is:
                alpha_pairs_changed += inner_L(i, oS)
                print("非邊界遍歷:第%d次迭代 樣本:%d, alpha優(yōu)化次數:%d" % (iter, i, alpha_pairs_changed))
            iter += 1
        # 遍歷一次后改為非邊界遍歷
        if entire_set:
            entire_set = False
        # 如果alpha沒有更新,計算全樣本遍歷
        elif (alpha_pairs_changed == 0):
            entire_set = True
        print("迭代次數: %d" % iter)
    return oS.b, oS.alphas


def calc_ws(alphas, data_arr, class_labels):
    """
        Function:
            計算w
        Parameters:
            data_arr - 數據矩陣
            class_labels - 數據標簽
            alphas - alphas值
        Returns:
            w - 計算得到的w
        Modify:
            2018-09-22
    """
    X = np.mat(data_arr)
    label_mat = np.mat(class_labels).transpose()
    m, n = np.shape(X)
    w = np.zeros((n, 1))
    for i in range(m):
        w += np.multiply(alphas[i] * label_mat[i], X[i, :].T)
    return w


def show_data_set(data_mat, label_mat):
    """
        Function:
            可視化數據集
        Parameters:
            data_mat - 數據矩陣
            label_mat - 數據標簽
        Returns:
            無
        Modify:
            2018-09-23
    """
    data_plus = []
    data_minus = []
    for i in range(len(data_mat)):
        if label_mat[i] > 0:
            data_plus.append(data_mat[i])
        else:
            data_minus.append(data_mat[i])
    # 轉換為numpy矩陣
    data_plus_np = np.array(data_plus)
    data_minus_np = np.array(data_minus)
    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()


if __name__ == '__main__':
    data_arr, class_labels = loadDataSet('./machinelearninginaction/Ch06/testSet.txt')
    # b, alphas = smoSimple(data_arr, class_labels, 0.6, 0.001, 40)
    # print(b)
    # print(alphas[alphas > 0])
    #
    # for i in range(100):
    #     if alphas[i] > 0.0:
    #         print(data_arr[i], class_labels[i])
    #
    # w = calc_ws(alphas, data_arr, class_labels)
    # show_classifier(data_arr, class_labels, w, b, alphas)
    #
    b, alphas = smoP(data_arr, class_labels, 0.6, 0.001, 40)
    w = calc_ws(alphas, data_arr, class_labels)
    show_classifier(data_arr ,class_labels,w, b,alphas)

1、簡化版SMO

簡化版SMO運行結果

分類可視化結果

在幾百個點組成的小規(guī)模數據集上峰锁，簡化版SMO算法的運行是沒有什么問題的萎馅，但是在更大的數據集上的運行速度就會變慢。簡化版SMO算法的第二個α的選擇是隨機的虹蒋，針對這一問題糜芳，我們可以使用啟發(fā)式選擇第二個α值，來達到優(yōu)化效果千诬。

2耍目、完整 Platt SMO 算法

Platt SMO算法是通過一個外循環(huán)來選擇第一個alpha值的，并且其選擇過程會在兩種方式之間進行交替：一種方式是在所有數據集上進行單遍掃描徐绑，另一種方式則是在非邊界alpha中實現(xiàn)單遍掃描邪驮。而所謂非邊界alpha指的就是那些不等于邊界0或C的alpha值。對整個數據集的掃描相當容易傲茄，而實現(xiàn)非邊界alpha值的掃描時毅访，首先需要建立這些alpha值的列表，然后再對這個表進行遍歷盘榨。同時喻粹，該步驟會跳過那些已知的不會改變的alpha值。

在選擇第一個alpha值后草巡，算法會通過一個內循環(huán)來選擇第二個alpha值守呜。在優(yōu)化過程中，會通過最大化步長的方式來獲得第二個alpha值山憨。在簡化版SMO算法中查乒，我們會在選擇 j 之后計算錯誤率 Ej 。但在這里郁竟，我們會建立一個全局的緩存用于保存誤差值玛迄，并從中選擇使得步長或者說Ei-Ej 最大的alpha值。

完整SMO分類結果

完整版的SMO算法計算速度棚亩，上面這個運行結果蓖议，1~2秒就可以跑完，而簡化版的SMO都是在10秒以上跑完讥蟆，完整版的SMO算法計算速度得到大大提升勒虾。

給定C的設置，圖中畫紅圈的樣本點為支持向量上的點攻询，是滿足算法的一種解从撼。如果數據集非線性可分，就會發(fā)現(xiàn)支持向量會在超平面附近聚集成團〉土悖可以看出婆翔，完整版SMO算法選出的支持向量樣點更多，更接近理想的分隔超平面掏婶。

三啃奴、非線性SVM

貼上第三、第四節(jié)全部代碼

from numpy import *
import numpy as np
import matplotlib.pyplot as plt


def loadDataSet(filename):
    """
        Function:
            讀取數據
        Parameters:
            fileName - 文件名
        Returns:
            dataMat - 數據矩陣
            labelMat - 數據標簽
        Modify:
            2018-09-21
    """
    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[-1]))
    return dataMat, labelMat



def selectJrand(i, m):
    """
        Function:
            隨機選擇alphaJ
        Parameters:
            i - alpha下標
            m - alpha參數個數
        Returns:
            j - 隨機選擇alpha另一個下標
        Modify:
            2018-09-21
    """
    j = i
    while (j == i):
        j = int(random.uniform(0, m))
    return j


def clip_alpha(aj, H, L):
    """
        Function:
            修剪alpha
        Parameters:
            aj - alpha值
            H - alpha上限
            L - alpha下限
        Returns:
            aj - alpah值
        Modify:
            2018-09-21
    """
    if aj > H:
        aj = H
    if L > aj:
        aj = L
    return aj


def show_data_set(data_mat, label_mat):
    """
        Function:
            可視化數據集
        Parameters:
            data_mat - 數據矩陣
            label_mat - 數據標簽
        Returns:
            無
        Modify:
            2018-09-23
    """
    data_plus = []
    data_minus = []
    for i in range(len(data_mat)):
        if label_mat[i] > 0:
            data_plus.append(data_mat[i])
        else:
            data_minus.append(data_mat[i])
    # 轉換為numpy矩陣
    data_plus_np = np.array(data_plus)
    data_minus_np = np.array(data_minus)
    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()


def calc_Ek_kernel(oS, k):
    """
        Function:
            計算誤差
        Parameters:
            oS - 數據結構
            k - 下標為k的數據
        Returns:
             Ek - 下標為k的數據誤差
        Modify:
            2018-09-22
    """
    fXk = float(np.multiply(oS.alphas, oS.label_mat).T * (oS.K[:, k]) + oS.b)
    Ek = fXk - float(oS.label_mat[k])
    return Ek


def select_J_kernel(i, oS, Ei):
    """
        Function:
            內循環(huán)啟發(fā)方式
        Parameters:
            i - 下標為i的數據的索引值
            oS - 數據結構
            Ei - 下標為i的數據誤差
        Returns:
            j, maxK - 下標為j或maxK的數據的索引值
            Ej - 下標為j的數據誤差
        Modify:
            2018-09-22
    """
    maxK = -1
    max_deltaE = 0
    Ej = 0
    oS.eCache[i] = [1, Ei]
    # 返回Ei不為0的索引值數組
    valid_eCache_list = nonzero(oS.eCache[:, 0].A)[0]
    # 有不為0的Ei
    if (len(valid_eCache_list)) > 1:
        for k in valid_eCache_list:
            if k == i: continue
            Ek = calc_Ek_kernel(oS, k)
            deltaE = abs(Ei - Ek)
            if (deltaE > max_deltaE):
                maxK = k
                max_deltaE = deltaE
                Ej = Ek
        return maxK, Ej
    # 沒有不為0的Ei
    else:
        j = selectJrand(i, oS.m)
        Ej = calc_Ek_kernel(oS, j)
    return j, Ej


def update_Ek_kernel(oS, k):
    """
        Function:
            計算Ek,并更新誤差緩存
        Parameters:
            oS - 數據結構
            k - 下標為k的數據的索引值
        Returns:
            無
        Modify:
            2018-09-22
    """
    Ek = calc_Ek_kernel(oS, k)
    # 更新Ei緩存
    oS.eCache[k] = [1, Ek]


def inner_l_kernel(i, oS):
    """
        Function:
            優(yōu)化的SMO算法
        Parameters:
            i - 下標為i的數據的索引值
            oS - 數據結構
        Returns:
            1 - 有任意一對alpha值發(fā)生變化
            0 - 沒有任意一對alpha值發(fā)生變化或變化太小
        Modify:
            2018-09-22
    """
    # 步驟1：計算誤差Ei
    Ei = calc_Ek_kernel(oS, i)
    if ((oS.label_mat[i] * Ei < -oS.tol) and (oS.alphas[i] < oS.C)) or ((oS.label_mat[i] * Ei > oS.tol) and (oS.alphas[i] > 0)):
        # 使用內循環(huán)啟發(fā)方式2選擇alpha_j,并計算Ej
        j, Ej = select_J_kernel(i, oS, Ei)
        # 保存更新前的aplpha值雄妥，使用深拷貝
        alpha_i_old = oS.alphas[i].copy()
        alpha_j_old = oS.alphas[j].copy()
        # 步驟2：計算上下界L和H
        if (oS.label_mat[i] != oS.label_mat[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 0
            # 步驟3：計算eta
        eta = 2.0 * oS.K[i, j] - oS.K[i, i] - oS.K[j, j]
        if eta >= 0:
            print("eta>=0")
            return 0
        # 步驟4：更新alpha_j
        oS.alphas[j] -= oS.label_mat[j] * (Ei - Ej) / eta
        # 步驟5：修剪alpha_j
        oS.alphas[j] = clip_alpha(oS.alphas[j], H, L)
        # 更新Ei至誤差緩存
        update_Ek_kernel(oS, j)
        if (abs(oS.alphas[j] - alpha_j_old) < 0.00001):
            print('alpha_j變化太小')
            return 0
        # 步驟6：更新alpha_i
        oS.alphas[i] += oS.label_mat[j] * oS.label_mat[i] * (alpha_j_old - oS.alphas[j])
        # 步驟7：更新b1和b2
        b1 = oS.b - Ei - oS.label_mat[i] * (oS.alphas[i] - alpha_i_old) * oS.K[i, i] - \
             oS.label_mat[j] * (oS.alphas[j] - alpha_j_old) * oS.K[i, j]
        b2 = oS.b - Ej - oS.label_mat[i] * (oS.alphas[i] - alpha_i_old) * oS.K[i, j] - \
             oS.label_mat[j] * (oS.alphas[j] - alpha_j_old) * oS.K[j, j]
        # 步驟8：根據b1和b2更新b
        if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]):
            oS.b = b1
        elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]):
            oS.b = b2
        else:
            oS.b = (b1 + b2) / 2.0
        return 1
    else:
        return 0


def smoP_kernel(data_arr, class_labels, C, toler, max_iter, k_tup=('lin', 0)):
    """
    Function:
        完整的線性SMO算法
    Parameters：
        dataMatIn - 數據矩陣
        classLabels - 數據標簽
        C - 松弛變量
        toler - 容錯率
        maxIter - 最大迭代次數
        kTup - 包含核函數信息的元組
    Returns:
        oS.b - SMO算法計算的b
        oS.alphas - SMO算法計算的alphas
    Modify:
        2018-09-23
    """
    oS = OptStructKernel(np.mat(data_arr), np.mat(class_labels).transpose(), C, toler, k_tup)
    iter = 0
    entire_set = True
    alpha_pairs_changed = 0
    # 遍歷整個數據集alpha也都沒有更新或者超過最大迭代次數則退出循環(huán)
    while (iter < max_iter) and ((alpha_pairs_changed > 0) or (entire_set)):
        alpha_pairs_changed = 0
        # 遍歷整個數據集
        if entire_set:
            for i in range(oS.m):
                # 使用優(yōu)化的SMO算法
                alpha_pairs_changed += inner_l_kernel(i, oS)
                print("全樣本遍歷:第%d次迭代 樣本:%d, alpha優(yōu)化次數:%d" % (iter, i, alpha_pairs_changed))
            iter += 1
        # 遍歷非邊界值
        else:
            # 遍歷不在邊界0和C的alpha
            non_bound_is = np.nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
            for i in non_bound_is:
                alpha_pairs_changed += inner_l_kernel(i, oS)
                print("非邊界遍歷:第%d次迭代 樣本:%d, alpha優(yōu)化次數:%d" % (iter, i, alpha_pairs_changed))
            iter += 1
        # 遍歷一次后改為非邊界遍歷
        if entire_set:
            entire_set = False
        # 如果alpha沒有更新,計算全樣本遍歷
        elif (alpha_pairs_changed == 0):
            entire_set = True
        print("迭代次數: %d" % iter)
    return oS.b, oS.alphas


def kernel_trans(X, A, k_tup):
    m, n = np.shape(X)
    K = np.mat(np.zeros((m, 1)))
    # 線性核函數,只進行內積
    if k_tup[0] == 'lin': K = X * A.T
    elif k_tup[0] == 'rbf':
        for j in range(m):
            delta_row = X[j, :] - A
            K[j] = delta_row * delta_row.T
        K = np.exp(K / (-1 * k_tup[1] ** 2))
    else: raise NameError('核函數無法識別')
    return K



class OptStructKernel(object):
    def __init__(self, data_mat_in, class_labels, C, toler, k_tup):
        self.X = data_mat_in
        self.label_mat = class_labels
        self.C = C
        self.tol = toler
        # 數據矩陣行數
        self.m = np.shape(data_mat_in)[0]
        # 根據矩陣行數初始化alpha參數為0
        self.alphas = np.mat(np.zeros((self.m, 1)))
        self.b = 0
        # 根據矩陣行數初始化誤差緩存院仿，第一列為是否有效的標志位贿堰，第二列為實際的誤差E的值
        self.eCache = np.mat(np.zeros((self.m, 2)))
        # 初始化核K
        self.K = np.mat(np.zeros((self.m, self.m)))
        # 計算所有數據的核K
        for i in range(self.m):
            self.K[:, i] = kernel_trans(self.X, self.X[i, :], k_tup)


def calc_ws(alphas, data_arr, class_labels):
    """
        Function:
            計算w
        Parameters:
            data_arr - 數據矩陣
            class_labels - 數據標簽
            alphas - alphas值
        Returns:
            w - 計算得到的w
        Modify:
            2018-09-22
    """
    X = np.mat(data_arr)
    label_mat = np.mat(class_labels).transpose()
    m, n = np.shape(X)
    w = np.zeros((n, 1))
    for i in range(m):
        w += np.multiply(alphas[i] * label_mat[i], X[i, :].T)
    return w


def show_classifier(dataMat, class_labels,w, b, alphas):
    """
        Function:
            可視化分類結果
        Parameters:
            dataMat - 數據矩陣
            w - 直線法向量
            b - 直線截距
            alphas - alphas值
        Returns:
            無
        Modify:
            2018-09-22
    """
    # 繪制樣本點
    data_plus = []
    data_minus = []
    for i in range(len(dataMat)):
        if class_labels[i] > 0:
            data_plus.append(dataMat[i])
        else:
            data_minus.append(dataMat[i])
    data_plus_np = np.array(data_plus)
    data_minus_np = np.array(data_minus)
    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 alpha > 0:
            x, y = dataMat[i]
            plt.scatter([x], [y], s=150, c='none', alpha=0.7, linewidth=1.5, edgecolor='red')
    plt.show()

def test_rbf(k1=1.3):
    """
        Function:
            利用核函數進行分類的徑向基測試函數
        Parameters:
            k1 - 使用高斯核函數的時候表示到達率
        Returns:
            無
        Modify:
            2018-09-26
    """
    data_arr, class_labels = loadDataSet('machinelearninginaction/Ch06/testSetRBF.txt')
    b, alphas = smoP_kernel(data_arr, class_labels, 200, 0.0001, 100, ('lin', k1))
    data_mat = mat(data_arr)
    labels_mat = mat(class_labels).transpose()
    sv_ind = nonzero(alphas.A > 0)[0]
    s_v_s = data_mat[sv_ind]
    label_sv = labels_mat[sv_ind]
    print('支持向量個數:%d' % shape(s_v_s)[0])
    m, n = shape(data_mat)
    error_count = 0
    for i in range(m):
        kernel_eval = kernel_trans(s_v_s, data_mat[i, :], ('lin', k1))
        predict = kernel_eval.T * multiply(label_sv, alphas[sv_ind]) + b
        if sign(predict) != sign(class_labels[i]): error_count += 1
    print('訓練集錯誤率: %.2f%%' % ((float(error_count) / m) * 100))

    data_arr, class_labels = loadDataSet('machinelearninginaction/Ch06/testSetRBF2.txt')
    data_mat = mat(data_arr)
    labels_mat = mat(class_labels).transpose()
    m, n = shape(data_mat)
    error_count = 0
    for i in range(m):
        kernel_eval = kernel_trans(s_v_s, data_mat[i, :], ('lin', k1))
        predict = kernel_eval.T * multiply(label_sv, alphas[sv_ind]) + b
        if sign(predict) != sign(class_labels[i]): error_count += 1
    print('測試集錯誤率: %.2f%%' % ((float(error_count) / m) * 100))

    # data_arr, class_labels = loadDataSet('machinelearninginaction/Ch06/testSet.txt')
    # b, alphas = smoP_kernel(data_arr, class_labels, 200, 0.0001, 100, ('lin', k1))
    # data_mat = mat(data_arr)
    # labels_mat = mat(class_labels).transpose()
    # sv_ind = nonzero(alphas.A > 0)[0]
    # s_v_s = data_mat[sv_ind]
    # label_sv = labels_mat[sv_ind]
    # print('支持向量個數:%d' % shape(s_v_s)[0])
    # m, n = shape(data_mat)
    # error_count = 0
    # for i in range(m):
    #     kernel_eval = kernel_trans(s_v_s, data_mat[i, :], ('lin', k1))
    #     predict = kernel_eval.T * multiply(label_sv,alphas[sv_ind]) + b
    #     if sign(predict) != sign(class_labels[i]): error_count += 1
    # print('訓練集錯誤率: %.2f%%' % ((float(error_count)/m)*100))
    # w = calc_ws(alphas, data_arr, class_labels)
    # show_classifier(data_arr, class_labels, w, b, alphas)


def img_vector(file_name):
    """
        Function:
            將32x32的二進制圖像轉換為1x1024向量
        Parameters:
            file_name - 文件名
        Returns:
            return_vect - 二進制圖像的1x1024向量
        Modify:
            2018-09-26
    """
    return_vect = zeros((1, 1024))
    f = open(file_name, 'r')
    for i in range(32):
        line_str = f.readline()
        for j in range(32):
            return_vect[0,32 * i + j] = int(line_str[j])
    return return_vect


def load_image(dir_name):
    """
        Function:
            加載圖片
        Parameters:
            dir_name - 圖片文件夾路徑
        Returns:
            training_mat - 數據矩陣
            hw_labels - 數據標簽
        Modify:
            2018-09-26
    """
    from os import listdir
    hw_labels = []
    training_file_list = listdir(dir_name)
    m = len(training_file_list)
    training_mat = zeros((m, 1024))
    for i in range(m):
        file_name_str = training_file_list[i]
        file_str = file_name_str.split('.')[0]
        class_num_str = int(file_str.split('_')[0])
        if class_num_str == 9: hw_labels.append(-1)
        else: hw_labels.append(1)
        training_mat[i,:] = img_vector('%s/%s' % (dir_name, file_name_str))
    return training_mat, hw_labels


def test_digits(k_tup=('rbf', 10)):
    """
        Function:
            基于SVM的手寫數字識別測試函數
        Parameters:
            k_tup - 包含核函數信息的元組
        Returns:
            無
        Modify:
            2018-09-26
    """
    data_arr, label_arr = load_image('./machinelearninginaction/Ch02/digits/trainingDigits')
    b, alphas = smoP_kernel(data_arr, label_arr, 200, 0.0001, 10, k_tup)
    data_mat = mat(data_arr)
    label_mat = mat(label_arr).transpose()
    sv_ind = nonzero(alphas.A > 0)[0]
    s_v_s = data_mat[sv_ind]
    label__sv = label_mat[sv_ind]
    print('支持向量個數:%d' % np.shape(s_v_s)[0])
    m, n = shape(data_mat)
    error_count = 0
    for i in range(m):
        kernel_eval = kernel_trans(s_v_s, data_mat[i, :], k_tup)
        predict = kernel_eval.T * multiply(label__sv, alphas[sv_ind]) + b
        if sign(predict) != sign(label_arr[i]): error_count += 1
    print('訓練集錯誤率: %.2f%%' % (float(error_count) / m))

    data_arr, label_arr = load_image('./machinelearninginaction/Ch02/digits/testDigits')
    data_mat = mat(data_arr)
    label_mat = mat(label_arr).transpose()
    m, n = shape(data_mat)
    error_count = 0
    for i in range(m):
        kernel_eval = kernel_trans(s_v_s, data_mat[i, :], k_tup)
        predict = kernel_eval.T * multiply(label__sv, alphas[sv_ind]) + b
        if sign(predict) != sign(label_arr[i]): error_count += 1
    print('測試集錯誤率: %.2f%%' % (float(error_count) / m))

if __name__ == '__main__':
    # data_arr, class_labels = loadDataSet('./machinelearninginaction/Ch06/testSetRBF.txt')
    # show_data_set(data_arr, class_labels)

    # test_rbf()

    test_digits(('rbf', 20))

非線性數據集可視化

運用徑向基核函數分類結果，k1=0.1

共有 100 個數據點，88 個為支持向量鹿驼，必須使用這些支持向量才能對數據進行正確分類鸣剪。

運用徑向基核函數分類結果零蓉，k1=1.3

運用徑向基核函數分類結果卿嘲，k1=2

該數據集在這個設置的某處存在著最優(yōu)值。如果降低 σ 淀弹，那么訓練錯誤率就會降低丹壕，但是測試錯誤率卻會上升，σ 過小會出現(xiàn)欠擬合現(xiàn)象薇溃；如果增大σ 菌赖，那么訓練錯誤率就會上升，但是測試錯誤率卻會降低沐序，σ 過大會出現(xiàn)出現(xiàn)過擬合現(xiàn)象琉用。

運用線性核函數分類結果，k1=1.3

使用線性核函數對非線性數據分類效果較差策幼。

對于第2節(jié)線性數據集辕羽，看看運用線性核函數分類結果也較為理想，見下圖垄惧。

運用線性核函數對線性數據分類結果，k1=1.3

四绰寞、基于SVM的手寫數字識別

之前我們使用了kNN算法實現(xiàn)了手寫數字識別到逊，但是它的缺點是存儲空間大，因為要保留所有的訓練樣本滤钱，而對于支持向量機而言觉壶，其需要保留的樣本少了很多（即只保留支持向量），但是能獲得可比的效果件缸，效率也更高铜靶。

現(xiàn)在看看基于SVM的手寫數字識別分類的效果，在 kNN.py 中代碼直接應用類別標簽他炊，而同支持向量機一起使用時争剿，類別標簽為 -1 或者 +1 已艰。因此，一旦碰到數字 9 蚕苇，則輸出類別標簽 -1 哩掺，否則輸出 +1 。

基于SVM的手寫數字識別分類的運行結果

可以看到訓練集錯誤率為0涩笤，測試集錯誤率也僅為0.01%嚼吞。

五、示例：應用scikit-learn構建SVM分類器

sklearn.svm模塊提供了很多模型供我們使用蹬碧，其中是基于libsvm實現(xiàn)的舱禽。下面我們使用svm.SVC對手寫數字進行識別分類。

可以看到恩沽，訓練和測試的時間總共加起來才4.9s誊稚。而且，測試集的錯誤率僅為1.37%飒筑∑酰可以說識別效果是挺不錯的。

六协屡、小結

支持向量機是一種分類器俏脊。之所以稱為“機”是因為它會產生一個二值決策結果，即它是一種決策“機”肤晓。它可用于線性/非線性分類爷贫，也可以用于回歸，支持向量機的泛化錯誤率較低补憾，也就是說它具有良好的學習能力漫萄，且學到的結果具有很好的推廣性。

核方法或者說核技巧會將數據（有時是非線性數據）從一個低維空間映射到一個高維空間盈匾，可以將一個在低維空間中的非線性問題轉換成高維空間下的線性問題來求解腾务。核方法不止在 SVM中適用，還可以用于其他算法中削饵。而其中的徑向基函數是一個常用的度量兩個向量距離的核函數岩瘦。

支持向量機是一個二類分類器。當用其解決多類問題時窿撬，則需要額外的方法對其進行擴展启昧。SVM 的效果也對優(yōu)化參數和所用核函數中的參數敏感。