支持向量机理论详解及使用二次规划包的底层实现

支持向量机

CSDN的Latex公式显示有误，可以自行复制编译查看

基础

支持向量机中涉及超平面以及点到超平面距离等概念，现在首先对相关概念进行介绍。

超平面

在$p$维空间中，超平面是$p-1$维的仿射子空间，其公式可以表达为
$${\beta }_0+{\beta }_1{x}_1+{\beta }_2{x}_2+\cdots +{\beta }_p{x}_p=\sum _{i=1}^p{\beta }_i{x}_i+{\beta }_0=0,$$
其法向量为
$$\omega =({\beta }_1,{\beta }_2,\cdots ,{\beta }_p),$$
当${\sum }_{i=1}^p{\beta }_i=0$时，其法向量为单位向量。

点到超平面距离

设$p$维空间中一点为
$$X=(x_1,x_2,\cdots ,x_n),$$
则该点到超平面的距离为
$$d=\frac{|\omega X+{\beta }_0|}{|\omega |},$$
该公式可由超平面法向量的概念推导得出，当${\sum }_{i=1}^p{\beta }_i=0$时，该公式可直接简化为
$$d=|\omega X+{\beta }_0|.$$

分割超平面

下面证明超平面可以用于对特征空间中的点进行二分类，这是最大间隔分类器的理论基础。

定义和性质

现在设有两类共$n$样本，每个样本均为$p$维空间中的一个点，其取值分别为
$$
x_1=\begin{Bmatrix}
x_{11}\x_{12}
\\vdots
\x_{1p}

\end{Bmatrix}
,x_2=\begin{Bmatrix}
x_{21}\x_{22}
\\vdots
\x_{2p}

\end{Bmatrix}
,\dotsb,
x_n=\begin{Bmatrix}
x_{n1}\x_{n2}
\\vdots
\x_{np}

\end{Bmatrix},
KaTeX parse error: Can't use function '$' in math mode at position 8: 其样本标签为$̲y_1,y_2,\dotsb,…
\beta_0 + \beta_1x_1 + \beta_2x_2+\dotsb+\beta_{p}x_p=0,
$$满足$$
\begin{cases}
\beta_0 + \beta_1x_{i1} + \beta_2x_{i2}+\dotsb+\beta_{p}x_{ip}>0,& \text{ if } y_i = 1, \
\beta_0 + \beta_1x_{i1} + \beta_2x_{i2}+\dotsb+\beta_{p}x_{ip}<0,& \text{ if } y_i=-1,
\end{cases}
$$
则称该超平面为分割超平面。

超平面分离定理

超平面可以用于二分类的原理即超平面分离定理，也称为凸集分离定理。

引理

设$S\subseteq R^n$为非空凸集，$y\subseteq R^n\setminus S$，则存在唯一的点$\bar{x}$，使得该点与$y$的距离最小，即有
∥ x ˉ − y ∥ = inf ⁡ { ∥ x − y ∥ ∣ x ∈ S } > 0 , \|\bar{x}-y\|=\inf \{\|x-y\| \mid x \in S\}>0, ∥xˉ−y∥=inf{∥x−y∥∣x∈S}>0,
同时由范数的等价性可以得知，这里的范数可以是任意一种。

引理证明

先证明其存在性，考虑单位超球
B = { z ∈ R n ∣ ∥ z ∥ ≤ 1 } , B=\left\{z \in R^{n} \mid\|z\| \leq 1\right\}, B={z∈Rn∣∥z∥≤1},
取足够大的正数$\beta$，使得
D = S ∩ ( { y } + β B ) ≠ ∅ , D=S \cap(\{y\}+\beta B) \neq \emptyset, D=S∩({y}+βB)=∅,
由$S$为闭集， { y } + β B \{y\}+\beta B {y}+βB为有界闭集可知，$D$是非空有界闭集，故$S$可以在$D$上取得一个点$\bar{x}$使其满足到带你$y$的距离最小，存在性得证。

再证唯一性，记上述最短距离为$r$，设存在另外一点$\overline{x^{\prime }}\in S$满足
$$\Vert \bar{x}-y\Vert =\Vert \overline{x^{\prime }}-y\Vert =r,$$
设$x^{\prime \prime }=\frac{1}{2}(\bar{x}+\overline{x^{\prime }})$，由$x^{\prime \prime }-y=\frac{1}{2}x-\frac{1}{2}y+\frac{1}{2}{x}^{\prime }-\frac{1}{2}y$可知，两边取范数有
$$\Vert x^{\prime \prime }-y\Vert \le \frac{1}{2}\Vert x-y\Vert +\frac{1}{2}\Vert x^{\prime }-y\Vert =r,$$
又由$S$是凸集可知，$x^{\prime \prime }$是$\bar{x}\in S$与$\overline{x^{\prime }}\in S$的一个凸组合，即$x^{\prime \prime }\in S$。

由$r$是$y$到$S$的最短距离知道
$$\Vert x^{\prime \prime }-y\Vert =r,$$
由平行四边形定律可知
$${\Vert \bar{x}-\overline{x^{\prime }}\Vert }^2+4{\Vert x^{\prime \prime }-y\Vert }^2=2\Vert \bar{x}-y{\Vert }^2+2{\Vert \overline{x^{\prime }}-y\Vert }^2,$$
将$r$代入对应值得到
$${\Vert \bar{x}-\overline{x^{\prime }}\Vert }^2=2r^2+2r^2-4r^2=0,$$
即$\overline{x^{\prime }}=\bar{x}$，唯一性得证。

定理

设$S_1,S_2\subseteq R^n$，若存在非零向量$p\in R^n$和$\alpha \in R$使得
$$
\begin{array}{l}
S_{1} \subseteq H^{-}=\left{x \in R^{n} \mid p^{T} x \leq \alpha\right}, \
S_{2} \subseteq H^{+}=\left{x \in R^{n} \mid p^{T} x \geq \alpha\right},
\end{array}

$$
则称超平面$H=\left\{x\in R^n\mid p^{T}x=\alpha \right\}$分离了集合$S_1$和$S_2$。

定理证明

由$S$是非空集合可知，对$S$外任意一点$y$，存在$x^{\prime }\in S$使得
∥ x ′ − y ∥ = inf ⁡ { ∥ x − y ∥ x ∈ S } > 0 , \left\|x^{\prime}-y\right\|=\inf \{\|x-y\| x \in S\}>0, ∥x′−y∥=inf{∥x−y∥x∈S}>0,
设$x\in S$，由$S$是凸集可知，对任意的$\lambda \in (0,1)$，有
$$z=\lambda x+(1-\lambda )x^{\prime }\in S,$$
故有
KaTeX parse error: Expected 'EOF', got '&' at position 34: …-y\right\|^{2} &̲ \leq\|z-y\|^{2…
去掉上式两边的$\Vert x^{\prime }-y{\Vert }^2$可以得到
$$\lambda {\Vert x-x^{\prime }\Vert }^2+2{\left(x^{\prime }-y\right)}^T\left(x-x^{\prime }\right)\ge 0,\forall \lambda \in (0,1),\forall x\in S,$$
令$\lambda \to 0^{+}$，有
$${\left(x^{\prime }-y\right)}^T\left(x-x^{\prime }\right)\ge 0,\forall x\in S,$$
记$p=y-x^{\prime }\ne 0$，有
$$p^T(x-x^{\prime })\le 0,$$
又记$\alpha =p^T{x}^{\prime }$，此时有
$$p^T\le \alpha ,$$
同时又由
$$p^{T}y-\alpha =p^T\left(y-x^{\prime }\right)={\left(y-x^{\prime }\right)}^T\left(y-x^{\prime }\right)={\Vert y-x^{\prime }\Vert }^2>0$$
可知
$$p^{T}x\le \alpha \le p^{T}y,$$
即将该集合分成了两部分。

最大间隔分类器

原理

如果数据可以被超平面分开，那在数据的特征空间上就存在无限个可以将数据分开的超平面，这时一般考虑最大间隔超平面，即使得样本到超平面的最小间隔最大化，同时，用于计算样本集合到超平面的最小间隔的所有样本点被称为支持向量（如果有多个点，那么这些点到超平面的距离相同）。

可以将最大间隔分类器中的超平面参数视为规划问题中的待求值，而所有样本点都是约束条件（一个样本点是一个约束条件），可以列出规划问题如下

$$
\operatorname{maximize}{\beta{0}, \beta_{1}, \cdots, \beta_{p}} M\

\left{\begin{matrix}
\sum_{j=1}^{p} \beta_{j}^{2}=1,\
y_{i}\left(\beta_{0}+\beta_{1} x_{i \mathrm{t}}+\beta_{2} x_{i z}+\cdots+\beta_{p} x_{i p}\right) \geqslant M, \quad i=1, \cdots, n\\beta_j无约束,\quad j=1,2,3,\dotsb,p,
\end{matrix}\right.
$$

其中目标函数${\mathrm{maximize}}_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p}M$，其实可以直接写成${\max }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p}M$，这样得到的结果其实是和原来一样的。

另外，第一个约束$ \sum_{j=1}^{p} \beta_{j}^{{2}=1$用于保证求出的结果中$M$就是最大间隔分类器中样本点与超平面的最小距离，如果将该约束去掉，对上述规划问题求出的结果而言，最小距离则化为$\frac{M}{\sum_{j=1}}{p}\beta_j^2}$；如果在去掉该约束的同时，将第二个约束修改为
$$\begin{gathered} \frac{1}{\Vert \omega \Vert }{y}_i({\beta }_0+{\beta }_1{x}_{i\mathrm{t}}+{\beta }_2{x}_{iz}+\cdots +{\beta }_p{x}_{ip})\ge M, \\ 或者 \\ {y}_i({\beta }_0+{\beta }_1{x}_{i\mathrm{t}}+{\beta }_2{x}_{iz}+\cdots +{\beta }_p{x}_{ip})\ge M\Vert \omega \Vert , \end{gathered}$$
此时$M$就仍然表示样本点集合到超平面的最小距离，若令$M=\frac{1}{\Vert \omega \Vert }$，则目标函数可以进一步写成${\max }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p}\frac{1}{\Vert \omega \Vert }$，即
$${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p}\Vert \omega \Vert .$$
第二个约束则是保证样本点全部落在超平面正确的一侧。

由于每个样本点都会影响最大间隔分类器中超平面的参数，所以两类样本中一旦出现距离较小的点，这些点将对超平面产生巨大的影响，这也容易因为某些点的出现而出现过拟合，反而会造成对测试集划分效果不好的情况，所以这就引入了对样本点划分没有那么严格的支持向量分类器（也称为软间隔分类器）。

代码说明

导入库

import numpy as np
from scipy.optimize import linprog

不等式约束矩阵

def GetA(x, y):
    x = np.array(x)
    y = np.array(y)
    height, width = x.shape
    y = -y.reshape((height, 1))
    tmp = x * y  # 计算得到超平面样本点的系数
    tmp = np.concatenate((tmp, y), axis=1)  # 增加常量的正负系数
    tmp = np.concatenate((tmp, np.ones((height, 1))), axis=1)  # 增加M的系数
    return tmp

使用最开始的规划问题进行编程，其中
$$y_i\left({\beta }_0+{\beta }_1{x}_{i\mathrm{t}}+{\beta }_2{x}_{iz}+\cdots +{\beta }_p{x}_{ip}\right)⩾M,i=1,\cdots ,n$$
改写为
$$y_i\left({\beta }_0+{\beta }_1{x}_{i\mathrm{t}}+{\beta }_2{x}_{iz}+\cdots +{\beta }_p{x}_{ip}\right)-M⩾0,i=1,\cdots ,n$$
由于最终使用的scipy.optimize.linprog中要求的规 划问题格式为：不等式的约束条件中符号均为小于等于，目标函数为求最小值，所以继续改写为
$$-y_i\left({\beta }_0+{\beta }_1{x}_{i\mathrm{t}}+{\beta }_2{x}_{iz}+\cdots +{\beta }_p{x}_{ip}\right)+M⩽0,i=1,\cdots ,n$$
若不等式约束设为$AX\le b$，则根据上述公式可以得到不等式约束矩阵为
$$A=\left[\begin{array}{ccccc}-x_{11}{y}_1& \cdots & -x_{1p}{y}_1& -y_1& 1\\ ⋮& \ddots & ⋮& ⋮& ⋮\\ -x_{n1}{y}_n& \cdots & -x_{np}{y}_n& -y_n& 1\end{array}\right],$$
其中待求变量为
$$X={\left[\begin{array}{cccccc}{\beta }_1& {\beta }_2& \cdots & {\beta }_p& {\beta }_0& M\end{array}\right]}^T,$$
不等式约束变量为
$$b={\left[\begin{array}{ccccc}0& 0& 0& \cdots & 0\end{array}\right]}^T,$$
目标函数为
$$\max M=\min (-M).$$

不等式约束向量

def GetB(x):
    x = np.array(x)
    height, width = x.shape
    return np.zeros((height,))  # 不等式约束条件右边的b全部为0，约束条件个数即为样本点个数

目标函数系数

def GetC(x):
    x = np.array(x)
    height, width = x.shape
    c_num = width + 2  # 需要求解的参数的个数 = 特征空间维度 + 超平面方程中的一个常量 + 所求大M
    tmp = np.zeros((c_num,))
    tmp[-1] = -1
    return tmp

决策变量的界

def GetBound(x):
    x = np.array(x)
    height, width = x.shape
    get = []
    for i in range(width + 1):
        get.append((None, None))
    get.append((0, None))  # 只有M需要约束其界大于等于0
    return tuple(get)

求解

def solve(x, y):
    A = GetA(x, y)
    B = GetB(x)
    C = GetC(x)
    Bound = GetBound(x)
    r = linprog(c=C, A_ub=A, b_ub=B, bounds=Bound)
    return r

值得注意的是，该求解最终得到的结果中的$M$并不是样本结合到超平面的最小间隔，需要进一步除以超平面法向量的模来求得真实间隔。

其中返回的结果结构如下

     con: array([], dtype=float64)
     fun: -54452226.63349579
 message: 'The algorithm terminated successfully and determined that the problem is unbounded.'
     nit: 4
   slack: array([4.72635938e+07, 8.77588747e+07, 7.34542083e+07, 3.05744090e+08,
       2.61286827e+08, 3.74980180e+08, 4.90151060e+08, 4.86659674e+08,
       4.48112515e+08, 4.65908568e+08, 1.98497176e+08, 6.13494814e+08,
       2.90681402e+08, 4.94361952e+08, 1.13727802e+08, 4.84489864e+08,
       1.39447820e+08, 7.81356956e+07, 2.64020192e+08, 1.25613750e+08,
       1.84441418e+08, 2.09625102e+08, 2.15064611e+08])
  status: 3
 success: False
       x: array([ 1.77960531e+08, -1.23565441e+08, -5.52446963e+07,  5.44522266e+07])

该结果中的$x$前三个数分别代表超平面
$$Ax+By+C=0$$
中的$A、B、C$，可根据该参数得到超平面方程。

结果

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-78sQY5yu-1678444236210)(D:\DESKTOP\我的坚果云\typora图片\image-20221116171624844.png)]

二次规划问题

根据上文中的推导，当损失函数设置为${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p,{\epsilon }_1,{\epsilon }_2,\cdots ,{\epsilon }_n}\Vert \omega \Vert$时，规划问题如下
$$
\operatorname{min}{\beta{0}, \beta_{1}, \cdots, \beta_{p}} |\omega|\

s.t.\left{\begin{matrix}
y_{i}\left(\beta_{0}+\beta_{1} x_{i \mathrm{t}}+\beta_{2} x_{i z}+\cdots+\beta_{p} x_{i p}\right) \geqslant 1, \quad i=1, \cdots, n\
\beta_j无约束,\quad j=0,1,2,3,\dotsb,p,
\end{matrix}\right.
KaTeX parse error: Can't use function '$' in math mode at position 26: …示的意义，可以将目标函数改写为$̲\operatorname{m…
\operatorname{min}{\beta{0}, \beta_{1}, \cdots, \beta_{p}} \frac{1}{2}\omega^T \omega\

s.t.\left{\begin{matrix}
y_{i}\left(\omega^Tx_i+b\right) \geqslant 1, \quad i=1, \cdots, n\
\omega无约束
\end{matrix}\right.
$$可以写出拉格朗日函数：$$
L(\omega,b,\lambda)=\frac{1}{2}\omega^{T\omega+\sum_{i=1}}{N}\lambda_i[1-y_{i}\left(\omega^Tx_i+b\right)],
$$
若求$\max L(\omega ,b,\lambda )$时，若所有样本点$x_i$均满足原规划问题中的约束，有$1-y_i\left({\omega }^T{x}_i+b\right)⩽0$，此时所有的${\lambda }_i$均取$0$时$L(\omega ,b,\lambda )$才能取到最大值，且最大值为$L(\omega ,b,\lambda )=\frac{1}{2}{\omega }^T\omega$；反之，若存在样本点$x_i$不满足原规划问题中的约束，对应位置的$1-y_i\left({\omega }^T{x}_i+b\right)⩾0$，此时若要使得$L(\omega ,b,\lambda )$取最大值，对应的${\lambda }_i$取正无穷时可以取得最大值即为正无穷。

综上所述，有
$$\underset{\omega ,b}{\min }\underset{\lambda }{\max }L(\omega ,b,\lambda )=\underset{\omega ,b}{\min }(+\infty ,\frac{1}{2}{\omega }^T\omega )=\underset{\omega ,b}{\min }\frac{1}{2}{\omega }^T\omega ,$$
此时原规划问题化为
$$
\min_{\omega,b} \max_{\lambda} L(\omega,b,\lambda)\

s.t.\quad \lambda_i\geqslant0, \quad i=1, \cdots, n\
$$其对偶问题为$$
\max_{\lambda} \min_{\omega,b} L(\omega,b,\lambda)\

s.t.\quad \lambda_i\geqslant0, \quad i=1, \cdots, n\
KaTeX parse error: Can't use function '$' in math mode at position 10: 求解对偶问题，令$̲\frac{\partial …
\begin{align}L(\omega,b,\lambda) & =\frac{1}{2}\omega^{T\omega+\sum_{i=1}}{n}\lambda_i-\sum_{i=1}^{n}\lambda_i y_i \omega^T x_i-\sum_{i=1}^{n}\lambda_iy_ib\
&=\frac{1}{2}\omega^{T\omega+\sum_{i=1}}{n}\lambda_i-\sum_{i=1}^{n}\lambda_i y_i \omega^T x_i,
\end{align}
KaTeX parse error: Can't use function '$' in math mode at position 4: 又令$̲\frac{\partial …
\begin{align}L(\omega,b,\lambda) & =\frac{1}{2}(\sum_{i=1}^{{n}\lambda_iy_ix_i)}T(\sum_{i=1}^{{n}\lambda_iy_ix_i)+\sum_{i=1}}{n}\lambda_i-\sum_{i=1}^{n}\lambda_i y_i (\sum_{i=1}^{{n}\lambda_iy_ix_i)}T x_i\
&=-\frac{1}{2}\sum_{i=1}^{{n}\sum_{j=1}}{n}\lambda_i\lambda_jy_iy_jx_i^{Tx_j+\sum_{i=1}}{n}\lambda_i,
\end{align}
$$此时原对偶问题化为$$
\max_{\lambda} -\frac{1}{2}\sum_{i=1}^{{n}\sum_{j=1}}{n}\lambda_i\lambda_jy_iy_jx_i^{Tx_j+\sum_{i=1}}{n}\lambda_i\

s.t.\left{\begin{matrix}
\lambda_i \geqslant 1, \quad i=1, \cdots, n\
\sum_{i=1}^{n}\lambda_i y_i=0
\end{matrix}\right.
KaTeX parse error: Can't use function '$' in math mode at position 23: …多数凸优化求解包中目标函数都是$̲\min$的要求，在目标函数上…
\min_{\lambda} \frac{1}{2}\sum_{i=1}^{{n}\sum_{j=1}}{n}\lambda_i\lambda_jy_iy_jx_i^{Tx_j-\sum_{i=1}}{n}\lambda_i\

s.t.\left{\begin{matrix}
\lambda_i \geqslant 1, \quad i=1, \cdots, n\
\sum_{i=1}^{n}\lambda_i y_i=0
\end{matrix}\right.
$$
由于上述推导过程中已经满足KKT条件，所以该二次规划问题的最优解即为原问题的最优解。

支持向量分类器

原理

支持向量分类器并不寻求最大的间隔，而是希望对大部分样本点划分正确，这样可以使得模型的泛化能力更强。

可以通过修改最大间隔分类器的规划问题满足新的要求，具体规划问题为
$$
\operatorname{maximize}{\beta{0}, \beta_{1}, \cdots, \beta_{p}, \varepsilon_1,\varepsilon_2,\cdots,\varepsilon_n} M\

s.t.\left{\begin{matrix}
\sum_{j=1}^{p} \beta_{j}^{2}=1,\
y_{i}\left(\beta_{0}+\beta_{1} x_{i \mathrm{t}}+\beta_{2} x_{i z}+\cdots+\beta_{p} x_{i p}\right) \geqslant M(1-\varepsilon_i), \quad i=1, \cdots, n\
\varepsilon_i\ge0,\quad \sum_{i=1}^{n}\varepsilon_i\le C,\
\beta_j无约束,\quad j=1,2,3,\dotsb,p,
\end{matrix}\right.
$$
其中${\epsilon }_i$的设置是允许一部分样本点被划分在间隔内甚至是被划分到超平面的错误一侧，又由于${\epsilon }_i$大于1时表示得到该约束条件的样本点被错误分类到了超平面的另一侧，所以$C$（需要由自己设置）代表该支持向量分类器分类错误（被划分到超平面的错误一侧）的最大样本点个数，且当$C=0$时支持向量分类器恢复为最大间隔分类器（因为所有的${\epsilon }_i$都将为0）。

若根据上文使用${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p}\Vert \omega \Vert$作为目标函数，此时该规划问题变为一个凸优化问题，其形式如下
$$
\operatorname{min}{\beta{0}, \beta_{1}, \cdots, \beta_{p}, \varepsilon_1,\varepsilon_2,\cdots,\varepsilon_n} |\omega|\

s.t.\left{\begin{matrix}
y_{i}\left(\beta_{0}+\beta_{1} x_{i \mathrm{t}}+\beta_{2} x_{i z}+\cdots+\beta_{p} x_{i p}\right) \geqslant 1-\varepsilon_i, \quad i=1, \cdots, n\
\varepsilon_i\ge0,\quad \sum_{i=1}^{n}\varepsilon_i\le C,\
\beta_j无约束,\quad j=1,2,3,\dotsb,p,
\end{matrix}\right.
$$
同时，在支持向量分类器中，刚好落在间隔上和落在间隔错误一侧的样本点被称为支持向量（比最大间隔分类器中的支持向量增加了落在间隔错误一侧的样本点），且只有支持向量会对支持向量分类器的参数产生影响。

但是支持向量分类器只对线性可分的样本点起作用，如果样本点集合不是线性可分的，就需要引入支持向量机。

代码编写推导

对上述目标函数为${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p}\Vert \omega \Vert$的规划问题，可以将第一行约束条件化为
$$y_i\left({\beta }_0+{\beta }_1{x}_{i\mathrm{t}}+{\beta }_2{x}_{iz}+\cdots +{\beta }_p{x}_{ip}\right)+{\epsilon }_i⩾1,i=1,\cdots ,n$$
可以得到其不等式约束矩阵为
$$A_1=\left[\begin{array}{cccccccc}{x}_{11}{y}_1& \cdots & x_{1p}{y}_1& y_1& 1& 1& \cdots & 0\\ ⋮& \ddots & ⋮& ⋮& ⋮& ⋮& \ddots & ⋮\\ x_{n1}{y}_n& \cdots & x_{np}{y}_n& y_n& 1& 0& \cdots & 1\end{array}\right],$$
待求变量为
$$X={\left[\begin{array}{ccccccccc}{\beta }_1& {\beta }_2& \cdots & {\beta }_p& {\beta }_0& {\epsilon }_1& {\epsilon }_2& \cdots & {\epsilon }_n\end{array}\right]}^T,$$
为了保持与上面约束条件的符号方向相同，第二行的约束条件可以化为
$$-{\epsilon }_1-{\epsilon }_2-\cdots -{\epsilon }_n⩾-C,$$
除此之外，还有约束条件${\epsilon }_i\ge 0$需要满足，将其约束矩阵均添加到总的约束矩阵中可以得到
$$A=\left[\begin{array}{cccccccc}{x}_{11}{y}_1& \cdots & x_{1p}{y}_1& y_1& 1& 1& \cdots & 0\\ ⋮& \ddots & ⋮& ⋮& ⋮& ⋮& \ddots & ⋮\\ x_{n1}{y}_n& \cdots & x_{np}{y}_n& y_n& 1& 0& \cdots & 1\\ 0& \cdots & 0& 0& -1& -1& \cdots & -1\\ 0& \cdots & 0& 0& 0& 1& \cdots & 0\\ ⋮& \ddots & ⋮& ⋮& ⋮& ⋮& \ddots & 0\\ 0& \cdots & 0& 0& 0& 1& \cdots & 1\end{array}\right],$$
同时可以写出约束向量为
$$b={\left[\begin{array}{ccccccccc}1& 1& \cdots & 1& -C& 0& 0& \cdots & 0\end{array}\right]}^T,$$
其中共有$n$个1、$n$个0和1个$C$，可以根据上述推导写出约束条件被化为
$$AX⩾b.$$

代码实现

导入库

import numpy as np
import matplotlib.pyplot as plt
import matplotlib
import cvxpy as cp
from cvxopt import matrix, solvers

系数矩阵获取

def GetA(x, y):
    x = np.array(x)
    y = np.array(y)
    height, width = x.shape
    y = y.reshape((height, 1))
    tmp = x * y  # 计算得到超平面样本点的系数
    tmp = np.concatenate((tmp, y), axis=1)  # 增加常量的正负系数
    tmp = np.concatenate((tmp, np.eye(height)), axis=1)  # 增加beta0的系数
    # 使得所有epsilon之和小于等于C
    bottom = np.concatenate((np.zeros((1, width + 1)), -np.ones((1, height))), axis=1)  
    tmp = np.concatenate((tmp, bottom), axis=0)
    # 使得所有epsilon的取值大于等于0
    bottom = np.concatenate((np.zeros((height, width + 1)), np.eye(height)), axis=1)
    tmp = np.concatenate((tmp, bottom), axis=0)
    return tmp

约束向量获取

def GetB(x, C):
    x = np.array(x)
    height, width = x.shape
    tmp = np.ones((height,))  # 前面由样本点得到的约束条件全是大于等于1
    tmp = np.append(tmp, -C)  # 所有epsilon之和应该小于等于C，这里因为用的是大于等于矩阵，所以使用-C
    tmp = np.concatenate((tmp,np.zeros(height,)))
    return tmp

求解实现

def Solve(x, y, C):
    x = np.array(x)
    y = np.array(y)
    height, width = x.shape  # n = height, p = width
    A = GetA(x, y)
    B = GetB(x, C)
    # print(A.shape)
    # 这里是目标函数，由于cvxopt库会自动在二次型前加1/2，所以对结果值×2在开方才是法向量的模
    P = matrix(np.diag(np.concatenate((np.ones((width,)),np.zeros((height + 1,))))))
    # 没有使用到后面不同参数相乘的部分
    q = matrix(np.zeros((height+width+1)))
    # 由于cvxopt库中默认是小于等于，所以需要将推导中的大于等于修改为小于等于
    G = matrix(-A*1.0)
    # print(G.size)
    h = matrix(-B)
    sol = solvers.qp(P, q, G, h)
    print("primal objective:      ",sol['primal objective'])
    print("primal objective（乘2并开方）:      ", (sol['primal objective'])**(1/2))
    return np.array(sol['x']), sol['primal objective']

实验

x1 = [-0.4, 0.9, -0.2, 0.3, -0.7, -1.2, -1.5, -0.3, -0.5, -0.6, 0]
y1 = [-0.2, 2, 0.3, 2.9, 1.1, 1.3, 1.8, 3.5, 2.9, 2.9, 1.6]
x2 = [3.3, 1, 2.7, 0.7, 3.2, 1.4, 0.5, 2.1, 1.6, 2, 1.1, 1.2]
y2 = [-1.1, -1.8, -1, -0.8, -0.2, 0, -0.8, 0, 0.4, 0.5, -1, -0.9]
data = []
target = []
for i in range(len(x1)):
    data.append([x1[i], y1[i]])
    target.append(-1)
for j in range(len(x2)):
    data.append([x2[j], y2[j]])
    target.append(1)
data = np.array(data)

x_get,M_get = Solve(data,target,2)
A,B,C = x_get[:3]
cos = np.sqrt(B**2/(A**2+B**2))

画图

d = float((1/np.sqrt(2*M_get))/cos)
xx1 = np.linspace(-5, 5)
k = -A / B
t = -C / B
yy1 = k * xx1 + t
yy_up1 = k * xx1 + (t + d)
# yy_down1 = k * xx1 + (t - 0.8 - (k * 0.5 + t))
yy_down1 = k * xx1 + (t - d)

# 画出超平面
plt.plot(xx1, yy1, 'k-')
plt.plot(xx1, yy_down1, 'k--')
plt.plot(xx1, yy_up1, 'k--')

plt.fill_between(xx1, yy_down1, yy_up1, color='green', alpha=.25)  

cm_dark = matplotlib.colors.ListedColormap(['r', 'b'])
plt.scatter(data[:, 0], data[:, 1],
            s=80, facecolors='none')
plt.scatter(data[:, 0], data[:, 1], c=target, cmap=cm_dark)

plt.xlim(xmin=-2, xmax=4)
plt.ylim(ymin=-2, ymax=4)

plt.show()

结果

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-UwpW5UZZ-1678444236212)(D:\DESKTOP\我的坚果云\typora图片\image-20221116171724324.png)]

可以看到，在设置$C$后，有部分点被划分在了间隔的错误一侧，但没有点被划分在超平面错误的一侧。

支持向量机

核函数与解函数

支持向量机的核心就是通过使用某些函数将当前的样本点映射到另一个特征空间再使用超平面进行划分，这种函数称为核函数。

通过对上述规划问题（二次规划）的求解，最后可以求得一个解函数为
$$f(x)=\sum _{i=1}^N{\alpha }_i{y}_{i}K\left(x,x_i\right)+{\beta }_0,$$
即，将样本点输入该函数就可以得到其在超平面的哪一侧，并且解函数中只涉及核函数的内积，并不需要清晰知道核函数的计算，知道其变换后内积的计算就足够了，常用的核函数（内积计算）如下

名称	核函数内积
线性核函数	$K\left(x,x^{\prime }\right)=⟨x,x_i⟩$
$d$阶多项式核函数	$K\left(x,x^{\prime }\right)={\left(1+⟨x,x^{\prime }⟩\right)}^d$
高斯核函数	$K\left(x,x^{\prime }\right)=\exp \left(-\gamma {\left|x-x^{\prime }\right|}^2\right)$

二次规划形式

现在以最大间隔分类器二次规划形式的推导过程为基础，对软间隔分类器的二次规划形式进行推导。

当损失函数设置为${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p,{\epsilon }_1,{\epsilon }_2,\cdots ,{\epsilon }_n}\Vert \omega \Vert$时,软间隔分类器的规划形式为
$${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p,{\epsilon }_1,{\epsilon }_2,\cdots ,{\epsilon }_n}\Vert \omega \Vert$$
$$s.t.\{\begin{array}{c}{y}_i\left({\omega }^T{x}_i+b\right)⩾1-{\epsilon }_i,i=1,\cdots ,n\\ {\epsilon }_i\ge 0,{\sum }_{i=1}^n{\epsilon }_i\le C,\end{array}$$
与最大间隔分类器类似的，其目标函数可以改写为${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p,{\epsilon }_1,{\epsilon }_2,\cdots ,{\epsilon }_n}\frac{1}{2}{\omega }^T\omega$，此时规划问题为
$${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p,{\epsilon }_1,{\epsilon }_2,\cdots ,{\epsilon }_n}\frac{1}{2}{\omega }^T\omega$$
$$s.t.\{\begin{array}{c}{y}_i\left({\omega }^T{x}_i+b\right)⩾1-{\epsilon }_i,i=1,\cdots ,n\\ {\epsilon }_i\ge 0,{\sum }_{i=1}^n{\epsilon }_i\le C,\end{array}$$

为了方便后续的求解，可以将约束条件中与$C$有关的部分去掉，并将目标函数修改为${\min }_{{\beta }_0,{\beta }_1,\cdots ,{\beta }_p,{\epsilon }_1,{\epsilon }_2,\cdots ,{\epsilon }_n}\frac{1}{2}{\omega }^T\omega +C{\sum }_{i=1}^n{\epsilon }_i$，即得到如下形式

KaTeX parse error: Undefined control sequence: \label at position 379: …{array}\right. \̲l̲a̲b̲e̲l̲{initial}

可以写出对应的拉格朗日函数：
$$L(\omega ,b,\epsilon ,\lambda ,\mu )=\frac{1}{2}{\omega }^T\omega +C\sum _{i=1}^n{\epsilon }_i-\sum _{i=1}^n{\lambda }_i\left[y_i\left({\omega }^T\cdot x_i+b\right)-1+{\epsilon }_i\right]-\sum _{i=1}^n{\mu }_i{\epsilon }_i,$$
其中${\lambda }_i,{\mu }_i⩾0$。

经过与最大间隔分类器中类似的分析，可以知道规划问题\ref{initial}可以写成如下形式
$$\begin{gathered} \underset{\omega ,b,\epsilon }{\min }\underset{\lambda ,\mu }{\max }L(\omega ,b,\epsilon ,\lambda ,\mu ) \\ \text{ s.t. }\{\begin{array}{c}{\lambda }_i⩾0,i=1,\cdots ,n\\ {\epsilon }_i⩾0,i=1,\cdots ,n,\end{array} \end{gathered}$$
其对偶问题为
$$\begin{gathered} \underset{\lambda ,\mu }{\max }\underset{\omega ,b,\epsilon }{\min }L(\omega ,b,\epsilon ,\lambda ,\mu ) \\ \text{ s.t. }\{\begin{array}{c}{\lambda }_i⩾0,i=1,\cdots ,n\\ {\epsilon }_i⩾0,i=1,\cdots ,n,\end{array} \end{gathered}$$
对该对偶问题进行分析，求${\min }_{\omega ,b,\epsilon }L(\omega ,b,\epsilon ,\lambda ,\mu )$，令其满足
$$\begin{array}{l}{\nabla }_{w}L(w,b,\epsilon ,\lambda ,\mu )=w-{\sum }_{i=1}^N{\lambda }_i{y}_i{x}_i=0,\\ {\nabla }_{b}L(w,b,\epsilon ,\lambda ,\mu )=-{\sum }_{i=1}^N{\lambda }_i{y}_i=0,\\ {\nabla }_{{\epsilon }_i}L(w,b,\epsilon ,\lambda ,\mu )=C-{\lambda }_i-{\mu }_i=0,\end{array}$$
可以求得
$$\begin{gathered} \omega =\sum _{i=1}^n{\lambda }_i{y}_i{x}_i, \\ \sum _{i=1}^n{\lambda }_i{y}_i=0, \\ C-{\lambda }_i-{\mu }_i=0, \end{gathered}$$
将其代入$L(\omega ,b,\epsilon ,\lambda ,\mu )$可得其最小值为
$$\underset{\omega ,b,\epsilon }{\min }L(w,b,\epsilon ,\alpha ,\mu )=-\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_j{x}_i^T{x}_j+\sum _{i=1}^n{\lambda }_i,$$
此时得到的对偶问题可以写成
$$\begin{gathered} \underset{\lambda ,\mu }{\max }-\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_j{x}_i^T{x}_j+\sum _{i=1}^n{\lambda }_i \\ \text{ s.t. }\{\begin{array}{c}{\sum }_{i=1}^n{\lambda }_i{y}_i=0\\ C-{\lambda }_i-{\mu }_i=0,\\ {\lambda }_i⩾0,i=1,\cdots ,n\\ {\epsilon }_i⩾0,i=1,\cdots ,n,\end{array} \end{gathered}$$
又由式$C-{\lambda }_i-{\mu }_i=0$可知，可将约束条件${\epsilon }_i⩾0$改写为$0⩽{\lambda }_i⩽C$，即得到
$$\begin{gathered} \underset{\lambda }{\max }-\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_j{x}_i^T{x}_j+\sum _{i=1}^n{\lambda }_i \\ \text{ s.t. }\{\begin{array}{c}{\sum }_{i=1}^n{\lambda }_i{y}_i=0\\ 0⩽{\lambda }_i⩽C,i=1,\cdots ,n\end{array} \end{gathered}$$

可以从上面的推导过程中得出，当原问题与对偶问题的解相等时，需要该问题符合以下KKT条件
$$\begin{gathered} {\nabla }_{w}L\left(\omega ,b,\epsilon ,\lambda ,\mu \right)=\omega -\sum _{i=1}^n{\omega }_i{y}_i{x}_i=0, \\ {\nabla }_{b}L\left(\omega ,b,\epsilon ,\lambda ,\mu \right)=-\sum _{i=1}^n{\lambda }_i{y}_i=0, \\ {\nabla }_{\xi }L\left(\omega ,b,\epsilon ,\lambda ,\mu \right)=C-\lambda -\mu =0, \\ {\lambda }_i\left[y_i\left({\omega }^T\cdot x_i+b\right)-1+{\epsilon }_i\right]=0, \\ {\mu }_i{\epsilon }_i=0, \end{gathered}$$
其中前四条与最大间隔分类器二次规划形式的推导中类似，不再赘述。针对最后一条式子，有

• 当${\epsilon }_i\ne 0$时，对应样本点在间隔的错误一侧，为使式(\ref{lagrange2})最大化，有${\mu }_i=0$，此时即${\mu }_i{\epsilon }_i=0$，且${\lambda }_i=C$;
• 当${\epsilon }_i=0$时，对应样本点在间隔的错误一侧，可有${\mu }_i⩾0$，此时即${\mu }_i{\epsilon }_i=0$，${\lambda }_i⩽C$;

可以知道${\mu }_i{\epsilon }_i=0$成立。

综上所述，对问题原问题的求解可以化为以下问题的求解：
$$\begin{gathered} \underset{\lambda }{\min }\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_j{x}_i^T{x}_j-\sum _{i=1}^n{\lambda }_i \\ \text{ s.t. }\{\begin{array}{c}{\sum }_{i=1}^n{\lambda }_i{y}_i=0\\ 0⩽{\lambda }_i⩽C,i=1,\cdots ,n\end{array} \end{gathered}$$
与最大间隔分类器二次规划形式的推导中类似，存在一个样本点满足其对应的$0<{\lambda }_i⩽C$，可由该样本点得到$y_j\left({\omega }^T\cdot x_i+b\right)-1=0$，即可以得到超平面中的参数存在其所求超平面的参数为
KaTeX parse error: Undefined control sequence: \label at position 40: …ambda_i y_ix_i,\̲l̲a̲b̲e̲l̲{finalomega}\\ …
若使用了核函数，则该二次规划问题写为
$$\begin{gathered} \underset{\lambda }{\min }\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_{j}K(x_i,x_j)-\sum _{i=1}^n{\lambda }_i \\ \text{ s.t. }\{\begin{array}{c}{\sum }_{i=1}^n{\lambda }_i{y}_i=0\\ 0⩽{\lambda }_i⩽C,i=1,\cdots ,n\end{array} \end{gathered}$$
此时类别判别函数为
$$f(x)=\mathrm{sgn}(\sum _{i=1}^N{\alpha }_i{y}_{i}K\left(x,x_i\right)+b).$$

代码编写

根据二次规划形式的推导可以知道，对支持向量机中超平面的求解可以化为对以下问题的求解
$$\begin{gathered} \underset{\lambda }{\min }\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_{j}K(x_i,x_j)-\sum _{i=1}^n{\lambda }_i \\ \text{ s.t. }\{\begin{array}{c}{\sum }_{i=1}^n{\lambda }_i{y}_i=0\\ 0⩽{\lambda }_i⩽C,i=1,\cdots ,n\end{array} \end{gathered}$$
若令
$$
P=\begin{bmatrix}
y_1y_1K(x_1,x_1)&\cdots &y_1y_nK(x_1,x_n)\
\vdots & \ddots &\vdots \
y_ny_1K(x_n,x_1)&\cdots &y_ny_nK(x_n,x_n)
\end{bmatrix},\quad \lambda = \begin{pmatrix}
\lambda_1\
\vdots \\lambda_n

\end{pmatrix},
$$则有$$
\sum_{i=1}^{n} \sum_{j=1}^{n} \lambda_{i} \lambda_{j} y_{i} y_{j}K(x_{i},x_{j})=\lambda^TP\lambda.
KaTeX parse error: Can't use function '$' in math mode at position 12: 同样的，目标函数中的$̲-\sum_{i=1}^{n}…
-\sum_{i=1}^{n} \lambda_{i}=\begin{bmatrix}
-1 & & \
& \ddots & \
& & -1
\end{bmatrix}\lambda = Q\lambda,
$$此时，目标函数为$$
\frac{1}{2} \sum_{i=1}^{n} \sum_{j=1}^{n} \lambda_{i} \lambda_{j} y_{i} y_{j}K(x_{i}, x_{j})-\sum_{i=1}^{n} \lambda_{i}=\frac{1}{2}\lambda^TP\lambda+Q\lambda.
$$若令$$
G=\begin{bmatrix}1&&\&\ddots& \&&1 \ -1&&\&\ddots& \&&-1\end{bmatrix},\quad h=\begin{bmatrix}C\\vdots \C \ 0\\vdots\0 \end{bmatrix},
KaTeX parse error: Can't use function '$' in math mode at position 7: 其中，矩阵$̲G$中前$n$行为单位矩阵，后…
G\lambda\leqslant h,
$$同样的，令$$
A=\begin{bmatrix}y_1,y_2,\cdots,y_n\end{bmatrix},\
b = 0,
KaTeX parse error: Can't use function '$' in math mode at position 7: 则约束条件$̲\sum_{i=1}^{n} …
A\lambda=b.
$$综上所述，原二次规划问题可以化为以下形式$$
\min_{\lambda}\frac{1}{2} \lambda^TP\lambda+Q\lambda \
\text { s.t. }\left{\begin{array}{c}
A\lambda=b,\
G\lambda\leqslant h,
\end{array}\right.
$$a,$$
此时，目标函数为
$$\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n{\lambda }_i{\lambda }_j{y}_i{y}_{j}K(x_i,x_j)-\sum _{i=1}^n{\lambda }_i=\frac{1}{2}{\lambda }^{T}P\lambda +Q\lambda .$$
若令
$$G=\left[\begin{array}{ccc}1& & \\ & \ddots & \\ & & 1\\ -1& & \\ & \ddots & \\ & & -1\end{array}\right],h=\left[\begin{array}{c}C\\ ⋮\\ C\\ 0\\ ⋮\\ 0\end{array}\right],$$
其中，矩阵$G$中前$n$行为单位矩阵，后$n$行为负单位矩阵；向量$h$中前$n$行均为$C$，后$n$行均为$0$，此时原二次规划的约束条件$0⩽{\lambda }_i⩽C,i=1,\cdots ,n$可以化为
$$G\lambda ⩽h,$$
同样的，令
$$\begin{gathered} A=\left[\begin{array}{c}{y}_1,y_2,\cdots ,y_n\end{array}\right], \\ b=0, \end{gathered}$$
则约束条件${\sum }_{i=1}^n{\lambda }_i{y}_i=0$可以化为
$$A\lambda =b.$$
综上所述，原二次规划问题可以化为以下形式
$$\begin{gathered} \underset{\lambda }{\min }\frac{1}{2}{\lambda }^{T}P\lambda +Q\lambda \\ \text{ s.t. }\{\begin{array}{c}A\lambda =b,\\ G\lambda ⩽h,\end{array} \end{gathered}$$