手把手教你实现SVM算法

什么是机器学习 （Machine Learning）

机器学习是研究计算机怎样模拟或实现人类的学习行为，以获取新的知识或技能，重新组织已有的知识结构使之不断改善自身的性能。它是人工智能的核心，是使计算机具有智能的根本途径，其应用遍及人工智能的各个领域。

机器学习的大致分类：

1）分类（模式识别）：要求系统依据已知的分类知识对输入的未知模式（该模式的描述）作分析，以确定输入模式的类属，例如手写识别（识别是不是这个数）。

2）问题求解：要求对于给定的目标状态,寻找一个将当前状态转换为目标状态的动作序列。

SVM一般是用来分类的（一般先分为两类，再向多类推广一生二，二生三，三生万物哈）

问题的描述

向量表示:假设一个样本有n个变量(特征)：Ⅹ= (X1,X2,…,Xn)T

样本表示方法：

SVM线性分类器

SVM从线性可分情况下的最优分类面发展而来。最优分类面就是要求分类线不但能将两类正确分开(训练错误率为0),且使分类间隔最大。SVM考虑寻找一个满足分类要求的超平面,并且使训练集中的点距离分类面尽可能的远,也就是寻找一个分类面使它两侧的空白区域(margin)最大。

过两类样本中离分类面最近的点且平行于最优分类面的超平面上H1,H2的训练样本就叫做支持向量。

图例：

问题描述：

假定训练数据 ：

可以被分为一个超平面：

进行归一化：

此时分类间隔等于：

即使得：最大间隔最大等价于使

最小

下面这两张图可以看一下，有个感性的认识。那个好？

看下面这张图：

下面我们要开始优化上面的式子，因为推导要用到拉格朗日定理和KKT条件，所以我们先了解一下相关知识。在求取有约束条件的优化问题时，拉格朗日乘子法（Lagrange Multiplier) 和KKT条件是非常重要的两个求取方法，对于等式约束的优化问题，可以应用拉格朗日乘子法去求取最优值；如果含有不等式约束，可以应用KKT条件去求取。当然，这两个方法求得的结果只是必要条件，只有当是凸函数的情况下，才能保证是充分必要条件。KKT条件是拉格朗日乘子法的泛化。之前学习的时候，只知道直接应用两个方法，但是却不知道为什么拉格朗日乘子法（Lagrange Multiplier) 和KKT条件能够起作用，为什么要这样去求取最优值呢？

拉格朗日乘子法和KKT条件

定义：给定一个最优化问题：

最小化目标函数：

制约条件：

定义拉格朗日函数为：

求偏倒方程

可以求得

的值。这个就是神器拉格朗日乘子法。

上面的拉格朗日乘子法还不足以帮我们解决所有的问题，下面引入不等式约束

最小化目标函数：

制约条件变为：

定义拉格朗日函数为：

可以列出方程：

新增加的条件被称为KKT条件

KKT条件详解

对于含有不等式约束的优化问题，如何求取最优值呢？常用的方法是KKT条件，同样地，把所有的不等式约束、等式约束和目标函数全部写为一个式子L(a, b, x)= f(x) + a*g(x)+b*h(x)，KKT条件是说最优值必须满足以下条件：

1. L(a, b, x)对x求导为零；

2. h(x) =0;

3. a*g(x) = 0;

求取这三个等式之后就能得到候选最优值。其中第三个式子非常有趣，因为g(x)<=0，如果要满足这个等式，必须a=0或者g(x)=0. 这是SVM的很多重要性质的来源，如支持向量的概念。

二. 为什么拉格朗日乘子法（Lagrange Multiplier) 和KKT条件能够得到最优值？

为什么要这么求能得到最优值？先说拉格朗日乘子法，设想我们的目标函数z = f(x), x是向量, z取不同的值，相当于可以投影在x构成的平面（曲面）上，即成为等高线，如下图，目标函数是f(x, y)，这里x是标量，虚线是等高线，现在假设我们的约束g(x)=0，x是向量，在x构成的平面或者曲面上是一条曲线，假设g(x)与等高线相交，交点就是同时满足等式约束条件和目标函数的可行域的值，但肯定不是最优值，因为相交意味着肯定还存在其它的等高线在该条等高线的内部或者外部，使得新的等高线与目标函数的交点的值更大或者更小，只有到等高线与目标函数的曲线相切的时候，可能取得最优值，如下图所示，即等高线和目标函数的曲线在该点的法向量必须有相同方向，所以最优值必须满足：f(x)的梯度 = a* g(x)的梯度，a是常数，表示左右两边同向。这个等式就是L(a,x)对参数求导的结果。（上述描述，我不知道描述清楚没，如果与我物理位置很近的话，直接找我，我当面讲好理解一些，注：下图来自wiki）。

而KKT条件是满足强对偶条件的优化问题的必要条件，可以这样理解：我们要求min f(x), L(a, b, x) = f(x) + a*g(x) + b*h(x)，a>=0，我们可以把f(x)写为：max_{a,b} L(a,b,x)，为什么呢？因为h(x)=0, g(x)<=0，现在是取L(a,b,x)的最大值，a*g(x)是<=0，所以L(a,b,x)只有在a*g(x) = 0的情况下才能取得最大值，否则，就不满足约束条件，因此max_{a,b} L(a,b,x)在满足约束条件的情况下就是f(x)，因此我们的目标函数可以写为 min_x max_{a,b} L(a,b,x)。如果用对偶表达式： max_{a,b} min_x L(a,b,x)，由于我们的优化是满足强对偶的（强对偶就是说对偶式子的最优值是等于原问题的最优值的），所以在取得最优值x0的条件下，它满足 f(x0) = max_{a,b} min_x L(a,b,x) = min_x max_{a,b} L(a,b,x) =f(x0)，我们来看看中间两个式子发生了什么事情：

f(x0) = max_{a,b} min_x L(a,b,x) = max_{a,b} min_x f(x) + a*g(x) + b*h(x) = max_{a,b} f(x0)+a*g(x0)+b*h(x0) = f(x0)

可以看到上述加黑的地方本质上是说 min_x f(x) + a*g(x) + b*h(x) 在x0取得了最小值，用Fermat定理，即是说对于函数 f(x) + a*g(x) + b*h(x)，求取导数要等于零，即

f(x)的梯度+a*g(x)的梯度+ b*h(x)的梯度 = 0

这就是KKT条件中第一个条件：L(a, b, x)对x求导为零。

而之前说明过，a*g(x) = 0，这时KKT条件的第3个条件，当然已知的条件h(x)=0必须被满足，所有上述说明，满足强对偶条件的优化问题的最优值都必须满足KKT条件，即上述说明的三个条件。可以把KKT条件视为是拉格朗日乘子法的泛化。

上面跑题了，下面我继续我们的SVM之旅。

经过拉格朗日乘子法和KKT条件推导之后

最终问题转化为：

最大化：

条件：

这个是著名的QP问题。决策面：

其中

为问题的优化解。

松弛变量（slack vaviable）

由于在采集数据的过程中，也可能有误差（如图）

所以我们引入松弛变量对问题进行优化。

式子就变为

最终转化为下面的优化问题：

其中的C是惩罚因子，是一个由用户去指定的系数，表示对分错的点加入多少的惩罚，当C很大的时候，分错的点就会更少，但是过拟合的情况可能会比较严重，当C很小的时候，分错的点可能会很多，不过可能由此得到的模型也会不太正确。

上面那个个式子看似复杂，现在我带大家一起推倒一下

……

…(草稿纸上，敲公式太烦人了)

最终得到：

最大化：

条件：

呵呵，是不是感觉和前面的式子没啥区别内，亲，数学就是这么美妙啊。

这个式子看起来beautiful，但是多数情况下只能解决线性可分的情况，只可以对线性可分的样本做处理。如果提供的样本线性不可分，结果很简单，线性分类器的求解程序会无限循环，永远也解不出来。但是不怕不怕。我们有杀手锏还没有出呢。接着咱要延伸到一个新的领域：核函数。嘻嘻，相信大家都应该听过这厮的大名，这个东东在60年代就提出来，可是直到90年代才开始火起来（第二春哈），主要是被Vapnik大大翻出来了。这也说明计算机也要多研读经典哈，不是说过时了就不看的，有些大师的论文还是有启发意义的。废话不多说，又跑题了。

核函数

那到底神马是核函数呢？

介个咱得先介绍一下VC维的概念。

为了研究经验风险最小化函数集的学习一致收敛速度和推广性，SLT定义了一些指标来衡量函数集的性能，其中最重要的就是VC维(Vapnik-Chervonenkis Dimension)。

VC维定义：对于一个指示函数（即只有0和1两种取值的函数）集，如果存在h个样本能够被函数集里的函数按照所有可能的2h种形式分开，则称函数集能够把h个样本打散，函数集的VC维就是能够打散的最大样本数目。

如果对任意的样本数，总有函数能打散它们，则函数集的VC维就是无穷大。

看图比较方便（三个点分类，线性都可分的）。

如果四个点呢？哈哈，右边的四个点要分为两个类，可能就分不啦。

如果四个点，一条线可能就分不过来啦。

一般而言,VC维越大, 学习能力就越强,但学习机器也越复杂。

目前还没有通用的关于计算任意函数集的VC维的理论,只有对一些特殊函数集的VC维可以准确知道。

N维实数空间中线性分类器和线性实函数的VC维是n+1。

Sin(ax)的VC维为无穷大。

对于给定的学习函数集,如何计算其VC维是当前统计学习理论研究中有待解决的一个难点问题，各位童鞋有兴趣可以去研究研究。

咱们接着要说说为啥要映射。

例子是下面这张图：

我们把横轴上端点a和b之间红色部分里的所有点定为正类，两边的黑色部分里的点定为负类。试问能找到一个线性函数把两类正确分开么？不能，因为二维空间里的线性函数就是指直线，显然找不到符合条件的直线。

但我们可以找到一条曲线，例如下面这一条：

显然通过点在这条曲线的上方还是下方就可以判断点所属的类别（你在横轴上随便找一点，算算这一点的函数值，会发现负类的点函数值一定比0大，而正类的一定比0小）。这条曲线就是我们熟知的二次曲线，它的函数表达式可以写为：

问题只是它不是一个线性函数，但是，下面要注意看了，新建一个向量y和a：

这样g(x)就可以转化为f(y)=<a,y>，你可以把y和a分别回带一下，看看等不等于原来的g(x)。用内积的形式写你可能看不太清楚，实际上f(y)的形式就是：

g(x)=f(y)=ay

在任意维度的空间中，这种形式的函数都是一个线性函数（只不过其中的a和y都是多维向量罢了），因为自变量y的次数不大于1。

看出妙在哪了么？原来在二维空间中一个线性不可分的问题，映射到四维空间后，变成了线性可分的！因此这也形成了我们最初想解决线性不可分问题的基本思路——向高维空间转化，使其变得线性可分。

而转化最关键的部分就在于找到x到y的映射方法。遗憾的是，如何找到这个映射，没有系统性的方法（也就是说，纯靠猜和凑）。具体到我们的文本分类问题，文本被表示为上千维的向量，即使维数已经如此之高，也常常是线性不可分的，还要向更高的空间转化。其中的难度可想而知。

为什么说f(y)=ay是四维空间里的函数?

大家可能一时没看明白。回想一下我们二维空间里的函数定义 g(x)=ax+b 变量x是一维的，为什么说它是二维空间里的函数呢？因为还有一个变量我们没写出来，它的完整形式其实是 y=g(x)=ax+b 即 y=ax+b 看看，有几个变量？两个。那是几维空间的函数？ 再看看 f(y)=ay 里面的y是三维的变量，那f(y)是几维空间里的函数？

用一个具体文本分类的例子来看看这种向高维空间映射从而分类的方法如何运作，想象一下，我们文本分类问题的原始空间是1000维的（即每个要被分类的文档被表示为一个1000维的向量），在这个维度上问题是线性不可分的。现在我们有一个2000维空间里的线性函数

f(x’)=<w’,x’>+b

注意向量的右上角有个 ’哦。它能够将原问题变得可分。式中的 w’和x’都是2000维的向量，只不过w’是定值，而x’是变量（好吧,严格说来这个函数是2001维的,哈哈），现在我们的输入呢，是一个1000维的向量x，分类的过程是先把x变换为2000维的向量x’，然后求这个变换后的向量x’与向量w’的内积，再把这个内积的值和b相加，就得到了结果，看结果大于阈值还是小于阈值就得到了分类结果。

你发现了什么？我们其实只关心那个高维空间里内积的值，那个值算出来了，分类结果就算出来了。而从理论上说， x’是经由x变换来的，因此广义上可以把它叫做x的函数（有一个x，就确定了一个x’，对吧，确定不出第二个），而w’是常量，它是一个低维空间里的常量w经过变换得到的，所以给了一个w 和x的值，就有一个确定的f(x’)值与其对应。这让我们幻想，是否能有这样一种函数K(w,x),他接受低维空间的输入值，却能算出高维空间的内积值<w’,x’>？

如果有这样的函数，那么当给了一个低维空间的输入x以后，

g(x)=K(w,x)+b

f(x’)=<w’,x’>+b

这两个函数的计算结果就完全一样，我们也就用不着费力找那个映射关系，直接拿低维的输入往g(x)里面代就可以了（再次提醒，这回的g(x)就不是线性函数啦，因为你不能保证K(w,x)这个表达式里的x次数不高于1哦）。

万幸的是，这样的K(w,x)确实存在（发现凡是我们人类能解决的问题，大都是巧得不能再巧，特殊得不能再特殊的问题，总是恰好有些能投机取巧的地方才能解决，由此感到人类的渺小），它被称作核函数（核，kernel），而且还不止一个，事实上，只要是满足了Mercer条件的函数，都可以作为核函数。核函数的基本作用就是接受两个低维空间里的向量，能够计算出经过某个变换后在高维空间里的向量内积值。几个比较常用的核函数，俄，教课书里都列过，我就不敲了（懒！）。

回想我们上节说的求一个线性分类器，它的形式应该是：

现在这个就是高维空间里的线性函数（为了区别低维和高维空间里的函数和向量，我改了函数的名字，并且给w和x都加上了 ’），我们就可以用一个低维空间里的函数（再一次的，这个低维空间里的函数就不再是线性的啦）来代替，

又发现什么了？f(x’) 和g(x)里的α，y，b全都是一样一样的！这就是说，尽管给的问题是线性不可分的，但是我们就硬当它是线性问题来求解，只不过求解过程中，凡是要求内积的时候就用你选定的核函数来算。这样求出来的α再和你选定的核函数一组合，就得到分类器啦！

明白了以上这些，会自然的问接下来两个问题：

1． 既然有很多的核函数，针对具体问题该怎么选择？

2． 如果使用核函数向高维空间映射后，问题仍然是线性不可分的，那怎么办？

第一个问题现在就可以回答你：对核函数的选择，现在还缺乏指导原则！各种实验的观察结果（不光是文本分类）的确表明，某些问题用某些核函数效果很好，用另一些就很差，但是一般来讲，径向基核函数是不会出太大偏差的一种，首选。（我做文本分类系统的时候，使用径向基核函数，没有参数调优的情况下，绝大部分类别的准确和召回都在85%以上。

感性理解，映射图：

常用的两个Kernel函数：

多项式核函数：

高斯核函数：

定义：

将核函数带入，问题又转化为线性问题啦，如下：

求

，其中

式子是有了，但是如何求结果呢？接下来一步一步的解决这个问题，并且通过动手编程使大家对这个有个问题有个直观的认识。（PS：大家都对LIBSVM太依赖了，这样无助于深入的研究与理解，而且我觉得自己动手实现的话会比较有成就感）

一．SMO算法的原理

SMO算法和以往的一些SVM改进算法一样，是把整个二次规划问题分解为很多较易处理的小问题，所不同的是，只有SMO算法把问题分解到可能达到的最小规模：每次优化只处理两个样本的优化问题，并且用解析的方法进行处理。我们将会看到，这种与众不同的方法带来了一系列不可比拟的优势。

对SVM来说，一次至少要同时对两个样本进行优化（就是优化它们对应的Lagrange乘子），这是因为等式约束的存在使得我们不可能单独优化一个变量。

所谓“最小优化”的最大好处就是使得我们可以用解析的方法求解每一个最小规模的优化问题，从而完全避免了迭代算法。

当然，这样一次“最小优化”不可能保证其结果就是所优化的Lagrange乘子的最终结果，但会使目标函数向极小值迈进一步。我们再对其它Lagrange乘子做最小优化，直到所有乘子都符合KKT条件时，目标函数达到最小，算法结束。

这样，SMO算法要解决两个问题：一是怎样解决两个变量的优化问题，二是怎样决定先对哪些Lagrange乘子进行优化。

二．两个Lagrange乘子的优化问题（子程序takeStep）

我们在这里不妨设正在优化的两个Lagrange乘子对应的样本正是第一个和第二个，对两个Lagrange乘子α1和α2，在不改变其他乘子的情况下，它们的约束条件应表达为正方形内的一条线段。如图所示：

在这条线段上求一个函数的极值，相当于一个一维的极值问题。我们可以把α1用α2表示，对α2求无条件极值，如果目标函数是严格上凹的，最小值就一定在这一极值点（极值点在区间内）或在区间端点（极值点在区间外）。α2确定后，α1也就确定下来了。因此我们先找到α2优化区间的上下限制，再在这个区间中对α2求最小值。

由图1我们容易得到α2的上下限应为：

L=max(0,α2－α1)，H=min(C,C+α2 –α1) , 若y1与y2异号；

L=max(0,α2 +α1－C), H=min(C, α2 +α1) ,若y1与y2同号；

令s=y1y2标志这两个样本是否同类，则有

L=max(0, α2+sα1－ 1/2 (s+1)C), H=min(C, α2 +sα1 –1/2 (s－1)C)

而α1和α2在本次优化中所服从的等式约束为：

α1＋sα2=α01+sα02=d

下面我们推导求最小值点α2的公式：由于只有α1，α2两个变量需要考虑，目标函数可以写成

Wolfe(α1,α2)=1/2 K11α21+1/2 K22α22+ sK12α1α2 + y1α1v1 +y2α2v2 －α1 －α2+常数

其中Kij=K(xi,xj) , vi=y3α03Ki3+…+ylα0lKil = ui+b0－ y1α01K1i – y2α01K2i

上标为0的量表示是本次优化之前Lagrange乘子的原值。

将α2用α1表示并代入目标函数：

Wolfe（α2）=1/2 K11(d-sα2)2+1/2 K22α22+sK12(d-sα2) α2

+y1(d-sα2)v1 – d+sα2+y2α2v2－α2+常数

对α2求导：

dWolfe(α2)/dα2

=－sK11(d－sα2)+K22α2－K12α2+sK12(d－sα2)－y2v2+s+y2v2－1 ＝0

如果Wolfe函数总是严格上凹的，即二阶导数K11+K22－2K12>0, 那么驻点必为极小值点，无条件的极值点就为

α2＝[s(K11－K12)d+y2(v1－v2)＋1－s] / (K11＋K22－2K12)

将d,v与α0,u之间关系代入，就得到用上一步的α02，u1,u2表示的α2的无条件最优点：

α2＝[α02(K11＋K22－2K12) +y2(u1－u2＋y2－y1)] / (K11＋K22－2K12)

令η=K11＋K22－2K12为目标函数的二阶导数，Ei=ui－yi为第i个训练样本的“误差”，这个式子又可以写为

α2＝α02＋y2(E1－E2)/η

除非核函数K不满足Mercer条件（也就是说不能作为核函数），η不会出现负值。但η=0是可以出现的情况。这时我们计算目标函数在线段两个端点上的取值，并将Lagrange乘子修正到目标函数较小的端点上：

f1=y1(E1+b)－α1K(x1,x1)­－sα2K(x1,x1)

f2=y2(E2+b)－sα1K(x1,x2)­－α2K(x2,x2)

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

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

WolfeL=L1f1+Lf2+1/2 L21K(x1,x1)+1/2 L2K(x2,x2)+sLL1K(x1,x2)

WolfeH=H1f1+Hf2+1/2 H21K(x1,x1)+1/2 H2K(x2,x2)+sHH1K(x1,x2)

当两个端点上取得相同的目标函数值时，目标函数在整条线段上的取值都会是一样的（因为它是上凹的），这时不必对α1，α2作出修正。

α2的无条件极值确定后，再考虑上下限的限制，最终的α2为

最后，由等式约束确定α1：

α1*=α1+s(α2－α2*)

三．选择待优化Lagrange乘子的试探找点法

事实上即使我们不采用任何找点法，只是按顺序抽取αi，αj的所有组合进行优化，目标函数也会不断下降，直到任一对αi，αj都不能继续优化，目标函数就会收敛到极小值。我们采取某种找点方法只是为了使算法收敛得更快。

这种试探法先选择最有可能需要优化的α2，再针对这样的α2选择最有可能取得较大修正步长的α1。这样，我们在程序中使用两个层次的循环：

内层循环（子程序examineExample）针对违反KKT条件的样本选择另一个样本与它配对优化（指优化它们的Lagrange乘子），选择的依据是尽量使这样一对样本能取得最大优化步长。对其中一个Lagrange乘子α2来说优化步长为|(E1－E2)/η|,但由于核函数估算耗时较大，我们只用|E1－E2|来大致估计有可能取得的步长大小。也就是说，选出使得|E1－E2|最大的样本作为第二个样本。需要注意的是，这样的步长估计是比较粗略的，选择出来的一对样本有时非但不能“一劳永逸”地“一步到位”，反而不能作出进一步调整，（例如η=0的情况，最小优化问题的二次型只是半正定的）。这时我们遍历所有非边界样本（非边界样本就是Lagrange乘子不在边界0或C上的样本），继续寻找能与α2配对优化的α1，如果这样的样本在非边界样本中找不到，再遍历所有样本。这两次遍历都是从随机位置开始的，以免算法总是在一开始遍历就向固定的方向偏差。在极端退化的情形，找不到与α2配对能作出进一步调整的α1，这时我们放弃第一个样本。

外层循环（主程序smo）遍历非边界样本或所有样本：优先选择遍历非边界样本，因为非边界样本更有可能需要调整，而边界样本常常不能得到进一步调整而留在边界上（可以想象大部分样本都很明显不可能是支持向量，它们的Lagrange乘子一旦取得零值就无需再调整）。循环遍历非边界样本并选出它们当中违反KKT条件的样本进行调整，直到非边界样本全部满足KKT条件为止。当某一次遍历发现没有非边界样本得到调整时，就遍历所有样本，以检验是否整个集合也都满足KKT条件。如果在整个集合的检验中又有样本被进一步优化，就有必要再遍历非边界样本。这样，外层循环不停地在“遍历所有样本”和“遍历非边界样本”之间切换，直到整个训练集都满足KKT条件为止。

以上用KKT条件对样本所作检验都是达到一定精度就可以了，例如正侧的非边界样本的输出ui可以在1的一定公差范围之内，通常这一公差（tolerance）取0.001，如果要求十分精确的输出算法就不能很快收敛。

四．每次最小优化后的重置工作

每做完一次最小优化，必须更新每个样本的误差（Error Cache），以便用修正过的分类面对其它样本再做KKT检验，以及选择第二个配对优化样本时估计步长之用。

更新Error Cache首先要重置阈值b 。我们可直接利用刚刚被优化的两个样本的信息在原阈值b0基础上作简单修正，而不需要调用所有支持向量重新计算b 。最小优化后的α1*如果不在边界上，b的计算公式为：

b1=E1+y1(α1*－α10)K(x1,x1)+y2(α2*－α20)K(x1,x2)+b0

最小优化后的α2*如果不在边界上，b的计算公式为：

b2=E2+y1(α1*－α10)K(x1,x2)+y2(α2*－α20)K(x2,x2)+b0

α1*，α2*都不在边界上时，b1和b2是相等的。两个Lagrange乘子都在边界上时，b1和b2以及它们之间的数都可作为符合KKT条件的阈值。这时SMO算法选择最安全的b1 ，b2之中点作为阈值。

非线性的情况，误差的计算要用到所有已找到的支持向量及它们的Lagrange乘子：

线性的情况则是先重置分类超平面的法向量w，再根据uj=w’xj-b计算输出uj和误差Ej=uj－yj 。同阈值的重置一样，法向量的重置也不需要调用所有的支持向量，只需在原来的法向量基础上作改动：

w*=w+y1(α1*－α1)x1+y2(α2*－α2)x2

大部分重置工作都是以简单的非循环计算来完成的，这使得需要做很多次最小优化的SMO算法不必在每次优化后的重置中花费太多时间。但是我们也看到，非线性的情况误差的重置必须与所有支持向量逐个计算核函数，而且核函数的计算本身就比点积复杂，于是非线性的情况误差的重置将成为算法速度的瓶颈。

四、算法的流程图

五、Platt大神的逻辑代码