机器学习入门笔记：（4.3）SMO算法

前言

之前的博客中，已经介绍了SVM的原理:

机器学习入门学习笔记：（4.1）SVM算法

机器学习入门学习笔记：（4.2）核函数和软间隔

最后我们得到的优化问题如下：

maxα∑i=1mαi−12∑i=1m∑j=1mαiαjyiyjxixj,s.t.∑i=1mαiyi=00≤αi≤C,i=1,2,...,m

这个问题的解就是一系列的α，这些α会使得上面的式子有最大值。

这个式子是引入了软间隔后的支持向量机的问题，再进一步，用上核函数，就可以表示为：

maxα∑i=1mαi−12∑i=1m∑j=1mαiαjyiyjK(xi,xj),s.t.∑i=1mαiyi=00≤αi≤C,i=1,2,...,m

一般来说，我们都不是很喜欢求最大值的问题，而是求最小值，所以将上面的问题换成求最小值的形式：

minα12∑i=1m∑j=1mαiαjyiyjK(xi,xj)−∑i=1mαi,s.t.∑i=1mαiyi=00≤αi≤C,i=1,2,...,m

我们之后的讨论都会围绕着这个问题进行。

SMO算法

为了解决二次规划问题，人们提出许多高效的算法。其中比较典型的一个就是SMO(Sequential Minimal Optimization)算法。SMO算法由Microsoft Research的John C. Platt在1998年提出，并成为最快的二次规划优化算法，特别针对线性SVM和数据稀疏时性能更优。关于SMO最好的资料就是他本人写的《Sequential Minimal Optimization A Fast Algorithm for Training Support Vector Machines》了。

（摘自支持向量机（五）SMO算法）

SMO 概念

SMO的基本思路就是：先固定αi之外的所有参数，然后求αi的极值。但是问题中存在约束条件：∑mi=1αiyi=0。如果固定了αi之外的其他变量，则αi也会被固定，可以由其他的变量导出。于是，一次只留一个参数，固定其余参数的方法在这里是不适用的，但是这个思想却给了我们不错的启发。那么，SMO可以每次选择两个变量αi和αj，并固定其他参数。这样，在参数初始化之后，SMO不断迭代重复下面的步骤，直至收敛：

• 选取一对新的αi和αj;
• 固定αi和αj之外的参数，求解前面的优化问题，获取更新后的αi和αj。

假设选取的两个变量为α1和α2，那么由于其余参数均被固定，目标函数最后也只与α1和α2有关。

由约束条件∑mi=1αiyi=0有：α1y1+α2y2=C，其中C为常数。

为简化表示，我们用Kij表示K(xi,xj)。

对原始问题进行化简：

minα12∑i=1m∑j=1mαiαjyiyjK(xi,xj)−∑i=1mαi

=minα∑i=1m[αiα1yiy1K(xi,x1)+αiα2yiy2K(xi,x2)+∑j=3mαiαjyiyjK(xi,xj)]−α1−α2−∑i=3mαi

=minα12[α21y21K(x1,x1)+α1α2y1y2K(x1,x2)+∑j=3mα1αjy1yjK(x1,xj)]+12[α2α1y2y1K(x2,x1)+α22y22K(x2,x2)+∑j=3mα2αjy2yjK(x2,xj)]+12∑i=3m[αiα1yiy1K(xi,x1)+αiα2yiy2K(xi,x2)+∑j=3mαiαjyiyjK(xi,xj)]−α1−α2−∑i=3mαi

=minα12[α21y21K11+α22y22K22+2α1α2y1y2k12+∑i=3m∑j=3mαiαjyiyjK(xi,xj)+2α1y1v1+2α2y2v2]−α1−α2−∑i=3mαi

（其中，因为核函数满足K12=K21，所以直接合在一起了；为了便于表示，令vi=∑mj=3yjαjkij）

由于除了α1和α2之外的αi,i=3,4,...,m都被固定了，所以都是常数。

=minα12[α21y21K11+α22y22K22+2α1α2y1y2k12+2α1y1v1+2α2y2v2]−α1−α2+C

其中C是任意常数。

解的范围

现在的结果是有关α1和α2的表达式，为了进一步化简我们还可以将α1用α2表示。但是在此之前，还有一个问题，那就是α1和α2的取值范围是多少？尽管有了约束条件∑mi=1αiyi=0，但是在y1和y2取不同值时，也会有不同的约束关系产生。

（摘自John C. Platt的论文）

前面推出的α1和α2的约束关系为：α1y1+α2y2=C。

我们还有些已知的条件：

0≤α≤C

y1和y2都是输出的标签，为±1，显然：y21=y22=1。

如上图所示，分两种情况讨论：

1. 如果y1≠y2，则y1和y2一定异号，那么约束关系变为：α1−α2=k，k是任意常数。k具体是多少我们并不关心，但是我们知道α1和α2的取值都落在途中的直线上。k无非就是一个截距，随着k的变化，这根直线在方框内会上下移动，交点也变，但是一定要在方框范围内，所以边界一定会落在方框与直线的交点上。假设几种可能的情况就不难推算出α 的范围了假设只考虑α2的范围，设L为α2可能的最小取值，H为α2可能的最大取值。

   L=max(0,−k),H=min(C,C−k)L=max(0,α2−α1),H=min(C,C+α2−α1)

2. 如果y1=y2，则y1和y2一定同号，那么约束关系变为：α1+α2=k，k是任意常数。与前面同样分析。只考虑α2的范围，设L为α2可能的最小取值，H为α2可能的最大取值。根据图中所示，可以得到如下关系：

   L=max(0,k−C),H=min(C,k)L=max(0,α1+α2−C),H=min(C,α1+α2)

通过上面的讨论，我们得到了α2的可能取值范围：L≤α2≤H

1. 当y1≠y2时，L=max(0,α2−α1),H=min(C,C+α2−α1)；
2. 当y1=y2时，L=max(0,α1+α2−C),H=min(C,α1+α2)。

同理，α1的范围与α2是一样的。这个范围先保留，后面再用。

求解优化问题

转化为一元函数求极值点

接下来，将α1用α2表示。将约束关系：α1y1+α2y2=C左右同时乘上y1得到：

α1=(ζ−α2y2)y1

，这里的ζ是常数，为了不跟原式子中的C混淆，换成ζ表示这个常数。

好了，我们已经将α1用α2表示出来了，可以代回到前面的问题中了：

minα12[α21y21K11+α22y22K22+2α1α2y1y2k12+2α1y1v1+2α2y2v2]−α1−α2+C其中vi=∑j=3myjαjkij

我们有y21=y22=1，以及α1=(ζ−α2y2)y1，代入化简得到：

minα12[(ζ−α2y2)2k11+k22α22+2(ζ−α2y2)α2y2k12+2(ζ−α2y2)v1+2α2y2v2]−(ζ−α2y2)y1−α2+C

常数可以去掉，不影响结果：

minα12[(ζ−α2y2)2k11+k22α22+2(ζ−α2y2)α2y2k12+2(ζ−α2y2)v1+2α2y2v2]−(ζ−α2y2)y1−α2

上面的问题已经化成了单变量的优化问题了，使用常规套路，求偏导取0，即可解出α2的值。

对目标函数求偏导数：

∂Φ∂α2=α2(K11+K22−2K12)−K11ζy2+K12ζy2−y2v1+y2v2+y1y2−1

令这个偏导数为0可以求出新的α2，利用

α1=(ζ−α2y2)y1

这个关系，又可以求出新的α1，这两新求出的值即为我们使用SMO算法优化之后的结果。为了与原始的α1和α2值区分，我们将这两个新的值标记为α∗1和α∗2。

修改后的几个条件如下：

∂Φ∂α∗2=α∗2(K11+K22−2K12)−K11ζy2+K12ζy2−y2v1+y2v2+y1y2−1=0

ζ=α1y1+α2y2=α∗1y1+α∗2y2

其中vi=∑j=3myjαjkij

为了后面表示简便，还要给出几个关系：

如果把vi直接代入，结果太复杂了。还需要对vi做一些变换：

前面的博客中，我们已经推导过SVM的数学模型的最终结果：f(xi)=ωTx+b=∑mj=1αjyjK(xi,xj)+b

展开f(x)看看：

f(x1)=α1y1K11+α2y2K12+∑j=3mαjyjK1j+b=α1y1K11+α2y2K12+v1+b

f(x2)=α1y1K12+α2y2K22+∑j=3mαjyjK2j+b=α1y1K12+α2y2K22+v2+b

所以：

v1=f(x1)−α1y1K11−α2y2K12−b

v2=f(x2)−α1y1K12−α2y2K22−b

好了，条件基本都得到了，再列一下我们要用的条件：

∂Φ∂α∗2=α∗2(K11+K22−2K12)−K11ζy2+K12ζy2−y2v1+y2v2+y1y2−1=0

ζ=α1y1+α2y2(=α∗1y1+α∗2y2)

v1=f(x1)−α1y1K11−α2y2K12−b

v2=f(x2)−α1y1K12−α2y2K22−b

联立，化简得：(注：ζ=α1y1+α2y2，用旧的α1和α2，因为我们最后要表示出新的α∗2）

α∗2(K11+K22−2K12)=(K11+K22−2K12)α2+y2[y2−y1+f(x1)−f(x2)]

设预测值与真实值之差为Ei：Ei=f(xi)−yi

继续化简：

α∗2=α2+y2E1−E2(K11+K22−2K12)

再记η=(K11+K22−2K12)：

α∗2=α2+y2E1−E2η,η=(K11+K22−2K12)

注意，α∗2是经过优化后求出的解，α2是之前的值。

当然还有α2的范围约束，前面我们已经推导了：

了α2的可能取值范围：L≤α2≤H

1. 当y1≠y2时，L=max(0,α2−α1),H=min(C,C+α2−α1)；
2. 当y1=y2时，L=max(0,α1+α2−C),H=min(C,C+α1+α2)。

所以求出了优化后的α∗2后，还需要经过一个范围的约束：

αnew2=⎧⎩⎨⎪⎪H,α∗2>Hα∗2,L≤α∗2≤HL,α∗2<L

求解α1

知道了α2求α1就很容易了：

由约束条件：α1y1+α2y2=αnew1y1+αnew2y2=ζ得到：

αnew1=α1+y1y2(α2−αnew2)

取临界情况

前面推导的结果：

α∗2=α2+y2E1−E2η,η=(K11+K22−2K12)

大部分情况下，都有η=(K11+K22−2K12)>0，但是在不满足这个条件时，α∗2需要取临界值。

1. η<0,当核函数K不满足Mercer定理时，矩阵K非正定；
2. η=0,样本x1与x2输入特征相同；

也可以换个方式来理解：

原问题：

minα12[(ζ−α2y2)2k11+k22α22+2(ζ−α2y2)α2y2k12+2(ζ−α2y2)v1+2α2y2v2]−(ζ−α2y2)y1−α2+C

其一阶偏导数为：

∂Φ∂α2=α2(K11+K22−2K12)−K11ζy2+K12ζy2−y2v1+y2v2+y1y2−1

二阶偏导数为：

∂2Φ∂α22=η=(K11+K22−2K12)

这个η就是原问题的二阶偏导数，根据函数的性质来看：

1. 当η<0时，目标函数为凸函数，没有极小值，最小值会在边界取得；
2. 当η=0时，目标函数为单调函数，很明显，最小值或者最大值都会在边界上取得。

所以，当η≤0时，把α∗2=L和α∗2=H分别代入α1y1+α2y2=αnew1y1+αnew2y2=ζ解出α∗1=L1和α∗1=H1，其中令s=y1y2：

L1=α1+s(α2−L)

H1=α1+s(α2−H)

代回到目标函数中可以求出对应的两个可能值ΨL和ΨH，最后取两者中更小的那个就是最小值了。

代入之前先看看目标函数：

Ψ=12[α21K11+α22K22+2α1α2y1y2k12+2α1y1v1+2α2y2v2]−α1−α2

因为v1和v2的存在，展开后还是有些不太好看的。vi=∑mj=3yjαjkij这东西不好化简，所以使用vi的另一种表示形式：

vi=f(xi)−α1y1Ki1−α2y2Ki2−b

代入Ψ中：(注：y21=y22=1)

Ψ=12α21K11+12α22K22+α1α2y1y2K12+α1y1(f(x1)−α1y1K11−α2y2K12−b)+α2y2(f(x2)−α1y1K12−α2y2K22−b)−α1y21−α2y22

Ψ=12α21K11+12α22K22+α1α2y1y2K12+α1(y1f(x1)−α1K11−α2y1y2K12−y1b−y21)+α2(y2f(x2)−α1y1y2K12−α2K22−by2−y22)

下面令Ψ=12α21K11+12α22K22+α1α2y1y2k12+α1f1+α2f2，则：

f1=y1f(x1)−α1K11−α2y1y2K12−y1b−y21=y1(E1−b)−α1K11−α2y1y2K12

f2=y2f(x2)−α1y1y2K12−α2K22−by2−y22=y2(E2−b)−α1y1y2K12−α2K22

使用上面的这些式子，将L1=α1+s(α2−L)、H1=α1+s(α2−H)以及s=y1y2代入可以得到如下结果：

这里是论文中的结果，我就偷懒不写步骤了。前面推导的式子联立就可以得到上面的6个式子。

计算ω和b

首先看ω，前面我们可以解出α，根据公式：ω=∑mi=1αiyixi就可以求出来ω。

因为除了α1和α2之外的αi都被固定了，所以优化前后都会有如下关系：

ω=α1y1x1+α2y2x2+∑i=3mαiyixi

ωnew=αnew1y1x1+αnew2y2x2+∑i=3mαiyixi

两式做差即可求出新的ω：

ωnew=ω+y1x1(αnew1−α1)+y2x2(αnew2−α2)

然后是b，我们没有直接的公式来计算，只能通过KKT条件间接求出来。

（摘自Platt的论文）

这是原优化问题的KKT条件：

1. 当αi=0时，分类是正确的；
2. 当0≤αi≤C时，这时的样本点是支持向量，处在边界上；
3. 当αi=C时，位于边界之间。

参考上面的KKT条件进行分类讨论：

1. 如果0<α1<C，则(x1,y1)为支持向量，满足yi(∑mi=1αiyiKi1+b1)=1：

   αnew1y1K11+αnew2y2K21+∑i=3mαiyiKi1+bnew1=y1

   因为：y1−∑mi=3αiyiKi1=y1−f(x1)+α1y1K11+α2y2K21+b

   所以：

   bnew1=y1−f(x1)+α1y1K11+α2y2K21+b−αnew1y1K11−αnew2y2K21
   bnew1=−E1−y1K11(αnew1−α1)−y2K21(αnew2−α2)+b

   其中α1,α2,b为旧的值，αnew1,αnew2,bnew为优化后的值。

2. 如果0<α2<C，则(x2,y2)为支持向量，同理可以得到：

   bnew2=−E2−y1K12(αnew1−α1)−y2K22(αnew2−α2)+b

3. 如果同时有0<α1<C，0<α2<C，那么bnew1=bnew2。

4. 如果均不满足0<αi<C，这时选择他们的中点：bnew=bnew1+bnew22

启发式选择变量

由Osuna定理知道，只要选取的αi和αj中有一个不满足KKT条件的，那么最初的那个目标函数就会在迭代之后减小。直观来看，如果KKT条件违背的越厉害，那么变量更新后可能导致目标函数减小的更厉害。于是，SMO可以首先取违背KKT条件程度最大的向量，第二个变量选取一个使目标函数减小最快的变量。

（摘自【机器学习详解】SMO算法剖析）

选取第一个变量

第一个变量的选择称为外循环.

首先遍历整个样本集，选择违反KKT条件的αi作为第一个变量，接着依据相关规则选择第二个变量(见下面分析),对这两个变量采用上述方法进行优化。

当遍历完整个样本集后，遍历非边界样本集(0<αi<C)中违反KKT的αi作为第一个变量，同样依据相关规则选择第二个变量，对此两个变量进行优化。

当遍历完非边界样本集后，再次回去，遍历整个样本集寻找，即在整个样本集与非边界样本集上来回切换，寻找违反KKT条件的αi作为第一个变量。直到遍历完整个样本集后，没有违反KKT条件αi，然后退出。
边界上的样本对应的ααi=0或者αi=C，在优化过程中很难变化，然而非边界样本0<αi<C会随着对其他变量的优化有大的变化。

第二个变量的选择

SMO称第二个变量的选择过程为内循环.

假设在外循环中找到的第一个变量记为α1，那么第二个变量的选择我希望能使α2有较大的变化。由于α2是依赖于|E1−E2|，当Ei为正时，那么选择最小的Ei作为E2 ,如果Ei为负，选择最大Ei作为E2，通常为每个样本的Ei保存在一个列表中，选择最大的|E1−E2|来近似最大化步长。
有时按照上述的启发式选择第二个变量，不能够使得函数值有足够的下降，这时按下述步骤:

首先在非边界集上选择能够使函数值足够下降的样本作为第二个变量，
如果非边界集上没有，则在整个样本集上选择第二个变量，
如果整个样本集依然不存在，则重新选择第一个变量。

后记

SVM从头到尾推一边，发现许多概念都更熟悉了。许多推导步骤，我也尽可能写在博客中了，工作量的确挺大的。下次再把SVM的代码仔细研究研究。

参考资料：

《机器学习》周志华

【机器学习详解】SMO算法剖析

[支持向量机（五）SMO算法](http://www.cnblogs.com/jerrylead/archive/2011/03/18/1988419.html)

解密SVM系列（三）：SMO算法原理与实战求解