支持向量机之SMO算法

算法的提出

众所周知,支持向量机的学习问题可以形式化为求解凸二次规划问题。这样的凸二次规划问题具有全局最优解,并且有许多优化算法可以用于这一问题的求解,但是当训练样本容量很大的时候,这些算法往往会变得非常低效。其实这个问题也蛮显然的,在SVM中,目标函数\(\min_{\alpha} {1\over 2}\sum_{i=1}^N\sum_{j=1}^N \alpha_i \alpha_j y_iy_jK(x_j , x_j) - \sum_{i=1}^N\alpha_i\) 中,变量是拉格朗日乘子,一个变量\(\alpha_i\)对应于一个样本点\((x_i , y_j)\),而变量的总数等于训练样本的容量N,如果训练样本容量N非常大的话,\(\alpha_i\)的求解会变得很慢,这也是为什么要提出SMO算法的原因

SMO优化算法

SMO算法的全名是序列最小最优化算法,由Microsoft Research的John C. Platt在1998年提出,并成为最快的二次规划优化算法,特别针对线性SVM和数据稀疏时性能更优。关于SMO最好的资料就是他本人写的《Sequential Minimal Optimization A Fast Algorithm for Training Support Vector Machines》了。

细节总结

首先这里先列出SVM目标函数: \[\max_{\alpha} W(\alpha) = \sum_{i=1}^m\alpha_i - {1\over 2}\sum_{i,j=1}^my^{(i)}y^{(j)}\alpha_i \alpha_jK(x^{(1)} , x^{(2)}) \\ s.t. 0\le \alpha_i \le C , i = 1,\ldots,m \\ \sum_{i=1}^m\alpha_iy^{(i)} = 0\]

要解决的是在参数\(\lbrace \alpha_1, \alpha2,\ldots,\alpha_n \rbrace\)上求最大值\(W\)的问题,至于\(x^{(i)}\)\(y^{(i)}\)都是已知数(都是训练集里的)。\(C\)由我们预先设定,也就是惩罚因子,也是已知数。

SMO就是用来求\(\alpha_i\)的,SMO是一种启发式算法,其基本思路是:如果所有变量的解都满足此最优化问题的KKT条件(Karush-Kuhn-Tucker conditions),那么这个最优化问题的解也就得到了。因为KKT条件本来就是该优化问题的充分必要条件。否则,选择两个变量,固定其他变量,针对这两个变量构建一个二次规划问题。这个二次规划问题关于这两个变量的解应该更接近原始二次规划问题的解,因为这会使得原始二次规划问题的目标函数值变得更小。重要的是,这时子问题可以通过解析的方法求解,这样就可以大大提高整个算法的计算速度。子问题有两个变量,一个是违反KKT条件最严重的哪一个,另一个由约束条件自动确定,如此,SMO算法将原问题不断分解为子问题并对子问题进行求解,进而达到求解原问题的目的。

按照坐标上升的思路,我们首先固定除\(\alpha_1\)以外的所有参数,然后在\(\alpha_1\)上求极值。等一下,这个思路有问题,因为如果固定\(\alpha_1\)以外的所有参数,那么\(\alpha_1\)将不再是变量(可以由其他值推出),因为问题中规定了:\[\alpha_1y^{(1)} = -\sum_{i=2}^m\alpha_iy^{(i)}\] 因此,我们需要一次选取两个参数做优化,比如\(\alpha_1\)\(\alpha_2\),此时\(\alpha_2\)可以由\(\alpha_1\)和其他参数表示出来。这样回代到\(W\)中,\(W\)就只是关于\(\alpha_1\)的函数了,可解。

这样,SMO的主要步骤如下: SMO的主要步骤

意思是,第一步选取一对\(\alpha_i\)\(\alpha_j\),选取方法使用启发式方法(后面讲)。第二步,固定除\(\alpha_i\)\(\alpha_j\)之外的其他参数,确定\(W\)极值条件下的\(\alpha_i\)\(\alpha_j\)\(\alpha_i\)表示。

SMO之所以高效就是因为在固定其他参数后,对一个参数优化过程很高效。

下面讨论具体方法:

假设我们选取了初始值\(\lbrace \alpha_1, \alpha2,\ldots,\alpha_n \rbrace\)满足了问题中的约束条件。接下来,我们固定\(\lbrace \alpha_3, \alpha_4,\ldots,\alpha_n \rbrace\),这样\(W\)就是\(\alpha_1\)\(\alpha_2\)的函数。并且\(\alpha_1\)\(\alpha_2\)满足条件: \[ \alpha_1y^{(1)} + \alpha_2y^{(2)} = -\sum_{i=3}^m\alpha_iy^{(i)} \] 由于\(\lbrace \alpha_3, \alpha_4,\ldots,\alpha_n \rbrace\)都是已知固定值,因此为了方面,可将等式右边标记成实数值\(\zeta\). \[ \alpha_1y^{(1)} + \alpha_2y^{(2)} = \zeta \]\(y^{(1)}\)\(y^{(2)}\)异号时,也就是一个为1,一个为-1时,他们可以表示成一条直线,斜率为1。如下图: \alpha_1 -\alpha_2 =\zeta

横轴是\(\alpha_1\),纵轴是\(\alpha_2\)\(\alpha_1\)\(\alpha_2\)既要在矩形方框内,也要在直线上,因为\(0\le \alpha_i \le C\) , 因此 \[ L = \max(0,\alpha_2 -\alpha_1) \quad H = \min(C, C+\alpha_2 -\alpha_1) \] 为什么呢?注意这里的\(\alpha_2 -\alpha_1 = -\zeta\),而\(-\zeta\)显然就是这条直线的截距,那么可想而知为什么要在这里这样写了。

同理,当\(y^{(1)}\)\(y^{(2)}\)同号的时候, \[ L = \max (0,\alpha_2 + \alpha_1 -C) \quad H = \min (C, \alpha_2 + \alpha1) \] 然后我们打算将\(\alpha_1\)\(\alpha_2\)表示: \[ \alpha_1 = (\zeta - \alpha_2y^{(2)})y^{(1)} \] 然后反代入\(W\)中,得: \[ W(\alpha_1 , \alpha_2 , \ldots , \alpha_m) = W((\zeta - \alpha_2y^{(2)})y^{(1)}, \alpha_2 , \ldots , \alpha_m) \] 展开后\(W\)可以表示成\(a\alpha_2^2 + b\alpha_2 + c\). 其中a,b,c是固定值。这样,通过对\(W\)进行求导可以得到\(\alpha_2\),然而要保证\(\alpha_2\)满足\(L\le \alpha_2 \le H\),我们使用\(\alpha_2^{new , unclipped}\)表示求导求出来的\(\alpha_2\),然而最后的\(\alpha_2\),要根据下面情况得到: \[ \alpha_2^{new} = \begin{cases} H, & \text{if $\alpha_2^{new,unclipped} \gt H$ } \\ \alpha_2^{new,unclipped} , & \text{if $L\le \alpha_2^{new,unclipped} \le H$ } \\ L, & \text{if $\alpha_2^{new,unclipped} \lt L$} \end{cases} \] 这样得到\(\alpha_2^{new}\)后,我们可以得到\(\alpha_1\)的新值\(\alpha_1^{new}\).

下面来找到启发式搜索的方法和求b值的公式。

这里首先定义特征到结果的输出函数为: \[ u = \vec w \cdot \vec x -b \tag{1} \] 原始的优化问题为: \[ \min_{w,b} {1\over 2}\left \| \vec w\right \|^2 \quad s.t. y_i(\vec w \cdot \vec x_i -b) \ge 1 , \forall i \tag{3} \] 求导得到: \[ \vec w = \sum_{i=1}^Ny_i\alpha_i \vec x_i , \quad b = \vec w \cdot \vec x_k - y_k \quad for \quad some \quad \alpha_k \gt 0. \tag{7} \] 经过对偶后为: \[ \min_{\alpha}\psi (\vec \alpha) = \min_{\alpha}{1\over 2}\sum_{i=1}^N \sum_{j=1}^Ny_iy_j(\vec x_i \cdot \vec x_j)\alpha_i \alpha_j - \sum_{i=1}^N\alpha_i , \\ s.t. \alpha_i \ge 0 , \forall i \\ \sum_{i=1}^Ny_i\alpha_i = 0. \]

这里与W函数是一样的,只是符号求反后,变成求最小值了。\(y_i \in \lbrace +1,-1 \rbrace\)

经过加入松弛变量\(\zeta\)后,模型修改为: \[ \min_{w,b,\zeta} {1\over 2}\left \| \vec w \right \|^2 + C \sum_{i=1}^N\xi_i \tag{8}\]

\[ 0 \le \alpha_i \le C , \forall i. \tag{9} \\ \]

由公式(7)代入(1)中可知, \[ u = \sum_{j=1}^Ny_j\alpha_j K(\vec x_j , \vec x)-b \tag{10} \] 这个过程和之前的对偶过程一样。

重新整理我们要求的问题为: \[ \min_{\alpha}\psi (\vec \alpha) = \min_{\alpha} {1\over 2}\sum_{i=1}^N \sum_{j=1}^Ny_iy_jK(\vec x_i, \vec x_j)\alpha_i \alpha_j - \sum_{i=1}^N\alpha_i , \\ 0 \le \alpha_i \le C , \forall i , \\ \sum_{i=1}^Ny_i\alpha_i = 0. \tag(11) \\ \] 与之对应的KKT条件为: \[ \alpha_i = 0 \leftrightarrow y_i u_i \ge 1. \\ 0 \lt \alpha_i \lt C \leftrightarrow y_iu_i = 1 \\ \alpha_i = C \leftrightarrow y_iu_i \le 1 \tag{12} \] 这个KKT条件说明,在两条间隔线外面的点,对应前面的系数\(\alpha_i\)为0,在两条间隔线里面的对应\(\alpha_i\)为C,在两条间隔线上的对应的系数\(\alpha_i\)在0和C之间。

将之前的\(L\)\(H\)重新拿过来: \[ L = \max(0,\alpha_2 -\alpha_1) \quad H = \min(C, C+\alpha_2 -\alpha_1) \tag{13} \]

\[ L = \max (0,\alpha_2 + \alpha_1 -C) \quad H = \min (C, \alpha_2 + \alpha1) \tag{14} \] 之前我们将问题进行到这里,然后说将\(\alpha_1\)\(\alpha_2\)表示后代入\(W\)中,这里将代入\(\psi\)中,得: \[ \psi = {1\over 2}K_{11}\alpha_1^2 + {1\over 2}K_{22}\alpha_2^2 + sK_{12}\alpha_1 \alpha_2 + y_1\alpha_1v_1 + y_2\alpha_2 v_2 - \alpha_1 - \alpha_2 + \psi_{constant} \tag{24} \] 其中: \[ K_{ij} = K(\vec x_i , \vec x_j) , \\ v_i = \sum_{j = 3}^Ny_j\alpha_j^*K_ij = u_i+b^* - y_1\alpha_1^*K_{1i} - y_2\alpha_2^*K_{2i} , \tag{25} \] 这里的\(\alpha_1^*\)\(\alpha_2^*\)代表某次迭代前的原始值,因此是常数,而\(\alpha_1\)\(\alpha_2\)是变量,待求。公式(24)中的最后一项是常数。

由于\(\alpha_1\)\(\alpha_2\)满足以下公式: \[ y_1\alpha_1^* + y_2\alpha_2^* = -\sum_{i=3}^n y_i\alpha_i^* = y_1\alpha_1 + y_2\alpha_2\] 因为\(\alpha_i^* (i \gt 2)\)的值是固定值,在迭代前后不会变。

那么用\(s\)表示\(y_1y_2\),上式两边乘以\(y_1\)时,变为: \[ \alpha_1 + s\alpha_2 = \alpha_1^* + s\alpha_2^* = w \tag{26} \] 其中: \[ w = -y_1\sum_{i=3}^ny_i\alpha_i^* \] 代入(24)中,得: \[ \psi = {1\over 2}K_{11}(w-s\alpha_2)^2 + {1\over 2}K_{22}\alpha_2^2 + sK_{12}(w-s\alpha_2) \alpha_2 + y_1(w-s\alpha_2)v_1 + y_2\alpha_2 v_2 - (w-s\alpha_2) - \alpha_2 + \psi_{constant} \tag{27} \] 这时候只有\(\alpha_2\)是变量了,求导: \[ {d\psi \over d\alpha_2} = -sK_{11}(w-s\alpha_2) + K_{22}\alpha_2 - K_{12}\alpha_2 + sK_{12}(w-s\alpha_2) - y_2v_1 + s + y_2v_2 -1 = 0 \tag{28} \] 如果\(\psi\)的二阶导数大于0(凹函数),那么一阶导数为0时,就是极小值了。

假设其二阶导数为0(一般成立),那么上式化简为: \[ \alpha_2 (K_{11}+ K{22}-2K{12}) = s(K_{11} - K_{12})w + y_2(v_1 - v_2) + 1 -s \tag{29} \]\(w\)和\(v\)代入后,继续化简推导,得(推导了六七行推出来了) \[ \alpha_2 (K_{11}+ K{22}-2K{12}) = \alpha_2^*(K_{11}+ K{22}-2K{12})+ y_2(u_1-u_2+y_2-y_1) \tag{30} \] 我们使用\(\eta\)来表示: \[ \eta = K_{11}+ K{22}-2K{12} \] 通常情况下目标函数是正定的,也就是说,能够在直线约束方向上求得最小值,并且\(\eta \gt 0\)

那么我们在(30)两边都除以\(\eta\)可以得到: \[ \alpha_2^{new} = \alpha_2 + {y_2(E_1-E_2) \over \eta} \tag{16} \] 这里我们使用\(\alpha_2^{new}\)表示优化后的值,\(\alpha_2\)是迭代前的值,\(E_i = u_i-y_i\)

与之前提到的一样\(\alpha_2^{new}\)不是最终迭代后的值,需要进行约束: \[ \alpha_2^{new,clipped} = \begin{cases} H, & \text{if $\alpha_2^{new} \gt H$ } \\ \alpha_2^{new} , & \text{if $L\le \alpha_2^{new} \le H$ } \\ L, & \text{if $\alpha_2^{new} \lt L$} \end{cases} \tag{17} \] 那么: \[ \alpha_1^{new} = \alpha_1 + s(\alpha_2 - \alpha_2^{new,clipped}) \tag{18} \] 特殊情况下,\(\eta\)可能不为正,如果核函数\(K\)不满足Mercer定理,那么目标函数可能变得非正定,\(\eta\)可能出现负值。即使\(K\)是有效的核函数,如果训练样本中出现相同的特征\(x\),那么\(\eta\)仍有可能为0。SMO算法在\(\eta\)不为正值的情况下仍有效。为保证有效性,我们可以推导出\(\eta\)就是\(\psi\)的二阶导数,\(\eta \lt 0\)\(\psi\)没有极小值,最小值在边缘处取到(类比\(y=-x^2\)),\(\eta = 0\)时更是单调函数了,最小值也在边缘处取得,而\(\alpha_2\)的边缘就是\(L\)\(H\)。这样将\(\alpha_2 =L\)\(\alpha_2 = H\)分别代入\(\psi\)中即可求得\(\psi\)的最小值,相应的\(\alpha_2 =L\)还是\(\alpha_2 = H\)也可以知道了。具体计算公式如下: 具体计算公式

至此,迭代关系式出了\(b\)的推导式以外,都已经推出。

\(b\)每一步都要更新,因为前面的KKT条件指出了\(\alpha_i\)\(y_iu_i\)的关系,而\(u_i\)\(b\)有关,在每一步计算出\(\alpha_i\)后,根据KKT条件来调整\(b\)

b的更新有几种情况: b的更新情况–来自罗林开的ppt

这里的界内指\(0 \lt \alpha_i \lt C\),界上就是等于0或者C了。 前面两个的公式推导可以根据\(y_1\alpha_1^* + y_2\alpha_2^* = -\sum_{i=3}^ny_i\alpha_i^* = y_1\alpha_1 + y_2\alpha_2\) 和对于\(0 \lt \alpha_i \lt C\)\(y_iu_i = 1\)的KKT条件推出。

这样全部参数的更新公式都已经介绍完毕,附加一点,如果使用的是线性核函数,我们就可以继续使用\(w\)了,这样不用扫描整个样本库来作内积了。

w值的更新方法为: \[ \vec w^{new} = \vec w + y_1(\alpha_1^{new}- \alpha_1)\vec x_1 + y_2(\alpha_2^{new,clipped}-\alpha_2)\vec x_2 \tag{22} \] 根据前面的 \[ \vec w = \sum_{i=1}^Ny_i\alpha_i \vec x_i , \quad b = \vec w \cdot \vec x_k - y_k \quad for \quad some \quad \alpha_k \gt 0. \tag{7} \] 公式推导出。

SMO中拉格朗日乘子的启发式选择方法

终于到了最后一个问题了,所谓的启发式选择方法主要思想是每次选择拉格朗日乘子的时候,优先选择样本前面系数\(0 \lt \alpha_i \lt C\)\(\alpha_i\)作优化(论文中称为无界样例),因为在界上(\(\alpha_i\)为0或C)的样例对应的系数\(\alpha_i\)一般不会更改。

这条启发式搜索方法是选择第一个拉格朗日乘子用的,比如前面的\(\alpha_2\)。那么这样选择的话,是否最后会收敛。可幸的是Osuna定理告诉我们只要选择出来的两个\(\alpha_i\)中有一个违背了KKT条件,那么目标函数在一步迭代后值会减小。违背KKT条件不代表\(0 \lt \alpha_i \lt C\),在界上也有可能会违背。是的,因此在给定初始值\(\alpha_i =0\)后,先对所有样例进行循环,循环中碰到违背KKT条件的(不管界上还是界内)都进行迭代更新。等这轮过后,如果没有收敛,第二轮就只针对\(0 \lt \alpha_i \lt C\)的样例进行迭代更新。

在第一个乘子选择后,第二个乘子也使用启发式方法选择,第二个乘子的迭代步长大致正比于\(\left | E_1-E_2 \right |\),选择第二个乘子能够最大化\(\left | E_1-E_2 \right |\)。即当\(E_1\)为正时选择负的绝对值最大的\(E_2\),反之,选择正值最大的\(E_2\)

最后的收敛条件是在界内(\(0 \lt \alpha_i \lt C\))的样例都能够遵循KKT条件,且其对应的\(\alpha_i\)只在极小的范围内变动。

至于如何写具体的程序,请参考John C. Platt在论文中给出的伪代码。

总结

这份SVM的讲义重点概括了SVM的基本概念和基本推导,中规中矩却又让人醍醐灌顶。起初让我最头疼的是拉格朗日对偶和SMO,后来逐渐明白拉格朗日对偶的重要作用是将w的计算提前并消除w,使得优化函数变为拉格朗日乘子的单一参数优化问题。而SMO里面迭代公式的推导也着实让我花费了不少时间。

对比这么复杂的推导过程,SVM的思想确实那么简单。它不再像logistic回归一样企图去拟合样本点(中间加了一层sigmoid函数变换),而是就在样本中去找分隔线,为了评判哪条分界线更好,引入了几何间隔最大化的目标。

之后所有的推导都是去解决目标函数的最优化上了。在解决最优化的过程中,发现了w可以由特征向量内积来表示,进而发现了核函数,仅需要调整核函数就可以将特征进行低维到高维的变换,在低维上进行计算,实质结果表现在高维上。由于并不是所有的样本都可分,为了保证SVM的通用性,进行了软间隔的处理,导致的结果就是将优化问题变得更加复杂,然而惊奇的是松弛变量没有出现在最后的目标函数中。最后的优化求解问题,也被拉格朗日对偶和SMO算法化解,使SVM趋向于完美。

常见的问题

Q

违反kkt条件是什么意思,是不是类似:我遍历时选中了一个样本,计算其输出\(u_i\),然后计算\(y_i*u_i\),发现 \(y_i*u_i \lt 1\),但是该样本对应的\(\alpha\)值为:\(\alpha_i=0\)

违反kkt条件就是违反上面的公式(12)。真正实现的时候,会允许有一定的误差(在Platt论文中的伪代码中称为tolerance),给定一个样本点\(x_i\)后,首先使用\(u_i\)计算\(y_i^estimate\),然后计算\(E=y_i^estimate - y_i\),然后计算\(r = E * y_i\)。违反kkt条件的判断条件是\((r \lt -tolerance \land \alpha <C) \lor (r \gt tolerance \land \alpha \gt 0)\)

Q

在启发式的选取第一个Langrange乘子的过程是否可以理解为:首先遍历所有样本,如果某个样本违反kkt条件则选中它为第一个Langrange乘子,然后启发式的选第二个Langrange乘子,最后做alpha值的优化,如果所有样本都满足kkt条件,则遍历样本,选择界内 \(0 \lt \alpha \lt C\)的样本的\(\alpha\)值作为第一个Langrange乘子,然后启发式的选第二个Langrange乘子,最后做\(\alpha\)值的优化。

思想是你说的那样。在伪代码中,第一遍遍历所有样例,如果样例违反了kkt条件,先看是否存在其他alpha非0又非C的样例,如果存在使用第二个Langrange乘子启发式搜索算法来找第二个alpha做优化,不存在也没关系,接下来还有两个循环,第一个循环随便挑一个非0又非C的alpha来优化,第二个随便挑一个样例来优化。上面提到的三个优化中,如果任何一个对alpha值做了更新,那么会提前退出。我说的有点绕,你可以看看具体的伪代码,自己理解一下。总是就是优先选择不在界上的alpha做优化,因为不在界上的更容易违反kkt条件。

热评文章