【数据科学系统学习】机器学习算法 # 西瓜书学习记录 [5] 支持向量机实践

摘要：本篇内容为机器学习实战第章支持向量机部分程序清单。支持向量机优点泛化错误率低，计算开销不大，结果易解释。注以上给出的仅是简化版算法的实现，关于完整的算法加速优化并应用核函数，请参照机器学习实战第页。

本篇内容为《机器学习实战》第 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 页。

不足之处，欢迎指正。