参考: 支持向量机整理
SVM 硬边界的结果如下:
$$
min \quad \frac{1}{2} \sum_{i=1}^m\sum_{j=1}^m \alpha_i\alpha_jy_iy_j \vec x_i \vec x_j - \sum_{i=1}^m\alpha_i
\\
s.t. \quad \alpha_i\ge0 \quad i=1...m
\\
\quad \sum_{i=1}^m \alpha_i y_i=0
$$
一. 数据准备
测试数据如下所示, 前两个为 -1, 后面三个为 1, 如下图可以看到分割线即为:
$$
y = x + 1
$$
import numpy as np
import matplotlib.pyplot as plt X = np.array([[1,3],
[0,2],
[0,0],
[2,0],
[2,2]]) x = np.linspace(-2, 3, 100) y = np.array([-1,-1,1,1,1]) plt.figure()
plt.scatter(X[:2,0],X[:2,1])
plt.scatter(X[2:,0],X[2:,1])
plt.plot(x, x+1)
plt.show()
二. 获取 QP 的参数并计算α
将下面的结果带入到二次规划问题中分别求得 P/p/G/h/A/b 的值.
这个过程不是很容易, 虽然数据量这么少, 我反反复复弄了好几遍最终才做对.
$$
min \quad \frac{1}{2} \sum_{i=1}^m\sum_{j=1}^m \alpha_i\alpha_jy_iy_j \vec x_i \vec x_j - \sum_{i=1}^m\alpha_i
\\
s.t. \quad \alpha_i\ge0 \quad i=1...m
\\
\quad \sum_{i=1}^m \alpha_i y_i=0
$$
按照下面的形式进行获取参数.
$$
min \quad \frac{1}{2}x^TPx + q^Tx\\
s.t. \quad Gx \le h\\
\quad\quad Ax = b
$$
# 需要将 X 中的数据彼此相乘, 得到一个 5*5 的矩阵
# 同时需要注意 y 的符号会影响
P = matrix([[10.0,6.0,0.0,-2.0,-8.0],
[6.0,4.0,0.0,0.0,-4.0],
[0.0,0.0,0.0,0.0,0.0],
[-2.0,0.0,0.0,4.0,4.0],
[-8.0,-4.0,0.0,4.0,8.0]]) # 为了得到一个常数, q 为 5*1 的矩阵, 转置后正好可以用
q = matrix(-1.0, (5,1)) # 首先将 ≥ 调整为 ≤, 然后按照向量的形式表示
# 结果 h 为 5*1 的矩阵
# 因此 G 为 5*5 的矩阵(α 是 5*1 矩阵)
G = matrix([[-1.0,0.0,0.0,0.0,0.0],
[0.0,-1.0,0.0,0.0,0.0],
[0.0,0.0,-1.0,0.0,0.0],
[0.0,0.0,0.0,-1.0,0.0],
[0.0,0.0,0.0,0.0,-1.0]]) h = matrix(0.0, (5,1)) # 结果为常数的形式, 因此 A 是一个 1*5 矩阵
A = matrix([1.0,1.0,-1.0,-1.0,-1.0]).T b = matrix(0.0, (1,1))
将上面的内容带入到二次规划的函数中进行求解.
sol = solvers.qp(P,q,G,h,A,b)
alpha = sol['x']
print(alpha)
pcost dcost gap pres dres
0: -1.4151e+00 -3.0463e+00 1e+01 3e+00 2e+00
1: -3.4780e-01 -2.4147e+00 2e+00 7e-16 8e-16
2: -9.2882e-01 -1.0856e+00 2e-01 3e-16 5e-16
3: -9.9882e-01 -1.0010e+00 2e-03 2e-16 3e-16
4: -9.9999e-01 -1.0000e+00 2e-05 1e-16 3e-16
5: -1.0000e+00 -1.0000e+00 2e-07 2e-16 2e-16
Optimal solution found.
[ 4.31e-01]
[ 5.69e-01]
[ 2.84e-01]
[ 5.88e-08]
[ 7.16e-01]
三. 根据α来计算w
目前已经求出了所有的$\alpha$, 根据下面的公式将所有的样本点数据带入求得$\vec w$. 根据$\alpha$的结果可以判断哪些是支持向量, 包括 index = 0, 1, 2, 4 都满足.
$$
\vec w=\sum_{i=1}^m \alpha_i y_i \vec x_i
$$
X0 = X[:,0].flatten()
X1 = X[:,1].flatten() w1 = (w*y*X0).sum()
w2 = (w*y*X1).sum() W = np.array([w1,w2]) print("w1=", w1)
print("w2=", w2)
w1= 1.0000000446896518
w2= -1.000000054585139
四. 根据w来求b
$\vec w$已经求出了, 这时候只要带入任何一个支持向量里面即可, 公式如下:
$$
y_i(\vec w ^T\vec x_i+b) = 1
$$
化简后得到:
$$
b = y_i - \vec w^T\vec x_i
$$
由上面计算可知, 第一个点在支持向量上面, 因此可以计算获得b值.
$$
b = y_1 - \vec w^T\vec x_1
$$
W = np.mat(W)
xx = np.mat(X[1,:].T)
xx = xx.T b = int(-1 - W*xx)
print("b=", b)
b= 1
所以最终的结果就是:
$$
x_1 - x_2 + 1 = 0
$$
将$x_1$换成$x$, 将$x_2$换成$y$, 则得到:
$$
y = x + 1
$$