机器学习-SVM（支持向量机）

推荐课程：【机器学习实战】第5期支持向量机 |数据分析|机器学习|算法|菊安酱_哔哩哔哩_bilibili

赞美菊神ヾ ( ゜ⅴ゜)ﾉ

一、什么是支持向量机？

支持向量机（Support Vector Machine, SVM）是一类按监督学习（supervised learning）方式对数据进行二元分类的广义线性分类器（generalized linear classifier），其决策边界是对学习样本求解的最大边距超平面（maximum-margin hyperplane）。

1.1 举个例子

在桌子上似乎有规律地放了两种颜色的球，要求你用一根棍子分离开他们，并且尽量再放更多的球之后，仍然适用。

SVM就是试图把棍放在最佳位置，好让在棍的两边有尽可能大的间隙（这被认为是最佳求解）。

并且，现在即使再放入更多的球，棍子仍然是一个很好的分界线。

但是，现在增加难度，如果将球散乱地放在桌子上，又该怎样进行划分呢？很明显，此时在二维平面中，这变成了一个线性不可分的问题，我们没有方法再用一根棍子将这些球分开了，那么怎么解决这样一个问题呢？

解决方法也很简单，我们可以使用一个核函数，将二维平面中的 '小球' 投影到三维空间，也许就可以在三维空间中，有可能找到这样一个平面将其分隔开来。（可以想象一下，用力拍向桌子，然后桌子上的球就被震到空中，瞬间抓起一张纸，插到两种球的中间。）

（话说，如果3维空间依旧找不到这样一个平面呢？没关系，我继续投四维呀╮(๑•́ ₃•̀๑)╭，总能找到一个维度解决线性不可分的问题。）

可以通过视频，更为直观地感受一下这个过程：支持向量机 SVM在线性不可分情况下进行分类可视化直观展示高维空间映射_哔哩哔哩_bilibili

而在二维平面的角度看这些球，这些球像是被一条曲线分开了。

更为规范的，我们把这些球叫做「data」，把棍子叫做「classififier」, 最大间隙 trick 叫做「optimization」，拍桌子叫做「kernelling」, 那张纸叫做「hyperplane」。

1.2 概述一下

当一个分类问题，数据是线性可分（linearly separable）的，也就是用一根棍就可以将两种小球分开的时候，我们只要将棍的位置放在让小球距离棍的距离最大化的位置即可，寻找这个最大间隔的过程，就叫做最优化。但是，一般的数据是线性不可分的，也就是找不到一个棍将两种小球很好的类。这个时候，我们就需要使用核函数 （kernel）将小球投影到多维空间（想象一下，将小球拍飞到空中，瞬间抓起一张纸，插到两种球的中间），而在多维空间中切分小球的平面，就是超平面（hyperplane）。如果数据集是N维的，那么超平面就是N-1维的。

Q：什么是支持向量？

A：把一个数据集正确分开的超平面可能有多个，而那个具有“最大间隔”的超平面就是SVM要寻找的最优解。而这个真正的最优解对应的两侧虚线所穿过的样本点，就是SVM中的支持样本点，称为“支持向量（support vector）”。支持向量到超平面的距离被称为间隔（margin）。

二、线性支持向量机

一个最优化问题通常有两个最基本的因素：

1）目标函数，也就是你希望什么东西的什么指标达到最好；

2）优化对象，你期望通过改变哪些因素来使你的目标函数达到最优。

在线性SVM算法中，目标函数显然就是那个“间隔”，而优化对象则是超平面。

2.1 超平面方程

在线性可分的二分类问题中，超平面其实就是一条直线，公式如下：

将原来的x轴看作x1轴，y轴看作x2轴，于是公式(1)转换为：

向量形式可以写成：

进一步可表示为：

其中，，b为截距。

2.2 间隔的计算公式

间隔的计算公式 = 点到直线距离推导公式，向量表示法：

其中， $\left \| \omega \right \|$ 是向量 $\omega$ 的模，假如，则，表示在空间中向量的长度。

我们评价分类器的好坏依据是分类间隔的的大小，即分类间隔W越大，我们认为这个超平面的分类效果越好。而追求分类间隔W的最大化也就是寻找 d的最大化。

2.3 约束条件

对于SVM的线性可分的二分类问题，我们希望：

1）超平面够将所有的样本点都正确分类。
2）超平面的位置应该是在间隔区域的中轴线上。
3）对于一个给定的超平面，我们能找到对应的支持向量，来计算距离d。

对于上述的三个约束条件，我们假设在平面空间中有红蓝两种点，对其分别标记为：

红色为正样本，标记为+1；
蓝色为负样本，标记为-1。

对每个样本点 $x_i$ 加上类别标签 $y_i$ ，则有

如果我们的超平面能够完全将红蓝两种样本点分离开（对应第一个约束条件），那么则有

如果要求在再高一点，假设超平面正好处于间隔区域的中轴线上（对应第二个约束条件），并且相应支持向量到超平面的距离为d（对应第三个约束条件），则公式可进一步写为：

对公式两边同时除以d，可得：

其中，。

简化方程，使用 $\omega$ 和 $d$ 代替 $\frac{\omega^T _d}{\left \| \omega _d \right \|}$ ， $\frac{b_d}{\left \| \omega _d \right \|}$ ，可得：

当等于 1或 -1时，对应的 $x_i$ 坐标就是支持向量的坐标。

以将上述方程糅合成一个约束方程：

这个方程就是SVM最优化问题的约束条件。

2.4 线性SVM最优化问题

当等于 1或 -1时，对应的 $x_i$ 坐标就是支持向量的坐标，相应的根据约束条件，我们可以知道只有超平面最大化时，才会等于 1或 -1。无论是等于1还是-1，即 $x_i$ 为支持向量坐标时，对于公式10来说，都有。

所以对于支持向量来说：

我们原来的任务是找到一组参数 $\omega,b$ 使得分类间隔 $W=2d$ 最大化，根据公式12 就可以转变为 $\left \| \omega \right \|$ 的最小化问题，也等效于 $\frac{1}{2}\left \| \omega \right \|^2$ 的最小化问题。我们之所以要在 $\left \| \omega \right \|$ 上加上平方和1/2的系数，是为了以后进行最优化的过程中对目标函数求导时比较方便，但这绝不影响最优化问题最后的解。

所以，线性SVM最优化问题的数学描述就是：

公式13，描述的是一个典型的不等式约束条件下的二次型函数优化问题，同时也是支持向量机的基本数学模型。

2.5 求解线性SVM最优化问题

有不等式约束条件的原始目标函数求解十分困难，我们需要先转化没有约束条件的新目标函数，再进行求导优化，简单来讲过程就是：构造拉格朗日函数 + 求解拉格朗日对偶函数。

最终得到新的目标函数和约束条件，方便使用高效优化算法SMO算法：

三、SMO算法

SMO算法就是序列最小优化（Sequential Minimal Optimization），它是由 John Platt 于1996年发布的专门用于训练SVM的一个强大算法（如同，梯度下降算法专门用于训练深度学习模型一般）。

简化版SMO算法，如下：

def smoSimple(dataMatIn, classLabels, C, toler, maxIter):"""smoSimpleArgs:dataMatIn    数据集classLabels  类别标签C   松弛变量(常量值)，允许有些数据点可以处于分隔面的错误一侧。控制最大化间隔和保证大部分的函数间隔小于1.0这两个目标的权重。可以通过调节该参数达到不同的结果。toler   容错率（是指在某个体系中能减小一些因素或选择对某个系统产生不稳定的概率。）maxIter 退出前最大的循环次数Returns:b       模型的常量值alphas  拉格朗日乘子"""dataMatrix = mat(dataMatIn)# 矩阵转置 和 .T 一样的功能labelMat = mat(classLabels).transpose()m, n = shape(dataMatrix)# 初始化 b和alphas(alpha有点类似权重值。)b = 0alphas = mat(zeros((m, 1)))# 没有任何alpha改变的情况下遍历数据的次数iter = 0while (iter < maxIter):# w = calcWs(alphas, dataMatIn, classLabels)# print("w:", w)# 记录alpha是否已经进行优化，每次循环时设为0，然后再对整个集合顺序遍历alphaPairsChanged = 0for i in range(m):# print 'alphas=', alphas# print 'labelMat=', labelMat# print 'multiply(alphas, labelMat)=', multiply(alphas, labelMat)# 我们预测的类别 y = w^Tx[i]+b; 其中因为 w = Σ(1~n) a[n]*lable[n]*x[n]fXi = float(multiply(alphas, labelMat).T*(dataMatrix*dataMatrix[i, :].T)) + b# 预测结果与真实结果比对，计算误差EiEi = fXi - float(labelMat[i])# 约束条件 (KKT条件是解决最优化问题的时用到的一种方法。我们这里提到的最优化问题通常是指对于给定的某一函数，求其在指定作用域上的全局最小值)# 0<=alphas[i]<=C，但由于0和C是边界值，我们无法进行优化，因为需要增加一个alphas和降低一个alphas。# 表示发生错误的概率：labelMat[i]*Ei 如果超出了 toler， 才需要优化。至于正负号，我们考虑绝对值就对了。'''# 检验训练样本(xi, yi)是否满足KKT条件yi*f(i) >= 1 and alpha = 0 (outside the boundary)yi*f(i) == 1 and 0<alpha< C (on the boundary)yi*f(i) <= 1 and alpha = C (between the boundary)'''if ((labelMat[i]*Ei < -toler) and (alphas[i] < C)) or ((labelMat[i]*Ei > toler) and (alphas[i] > 0)):# 如果满足优化的条件，我们就随机选取非i的一个点，进行优化比较j = selectJrand(i, m)# 预测j的结果fXj = float(multiply(alphas, labelMat).T*(dataMatrix*dataMatrix[j, :].T)) + bEj = fXj - float(labelMat[j])alphaIold = alphas[i].copy()alphaJold = alphas[j].copy()# L和H用于将alphas[j]调整到0-C之间。如果L==H，就不做任何改变，直接执行continue语句# labelMat[i] != labelMat[j] 表示异侧，就相减，否则是同侧，就相加。if (labelMat[i] != labelMat[j]):L = max(0, alphas[j] - alphas[i])H = min(C, C + alphas[j] - alphas[i])else:L = max(0, alphas[j] + alphas[i] - C)H = min(C, alphas[j] + alphas[i])# 如果相同，就没发优化了if L == H:print("L==H")continue# eta是alphas[j]的最优修改量，如果eta==0，需要退出for循环的当前迭代过程# 参考《统计学习方法》李航-P125~P128<序列最小最优化算法>eta = 2.0 * dataMatrix[i, :]*dataMatrix[j, :].T - dataMatrix[i, :]*dataMatrix[i, :].T - dataMatrix[j, :]*dataMatrix[j, :].Tif eta >= 0:print("eta>=0")continue# 计算出一个新的alphas[j]值alphas[j] -= labelMat[j]*(Ei - Ej)/eta# 并使用辅助函数，以及L和H对其进行调整alphas[j] = clipAlpha(alphas[j], H, L)# 检查alpha[j]是否只是轻微的改变，如果是的话，就退出for循环。if (abs(alphas[j] - alphaJold) < 0.00001):print("j not moving enough")continue# 然后alphas[i]和alphas[j]同样进行改变，虽然改变的大小一样，但是改变的方向正好相反alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])# 在对alpha[i], alpha[j] 进行优化之后，给这两个alpha值设置一个常数b。# w= Σ[1~n] ai*yi*xi => b = yj- Σ[1~n] ai*yi(xi*xj)# 所以：  b1 - b = (y1-y) - Σ[1~n] yi*(a1-a)*(xi*x1)# 为什么减2遍？ 因为是 减去Σ[1~n]，正好2个变量i和j，所以减2遍b1 = b - Ei- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i, :]*dataMatrix[i, :].T - labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[i, :]*dataMatrix[j, :].Tb2 = b - Ej- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i, :]*dataMatrix[j, :].T - labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[j, :]*dataMatrix[j, :].Tif (0 < alphas[i]) and (C > alphas[i]):b = b1elif (0 < alphas[j]) and (C > alphas[j]):b = b2else:b = (b1 + b2)/2.0alphaPairsChanged += 1print("iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))# 在for循环外，检查alpha值是否做了更新，如果在更新则将iter设为0后继续运行程序# 知道更新完毕后，iter次循环无变化，才推出循环。if (alphaPairsChanged == 0):iter += 1else:iter = 0print("iteration number: %d" % iter)return b, alphas

查看代码运行时间及结果：

# 惩罚因子C = 0.6
# 容错率toler = 0.001
# 退出前最大循环次数maxIter=40
%time b, alphas = smoSimple(dataArr, labelArr, 0.6, 0.001, 40)