浏览本文须要的背景知识点:硬距离反对向量机、松弛变量、一丢丢编程常识
一、引言
后面一节咱们介绍了一种最根底的反对向量机模型——硬距离反对向量机,该模型能对线性可分的数据集进行分类,但事实中的数据集往往是线性不可分的,那么这一节咱们来介绍反对向量机中的第二种——软距离反对向量机1(Soft-margin Support Vector Machine),来解决下面说的数据集线性不可分的问题。
二、模型介绍
原始模型
咱们先来看下硬距离反对向量机的原始模型如下:
$$\begin{array}{}\underset{w, b}{\max} \frac{1}{2} w^Tw \\\text { s.t. } \quad y_{i}\left(w^{T} x_{i}+b\right) \geq 1 \quad i=1,2, \cdots, N\end{array}$$
<center>式2-1</center>
约束条件 y(wx + b) 大于等于 1 是为了使得所有样本点都在正确的分类下,这也是为什么称为硬距离的起因。而当初数据集无奈用一个超平面齐全离开,这时就须要容许局部数据集不满足上述约束条件。
批改下面的代价函数,加上分类谬误的样本点的惩办项,其中 1(x) 在后面章节中提到过,即批示函数2(indicator function),当满足上面不等式时函数返回 1,不满足时函数返回 0。同时通过常数 C 来管制惩办力度,其中 C > 0,能够看到当 C 取无限大时,会迫使每一个样本点分类正确,等同于硬距离反对向量机。
$$\begin{array}{}\underset{w, b}{\min} \frac{1}{2} w^{T} w+C \sum_{i=1}^{N} 1_{x_{i}}\left(y_{i}\left(w^{T} x_{i}+b\right) \lt 1\right) \\\text { s.t. } \quad C>0\end{array}$$
<center>式2-2</center>
式2-2中的批示函数既不是凸函数又不是连续函数,解决起来比拟麻烦,这时能够用max函数来替换,模式如下:
$$\begin{array}{}\underset{w, b}{\min} \frac{1}{2} w^{T} w+C \sum_{i=1}^{N} \max\left(0, 1 - y_{i}\left(w^{T} x_{i}+b\right)\right) \\\text { s.t. } \quad C>0\end{array}$$
<center>式2-3</center>
这时再引入松弛变量3(slack variable),同时加上如下的条件,进一步将 max 函数去掉,式2-3就是软距离反对向量机的原始模型:
$$\begin{array}{}\underset{w, b, \xi}{\min} \frac{1}{2} w^{T} w+C \sum_{i=1}^{N} \xi_{i} \\\text { s.t. } \quad C>0 \quad \xi_{i} \geq 0 \quad \xi_{i} \geq 1-y_{i}\left(w^{T} x_{i}+b\right) \quad i=1,2, \cdots, N\end{array}$$
<center>式2-4</center>
对偶模型
与硬距离反对向量机的解决形式一样,同样对原始模型用拉格朗日乘子法进行转换:
(1)应用拉格朗日乘子法,引入两类拉格朗日参数 、,失去拉格朗日函数
(2)拉格朗日函数对 w 求偏导数并令其为零
(3)同硬距离反对向量机一样失去 w 的解析解
(4)拉格朗日函数对 b 求偏导数并令其为零
(5)拉格朗日函数对 求偏导数并令其为零
(6)失去 C = +
$$\begin{aligned}L(w, b, \xi, \lambda, \mu) &=\frac{1}{2} w^{T} w+C \sum_{i=1}^{N} \xi_{i}+\sum_{i=1}^{N} \lambda_{i}\left(1-\xi_{i}-y_{i}\left(w^{T} x_{i}+b\right)\right)-\sum_{i=1}^{N} \mu_{i} \xi_{i} & (1)\\\frac{\partial L}{\partial w} &=w-\sum_{i=1}^{N} \lambda_{i} y_{i} x_{i}=0 & (2)\\w &=\sum_{i=1}^{N} \lambda_{i} y_{i} x_{i} & (3)\\\frac{\partial L}{\partial b} &=\sum_{i=1}^{N} \lambda_{i} y_{i}=0 & (4)\\\frac{\partial L}{\partial \xi_{i}} &=C-\lambda_{i}-\mu_{i}=0 & (5)\\C &=\lambda_{i}+\mu_{i} & (6)\end{aligned}$$
<center>式2-5</center>
将式2-5失去的后果带回拉格朗日函数中:
(1)拉格朗日函数
(2)带入后开展括号
(3)能够看到(2)式中第2、5项互相对消,第3、8项互相对消,第7项为零,合并后得
$$\begin{aligned}L(w, b, \xi, \lambda, \mu) &=\frac{1}{2} w^{T} w+C \sum_{i=1}^{N} \xi_{i}+\sum_{i=1}^{N} \lambda_{i}\left(1-\xi_{i}-y_{i}\left(w^{T} x_{i}+b\right)\right)-\sum_{i=1}^{N} \mu_{i} \xi_{i} & (1)\\&=\frac{1}{2} \sum_{i=1}^{N} \sum_{j=1}^{N} \lambda_{i} \lambda_{j} y_{i} y_{j} x_{i}^{T} x_{i}+\sum_{i=1}^{N} \lambda_{i} \xi_{i}+\sum_{i=1}^{N} \mu_{i} \xi_{i}+\sum_{i=1}^{N} \lambda_{i}-\sum_{i=1}^{N} \lambda_{i} \xi_{i}-\sum_{i=1}^{N} \sum_{j=1}^{N} \lambda_{i} \lambda_{j} y_{i} y_{j} x_{i}^{T} x_{i}-\sum_{i=1}^{N} \lambda_{i} y_{i} b-\sum_{i=1}^{N} \mu_{i} \xi_{i} & (2)\\&=\sum_{i=1}^{N} \lambda_{i}-\frac{1}{2} \sum_{i=1}^{N} \sum_{j=1}^{N} \lambda_{i} \lambda_{j} y_{i} y_{j} x_{i}^{T} x_{i} & (3)\end{aligned}$$
<center>式2-6</center>
在来看下拉格朗日乘子参数的条件:
(1)拉格朗日乘子 的限度条件
(2)两边同时加上拉格朗日乘子
(3)联合式2-5(6)的后果失去
$$\begin{aligned}\mu_i &\ge 0 & (1) \\\lambda_i + \mu_i &\ge \lambda_i & (2) \\C &\ge \lambda_i & (3)\end{aligned}$$
<center>式2-7</center>
同硬距离反对向量机一样,引入 KKT 条件如下:
$$\left\{\begin{aligned}\nabla_{w, b, \xi} L(w, b, \xi, \lambda, \mu) &=0 \\\lambda_{i} & \geq 0 \\y_{i}\left(w^{T} x_{i}+b\right)-1+\xi_{i} & \geq 0 \\\lambda_{i}\left(y_{i}\left(w^{T} x_{i}+b\right)-1+\xi_{i}\right) &=0 \\\mu_{i} & \geq 0 \\\xi_{i} & \geq 0 \\\mu_{i} \xi_{i} &=0\end{aligned}\right.$$
<center>式2-8</center>
通过使用拉格朗日乘子法后,失去了原模型的拉格朗日对偶模型并且须要满足如上所示的 KKT 条件:
$$\begin{array}{}\underset{\lambda}{\max} \sum_{i=1}^{N} \lambda_{i}-\frac{1}{2} \sum_{i=1}^{N} \sum_{j=1}^{N} \lambda_{i} \lambda_{j} y_{i} y_{j} x_{i}^{T} x_{j} \\\text { s.t. } \quad \sum_{i=1}^{N} \lambda_{i} y_{i}=0 \quad 0 \leq \lambda_{i} \leq C \quad i=1,2, \cdots, N\end{array}$$
<center>式2-9</center>
三、算法步骤
将下面式2-9的对偶模型与硬距离反对向量机的对偶模型互相比拟后,会发现两者的差异就在于拉格朗日乘子的约束条件不同,前者只多了一个惩办因子 C 的上界,所以其求解的算法与硬距离反对向量机的大致相同,只有如下两个中央不一样。
算法中一个有改变的中央是新 _b 的上下界的计算:
(1)、(2)所有 都须要大于等于零并且小于等于 C
(3)分状况探讨 _b 的限度条件
(4)综合(1)式失去最终变量 _b 的限度条件
$$\begin{aligned}&0 \leq \lambda_{b}^{\text {new }} \leq C & (1)\\&0 \leq \lambda_{a}^{\text {old }}+y_{a} y_{b}\left(\lambda_{b}^{\text {old }}-\lambda_{b}^{\text {new }}\right) \leq C & (2)\\\Rightarrow&\left\{\begin{array}{}\lambda_{a}^{\text {old }}+\lambda_{b}^{\text {old }}-C \leq \lambda_{b}^{\text {new }} \leq \lambda_{a}^{\text {old }}+\lambda_{b}^{\text {old }}, \quad y_{a} y_{b}=+1 \\\lambda_{b}^{\text {old }}-\lambda_{a}^{\text {old }} \leq \lambda_{b}^{\text {new }} \leq C-\lambda_{a}^{\text {old }}+\lambda_{b}^{\text {old }}, \quad y_{a} y_{b}=-1\end{array}\right. & (3)\\\Rightarrow&\left\{\begin{array}{}\max \left(0, \lambda_{a}^{\text {old }}+\lambda_{b}^{\text {old }}-C\right) \leq \lambda_{b}^{\text {new }} \leq \min \left(C, \lambda_{a}^{\text {old }}+\lambda_{b}^{\text {old }}\right), \quad y_{a} y_{b}=+1 \\\max \left(0, \lambda_{b}^{\text {old }}-\lambda_{a}^{\text {old }}\right) \leq \lambda_{b}^{\text {new }} \leq \min \left(C, C-\lambda_{a}^{\text {old }}+\lambda_{b}^{\text {old }}\right), \quad y_{a} y_{b}=-1\end{array}\right. & (4)\end{aligned}$$
<center>式3-1</center>
算法中另一个有改变的中央是 KKT 条件的判断:
(1)思考 等于 0 的状况
(2)依据式2-5 中(6)的后果可知 等于 C
(3)因为 等于 C,依据式2-8的 KKT 条件可知 等于 0
(4)将 等于 0 带入 KKT 条件得
$$\begin{aligned}\lambda &=0 & (1)\\\mu &=C & (2)\\\xi &=0 & (3)\\y_{i}\left(w^{T} x_{i}+b\right)-1 & \geq 0 & (4)\end{aligned}$$
<center>式3-2</center>
(1)思考 等于 C 的状况
(2)依据式2-5 中(6)的后果可知 等于 0
(3)因为 不等于 0,依据式2-8的 KKT 条件可知
(4)因为 大于等于 0 联合(3)式可得
$$\begin{aligned}\lambda &=C & (1)\\\mu &=0 & (2)\\y_{i}\left(w^{T} x_{i}+b\right)-1 + \xi &=0 & (3)\\y_{i}\left(w^{T} x_{i}+b\right)-1 & \leq 0 & (4)\end{aligned}$$
<center>式3-3</center>
(1)思考 介于 0 与 C 之间的状况
(2)依据式2-5 中(6)的后果可知 也介于 0 与 C 之间
(3)依据式2-8的 KKT 条件可知
(4)因为 不等于 0,依据式2-8的 KKT 条件并联合(3)式可得
$$\begin{aligned}0 \lt \lambda &\lt C & (1)\\0 \lt \mu &\lt C & (2)\\\xi &=0 & (3)\\y_{i}\left(w^{T} x_{i}+b\right)-1 & = 0 & (4)\end{aligned}$$
<center>式3-4</center>
综合下面式3-2、式3-3、式3-4可得如下的新限度条件:
$$\left\{\begin{array}{}y_{i}\left(w^{T} x_{i}+b\right)-1 \geq 0, & \lambda=0 \\y_{i}\left(w^{T} x_{i}+b\right)-1=0, & 0 \lt \lambda \lt C \\y_{i}\left(w^{T} x_{i}+b\right)-1 \leq 0, & \lambda=C\end{array}\right.$$
<center>式3-5</center>
算法步骤请参考硬距离反对向量机的第三小节和上面的代码实现。
四、代码实现
应用 Python 实现
import numpy as npclass SMO: """ 软距离反对向量机 序列最小优化算法实现(Sequential minimal optimization/SMO) """ def __init__(self, X, y): # 训练样本特色矩阵(N * p) self.X = X # 训练样本标签向量(N * 1) self.y = y # 拉格朗日乘子向量(N * 1) self.alpha = np.zeros(X.shape[0]) # 误差向量,默认为负的标签向量(N * 1) self.errors = -y # 偏移量 self.b = 0 # 权重向量(p * 1) self.w = np.zeros(X.shape[1]) # 代价值 self.cost = -np.inf def fit(self, C = 1, tol = 1e-4): """ 算法来自 John C. Platt 的论文 https://www.microsoft.com/en-us/research/uploads/prod/1998/04/sequential-minimal-optimization.pdf """ # 更新变动次数 numChanged = 0 # 是否查看全副 examineAll = True while numChanged > 0 or examineAll: numChanged = 0 if examineAll: for idx in range(X.shape[0]): numChanged += self.update(idx, C) else: for idx in range(X.shape[0]): if self.alpha[idx] <= 0: continue numChanged += self.update(idx, C) if examineAll: examineAll = False elif numChanged == 0: examineAll = True # 计算代价值 cost = self.calcCost() if self.cost > 0: # 当代价值的变动小于答应的范畴时完结算法 rate = (cost - self.cost) / self.cost if rate <= tol: break self.cost = cost def update(self, idx, C = 1): """ 对下标为 idx 的拉格朗日乘子进行更新 """ X = self.X y = self.y alpha = self.alpha # 查看以后拉格朗日乘子是否满足KKT条件,满足条件则间接返回 0 if self.checkKKT(idx, C): return 0 if len(alpha[(alpha != 0)]) > 1: # 依照|E1 - E2|最大的准则寻找第二个待优化的拉格朗日乘子的下标 jdx = self.selectJdx(idx) # 对下标为 idx、jdx 的拉格朗日乘子进行更新,当胜利更新时间接返回 1 if self.updateAlpha(idx, jdx, C): return 1 # 当未更新胜利时,遍历不为零的拉格朗日乘子进行更新 for jdx in range(X.shape[0]): if alpha[jdx] == 0: continue # 对下标为 idx、jdx 的拉格朗日乘子进行更新,当胜利更新时间接返回 1 if self.updateAlpha(idx, jdx, C): return 1 # 当仍然没有没有更新胜利时,遍历为零的拉格朗日乘子进行更新 for jdx in range(X.shape[0]): if alpha[jdx] != 0: continue # 对下标为 idx、jdx 的拉格朗日乘子进行更新,当胜利更新时间接返回 1 if self.updateAlpha(idx, jdx, C): return 1 # 仍然没有更新时返回 0 return 0 def selectJdx(self, idx): """ 寻找第二个待优化的拉格朗日乘子的下标 """ errors = self.errors if errors[idx] > 0: # 当误差项大于零时,抉择误差向量中最小值的下标 return np.argmin(errors) elif errors[idx] < 0: # 当误差项小于零时,抉择误差向量中最大值的下标 return np.argmax(errors) else: # 当误差项等于零时,抉择误差向量中最大值和最小值的绝对值最大的下标 minJdx = np.argmin(errors) maxJdx = np.argmax(errors) if max(np.abs(errors[minJdx]), np.abs(errors[maxJdx])) - errors[minJdx]: return minJdx else: return maxJdx def calcB(self): """ 计算偏移量 别离计算每一个拉格朗日乘子不为零对应的偏移量后取其平均值 """ X = self.X y = self.y alpha = self.alpha alpha_gt = alpha[alpha > 0] # 拉格朗日乘子向量中不为零的数量 alpha_gt_len = len(alpha_gt) # 全副为零时间接返回 0 if alpha_gt_len == 0: return 0 # b = y - Wx,具体算法请参考文章中的阐明 X_gt = X[alpha > 0] y_gt = y[alpha > 0] alpha_gt_y = np.array(np.multiply(alpha_gt, y_gt)).reshape(-1, 1) s = np.sum(np.multiply(alpha_gt_y, X_gt), axis=0) return np.sum(y_gt - X_gt.dot(s)) / alpha_gt_len def calcCost(self): """ 计算代价值 依照文章中的算法计算即可 """ X = self.X y = self.y alpha = self.alpha cost = 0 for idx in range(X.shape[0]): for jdx in range(X.shape[0]): cost = cost + (y[idx] * y[jdx] * X[idx].dot(X[jdx]) * alpha[idx] * alpha[jdx]) return np.sum(alpha) - cost / 2 def checkKKT(self, idx, C = 1): """ 查看下标为 idx 的拉格朗日乘子是否满足 KKT 条件 1. alpha >= 0 2. alpha <= C 3. y * f(x) - 1 >= 0 4. alpha * (y * f(x) - 1) = 0 """ y = self.y errors = self.errors alpha = self.alpha r = errors[idx] * y[idx] if (alpha[idx] > 0 and alpha[idx] < C and r == 0) or (alpha[idx] == 0 and r > 0) or (alpha[idx] == C and r < 0): return True return False def calcE(self): """ 计算误差向量 E = f(x) - y """ X = self.X y = self.y alpha = self.alpha alpha_y = np.array(np.multiply(alpha, y)).reshape(-1, 1) errors = X.dot(X.T).dot(alpha_y).T[0] + self.b - y return errors def calcU(self, idx, jdx, C = 1): """ 计算拉格朗日乘子上界,使两个待优化的拉格朗日乘子同时大于等于0 依照文章中的算法计算即可 """ y = self.y alpha = self.alpha if y[idx] * y[jdx] == 1: return max(0, alpha[jdx] + alpha[idx] - C) else: return max(0.0, alpha[jdx] - alpha[idx]) def calcV(self, idx, jdx, C = 1): """ 计算拉格朗日乘子下界,使两个待优化的拉格朗日乘子同时大于等于0 依照文章中的算法计算即可 """ y = self.y alpha = self.alpha if y[idx] * y[jdx] == 1: return min(C, alpha[jdx] + alpha[idx]) else: return min(C, C + alpha[jdx] - alpha[idx]) def updateAlpha(self, idx, jdx, C = 1): """ 对下标为 idx、jdx 的拉格朗日乘子进行更新 依照文章中的算法计算即可 """ if idx == jdx: return False X = self.X y = self.y alpha = self.alpha errors = self.errors # idx 的误差项 Ei = errors[idx] # jdx 的误差项 Ej = errors[jdx] Kii = X[idx].dot(X[idx]) Kjj = X[jdx].dot(X[jdx]) Kij = X[idx].dot(X[jdx]) # 计算 K K = Kii + Kjj - 2 * Kij oldAlphaIdx = alpha[idx] oldAlphaJdx = alpha[jdx] # 计算 jdx 的新拉格朗日乘子 newAlphaJdx = oldAlphaJdx + y[jdx] * (Ei - Ej) / K U = self.calcU(idx, jdx, C) V = self.calcV(idx, jdx, C) if newAlphaJdx < U: # 当新值超过上界时,批改其为上界 newAlphaJdx = U if newAlphaJdx > V: # 当新值低于下界时,批改其为下界 newAlphaJdx = V if oldAlphaJdx == newAlphaJdx: # 当新值与旧值相等时,判断为未更新,间接返回 return False # 计算 idx 的新拉格朗日乘子 newAlphaIdx = oldAlphaIdx + y[idx] * y[jdx] * (oldAlphaJdx - newAlphaJdx) # 更新权重向量 self.w = self.w + y[idx] * (newAlphaIdx - oldAlphaIdx) * X[idx] + y[jdx] * (newAlphaJdx - oldAlphaJdx) * X[jdx] # 更新拉格朗日乘子向量 alpha[idx] = newAlphaIdx alpha[jdx] = newAlphaJdx oldB = self.b # 从新计算偏移量 self.b = self.calcB() # 从新计算误差向量 eps = np.finfo(np.float32).eps newErrors = [] for i in range(X.shape[0]): newError = errors[i] + y[idx] * (newAlphaIdx - oldAlphaIdx) * X[idx].dot(X[i]) + y[jdx] * (newAlphaJdx - oldAlphaJdx) * X[jdx].dot(X[i]) - oldB + self.b if np.abs(newError) < eps: newError = 0 newErrors.append(newError) self.errors = newErrors return True
五、第三方库实现
scikit-learn4 实现
from sklearn.svm import SVCsvc = SVC(kernel = "linear", C = 1)# 拟合数据svc.fit(X, y)# 权重系数w = svc.coef_# 截距b = svc.intercept_print("w", w, "b", b)
六、动画演示
下图展现了一个线性不可分的数据集,其中红色示意标签值为 -1 的样本点,蓝色代表标签值为 1 的样本点:
<center>图6-1</center>
下图展现了通过软距离反对向量机拟合数据后的后果,其中浅红色的区域为预测值为 -1 的局部,浅蓝色的区域则为预测值为 1 的局部:
<center>图6-2</center>
七、思维导图
<center>图7-1</center>
八、参考文献
- https://en.wikipedia.org/wiki...
- https://en.wikipedia.org/wiki...
- https://en.wikipedia.org/wiki...
- https://scikit-learn.org/stab...
残缺演示请点击这里
注:本文力求精确并通俗易懂,但因为笔者也是初学者,程度无限,如文中存在谬误或脱漏之处,恳请读者通过留言的形式批评指正
本文首发于——AI导图,欢送关注