机器学习：支持向量机（Support Vector Machine, SVM）

SVM是经典的监督型分类算法，广泛应用于机器学习、数据挖掘等领域。本质上SVM与线性回归方法类似，都是求一组权重系数。

感知机（perceptron）分类算法由于样本顺序和不同初值会导致解的多样性，解不唯一，亦不为最优，而SVM本质上是解最优目标方程的过程，解唯一，而且对于线性不可分的样本集会通过kernelling将样本映射到高维空间使其线性可分，再训练出分类边界，即超平面。因此SVM更具先进性。

将SVM的整个数学过程推到出来并不是一件简单的事，日后再完善orz。

一. 问题描述

如下图所示，SVM所解决的问题就是找出两组样本点间的一个最优分类边界，距离分类边界最近的那些点就是支持向量，在得到最终的分类边界时，支持向量到分类边界的距离也达到了最大。

设超平面方程为：$w^T \dot{} x_i+b=0$；
标签：$y_i=1~or~-1$，
其中$i$为样本的index

归一化：$y_i(w^T \dot{} x_i+b)\geqslant 1$

由点到面的距离公式可知此时支持向量到超平面的距离为：$\frac{1}{\| \|w\| \|}$

那么最大化$\frac{1}{\| \|w\| \|}$
可以转换为最小化$\| \|w\| \|^{2}$，此时的最优化问题如下：

$$\min \frac{1}{2}||w||^2 ~~~s.t.~~y_i(w^T \dot{} x_i+b)\geqslant1$$

二. 算法过程

首先讲上述优化问题通过拉格朗日乘子法和KKT条件转换为另一个更利于快速解决的优化问题，该优化方程为：

$$\max_{\alpha} \sum_{i} \alpha_{i}-\frac{1}{2}\sum_{i,j} \alpha_{i}\alpha_{j}y_iy_jx_i^Tx_j$$

$$s.t. ~~~ a_i\geqslant0,~\sum_{i}a_iy_i=0$$

这里，$a_i$是拉格朗日算子，接下来再用SMO算法解此最优化问题。

三. 相关代码

机器学习实战里面有实现SVM的python代码，如下：

<span style="font-size:24px;">from numpy import * 
import random 

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 


def selectJrand(i, m):  #生成一个随机数
    j=i 
    while(j==i): 
        j=int(random.uniform(0, m))  #生成一个[0, m]的随机数，int转换为整数。注意，需要import random
    return j 

def clipAlpha(aj, H, L):  #阈值函数
    if aj>H: 
        aj=H 
    if aj<L: 
        aj=L 
    return aj 

def smoSimple(dataMatIn, classLabels, C, toler, maxIter): 
    dataMatrix = mat(dataMatIn); labelMat = mat(classLabels).transpose()
    b = 0; m,n = shape(dataMatrix)
    alphas = mat(zeros((m,1)))
    iter = 0
    while(iter<maxIter):  #迭代次数
        alphaPairsChanged=0 
        for i in range(m):  #在数据集上遍历每一个alpha
            #print alphas 
            #print labelMat
            fXi = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + b
            #fXi=float(np.multiply(alphas, labelMat).T*dataMatrix*dataMatrix[i, :].T)+b  #.T也是转置
            Ei=fXi-float(labelMat[i]) 
            if((labelMat[i]*Ei<-toler) and (alphas[i]<C)) or ((labelMat[i]*Ei>toler) and (alphas[i]>0)): 
                j=selectJrand(i, m)  #从m中选择一个随机数，第2个alpha j
                fXj=float(multiply(alphas, labelMat).T*dataMatrix*dataMatrix[j, :].T)+b 
                Ej=fXj-float(labelMat[j]) 

                alphaIold=alphas[i].copy()  #复制下来，便于比较
                alphaJold=alphas[j].copy() 

                if(labelMat[i]!=labelMat[j]):  #开始计算L和H
                    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是alphas[j]的最优修改量，如果eta为零，退出for当前循环
                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 
                alphas[j]-=labelMat[j]*(Ei-Ej)/eta  #调整alphas[j] 
                alphas[j]=clipAlpha(alphas[j], H, L)  
                if(abs(alphas[j]-alphaJold)<0.00001):  #如果alphas[j]没有调整
                    print 'j not moving enough' 
                    continue 
                alphas[i]+=labelMat[j]*labelMat[i]*(alphaJold-alphas[j])  #调整alphas[i]
                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(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 'iter: %d i: %d, pairs changed %d' %(iter, i, alphaPairsChanged) 
        if(alphaPairsChanged==0): 
            iter+=1 
        else: 
            iter=0 
        print 'iteration number: %d' %iter
    return b, alphas 


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

    print b, alphas </span>

四. 参考

支持向量机通俗导论（理解SVM的三层境界）

手把手教你实现SVM算法（一）