机器学习--sklearn之支持向量机(SVM)

SVM是一个非常经典的监督学习算法。下面给出SVM对于二值分类的原理及推导过程。

1、问题转化

如下图所示:
在这里插入图片描述想要找一条直线 w x + b = 0 wx+b=0 wx+b=0将图中红蓝两类样本区分开,假设将这条线向左和向右分别平移,接触到样本点停止,得到两条新的直线,设它们分别为 w x + b = c    wx+b=c\; wx+b=c w x + b = − c wx+b=-c wx+b=c
w = w c , b = b c w=\frac{w}{c},b=\frac{b}{c} w=cw,b=cb,则这两条直线转化为 w x + b = 1    wx+b=1\; wx+b=1 w x + b = − 1 wx+b=-1 wx+b=1。此时,令 w x + b = 0 wx+b=0 wx+b=0也左右同时除以 c c c,仍得到 w x + b = 0 wx+b=0 wx+b=0
SVM的目的即找到直线 w x + b = 0 wx+b=0 wx+b=0,使得 w x + b = 1    wx+b=1\; wx+b=1 w x + b = − 1 wx+b=-1 wx+b=1之间的距离 d d d最大。
如果在决策边界上任意取两个点 x a x_a xa x b x_b xb,并带入决策边界的表达式,则有:
w x a + b = 0 wx_a+b=0 wxa+b=0
w x b + b = 0 wx_b+b=0 wxb+b=0
将两式相减,可以得到:
w ( x a − x b ) = 0 w(x_a-x_b)=0 w(xaxb)=0
所以参数向量 w w w的方向必然是垂直于我们的决策边界。
x p x_p xp为红色的点, x r x_r xr为蓝色的点,则:
w x p + b = 1 wx_p+b=1 wxp+b=1
w x r + b = − 1 wx_r+b=-1 wxr+b=1
两个式子相减,则有:
w ( x p − x r ) = 2 w(x_p-x_r)=2 w(xpxr)=2
由于:
在这里插入图片描述
在这里插入图片描述所以,只要保证:
m i n 1 2 ∣ ∣ w ∣ ∣ 2 min\quad \frac{1}{2}||w||^{2} min21w2
上式中的平方是因为L2范式中有根号。
我们的两条虚线表示的超平面,是数据边缘所在的点。所以对于任意样本 i i i,我们可以把决策函数写作:
在这里插入图片描述故有约束条件:
( w x i + b ) y i ≥ 1 (wx_{i}+b)y_{i}\geq1 (wxi+b)yi1
到这里,此问题已被转化为一有约束的优化问题:
{ m i n    1 2 ∣ ∣ w ∣ ∣ 2 ( w x i + b ) y i ≥ 1 \begin{cases} min\; \frac{1}{2}||w||^{2} \\ (wx_{i}+b)y_{i}\geq1 \end{cases} {min21w2(wxi+b)yi1
最小化函数是二次的,并且约束条件在参数w和b下是线性的,这样的最优化问题被称为“凸优化问题”。拉格朗日乘数法正好可以用来解决凸优化问题,这种方法也是业界常用的,用来解决带约束条件,尤其是带有不等式的约束条件的函数的数学方法。
关于凸优化理论请参考凸函数优化

2、拉格朗日对偶

此处引入拉格朗日对偶的概念。

1) 原始问题

假设 f ( x ) , c i ( x ) , h j ( x ) f_{(x)}, c_{i(x)}, h_{j(x)} f(x),ci(x),hj(x)是定义在 R n R^{n} Rn上的连续可微函数,考虑约束最优化问题:
m i n x ∈ R n f ( x ) \mathop{min}\limits_{x\in R^{n}}f_{(x)} xRnminf(x)
s . t . c i ( x ) ≤ 0    , i = 1 , 2 , … , k ; h j ( x ) = 0    , i = 1 , 2 , … , l s.t.\quad c_{i(x)}\leq0\;,i=1,2,…,k;\quad h_{j(x)}=0\;,i=1,2,…,l s.t.ci(x)0,i=1,2,,k;hj(x)=0,i=1,2,,l
称此约束最优化问题为原始最优化问题或原始问题。
首先,引进广义拉格朗日函数:
L ( x , α , β ) = f ( x ) + ∑ i = 1 k α i c i ( x ) + ∑ j = 1 l β j h j ( x ) L_{(x,\alpha,\beta)}=f_{(x)}+\sum_{i=1}^{k}\alpha_{i}c_{i(x)}+\sum_{j=1}^{l}\beta_{j}h_{j(x)} L(x,α,β)=f(x)+i=1kαici(x)+j=1lβjhj(x)
这里, x = ( x ( 1 ) , x ( 2 ) , … , x ( n ) ) T ∈ R n x=(x^{(1)},x^{(2)},…,x^{(n)})^{T}\in R^{n} x=(x(1),x(2),,x(n))TRn α i ,    β j \alpha_{i},\; \beta_{j} αi,βj是拉格朗日乘子, α i ≥ 0 \alpha_{i}\geq0 αi0。考虑 x x x的函数:
θ p ( x ) = m a x α ,    β ,    α i ≥ 0 L ( x , α , β ) \theta_{p(x)}=\mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}L_{(x,\alpha,\beta)} θp(x)=α,β,αi0maxL(x,α,β)
假设给定某个 x x x,如果 x x x违反原始问题的约束条件,即存在某个 i i i使得 c i ( x ) > 0    c_{i(x)}>0\; ci(x)>0或者存在某个 j j j使得 h j ( x ) ≠ 0 h_{j(x)}\ne0 hj(x)=0,那么就有:
θ p ( x ) = m a x α ,    β ,    α i ≥ 0 [ f ( x ) + ∑ i = 1 k α i c i ( x ) + ∑ j = 1 l β j h j ( x ) ] = + ∞ \theta_{p(x)}=\mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}[f_{(x)}+\sum_{i=1}^{k}\alpha_{i}c_{i(x)}+\sum_{j=1}^{l}\beta_{j}h_{j(x)}]=+\infty θp(x)=α,β,αi0max[f(x)+i=1kαici(x)+j=1lβjhj(x)]=+
因为若某个 i i i使得 c i ( x ) > 0 c_{i(x)}>0 ci(x)>0,则可以令 α i → + ∞ \alpha_{i}\rightarrow+\infty αi+,若某个 j j j使 h j ( x ) ≠ 0 h_{j(x)}\ne0 hj(x)=0,则可令 β j → + ∞ \beta_{j}\rightarrow+\infty βj+。以上两种情况均可以使    θ p ( x ) = + ∞ \;\theta_{p(x)}=+\infty θp(x)=+
相反地,若 x x x满足约束条件,则:
θ p ( x ) = m a x α ,    β ,    α i ≥ 0 [ f ( x ) + ∑ i = 1 k α i c i ( x ) ] \theta_{p(x)}=\mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}[f_{(x)}+\sum_{i=1}^{k}\alpha_{i}c_{i(x)}] θp(x)=α,β,αi0max[f(x)+i=1kαici(x)]
由于 α i ≥ 0 \alpha_{i}\geq0 αi0 c i ( x ) ≤ 0 c_{i(x)}\leq0 ci(x)0,故 ∑ i = 1 k α i c i ( x ) ≤ 0 \mathop{\sum}\limits_{i=1}^{k}\alpha_{i}c_{i(x)}\leq0 i=1kαici(x)0。所以:
θ p ( x ) = f ( x ) \theta_{p(x)}=f_{(x)} θp(x)=f(x)
即:
θ p ( x ) = { f ( x ) , x    满足原问题约束 + ∞ , 其他 \theta_{p(x)}=\begin{cases} f_{(x)}, & \text {x\;满足原问题约束} \\ +\infty, & \text {其他} \end{cases} θp(x)={f(x),+,x满足原问题约束其他
如果考虑极小化问题
m i n x θ p ( x ) = m i n x m a x α ,    β ,    α i ≥ 0 L ( x , α , β ) \mathop{min}\limits_{x}\theta_{p(x)}=\mathop{min}\limits_{x} \mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}L_{(x,\alpha,\beta)} xminθp(x)=xminα,β,αi0maxL(x,α,β)
θ p ( x ) \theta_{p(x)} θp(x)的表达式可知, m i n x θ p ( x ) \mathop{min}\limits_{x}\theta_{p(x)} xminθp(x)是与原始最优化问题等价的,也就是说它们的解相同。问题 m i n x m a x α ,    β ,    α i ≥ 0 L ( x , α , β ) \mathop{min}\limits_{x} \mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}L_{(x,\alpha,\beta)} xminα,β,αi0maxL(x,α,β)被称为广义拉格朗日函数的极小极大问题。

2)对偶问题

定义
θ D ( α , β ) = m i n x L ( x , α , β ) \theta_{D(\alpha,\beta)}=\mathop{min}\limits_{x}L_{(x,\alpha,\beta)} θD(α,β)=xminL(x,α,β)
再考虑极大化 θ D ( α , β ) \theta_{D(\alpha,\beta)} θD(α,β),即:
m a x α ,    β ,    α i ≥ 0 θ D ( α , β ) = m a x α ,    β ,    α i ≥ 0 m i n x L ( x , α , β ) \mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}\theta_{D(\alpha,\beta)}=\mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}\mathop{min}\limits_{x}L_{(x,\alpha,\beta)} α,β,αi0maxθD(α,β)=α,β,αi0maxxminL(x,α,β)
问题 m a x α ,    β ,    α i ≥ 0 m i n x L ( x , α , β )    \mathop{max}\limits_{\alpha,\; \beta,\; \alpha_{i}\geq0}\mathop{min}\limits_{x}L_{(x,\alpha,\beta)}\; α,β,αi0maxxminL(x,α,β)称为广义拉格朗日函数的极大极小问题。
若目标函数与所有约束函数皆为凸函数,则原始问题和它的对偶问题是等价的。
而对于 S V M SVM SVM,其所有函数皆为凸函数。

3、SVM的对偶问题

将函数从最初形态转换为拉格朗日乘数形态

{ m i n    1 2 w 2 ( w x i + b ) y i ≥ 1 \begin{cases} min\; \frac{1}{2}w^{2} \\ (wx_{i}+b)y_{i}\geq1 \end{cases} {min21w2(wxi+b)yi1
其广义拉格朗日函数为:
L ( w , b , α ) = 1 2 w 2 + ∑ i = 1 N α i [ 1 − ( w x i + b ) y i ] L_{(w,b,\alpha)}=\frac{1}{2}w^{2}+\sum_{i=1}^{N}\alpha_{i}[1-(wx_{i}+b)y_{i}] L(w,b,α)=21w2+i=1Nαi[1(wxi+b)yi]
所以损失函数能从最初形态转换为拉格朗日乘数形态:
m i n w ,    b m a x α i ≥ 0 L ( w , b , α ) \mathop{min}\limits_{w,\;b} \mathop{max}\limits_{\alpha_{i}\geq0}L_{(w,b,\alpha)} w,bminαi0maxL(w,b,α)
怎么理解这个式子呢?首先,我们第一步先执行 m a x max max,即最大化 L ( w , b , α ) L_{(w,b,\alpha)} L(w,b,α),那就有两种情况:

  • ( w x i + b ) y i > 1 (wx_{i}+b)y_{i}>1 (wxi+b)yi>1,函数的第二部分 ∑ i = 1 N α i [ 1 − ( w x i + b ) y i ] \sum_{i=1}^{N}\alpha_{i}[1-(wx_{i}+b)y_{i}] i=1Nαi[1(wxi+b)yi]就一定为负,式子 1 2 w 2 \frac{1}{2}w^{2} 21w2就要减去一个正数,此时若要最大化 L ( w , b , α ) L_{(w,b,\alpha)} L(w,b,α),则 α \alpha α必须取到0。
  • ( w x i + b ) y i < 1 (wx_{i}+b)y_{i}<1 (wxi+b)yi<1,函数的第二部分 ∑ i = 1 N α i [ 1 − ( w x i + b ) y i ] \sum_{i=1}^{N}\alpha_{i}[1-(wx_{i}+b)y_{i}] i=1Nαi[1(wxi+b)yi]就一定为正,式子 1 2 w 2 \frac{1}{2}w^{2} 21w2就要加上一个正数,此时若要最大化 L ( w , b , α ) L_{(w,b,\alpha)} L(w,b,α),则 α \alpha α必须取到正无穷。

若把函数第二部分当作一个惩罚项来看待,则 ( w x i + b ) y i (wx_{i}+b)y_{i} (wxi+b)yi大于1时函数没有受到惩罚,而 ( w x i + b ) y i (wx_{i}+b)y_{i} (wxi+b)yi小于1 时函数受到了极致的惩罚,即加上了一个正无穷项,函数整体永远不可能取到最小值。所以第二步,我们执行 m i n min min的命令,求解函数整体的最小值,我们就永远不能让 α \alpha α必须取到正无穷的状况出现,即是说永远不让 ( w x i + b ) y i < 1 (wx_{i}+b)y_{i}<1 (wxi+b)yi<1的状况出现,从而实现了求解最小值的同时让约束条件被满足。
现在, L ( w , b , α ) L_{(w,b,\alpha)} L(w,b,α)就是我们新的损失函数了,我们的目标是要通过先最大化,再最小化它来求解参数向量 w w w和截距 b b b的值。
接下来,我们需要将拉格朗日函数转换为拉格朗日对偶函数。

将拉格朗日函数转换为拉格朗日对偶函数

为什么要转换?

要求极值,最简单的方法还是对参数求导后让一阶导数等于0。我们先来试试看对拉格朗日函数求极值,在这里我 们对参数向量 w w w和截距 b b b分别求偏导并且让他们等于0。这个求导过程比较简单:
在这里插入图片描述由于两个求偏导结果中都带有未知的拉格朗日乘数 α i \alpha_i αi,因此我们还是无法求解出 w w w b b b,我们必须想出一种方法来求解拉格朗日乘数 α i \alpha_i αi。幸运地是,拉格朗日函数可以被转换成一种只带有 α i \alpha_i αi,而不带有 w w w b b b的形式,这种形式被称为拉格朗日对偶函数。在对偶函数下,我们就可以求解出拉格朗日乘数 α i \alpha_i αi,然后带入到上面推导出的(1)和(2)式中来求解 w w w b b b

为什么能够进行转换?

对于任何一个拉格朗日函数 L ( x , α ) = f ( x ) + ∑ i = 1 q α i h i ( x ) L(x,\alpha)=f(x)+\sum_{i=1}^{q}\alpha_ih_{i}(x) L(x,α)=f(x)+i=1qαihi(x),都存在一个与它对应的对偶函数 g ( α ) g(\alpha) g(α),只带有拉格朗日乘数 α \alpha α作为唯一的参数。如果 L ( x , α ) L(x,\alpha) L(x,α)的最优解存在并可以表示为 m i n x L ( x , α ) \mathop{min}\limits_{x}L(x,\alpha) xminL(x,α),并且对偶函数的最优解也存在并可以表示为 m i n α g ( α ) \mathop{min}\limits_{\alpha}g(\alpha) αming(α),则我们可以定义对偶差异,即拉格朗日函数的优解与其对偶函数的优解之间的差值:在这里插入图片描述如果 Δ = 0 \Delta=0 Δ=0,则称 L ( x , α ) L(x,\alpha) L(x,α)与其对偶函数之间存在强对偶关系,此时我们就可以通过求解其对偶函数的优解来替代求解原始函数的优解。那强对偶关系什么时候存在呢?则这个拉格朗日函数必须满足KKT条件:
在这里插入图片描述之前我们已经让拉格朗日函数上对参数 w w w b b b的求导为0,得到了式子:
在这里插入图片描述并且在我们的函数中,我们通过先求解最大值再求解最小值的方法使得函数天然满足:
在这里插入图片描述所以接下来,我们只需要再满足一个条件:
在这里插入图片描述这个条件其实很容易满足,能够让 ( w x i + b ) y i − 1 = 0 (wx_{i}+b)y_{i}-1=0 (wxi+b)yi1=0的就是落在虚线的超平面上的样本点,即我们的支持向量。所有不是支持向量的样本点则必须满足 α i = 0 \alpha_i=0 αi=0。满足这个式子说明了,我们求解的参数 w w w b b b以及求解的超平面的存在,只与支持向量相关,与其他样本点都无关。现在KKT的五个条件都得到了满足,我们就可以使用 L ( w , b , α ) L_{(w,b,\alpha)} L(w,b,α)的对偶函数来求解 α \alpha α了。

进行转换

原问题的对偶问题:
m a x α i ≥ 0 m i n w ,    b L ( w , b , α ) \mathop{max}\limits_{\alpha_{i}\geq0}\mathop{min}\limits_{w,\;b}L_{(w,b,\alpha)} αi0maxw,bminL(w,b,α)
先对 L ( w , b , α ) L_{(w,b,\alpha)} L(w,b,α)求极小值:
∂ L ∂ w = 0 ⇒ w = ∑ i = 1 N α i x i y i \frac{\partial L}{\partial w}=0\Rightarrow w=\sum_{i=1}^{N}\alpha_{i}x_{i}y_{i} wL=0w=i=1Nαixiyi
∂ L ∂ b = 0 ⇒ ∑ i = 1 N α i y i = 0 \frac{\partial L}{\partial b}=0\Rightarrow \sum_{i=1}^{N}\alpha_{i}y_{i}=0 bL=0i=1Nαiyi=0
将结果带入对偶问题,可将其转化为下式:
m a x α i ≥ 0    1 2 ( ∑ i = 1 N α i x i y i ) 2 + ∑ i = 1 N α i [ 1 − ( ( ∑ j = 1 N α j x j y j ) x i + b ) y i ] \mathop{max}\limits_{\alpha_{i}\geq0}\;\frac{1}{2}(\sum_{i=1}^{N}\alpha_{i}x_{i}y_{i})^{2}+\sum_{i=1}^{N}\alpha_{i}[1-((\sum_{j=1}^{N}\alpha_{j}x_{j}y_{j})x_{i}+b)y_{i}] αi0max21(i=1Nαixiyi)2+i=1Nαi[1((j=1Nαjxjyj)xi+b)yi]
= m a x α i ≥ 0    1 2 ( ∑ i = 1 N α i x i y i ) 2 + ∑ i = 1 N α i [ 1 − ( ∑ j = 1 N α j x j y j ) x i y i ] =\mathop{max}\limits_{\alpha_{i}\geq0}\;\frac{1}{2}(\sum_{i=1}^{N}\alpha_{i}x_{i}y_{i})^{2}+\sum_{i=1}^{N}\alpha_{i}[1-(\sum_{j=1}^{N}\alpha_{j}x_{j}y_{j})x_{i}y_{i}] =αi0max21(i=1Nαixiyi)2+i=1Nαi[1(j=1Nαjxjyj)xiyi]
= m a x α i ≥ 0    1 2 ( ∑ i = 1 N α i x i y i ) 2 + ∑ i = 1 N α i − ∑ i = 1 N ∑ j = 1 N α i x i y i α j x j y j =\mathop{max}\limits_{\alpha_{i}\geq0}\;\frac{1}{2}(\sum_{i=1}^{N}\alpha_{i}x_{i}y_{i})^{2}+\sum_{i=1}^{N}\alpha_{i}-\sum_{i=1}^{N}\sum_{j=1}^{N}\alpha_{i}x_{i}y_{i}\alpha_{j}x_{j}y_{j} =αi0max21(i=1Nαixiyi)2+i=1Nαii=1Nj=1Nαixiyiαjxjyj
= m a x α i ≥ 0    − 1 2 ∑ i = 1 N ∑ j = 1 N α i α j y i y j < x i x j > + ∑ i = 1 N α i =\mathop{max}\limits_{\alpha_{i}\geq0}\;-\frac{1}{2}\sum_{i=1}^{N}\sum_{j=1}^{N}\alpha_{i}\alpha_{j}y_{i}y_{j}<x_{i}x_{j}>+\sum_{i=1}^{N}\alpha_{i} =αi0max21i=1Nj=1Nαiαjyiyj<xixj>+i=1Nαi
约束条件为:
∑ i = 1 N α i y i = 0 α i ≥ 0 \sum_{i=1}^{N}\alpha_{i}y_{i}=0\quad \alpha_{i}\geq0 i=1Nαiyi=0αi0
K K T KKT KKT条件为:
α i [ 1 − ( w x i + b ) y i ] = 0 \alpha_{i}[1-(wx_{i}+b)y_{i}]=0 αi[1(wxi+b)yi]=0
故原始问题被转化为:
{ m a x α i ≥ 0    − 1 2 ∑ i = 1 N ∑ j = 1 N α i α j y i y j < x i x j > + ∑ i = 1 N α i ∑ i = 1 N α i y i = 0 α i ≥ 0 α i [ 1 − ( w x i + b ) y i = 0 \begin{cases} \mathop{max}\limits_{\alpha_{i}\geq0}\;-\frac{1}{2}\mathop{\sum}\limits_{i=1}^{N}\mathop{\sum}\limits_{j=1}^{N}\alpha_{i}\alpha_{j}y_{i}y_{j}<x_{i}x_{j}>+\mathop{\sum}\limits_{i=1}^{N}\alpha_{i} \\ \mathop{\sum}\limits_{i=1}^{N}\alpha_{i}y_{i}=0\quad \alpha_{i}\geq0 \\ \alpha_{i}[1-(wx_{i}+b)y_{i}=0 \end{cases} αi0max21i=1Nj=1Nαiαjyiyj<xixj>+i=1Nαii=1Nαiyi=0αi0αi[1(wxi+b)yi=0

4、SMO算法(序列最小优化算法)

优化此问题时可以采用SMO算法。
假设此问题中只有 α 1 \alpha_{1} α1 α 10 \alpha_{10} α10,则SMO算法的步骤如下。
1)先固定 α 3 \alpha_{3} α3 α 10 \alpha_{10} α10,则最优化函数变为:
m a x α i ≥ 0    − 1 2 ∑ i = 1 2 ∑ j = 1 2 α i α j y i y j < x i x j > + ∑ i = 1 2 α i \mathop{max}\limits_{\alpha_{i}\geq0}\;-\frac{1}{2}\mathop{\sum}\limits_{i=1}^{2}\mathop{\sum}\limits_{j=1}^{2}\alpha_{i}\alpha_{j}y_{i}y_{j}<x_{i}x_{j}>+\mathop{\sum}\limits_{i=1}^{2}\alpha_{i} αi0max21i=12j=12αiαjyiyj<xixj>+i=12αi
约束条件变为:
α 1 y 1 + α 2 y 2 = − ∑ i = 3 10 α i y i    为一常数 \alpha_{1}y_{1}+\alpha_{2}y_{2}=-\sum_{i=3}^{10}\alpha_{i}y_{i}\text{\;为一常数} α1y1+α2y2=i=310αiyi为一常数
则此优化问题变成了一个二维变量( α 1 , α 2 \alpha_{1},\alpha_{2} α1α2)的优化问题,使原来的十维优化问题大大简化。
这时,可以先算出最优的 α 1 \alpha_{1} α1 α 2 \alpha_{2} α2
2)固定 α 1 \alpha_{1} α1 α 4 \alpha_{4} α4 α 10 \alpha_{10} α10,重复以上操作并代入最优的 α 1 \alpha_{1} α1,得到最优的 α 2 \alpha_{2} α2 α 3 \alpha_{3} α3
3)固定 α 1 \alpha_{1} α1 α 2 \alpha_{2} α2 α 5 \alpha_{5} α5 α 10 \alpha_{10} α10,重复以上操作并代入最优的 α 1 \alpha_{1} α1 α 2 \alpha_{2} α2,得到最优的 α 3 \alpha_{3} α3 α 4 \alpha_{4} α4
……
4)得到所有的最优 α \alpha α,将其代入 w = ∑ i = 1 10 α i x i y i w=\mathop{\sum}\limits_{i=1}^{10}\alpha_{i}x_{i}y_{i} w=i=110αixiyi,得到最优的 w w w值。
5)若得到的 α i = 0 \alpha_{i}=0 αi=0,则表示第 i i i个约束条件不起作用;若得到的 α i ≠ 0 \alpha_{i}\ne0 αi=0,则它对应的 x i x_{i} xi y i y_{i} yi就是支持向量。
6)通过 y = ( ∑ i = 1 N α i x i y i ) T x = ∑ i = 1 N α i y i x i T x    y=(\mathop{\sum}\limits_{i=1}^{N}\alpha_{i}x_{i}y_{i})^{T}x=\mathop{\sum}\limits_{i=1}^{N}\alpha_{i}y_{i}x_{i}^{T}x\; y=(i=1Nαixiyi)Tx=i=1NαiyixiTx得到某样本点对应的标签。

5、SVM实现(sklearn)

线性

1、导入需要的库

from sklearn.datasets import make_blobs
from sklearn.svm import SVC
import matplotlib.pyplot as plt
import numpy as np

2、实例化数据集,可视化数据集

X,y = make_blobs(n_samples=50, centers=2, random_state=0,cluster_std=0.6)
plt.scatter(X[:,0],X[:,1],c=y,s=50,cmap="rainbow")#rainbow彩虹色
plt.xticks([])
plt.yticks([])
plt.show()

可见得到的数据集分布为:
在这里插入图片描述

3. 画决策边界

获取样本构成的平面,作为一个对象
#首先要有散点图
plt.scatter(X[:,0],X[:,1],c=y,s=50,cmap="rainbow")
ax = plt.gca() #获取当前的子图,如果不存在,则创建新的子图
用函数meshgrid制作网格
#获取平面上两条坐标轴的最大值和最小值
xlim = ax.get_xlim()
ylim = ax.get_ylim()
#在最大值和最小值之间形成30个规律的数据
axisx = np.linspace(xlim[0],xlim[1],30)
axisy = np.linspace(ylim[0],ylim[1],30)
 
axisy,axisx = np.meshgrid(axisy,axisx)
#我们将使用这里形成的二维数组作为我们contour函数中的X和Y
#使用meshgrid函数将两个一维向量转换为特征矩阵
#核心是将两个特征向量广播,以便获取y.shape * x.shape这么多个坐标点的横坐标和纵坐标
 
xy = np.vstack([axisx.ravel(), axisy.ravel()]).T
#其中ravel()是降维函数,vstack能够将多个结构一致的一维数组按行堆叠起来
#xy就是已经形成的网格,它是遍布在整个画布上的密集的点
plt.scatter(xy[:,0],xy[:,1],s=1,cmap="rainbow")

得到图像为:
在这里插入图片描述
在以上代码中,meshgrid和vstack的作用可以用下面这段代码说明:

#理解函数meshgrid和vstack的作用
a = np.array([1,2,3])
b = np.array([7,8])
#两两组合,会得到多少个坐标?
#答案是6个,分别是 (1,7),(2,7),(3,7),(1,8),(2,8),(3,8)
v1,v2 = np.meshgrid(a,b)
print('v1: ', v1)
print('v2: ', v2)
v = np.vstack([v1.ravel(), v2.ravel()]).T
print('v: ', v)

最终得到的v1,v2和v分别为:

v1:  [[1 2 3]
 [1 2 3]]
v2:  [[7 7 7]
 [8 8 8]]
v:  [[1 7]
 [2 7]
 [3 7]
 [1 8]
 [2 8]
 [3 8]]
建模,通过fit计算出对应的决策边界
clf = SVC(kernel = "linear").fit(X,y)#计算出对应的决策边界
Z = clf.decision_function(xy).reshape(axisx.shape)
print(Z.shape)
#重要接口decision_function,返回每个输入的样本所对应的到决策边界的距离
#然后再将这个距离转换为axisx的结构,这是由于画图的函数contour要求Z的结构必须与X和Y保持一致

所以Z的shape为(30, 30)。

用contour函数绘图
#首先要有散点图
plt.scatter(X[:,0],X[:,1],c=y,s=50,cmap="rainbow")
ax = plt.gca() #获取当前的子图,如果不存在,则创建新的子图
#画决策边界和平行于决策边界的超平面
ax.contour(axisx,axisy,Z
           ,colors="k"
           ,levels=[-1,0,1] #画三条等高线,分别是Z为-1,Z为0和Z为1的三条线
           ,alpha=0.5#透明度
           ,linestyles=["--","-","--"])
 
ax.set_xlim(xlim)#设置x轴取值
ax.set_ylim(ylim)

得到图像:
在这里插入图片描述

可以将整个绘图过程封装为一个函数
#将上述过程包装成函数:
def plot_svc_decision_function(model,ax=None):
    if ax is None:
        ax = plt.gca()
    xlim = ax.get_xlim()
    ylim = ax.get_ylim()
    
    x = np.linspace(xlim[0],xlim[1],30)
    y = np.linspace(ylim[0],ylim[1],30)
    Y,X = np.meshgrid(y,x) 
    xy = np.vstack([X.ravel(), Y.ravel()]).T
    P = model.decision_function(xy).reshape(X.shape)
    
    ax.contour(X, Y, P,colors="k",levels=[-1,0,1],alpha=0.5,linestyles=["--","-","--"]) 
    ax.set_xlim(xlim)
    ax.set_ylim(ylim)

4. 查看相关信息

clf.predict(X)
#根据决策边界,对X中的样本进行分类,返回的结构为n_samples
clf.score(X,y)
#返回给定测试数据和标签的平均准确度
clf.support_vectors_
#返回支持向量坐标
clf.n_support_#array([2, 1])
#返回每个类中支持向量的个数

环形

1、实例化数据集,可视化数据集

from sklearn.datasets import make_circles
X,y = make_circles(100, factor=0.1, noise=.1)
plt.scatter(X[:,0],X[:,1],c=y,s=50,cmap="rainbow")
plt.show()

可见数据集在图中的位置为:
在这里插入图片描述

2、 用线性SVM画决策边界并打分

clf = SVC(kernel = "linear").fit(X,y)
plt.scatter(X[:,0],X[:,1],c=y,s=50,cmap="rainbow")
plot_svc_decision_function(clf)
clf.score(X,y)

得到图像:
在这里插入图片描述
准确率为0.68。
明显,现在线性SVM已经不适合于我们的状况了,我们无法找出一条直线来划分我们的数据集,让直线的两边分别是两种类别。这个时候,如果我们能够在原本的X和y的基础上,添加一个维度r,变成三维,我们可视化这个数据,来看看添加维度让我们的数据如何变化。

3、为非线性数据增加维度并绘制3D图像

#定义一个由x计算出来的新维度r
r = np.exp(-(X**2).sum(1))

from mpl_toolkits import mplot3d
#定义一个绘制三维图像的函数
#elev表示上下旋转的角度
#azim表示平行旋转的角度
def plot_3D(elev=30,azim=30,X=X,y=y,r=r):
    ax = plt.subplot(projection="3d")  # 这一步需要导入mplot3d
    ax.scatter3D(X[:,0],X[:,1],r,c=y,s=50,cmap='rainbow')
    ax.view_init(elev=elev,azim=azim)
    ax.set_xlabel("x")
    ax.set_ylabel("y")
    ax.set_zlabel("r")
    plt.show()
    
plot_3D()

得到图像:
在这里插入图片描述
可以看见,此时我们的数据明显是线性可分的了:我们可以使用一个平面来将数据完全分开,并使平面的上方的所有数据点为一类,平面下方的所有数据点为另一类。
使用ipywidgets库中的函数可以通过调节3D图的视角更清晰地了解数据(限于jupyter notebook用户)。

from ipywidgets import interact,fixed
interact(plot_3D,elev=[0,30,60,90],azip=(-180,180),X=fixed(X),y=fixed(y))
plt.show()

在这里插入图片描述

4、用非线性SVM画决策边界

clf = SVC(kernel = "rbf").fit(X,y)
plt.scatter(X[:,0],X[:,1],c=y,s=50,cmap="rainbow")
plot_svc_decision_function(clf)

得到图像:
在这里插入图片描述
此时我们的数据在三维空间中,我们的超平面就是一个二维平面。明显我们可以用一个平面将两类数据隔开,这个平面就是我们的决策边界了。我们刚才做的,计算r,并将r作为数据的第三维度来将数据升维的过程,被称为“核变换”,即是将数据投影到高维空间中,以寻找能够将数据完美分割的超平面,即是说寻找能够让数据线性可分的高维空间。为了详细解释这个过程,我们需要引入SVM中的核心概念:核函数。

6、核函数

为了能够找出非线性数据的线性决策边界,我们需要将数据从原始的空间 x x x投射到新空间 Φ ( x ) \Phi(x) Φ(x)中。 Φ ( x ) \Phi(x) Φ(x)是一个映射函数,它代表了某种非线性的变换,如同我们之前所做过的使用r来升维一样,这种非线性变换看起来是一种非常有效的方式。使用这种变换,线性SVM的原理可以被很容易推广到非线性情况下,其推导过程和逻辑都与线性SVM一模一样,只不过在定义决策边界之前,我们必须先对数据进行升维度,即将原始的 x x x转换成 Φ ( x ) \Phi(x) Φ(x)
如此,非线性SVM的损失函数的初始形态为:
在这里插入图片描述同理,非线性SVM的拉格朗日函数和拉格朗日对偶函数为:
在这里插入图片描述使用同样的推导方式,让拉格朗日函数满足KKT条件,并在拉格朗日函数上对每个参数求导,经过和线性SVM相同的变换后,就可以得到拉格朗日对偶函数。同样使用梯度下降或SMO等方式对 α \alpha α进行求解,最后可以求得决策边界,并得到最终的决策函数:
在这里插入图片描述这种变换非常巧妙,但也带有一些实现问题。 首先,我们可能不清楚应该什么样的数据应该使用什么类型的映射函数来确保可以在变换空间中找出线性决策边界。极端情况下,数据可能会被映射到无限维度的空间中,这种高维空间可能不是那么友好,维度越多,推导和计算的难度都会随之暴增。其次,即使已知适当的映射函数,我们想要计算类似于 Φ ( x i ) ⋅ Φ ( x t e s t ) \Phi(x_i) \cdot \Phi(x_{test}) Φ(xi)Φ(xtest)这样的点积,计算量可能会无比巨大,要找出超平面所付出的代价是非常昂贵的。
而解决这些问题的数学方式,叫做“核技巧”(Kernel Trick),是一种能够使用数据原始空间中的向量计算来表示升维后的空间中的点积结果的数学方式。具体表现为, K ( u , v ) = Φ ( u ) ⋅ Φ ( v ) K(u,v)=\Phi(u) \cdot \Phi(v) K(u,v)=Φ(u)Φ(v)。而这个原始空间中的点积函数,就被叫做“核函数”(Kernel Function)。
选用不同的核函数,就可以解决不同数据分布下的寻找超平面问题。在SVC中,这个功能由参数“kernel”和一系列与核函数相关的参数来进行控制。参数“kernel"在sklearn中可选以下几种:
在这里插入图片描述

7、在不同数据集上对比各个核函数效果

导入需要的库

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from sklearn import svm  # 也可以用from sklearn.svm import SVC
from sklearn.datasets import make_circles, make_moons, make_blobs,make_classification

实例化数据集,可视化数据集

n_samples = 100
 
datasets = [
    make_moons(n_samples=n_samples, noise=0.2, random_state=0),
    make_circles(n_samples=n_samples, noise=0.2, factor=0.5, random_state=1),
    make_blobs(n_samples=n_samples, centers=2, random_state=5),#分簇的数据集
    make_classification(n_samples=n_samples,n_features = 2,n_informative=2,n_redundant=0, random_state=5)
                #n_features:特征数,n_informative:带信息的特征数,n_redundant:不带信息的特征数
    ]
 
Kernel = ["linear","poly","rbf","sigmoid"]
 
#四个数据集分别是什么样子呢?
for X,Y in datasets:
    plt.figure(figsize=(5,4))
    plt.scatter(X[:,0],X[:,1],c=Y,s=50,cmap="rainbow")

所以四个数据集在图上的分布为:
在这里插入图片描述在这里插入图片描述在这里插入图片描述在这里插入图片描述

创建画布以导入结果

nrows=len(datasets)
ncols=len(Kernel) + 1  # 多一个装原图像
 
fig, axes = plt.subplots(nrows, ncols,figsize=(20,16))

SVM决策边界

#第一层循环:在不同的数据集中循环
for ds_cnt, (X,Y) in enumerate(datasets):
    
    #在图像中的第一列,放置原数据的分布
    ax = axes[ds_cnt, 0]
    if ds_cnt == 0:
        ax.set_title("Input data")
    ax.scatter(X[:, 0], X[:, 1], c=Y, zorder=10, cmap=plt.cm.Paired,edgecolors='k')
    ax.set_xticks(())
    ax.set_yticks(())
    
    #第二层循环:在不同的核函数中循环
    #从图像的第二列开始,一个个填充分类结果
    for est_idx, kernel in enumerate(Kernel):
        
        #定义子图位置
        ax = axes[ds_cnt, est_idx + 1]
        
        #建模
        clf = svm.SVC(kernel=kernel, gamma=2).fit(X, Y)
        score = clf.score(X, Y)
        
        #绘制图像本身分布的散点图
        ax.scatter(X[:, 0], X[:, 1], c=Y
                   ,zorder=10
                   ,cmap=plt.cm.Paired,edgecolors='k')
        #绘制支持向量
        ax.scatter(clf.support_vectors_[:, 0], clf.support_vectors_[:, 1], s=50,
                    facecolors='none', zorder=10, edgecolors='k')# facecolors='none':透明的
        
        #绘制决策边界
        x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
        y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
        
        #np.mgrid,合并了我们之前使用的np.linspace和np.meshgrid的用法
        #一次性使用最大值和最小值来生成网格
        #表示为[起始值:结束值:步长]
        #如果步长是复数,则其整数部分就是起始值和结束值之间创建的点的数量,并且结束值被包含在内
        XX, YY = np.mgrid[x_min:x_max:200j, y_min:y_max:200j]
        #np.c_,类似于np.vstack的功能
        Z = clf.decision_function(np.c_[XX.ravel(), YY.ravel()]).reshape(XX.shape)
        #填充等高线不同区域的颜色
        ax.pcolormesh(XX, YY, Z > 0, cmap=plt.cm.Paired)
        #绘制等高线
        ax.contour(XX, YY, Z, colors=['k', 'k', 'k'], linestyles=['--', '-', '--'],
                    levels=[-1, 0, 1])
        
        #设定坐标轴为不显示
        ax.set_xticks(())
        ax.set_yticks(())
        
        #将标题放在第一行的顶上
        if ds_cnt == 0:
            ax.set_title(kernel)
            
        #为每张图添加分类的分数   
        ax.text(0.95, 0.06, ('%.2f' % score)
                , size=15
                , bbox=dict(boxstyle='round', alpha=0.8, facecolor='white')
                    #为分数添加一个白色的格子作为底色
                , transform=ax.transAxes #确定文字所对应的坐标轴,就是ax子图的坐标轴本身
                , horizontalalignment='right' #位于坐标轴的什么方向
               )
plt.show()

得到结果为:
在这里插入图片描述可以观察到,线性核函数和多项式核函数在非线性数据上表现会浮动,如果数据相对线性可分,则表现不错,如果是像环形数据那样彻底不可分的,则表现糟糕。在线性数据集上,线性核函数和多项式核函数即便有扰动项也可以表现不错,可见多项式核函数是虽然也可以处理非线性情况,但更偏向于线性的功能。
Sigmoid核函数就比较尴尬了,它在非线性数据上强于两个线性核函数,但效果明显不如rbf,它在线性数据上完全比不上线性的核函数们,对扰动项的抵抗也比较弱,所以它功能比较弱小,很少被用到。
高斯径向基核函数基本在任何数据集上都表现不错,属于比较万能的核函数。我个人的经验是,无论如何先试试看高斯径向基核函数,它适用于核转换到很高的空间的情况,在各种情况下往往效果都很不错,如果rbf效果不好,那我们再试试看其他的核函数。另外,多项式核函数多被用于图像处理之中。

8、注意事项

8.1 数据标准化

SVM严重受到数据量纲的影响,所以在进行分类前要先对数据进行标准化处理。使用乳腺癌数据集举例:

from sklearn.datasets import load_breast_cancer
from sklearn.svm import SVC
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import numpy as np
from time import time
import datetime
import pandas as pd

data = load_breast_cancer()
X = data.data
y = data.target
Xtrain, Xtest, Ytrain, Ytest = train_test_split(X,y,test_size=0.3,random_state=420)

data = pd.DataFrame(X)
data.describe([0.01,0.05,0.1,0.25,0.5,0.75,0.9,0.99]).T  # 描述性统计
#从mean列和std列可以看出严重的量纲不统一
#从1%的数据和最小值相对比,90%的数据和最大值相对比,查看是否是正态分布或偏态分布,如果差的太多就是偏态分布,谁大方向就偏向谁
#可以发现数据大的特征存在偏态问题
#这个时候就需要对数据进行标准化

from sklearn.preprocessing import StandardScaler
X = StandardScaler().fit_transform(X)#将数据转化为0,1正态分布
data = pd.DataFrame(X)
data.describe([0.01,0.05,0.1,0.25,0.5,0.75,0.9,0.99]).T  # 均值很接近,方差为1了

8.2 选取与核函数相关的参数

在这里插入图片描述在知道如何选取核函数后,我们还要观察一下除了kernel之外的核函数相关的参数。对于线性核函数,kernel是唯一能够影响它的参数,但是对于其他三种非线性核函数,他们还受到参数gamma,degree以及coef0的影响。参数gamma就是表达式中的 γ \gamma γ,degree就是多项式核函数的次数 d d d,参数coef0就是常数项 r r r。其中,高斯径向基核函数受到gamma的影响,而多项式核函数受到全部三个参数的影响。
在这里插入图片描述一般我们直接使用学习曲线或者网格搜索来帮助我们查找最佳的参数组合。
仍用上一个section中的乳腺癌数据集举例:

score = []
gamma_range = np.logspace(-10, 1, 50) #返回在对数刻度上均匀间隔的数字
for i in gamma_range:
    clf = SVC(kernel="rbf",gamma = i,cache_size=5000).fit(Xtrain,Ytrain)
    score.append(clf.score(Xtest,Ytest))
    
print(max(score), gamma_range[score.index(max(score))])
plt.plot(gamma_range,score)
plt.show()

得到学习曲线:
在这里插入图片描述
通过学习曲线,很容就找出了rbf的最佳gamma值。
但对于多项式核函数来说,一切就没有那么容易了,因为三个参数共同作用在一个数学公式上影响它的效果。因此,我们往往使用网格搜索来共同调整三个对多项式核函数有影响的参数。
依然使用乳腺癌数据集举例。

from sklearn.model_selection import StratifiedShuffleSplit#用于支持带交叉验证的网格搜索
from sklearn.model_selection import GridSearchCV#带交叉验证的网格搜索
 
time0 = time()
 
gamma_range = np.logspace(-10,1,20)
coef0_range = np.linspace(0,5,10)
 
param_grid = dict(gamma = gamma_range
                  ,coef0 = coef0_range)
cv = StratifiedShuffleSplit(n_splits=5, test_size=0.3, random_state=420)#将数据分为5份,5份数据中测试集占30%
grid = GridSearchCV(SVC(kernel = "poly",degree=1,cache_size=5000
                        ,param_grid=param_grid
                        ,cv=cv)
grid.fit(X, y)
 
print("The best parameters are %s with a score of %0.5f" % (grid.best_params_, 
grid.best_score_))

最终得到:

The best parameters are {'coef0': 0.0, 'gamma': 0.18329807108324375} with a score of 0.96959

9、硬间隔与软间隔:重要参数C

到这里,我们已经了解了线性SVM的基本原理,以及SVM如何被推广到非线性情况下,还了解了核函数的选择和应用。但实际上,我们依然没有完全了解sklearn当中的SVM用于二分类的全貌。我们之前在理论推导中使用的数据都有一个特点,那就是他们或是完全线性可分,或者是非线性的数据。在我们对比核函数时,实际上用到了一种不同的数据,那就是不完全线性可分的数据集。比如说如下数据集:
在这里插入图片描述
这个数据集和我们最开始介绍SVM如何工作的时候的数据集一模一样,除了多了P和Q两个点。我们注意到,虽然决策边界 B 1 B_1 B1的间隔已经非常宽了,然而点P和Q依然被分错了类别,相反,边际比较小的 B 2 B_2 B2却正确地分出了点P和Q的类别。这里并不是说 B 2 B_2 B2此时此刻就是一条更好的边界了,与之前的论述中一致,如果我们引入更多的训练数据,或引入测试数据, B 1 B_1 B1更加宽敞的边界可以帮助它有更好的表现。但是,和之前不一样,现在即便是边际最大的决策边界 B 1 B_1 B1的训练误差也不可能为0了。此时,我们就需要引入“软间隔”的概念:
当两组数据是完全线性可分,我们可以找出一个决策边界使得训练集上的分类误差为0,这两种数据就被称为是存在“硬间隔”的。当两组数据几乎是完全线性可分的,但决策边界在训练集上存在较小的训练误差,这两种数据就被称为是存在“软间隔”。
我们可以通过调整我们对决策边界的定义,将硬间隔时得出的数学结论推广到软间隔的情况上,让决策边界能够忍受一小部分训练误差。这个时候,我们的决策边界就不是单纯地寻求最大边际了,因为对于软间隔地数据来说,边
际越大被分错的样本也就会越多,因此我们需要找出一个“最大边际”与“被分错的样本数量”之间的平衡。
在这里插入图片描述在这里插入图片描述在这里插入图片描述其中C是用来控制惩罚项的惩罚力度的系数。
我们的拉格朗日函数为(其中 μ \mu μ是第二个拉格朗日乘数):
在这里插入图片描述需要满足的KKT条件为:
在这里插入图片描述拉格朗日对偶函数为:
在这里插入图片描述在这里插入图片描述参数C用于权衡“训练样本的正确分类”与“决策函数的边际最大化”两个不可同时完成的目标,希望找出一个平衡点来让模型的效果最佳。
在这里插入图片描述在实际使用中,C和核函数的相关参数(gamma,degree等等)们搭配,往往是SVM调参的重点。与gamma不同,C没有在对偶函数中出现,并且是明确了调参目标的,所以我们可以明确我们究竟是否需要训练集上的高精确度来调整C的方向。默认情况下C为1,通常来说这都是一个合理的参数。 如果我们的数据很嘈杂,那我们往往减小C。当然,我们也可以使用网格搜索或者学习曲线来调整C的值。

#调线性核函数
score = []
C_range = np.linspace(0.01,30,50)
for i in C_range:
    clf = SVC(kernel="linear",C=i,cache_size=5000).fit(Xtrain,Ytrain)
    score.append(clf.score(Xtest,Ytest))
print(max(score), C_range[score.index(max(score))])
plt.plot(C_range,score)
plt.show()
 
#换rbf
score = []
C_range = np.linspace(0.01,30,50)
for i in C_range:
    clf = SVC(kernel="rbf",C=i,gamma = 0.012742749857031322,cache_size=5000).fit(Xtrain,Ytrain)
    score.append(clf.score(Xtest,Ytest))
    
print(max(score), C_range[score.index(max(score))])
plt.plot(C_range,score)
plt.show()
 
#进一步细化
score = []
C_range = np.linspace(5,7,50)
for i in C_range:
    clf = SVC(kernel="rbf",C=i,gamma = 
0.012742749857031322,cache_size=5000).fit(Xtrain,Ytrain)
    score.append(clf.score(Xtest,Ytest))
    
print(max(score), C_range[score.index(max(score))])
plt.plot(C_range,score)
plt.show()
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

cofisher

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值