SMO方法的实现及证明

1．

问题的阐述

SVM是从线性可分情形下的最优分类面发展而来。基本思想可以用图(2-16)的二维情况说明。

图2-16 线性可分情况下的最优分类线

图中实心点和空心点代表两类样本，H为分类线H1,H2分别为过各类中离分界线最近的样本且平行于分类线的直线，它们之间的距离叫做分类间隔(margin)。所谓最优分类线就是要求分类线不仅能将两类正确分开(训练错误率为0)，而且使分类间隔最大。推广到一般线性可分情形，假设分类方程为+b=0，对其进行归一化，样本集

(xi,yi),i=1,2,K,n,x∈Rd,y∈{+1,−1}，满足

yi(+b)−1≥0

构造损失函数作为目标函数及约束条件，即：

minimize: W(α)=

w2

2

(1)

+C∑ξi

i

(2-a) (2-b) (2-c)

subject to yiwTxi+b≥1−ξi,∀i

()ξi≥0,∀i

经过拉格朗日变换以及KKT定理推导，式子变为：

1

minimize:W(α)=∑αi−∑αiαjyiyjxixj

2i,ji

subject to 0≤αi≤C

(3)

∑αy

ii

i

=0

引入核函数，最后的目标函数变为：

1nn

maximize :W(α)=∑αi−∑∑αiαjyiyjK(xi,xj)2i=1j=1i=1

n

subject to ∑αiyi=00≤αi≤C,∀i

i=1

n

(4)

改写为矩阵相乘的格式，得到：

1T

αQα−eTα

2

subject to yTα=00≤αi≤C,i=1,....,lmin

f(α)=

(5)

其中e为全1向量，C为所有变量的上界，Q为l×l的半正定矩阵。训练向量xi通过φ函数被映射到更高维的空间(可能为无穷维)， Qij≡yiyjK(xi,xj)，其中K(x,y)=φ(x)φ(y)为核函数。简言之，即通过训练样本寻找α，使得

T

1T

αQα−eTα最小。 2

2. 最小贯序列方法(SMO：Sequential Minimal Optimization)

SMO是改进的SVM训练算法。同其他SVM训练 算法一样，SMO将一个大的QP问题分解为一系列小的QP问题。与其他算法不一样的是，SMO处理的是规模最小的QP问题，因此能够快速解决并获得解析解，能够有效的改善空间和时间复杂度。

SMO的优点：

1. 即使很多拉格朗日乘子在界上，SMO仍然能够较好地处理； 2. SMO训练线性核SVM时性能优异，因为SMO的计算时间主要集中在SVM的计算上，而线性核SVM的计算可以表示为简单的内积，而非线性核的和；

3. SMO能够很好地处理输入为稀疏向量的情形，因为核计算时间的减少使SVM训练过程加速。因为块分解法花费的时间主要在跑QP代码上，因此无法有效利用SVM的线性特性或者输入向量的稀疏性。 4. SMO最大的优点在于能够有效处理大规模QP问题，因为它减小训练集规模优于块分解方法。

使用最小贯序列方法（SMO）每次仅选择两个元素共同优化，在其他参数固定的前提下，找到这两个参数的最优值，并更新相应的α向量，而这两个点乘子的优化可以获得解析解。尽管需要更多的迭代才收敛，每次迭代需要很少的操作，因此算法在整体上的速度上有数量级的提高。包括

收敛时间在内，算法其他特征没有矩阵操作，它不需要存储核矩阵，所需要的内存大小和训练数据集的大小线性增长，因此也有效降低了空间复杂度。

通常来说，获得精确的最优值是不太现实的，因此需要定义近似最优条件，因此后面提及优化一词时，将代表近似最优解。实现技术主要考虑两个步骤：如何选择工作集；如何更新两个点乘子。 2.1工作集选择

因为SMO方法每一步选择和更新两个拉格朗日乘子，并且至少有一个

拉格朗日乘子在更新前违背了KKT条件，因此根据Osuna原理每一步能够减小目标函数值。因此迭代收敛能够得到保证。为了能够快速收敛，SMO启发式地选择两个拉格朗日乘子进行优化。 2.1.1 Platt 方法[1]

根据KKT条件，可以得出：

αi=0⇔yif(xi)≥1 0<αil

（6a） （6b） （6c）

其中，

f(xi)=∑j=1yjαjK(xi⋅xj)+b

（7）

Platt方法通过式（6）检测αi是否满足KKT条件。它将启发式选择第一个和第二个拉格朗日乘子视为的过程。选择第一个乘子作为SMO算法的外层循环。任何违反KKT条件(6)的乘子（或样本）都是合法的第一个拉格朗日乘子。外层循环中分两种情况选择第一个乘子：

A. 在所有乘子中选择（第1次或者αi∈(0,C)没有违反KKT条件乘子的情

形）；

B. 在界内αi∈(0,C)中选择（其它情形，这是常规情况）。

外层循环分两种情况选择第一个乘子的理由：当SMO算法有进展时，界上的乘子（即等于0或者C的乘子）更新后很有可能仍然停留在界上，界内的乘子才发生改变，因此首先在界内乘子中选择，可以加快算法的运行时间。

界上乘子更新后很有可能仍然停留在界上的解释：

上图红点表示待更新的两个乘子的对应坐标，其中α2在界上α2=0，假

设核矩阵正定，根据

newold

=α2+α2

y2(E1−E2)η

η>0,因此若α2更新后值发生变化，必须y2(E1−E2)>0，而在后面的工作集选

择并不能够确保该式成立，若y2(E1−E2)≤0，则α2的值保持不变，因此更新失败，仍然留在界上。其他情形依次类推。

总言之，外循环首先对整个训练集进行迭代，确定是否每一个实例都

违背KKT条件。如果一个实例违背了KKT条件，它就作为待优化的候选乘子。遍历一次整个训练集后（最多在整个训练集寻找一次违背KKT条件的

α2），外循环对所有拉格朗日乘子既非0也非C(non-bound，注意non-bound

是针对0≤αi≤C这一约束条件而言的)的例子进行迭代，重复进行，每个例子检查是否满足KKT条件，违背该条件的实例作为优化的候选实例。外循环不断重复对这些非边界实例的遍历直到所有的实例在ε内满足KKT条件。外循环返回并对整个训练集进行迭代。

需要注意的事KKT条件是在精度ε内检查的。一般来说，ε设为10−3。

识别系统一般不要求在较高精确度上检查KKT条件：例如正界输出0.999到1.001之间的数是可接受的。如果要求很高的输出精度，SMO算法将无法快速收敛。

一旦第一个乘子被选中，SMO以最大化优化时的步长为目标来选择第二

个乘子。现在计算核函数K需要一定的时间花销，根据

newold

=α2+α2

y2(E1−E2)η

其中

Ei=∑j=1αjyjK(xj,xi)−yi

n

Fi=∑j=1αjyjK(xj,xi)+b−yi =f(x)i−yi=Ei+b

n

第二个乘子的迭代步长大致正比于E1−E2，因此应选择第二个乘子最大化

E1−E2，即当E1为正时最小化E2，E1当为负时最大化E2。

使用上述第二个乘子的启发式选择在一些特殊的情况，SMO无法有正的进

展（positive progress）。例如，第一个和第二个乘子指向同一个输入向量。在这种情况，SMO对第二个乘子的启发式选择使用分层的方法指导它寻找到一个拉格朗日乘子对能够实现正的进展，所谓正的进展通过乘子优化后两个乘子的增量非零。分层的启发式选择第二个乘子的方法描述如下，如果上述的启发式选择无法实现正的进展，SMO开始在非边界实例中迭代，寻找能够实现正的进展的实例，如果没有非边界实例能够实现正的进展，SMO开始在整个训练集中迭代指导寻找到能够实现正的进展的实例。为了避免SMO倾向于训练集中位置靠前的训练集，在非边界实例中和在整个训练集中寻找都从一个随机位置开始。在一些极其特殊的情形，没有实例能够满足要求，此时，第一个实例将被忽略，SMO继续寻找第一个实例。

b,f(xi),F的更新：

new

为使乘子α1或者α2的KKT条件成立。若α1new在界内，则有E1=0根据

（1） b的更新

F1new=∑j=1αjyjK1j+b1new−yi=0F

old1

n

=∑j=1αjyjK1j+b−yi

n

可以得到：

newold

b1new=b−F1old−y1(α1new−α1old)K(x1,x1)−y2(α2−α2)K(x2,x1)

new

同理，若α2在界内，则有

newnewold

b2=b−F2old−y1(α1new−α1old)K(x2,x1)−y2(α2−α2)K(x2,x2) new

newnewold

和α2都在界内，根据α2=α2+

如果α1

y2(F1−F2)κold

=α2+

y2(E1−E2)κ不难得出

new

b1new=b2，则有

newA. 若α1new和α2有一个在界内，b取对应binew（i=1,2）的，如α1new在界内，b取

b1new；

newnewB. 若α1new和α2都在界上，那么b1new和b2之间的任意数都满足KKT条件，都可作

new

b1new+b2=。

2

为b的更新值，一般取b

new

new

下面给出bnew=tb1new+(1−t)b2,t∈(0,1)满足KKT条件的证明。由于两个乘

new

子都在界上，因此α2经过了裁剪，令

newold

α2=α2+λy2(F1−F2)η,λ∈[0,1)

(*)

化简得

i.

F1new=(1−λ)(1−t)(F1−F2)F

new2

=−t(1−λ）(F1−F2)

λ∈[0,1),t∈(0,1)

new

若α1new=α2（同为C或者0）

new

=C，由于（*）式，所以y2(F1−F2)>0，1) y1≠y2。不妨令α1new=α2

亦可得到y1F1new<0，因此更新的两个乘子都满足KKT可得y2F2new<0，

new

=0的情形。 条件，同理可证明α1new=α2

new

2) y1=y2，α1=α2=α1new=α2（此时，λ=0），此时最多只有一个乘子

满足KKT条件，根据（16-b）计算L和H，可得到L=H，因此遇到该情形时，将不会进行后面的更新。

ii.

new若α1new≠α2（一个为0，另一个为C）

new

=C,α1new=0，由于（*）式，所以y2(F1−F2)>0，则1) y1=y2不妨令α2

有y2F2new<0,亦可得到y1F1new>0，因此更新的两个乘子都满足KKT条

new

=0,α1new=C的情形。 件，同理可证明α2

new

2) y1≠y2，α2=α2,α1new=α1,y1F1new与y2F2new同号，不可能同时满足KKT

条件。根据（16-a）计算L和H，可得到L=H，因此遇到该情形时，将不会进行后面的更新。

（2） f(x)i的更新

f(x)i

new

new

=∑j=1αnewyK(x,x)+b jjji

n

（3） F(x)i的更新

F(x)i

new

=f(x)i

new

−yi

procedure examineExample(i2) y2 = target[i2]

alph2 = Lagrange multiplier for i2

F2 = SVM output on point[i2] – y2 (check in error cache) r2 = F2*y2

if ((r2 < -tol && alph2 < C) || (r2 > tol && alph2 > 0)) {

if (number of non-zero & non-C alpha > 1) {

i1 = result of second choice heuristic (section 2.2) if takeStep(i1,i2) return 1 }

loop over all non-zero and non-C alpha, starting at a random point {

i1 = identity of current alpha if takeStep(i1,i2) return 1 }

loop over all possible i1, starting at a random point {

i1 = loop variable if (takeStep(i1,i2) return 1

} } return 0 endprocedure main routine: numChanged = 0; examineAll = 1;

while (numChanged > 0 | examineAll) {

numChanged = 0; if (examineAll)

loop I over all training examples numChanged += examineExample(I) else

loop I over examples where alpha is not 0 & not C numChanged += examineExample(I) if (examineAll == 1) examineAll = 0

else if (numChanged == 0) examineAll = 1 }

该工作集选择方法存在一些不足，导致算法效率较低。因为它计算和

使用一个单一阈值(b)的方式，可能会带来一些不必要的低效。在任何情况，SMO根据两个已经优化的乘子来修改b。但是，当算法检查剩下的实例是否违背KKT条件时，很有可能存在一个不同的b值能够工作得更好。所以在SMO算法中，很有可能存在这种情况，尽管α已经迭代得到了一个满足bup≥blow的值，SMO却因为当前的b值不合适而无法检测到该情况。在这种情况下，SMO为了进行更新步骤不得不付出大量且不必要的资源以搜索第二个乘子。

2.1.2 Keerthi 方法[2]

(1)步骤

把训练样例分为五个集合: 1)

I0={i:0<αi2) 计算bup和blow,α1和α2就是bup和blow对应的乘子：

(8)

bup=min{Ei:i∈I0∪I1∪I2}blow=max{Ei:i∈I0∪I3∪I4}其中 Ei=∑j=1αjyjK(xj,xi)−yi

3) 计算bup≥blow，则终止迭代，否则进行两个乘子的更新操作。 (2)证明

根据KKT定理有：

n

(9)

∇f(α)+by=λ−μ

(10-a) (10-b)

λiαi>0,μi(C−αi)=0,λi≥0,μi≥0

注意：∇f(α)=Qα−e，即为f(α)的梯度，结合（10-a）和（10-b）可以得到：

当αi>0时，λi=0,μi≥0，则有

∇f(α)i+byi≤0

（11-a）

当αi∇f(α)i+byi≥0

（11-b）

通过（9-a）和（9-b）可知，当αi>0,yi=−1或者αiyi∇f(α)i≥−b

（12-a）

同理，当αi>0,yi=1或者αiyi∇f(α)i≤−b

（12-b）

因此，若α*为可行解，设m(α*)=min(∇f(α*)i),i∈Iup以及M(α*)=max(∇f(α*)i),i∈Ilow，则有m(α*)≥M(α*)，实际上，已经证明α*为可行驻点当且仅当m(α*)≥M(α*),因此可作为终止条件。 不难得到

∇f(α)i=yi∑j=1αjyjK(xi,xj)−1

l

（13）

由（10-a）及（10-b）可以得到，当0<αi∇f(α)i+byi=0，

l

（14）

将式（13）代入式子（14），化简得到

yi(∑j=1αjyjK(xi,xj)+b)=1

（15）

因此式（12-a）中的b即为式（2-b）中的b。

注意该方法较Platt方法而言，不会出现界上乘子（即等于0或者C的乘子）更新后很有可能仍然停留在界上的情况，以下图为例，红点表示更新前乘子对的对应值，(α2=0，在界上)，α2对应的E值为blow，因此y2=−1，根据

y2(E1−E2)

old

+α2new=α2

η

η>0，易有y2(E1−E2)>0，因此α2的值必定能够更新成功，依次类推，所以该方法

不需要像Platt方法那样先考虑更新非边界值。

该方法一些优化技巧：

1） 假定在任何情况下，Ei对所有i可用。假设i_low是blow对应的索引，i_up是bup对应的索引，此时，寻找特定的用于优化的i非常方便。例如，假设

i∈I1∪I2，我们只需要检查Ei现了违反对{i_low,i}。因此，寻找违反对是同时进行的，这与原始的SMO算法是不同的。正如下面所看到的，我们通过高效的更新方法计算和使用（i_low，blow）和(i_up,bup)。

2) 为了达到高效，我们将在更新αi,i∈I0做出大量努力，缓存Ei，i∈I0被保持和更新。对所有的i∈I0都满足优化条件时，将检查所有的索引值进行优化。

3) 一些额外的步骤需要添加到TakeStep过程中。当使用（i2,i1）作为索引对成功更新后,令I=I0∪{i1,i2}。我们只计算不完全的（i_low，blow）和(i_up,bup)，也就是仅使用I（I表示被更新过的元素，因此索引值i对应的Ei有可能为0或者C）。注意这些额外的步骤只需要少量的花费，因为缓存Ei，i∈I0是可用的，并且能很快地更新Ei1和Ei2。因为（i_low，blow）和

(i_up,bup)可能从{i2,i1}中产生并且在后续的步骤中作为i1的选择，因此我们将Ei1和Ei2的值存储在缓存中。

改进：实验表明直接从活动集中寻找(i2,i1)比上述算法更加高效，作者考虑使用上述步骤可能是受Platt的影响，认为界上乘子很有可能更新失败，实际上，使用该方法更新边界乘子一定成功（已证明），因此很有可能得到更大的步长加速收敛。

4） 当仅作用在αi,i∈I0时，也即一个ExamineAll=0的循环。注意当bup≤blow−2τ在某个点成立，说明对应的乘子为最大违反对(maximal violating pair)，否则跳出ExamineAll=0的循环。

5) Keerthi提出了两种实现在I0集合中选择最大违反对的方法（ExamineAll=0）

方法1：遍历i2∈I0，对每个i2，检验是否违反KKT条件，如果违反，选择方法2：总是选择最坏的违反对；选择i2=i_low,i1=i_up。

需要注意的是，方法1和Platt方法除了判断是否满足优化条件外其余都是相合适的i1，例如，检验是否满足Ei2≤blow−2τ,如果满足，可选择i1=i_low。

同的。方法2可以认为是方法1的进一步优化，因为在方法2有效利用了缓存寻找违反对。

6）当所有i∈I0对应的乘子都满足最优条件时，即I0中找不到违反对，算法返回到(ExamineAll=1)，在整个训练集中检查最有条件。我们依次遍历所有索引。因为通过I0我们计算出部分的（i_low，blow）和(i_up,bup)，检查每一个索引i并对这些值进行更新。对于一个给定的i，先计算Ei，然后通过当前的（i_low，blow）和(i_up,bup)检查最优条件，如果不违反，Ei用于更新这些值。例如，若i∈I1∪I2并且Ei≤blow−2τ,表明(i,i_low)是一个违反对，对应两个乘子进行更新操作，否则Ei用于更新(i_up,bup)，也即Ei7） 当ExamineAll=1并且无违反对出现时，所有索引对应的乘子满足最优条件，优化结束。

2．2 更新两个乘子α`1和α`2[3]

对条件∑j=1αjyj=0需要在迭代中实现，意味着每步能优化的乘子的最小个

l

数为2；无论何时一个乘子被更新，至少需要调整另一个乘子来保持条件成立。因此首先需要对首先被更新的乘子进行上下界约束，使得两个乘子同时满足

0<αi1) 计算上下界L和H,若L=H时返回。

oldold

L=max0,α2−α1old,H=minC,C+α2−α1old，if y1≠y2 （16-a）

oldoldL=max(0,α2+α1old−C)，H=min(C,α2+α1old)，if y1=y2

()()（16-b）

2) 计算Ws的二阶导数

η=2φ(x1)Tφ(x2)−φ(x1)Tφ(x1)−φ(x2)Tφ(x2)

3) 更新α2

newold

α2=α2+

n

y2(E1−E2)η

Ei=∑j=1αjyjK(xj,xi)−yi

4) 计算修剪后的α2

new

⎧H，α2≥H⎪newnew=⎨α2，L≤α2≤H ⎪new

≤LL，α2⎩

new,clip

α2

5) 更新α1

oldnew,clip

α1new=α1old+y1y2(α2−α2)

6) 更新所有Ei

newoldEinew=Eiold+y1(α1new−α1old)K(x1,xi)+y2(α2−α2)K(x2,xi)

在一些特殊的场合，η将不为正。当核函数不满足Mercer条件时将会计算

得到负的η，会让目标函数变为不定的。即使核函数正确也可能会计算得到一个为0的η,假如多个实例具有相同的输入向量。在任何情况下，SMO都能够起作用，即使η不为正，此时应该在线段两端点处计算目标函数Ψ的值。

f1=y1E1−α1K11−sα2K12;f2=y2E2−sα1K12−α2K22;L1=α1+s(α2−L);

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

121

L1K11+L2K22+sLL1K12;22121

ψH=H1f1+Hf2+L1K11+H2K22+sHH1K12;

22

SMO将拉格朗日乘子移到目标函数最小对应的端点处。如果端点处对应的目标函数相同（在一个小的误差范围ε内），并且核满足Mercer条件，此时，更新乘子的步骤将不能够有正的进展，此时将重新选择第二个乘子。 2．3 计算b

当所有乘子都满足KKT条件后，所有的α值确定，也即经过核映射后的超平面法向量确定。只需确定b的值，分类平面便确定。根据原始式

ψL=L1f1+Lf2+

minimize: W(α)=

w

2

2

+C∑ξi

i

subject to yiwTxi+b≥1−ξi,∀i

()

ξi≥0,∀i

由于w随着α值确定而确定，所以转换为

minimize: W'(b)=∑ξi

i

subject to yiwTxi+b≥1−ξi,∀i

()

ξi≥0,∀i

yiwTxi+b≤1,∀i

将所有非支持向量舍弃，对于所有支持向量，由KKT条件，有

()ξi=1−yi(wTxi+b),∀i

所以最后优化问题转化为

maximize: W''(b)=b∑yi

i

subject to 1−yiwTxi+b≥0,∀i

()

约束式会得到b的上界upperBound和下界lowerBound，即一个取值区间。若∑yi为正，则取上界，若∑yi为负，则取下界，若∑yi为0，则取上下界的均值。

i

i

i

2.2.2 证明 定义:

vi=∑j=3αjyjK(xi,xj)=f(xi)−∑j=1αjyjK(xi,xj)−b,i=1,2

l

2

考虑到将α1和α2的函数作为目标：

11222

W(α1,α2)=(α1+α2)−α12φ(x1)−α2φ(x2)22 T

−α1α2y1y2φ(x1)φ(x2)−y1α1v1−y2α2v2+常数

l

l

注意约束∑j=1αjoldyj=∑j=1αjyj=0意味着条件：

old

α1+sα2=常数=α1old+sα2=γ

new

这里s=y1y2，利用这个方程可以用α2计算α1new，此约束下目标函数可写为：

W(α2)=γ−sα2+α2−

112

K11(γ−sα2)2−α2K22

22

−sK12(γ−sα2)α2−y1(γ−sα2)v1−y2α2v2+常数

对W(α2)求偏导，驻点应满足:

∂W(α2)

=1−s+sK11(γ−sα2)−α2K22

∂α2

由此得出：

+K12α2−sK12(γ−sα2)+y2v1−y2v2 =0

new

(K11+K12+2K22)=1−s+γs(K11−K12)+y2(v1−v2)α2

=y2(y2−y1+γy1(K11−K12)+v1−v2)

令κ=K11+K12+2K22，因此

newα2κy2=y2−y1+f(x1)−∑j=1αjyjK(x1,xj)+γy1K11

2

−f(x2)+∑j=1αjyjK(x2,xj)−γy2K12

2

=α2y2κ+(f(x1)−y1−b)−(f(x2)−y2−b)

给出： 得证。

newold需要注意的是，尽管有α2=α2−

newold

α2=α2+

y2(E1−E2)κ

y2(E1−E2)κold

的形式，α2实际上是凑出

old

来的，目的是为了方便计算，在等号右边的式子中，化简后最终将看不到α2的

身影（类似y=-x+(x-1)）。

这种算法得到了非常广泛的应用，然而f(αk)与f(αk+1)并无显式的对比，目

前尚并无收敛的充分证明，因此带有启发性质。是一种重要的思路。另一种思路是通过f(αk)与f(αk+1)进行显式的对比，在每一步中f(αk+1)newold

令αk+1=αk+d，其中d1=α1new−α1old,d2=α2，令s=y1y2，由于−α2

newold

，可以得到d1=−sd2，则有 α1old+sα2=α1new+sα2

f(αk+1)−f(αk)≈∇f(αk)d

=∇f(αk)1d1+∇f(αk)2d2=d2(∇f(α)2−s∇f(α)1)=y2d2(F2−F1)

y2(E1−E2)k

k

由于d2=λf(α

κk+1

,λ∈(0,1]，当且仅当被裁剪时λ∈(0,1)；

<0

)−f(α

k

2(E1−E2))=−λκ 因此，上述算法计算f(αk)能够保证收敛。

[1] Chih-Jen Lin.On the Convergence of the Decomposition Method for Support Vector Machines

[2] S.S. Keerthi Improvements to Platt’s SMO Algorithm for SVM Classifier Design [3] John C. Platt Sequential Minimal Optimization: A Fast Algorithm for Training Support Vector Machines

[4] Osuna, E. Improved Training Algorithm for Support Vector Machines