前言
- 在支持向量机中,理论逻辑很简单:最大化最小的几何间隔。但是实际编写代码过程中有一个小点需要注意。总是把二分类的类别分为
0
和1
,这样就导致我的目标函数跟算法描述的就不一样,所以求解结果就不正确。 - 同时还有第二个要注意的就是凸优化包
cvxpy
中各种运算的表示方法,比如凸优化中常见的二次方程的表示,变量的默认形状等,要查看官方文档才能熟悉。变量的默认形状为列向量。
参考:
分析
支持向量机算法中,我们的训练数据除了是两种类别以外,类别的编号也有要求,分别是1
、-1
,只有这样,我们才能求每个样本所对应的函数间隔
D
=
y
i
(
w
x
i
+
b
)
D=y_i(wx_i+b)
D=yi(wxi+b)和几何间隔
D
=
y
i
w
x
i
+
b
∥
w
∥
2
D=y_i\frac{wx_i+b}{\|w\|_2}
D=yi∥w∥2wxi+b,在这种类别标签的情况下,预测值
w
x
i
+
b
wx_i+b
wxi+b与真实值
y
i
y_i
yi之间的乘积才有意义。
预测与真实相同,乘积才会是大于0的;预测与真实相反,乘积就是小于0的。只有这样后面的目标函数最大化几何间隔才有意义。
重新表述一下可分数据集上支持向量机的目标函数和约束条件:
max
λ
s.t.
y
i
w
x
i
+
b
∥
w
∥
2
≥
λ
\begin{align*} &\textbf{max}&\lambda \\ &\textbf{s.t.} &y_i\frac{wx_i+b}{\|w\|_2}\ge\lambda \end{align*}
maxs.t.λyi∥w∥2wxi+b≥λ
如果我们使用类别标签为0
,1
,那么当错误分类时,几何间隔为0
,无法指导参数修改。所以必须要使用1
,-1
。
数据集线性可分情况下的支持向量机
此时有两种求法,一种是用原始算法,直接用cvxpy
函数包求解原始的这个凸优化问题,并把问题变为下述形式:
max
t
∥
w
∥
2
s.t.
y
i
(
w
x
i
+
b
)
≥
t
\begin{align*} &\textbf{max}&\frac{t}{\|w\|_2} \\ &\textbf{s.t.} &y_i{(wx_i+b)}\ge t \end{align*}
maxs.t.∥w∥2tyi(wxi+b)≥t
由于通过同比例放大w,b
可以实现条件中左边的乘积大小的任意变换,所以我们修改t
为
1
1
1。上述凸优化问题就变为:
min
∥
w
∥
2
s.t.
y
i
(
w
x
i
+
b
)
≥
1
\begin{align*} &\textbf{min}&{\|w\|_2} \\ &\textbf{s.t.} &y_i{(wx_i+b)}\ge 1 \end{align*}
mins.t.∥w∥2yi(wxi+b)≥1
第二种方法就是将原始问题使用拉格朗日乘子法变换为对偶问题,将加入条件和拉格朗日乘子的拉格朗日函数进行求导,并将求导得到的关系式带入拉格朗日函数,这样就可以得到对偶问题。
原始问题凸优化包解法
import numpy as np
import cvxpy as cp
#生成符合要求的样本数据
np.random.seed(3)
train_x=np.random.randn(2,2)
weight_x=np.random.randn(2)
bias_x=np.random.randn(1)
train_y=np.where(train_x@weight_x+bias_x<0,-1,1).reshape(2,-1)
print(train_x,train_y)
#求解对偶凸优化问题
w=cp.Variable(2)
b=cp.Variable(1)
obj=cp.Minimize(1/2*cp.sum_squares(w))
cons=[train_y[0]*(train_x@w+b)[0]>=1,train_y[1]*(train_x@w+b)[1]>=1]
prob=cp.Problem(obj,cons)
result=prob.solve()
#输出拉格朗日乘子的和最优化的目标函数值以及权重w
print(f'w.value,result,b.value:{w.value,result,b.value}')
结果如下:
[[ 1.78862847 0.43650985]
[ 0.09649747 -1.8634927 ]] [[-1]
[ 1]]
w.value,result,b.value:(array([-0.41507783, -0.56418804]), 0.24529887505030906, array([-0.01130633]))
对偶问题凸优化包解法
本例子中我们使用了小规模的数据,只有两个样本,所以这两个样本肯定都是支持向量,也就是对应的拉格朗日乘子都不为0,对于大规模样本数据的情况,如果不在分界面上,那么对应的拉格朗日乘子为0,也就不是支持向量。拉格朗日乘子不为0的肯定就是支持向量。
import numpy as np
import cvxpy as cp
#生成符合要求的样本数据
np.random.seed(3)
train_x=np.random.randn(2,2)
weight_x=np.random.randn(2)
bias_x=np.random.randn(1)
train_y=np.where(train_x@weight_x+bias_x<0,-1,1).reshape(2,-1)
print(train_x,train_y)
#求解对偶凸优化问题
alpha=cp.Variable(2)
obj=cp.Minimize(1/2*cp.quad_form(alpha,(train_x@train_x.T)*(train_y@train_y.T))-cp.sum(alpha))
cons=[alpha>=0,train_y.T@alpha>=0,train_y.T@alpha<=0]
prob=cp.Problem(obj,cons)
result=prob.solve()
#输出拉格朗日乘子的和最优化的目标函数值以及权重w
print(f'alpha.value,result,w:{alpha.value,result,np.array(alpha.value.reshape(2,-1)*train_y*train_x).sum(axis=0)}')
#检验支持向量机求出的分离面是否与这两个样本之间的连线垂直
w=np.array(alpha.value.reshape(2,-1)*train_y*train_x).sum(axis=0)
np.array(train_x[0,:]-train_x[1,:]).reshape(1,-1)@np.array([-w[1]/w[0],1]).reshape(2,-1)
结果如下:
[[ 1.78862847 0.43650985]
[ 0.09649747 -1.8634927 ]] [[-1]
[ 1]]
alpha.value,result,w:(array([0.24529888, 0.24529888]), -0.24529887505030898, array([-0.41507783, -0.56418804]))
array([[-4.4408921e-16]])
至于参数b,我们可以通过支持向量所对应的等式求出。
我们可以观察一下原始问题与对偶问题的解答是否一致。
数据集线性不可分情况下的线性支持向量机与软间隔最大化
软间隔顾名思义就是给原来的间隔留下一点宽容量,给那些不容易分正确的留一点余地。同时对于这些余地进行惩罚所得到的分割面,以上面线性可分的数据做演示。
import numpy as np
import cvxpy as cp
#生成符合要求的样本数据
np.random.seed(3)
train_x=np.random.randn(2,2)
weight_x=np.random.randn(2)
bias_x=np.random.randn(1)
train_y=np.where(train_x@weight_x+bias_x<0,-1,1).reshape(2,-1)
print(train_x,train_y)
#求解对偶凸优化问题
w=cp.Variable(2)
b=cp.Variable(1)
softgap=cp.Variable(2)
obj=cp.Minimize(1/2*cp.sum_squares(w)+100*cp.sum(softgap))
cons=[train_y[0]*(train_x@w+b)[0]>=1-softgap[0],train_y[1]*(train_x@w+b)[1]>=1-softgap[1],softgap>=0]
prob=cp.Problem(obj,cons)
result=prob.solve()
#输出拉格朗日乘子的和最优化的目标函数值以及权重w
print(f'w.value,result,b.value,softgap.value:{w.value,result,b.value,softgap.value}')
结果如下:
[[ 1.78862847 0.43650985]
[ 0.09649747 -1.8634927 ]] [[-1]
[ 1]]
w.value,result,b.value,softgap.value:(array([-0.41507783, -0.56418804]), 0.24529887505030906, array([-0.01130633]), array([-2.90746355e-22, 1.53881391e-22]))
可以看出在线性可分的的情况下软间隔不起作用。
那么我们制造一些线性不可分的数据,来测试一下。
import numpy as np
import cvxpy as cp
#生成符合要求的样本数据
np.random.seed(3)
train_x=np.array([[0,0],[1,0],[2,0]])
train_y=np.array([-1,1,-1]).reshape(3,-1)
print(train_x,train_y)
#求解对偶凸优化问题
w=cp.Variable(2)
b=cp.Variable(1)
softgap=cp.Variable(3)
obj=cp.Minimize(1/2*cp.sum_squares(w)+0.1*cp.sum(softgap))
cons=[train_y[0]*(train_x@w+b)[0]>=1-softgap[0],train_y[1]*(train_x@w+b)[1]>=1-softgap[1],train_y[2]*(train_x@w+b)[2]>=1-softgap[2],softgap>=0]
prob=cp.Problem(obj,cons)
result=prob.solve()
#输出拉格朗日乘子的和最优化的目标函数值以及权重w
print(f'w.value,result,b.value,softgap.value:{w.value,result,b.value,softgap.value}')
结果如下:
[[0 0]
[1 0]
[2 0]] [[-1]
[ 1]
[-1]]
w.value,result,b.value,softgap.value:(array([9.07653476e-18, 0.00000000e+00]), 0.2, array([-1.]), array([ 8.59013373e-23, 2.00000000e+00, -8.59013423e-23]))