分步实现机器学习 V - 支持向量机指南

Ryukkkk

5.00/5 (3投票s)

2019年5月14日

CPOL

4分钟阅读

6792

易于实现的机器学习

引言

支持向量机（SVM）是基于特征空间最大间隔的分类器。SVM的学习策略是最大化间隔，这可以转化为一个凸二次规划问题。在学习SVM算法之前，我先介绍一些术语。

函数间隔定义为：对于给定的训练集T和超平面(w,b)，超平面(w,b) 与样本(xi, yi)之间的函数间隔为

$\hat{\gamma}_{i}=y_{i}\left(w\cdot x+b\right)$

超平面(w,b)与训练集T之间的函数间隔是 $\hat{\gamma}_{i}$ 的最小值。

$\hat{\gamma}=\min\limits_{i=1,...,N} \hat{\gamma}_{i}$

函数间隔表示分类结果的置信度。如果超参数(w,b) 成比例地变化，例如，将(w,b)变为(2w,2b)。尽管超平面没有改变，但函数间隔会加倍。因此，我们对w 施加一些约束以使超平面明确化，例如归一化||w|| = 1。然后，该间隔被称为几何间隔：对于给定的训练集T和超平面(w,b)，超平面(w,b) 与样本(xi, yi)之间的函数间隔为

${\gamma}_{i}=y_{i}\left(\frac{w}{||w||}\cdot x_{i}+\frac{b}{||w||}\right)$

类似地，超平面(w,b)与训练集T之间的几何间隔是 ${\gamma}_{i}$ 的最小值。

${\gamma}=\min\limits_{i=1,...,N} {\gamma}_{i}$

现在，我们可以总结函数间隔和几何间隔之间的关系。

$\gamma_{i}=\frac{\hat{\gamma_{i}}}{||w||}$

$\gamma=\frac{\hat{\gamma}}{||w||}$

SVM模型

SVM模型由优化问题、优化算法和分类组成。

优化问题

当数据集线性可分时，支持向量是离超平面最近的样本，如下图所示。

H1 和H2 上的样本是支持向量。H1 和H2 之间的距离称为硬间隔。那么，SVM的优化问题是

$\min\limits_{w,b}\ \frac{1}{2}||w||^{2}\\ s.t.\ y_{i}\left(w\cdot x+b\right)-1\geq 0,\ i=1,2,...,N$

当数据集不可线性分割时，训练集中的一些样本不满足函数间隔大于或等于1的条件。为了解决这个问题，我们为每个样本 (xi, yi)添加一个松弛变量 $\xi_{i}\geq0$ 。然后，约束条件为

$y_{i}\left(w\cdot x+b\right)\geq 1-\xi_{i}$

同时，为每个松弛变量添加成本。目标函数为

$\frac{1}{2}||w||^{2}+C\sum_{i=1}^{N}\xi_{i}$

其中C 是惩罚系数。当C 较大时，误分类的惩罚会增加，而误分类的惩罚会减少。那么，SVM的优化问题是

$\min\limits_{w,b}\ \frac{1}{2}||w||^{2}+C\sum_{i=1}^{N}\xi_{i}\\ s.t.\ y_{i}\left(w\cdot x+b\right)\geq 1-\xi_{i},\ i=1,2,...,N\\ \xi_{i}\geq0,\ i=1,2,...,N$

支持向量位于间隔边界上或边界与超平面之间，如下图所示。在这种情况下，间隔称为软间隔。

当数据集不可线性分割时，需要应用核技巧将数据从原始空间转换为特征空间。用于核技巧的函数称为核函数，定义为

$K\left(x,z\right)=\phi\left(x\right)\cdot \phi\left(z\right)$

其中 $\phi\left(x\right)$ 是从输入空间到特征空间的映射，即：

$\phi\left(x\right):\chi\rightarrow\mathcal{H}$

核技巧的代码如下所示：

def kernelTransformation(self, data, sample, kernel):
        sample_num, feature_dim = np.shape(data)
        K = np.zeros([sample_num])
        if kernel == "linear":  # linear function
            K = np.dot(data, sample.T)
        elif kernel == "poly":  # polynomial function
            K = (np.dot(data, sample.T) + self.c) ** self.n
        elif kernel == "sigmoid":  # sigmoid function
            K = np.tanh(self.g * np.dot(data, sample.T) + self.c)
        elif kernel == "rbf":  # Gaussian function
            for i in range(sample_num):
                delta = data[i, :] - sample
                K[i] = np.dot(delta, delta.T)
            K = np.exp(-self.g * K)
        else:
            raise NameError('Unrecognized kernel function')
        return K

应用核技巧后，SVM的优化问题为

$\min\limits_{\alpha}\ \ \frac{1}{2}\sum_{i=1}^{N}\sum_{j=1}^{N}\alpha_{i}\alpha_{j}y_{i}y_{j}K\left(x_{i},y_{j}\right)-\sum_{i=1}^{N}\alpha_{i}\\s.t.\ \ \sum_{i=1}^{N}\alpha_{i}y_{i}=0\\ 0\leq\alpha_{i}\leq C,i=1,2,...,N$

其中 $\alpha_{i}$ 是拉格朗日乘子。

优化算法

使用传统的凸二次规划算法可以解决SVM优化问题。但是，当训练集很大时，算法会花费很长时间。为了解决这个问题，Platt提出了一种称为序列最小最优化（SMO）的高效算法。

SMO是一种启发式算法。当所有变量都满足KKT条件时，优化问题就解决了。否则，选择两个变量，固定其他变量，并构造一个只有这两个变量的凸二次规划问题。该问题有两个变量：一个选择违反KKT条件的alpha，另一个由约束决定。设 $\alpha_{1,}\alpha_{2}$ 为变量，固定 $\alpha_{3},\alpha_{4},...,\alpha_{N}$ 。因此， $\alpha_{1}$ 由下式计算：

$\alpha_{1}=-y_{1}\sum_{i=2}^{N}\alpha_{i}y_{i}$

如果 $\alpha_{2}$ 确定，则 $\alpha_{1}$ 也随之确定。在SMO的每次迭代中，都会更新两个变量，直到问题解决。然后，SVM的优化问题为：

$\min\limits_{\alpha_{1},\alpha_{2}} W\left(\alpha_{1},\alpha_{2}\right)=\frac{1}{2}K_{11}\alpha_{1}^{2}+\frac{1}{2}K_{22}\alpha_{2}^{2}+y_{1}y_{2}K_{12}\alpha_{1}\alpha_{2}-\\\left(\alpha_{1}+\alpha_{2}\right)+y_{1}\alpha_{1}\sum_{i=3}^{N}y_{i}\alpha_{i}K_{i1}+y_{2}\alpha_{2}\sum_{i=3}^{N}y_{i}\alpha_{i}K_{i2}\\s.t.\ \ \ \alpha_{1}y_{1}+\alpha_{2}y_{2}=-\sum_{i=3}^{N}y_{i}\alpha_{i}\\ 0\leq\alpha_{i}\leq C,i=1,2$

当只有一个变量时，这是一个简单的二次规划问题，如下所示：

因为约束是：

$\alpha_{1}y_{1}+\alpha_{2}y_{2}=-\sum_{i=3}^{N}y_{i}\alpha_{i} = k$

当 $y_{1}=y_{2}$ 时， $\alpha_{1}+\alpha_{2}=k$

当 $y_{1}\ne y_{2}$ 时，因为 $y_{1},y_{2}\in\left\{1，-1\right\}$ ，所以 $\alpha_{1}-\alpha_{2}=k$ 。

优化后的值 $\alpha_{2}^{new}$ 遵循：

$L\leq\alpha_{2}^{new}\leq H$

其中 L 和 H 是对角线的下界和上界，计算方法如下：

$L,H=\left\{ \begin{aligned} L= max\left(0, \alpha_{2}^{old}-\alpha_{1}^{old}\right),H= min\left(C,C+ \alpha_{2}^{old}-\alpha_{1}^{old}\right) if\ y_{1}\ne y_{2} \\ L= max\left(0, \alpha_{2}^{old}+\alpha_{1}^{old}-C\right),H= min\left(C,\alpha_{2}^{old}+\alpha_{1}^{old}\right) if\ y_{1}= y_{2} \\ \end{aligned} \right.$

未裁剪的优化值 $\alpha_{2}^{new,unc}$ 遵循：

$\alpha_{2}^{new,unc}=\alpha_{2}^{old}+\frac{y_{2}\left(E_{1}-E{2}\right)}{\eta}$

其中 E₁和 E₂是预测值 g(x) 与真实值之间的差值。g(x) 定义为：

$g\left(x\right)=\sum_{i=1}^{N}\alpha_{i}y_{i}K\left(x_{i},x\right)+b$

$\eta=K_{11}+K_{22}-2K_{12}$

至此，我们得到了 $\alpha_{1,}\alpha_{2}$ 的可行解。

$\alpha_{2}^{new}=\left\{\begin{aligned} &H,&\alpha_{2}^{new,unc}>H\\ &\alpha_{2}^{new,unc},&L\leq\alpha_{2}^{new,unc}\leq H\\ &L,&\alpha_{2}^{new,unc}<L\end{aligned}\right.$

$\alpha_{1}^{new}=\alpha_{1}^{old}+y_{1}y_{2}\left(\alpha_{2}^{old}-\alpha_{2}^{new}\right)$

SMO中有两个循环，即外循环和内循环。

在外循环中，选择违反KKT条件的 alpha，即：

$\alpha_{i}=0\Leftrightarrow y_{i}g\left(x_{i}\right)\geq1\\ 0<\alpha_{i}<C\Leftrightarrow y_{i}g\left(x_{i}\right)=1\\ \alpha_{i}=C\Leftrightarrow y_{i}g\left(x_{i}\right)\leq1$

首先，搜索满足 $0<\alpha_{i}<C$ 的样本。如果所有样本都满足条件，则搜索整个数据集。

在内循环中，首先搜索使得 $\left|E_{1}-E_{2}\right|$ 最大的 $\alpha_{2}$ 。如果选择的 $\alpha_{2}$ 没有足够减小，则首先搜索边界上的样本。如果选择的 $\alpha_{2}$ 减小得足够多，则停止搜索。否则，搜索整个数据集。如果在搜索整个数据集后找到一个可行的 $\alpha_{2}$ ，我们将选择一个新的 $\alpha_{1}$ 。

在每次迭代中，我们通过以下方式更新b ：

$b_{1}^{new}=-E_{1}-y_{1}K_{11}\left(\alpha_{1}^{new}-\alpha_{1}^{old}\right)-y_{2}K_{21}\left(\alpha_{2}^{new}-\alpha_{2}^{old}\right)+b^{old}$

$b_{2}^{new}=-E_{2}-y_{1}K_{12}\left(\alpha_{1}^{new}-\alpha_{1}^{old}\right)-y_{2}K_{22}\left(\alpha_{2}^{new}-\alpha_{2}^{old}\right)+b^{old}$

$b=\left\{\begin{aligned} &b_{1}^{new},&if\ 0<\alpha_{1}^{new}\\ &b_{2}^{new},&if\ 0<\alpha_{2}^{new}\\ &\frac{b_{1}^{new}+b_{2}^{new}}{2},&else\end{aligned}\right.$

为了方便起见，我们必须存储 $E_{i}$ 的值，计算方法如下：

$E_{i}^{new}=\sum_{s}y_{j}\alpha_{j}K\left(x_{i},y_{j}\right)+b^{new}-y_{i}$

搜索和更新的代码如下所示：

def innerLoop(self, i):
        Ei = self.calculateErrors(i)
        if self.checkKKT(i):

            j, Ej = self.selectAplha2(i, Ei)          # select alpha2 according to alpha1

            # copy alpha1 and alpha2
            old_alpha1 = self.alphas[i]
            old_alpha2 = self.alphas[j]

            # determine the range of alpha2 L and H      in page of 126
            # if y1 != y2    L = max(0, old_alpha2 - old_alpha1), 
                             H = min(C, C + old_alpha2 - old_alpha1)
            # if y1 == y2    L = max(0, old_alpha2 + old_alpha1 - C), 
                             H = min(C, old_alpha2 + old_alpha1)
            if self.train_label[i] != self.train_label[j]:
                L = max(0, old_alpha2 - old_alpha1)
                H = min(self.C, self.C + old_alpha2 - old_alpha1)
            else:
                L = max(0, old_alpha2 + old_alpha1 - self.C)
                H = min(self.C, old_alpha2 + old_alpha2)

            if L == H:
                # print("L == H")
                return 0

            # calculate eta in page of 127 Eq.(7.107)
            # eta = K11 + K22 - 2K12
            K11 = self.K[i, i]
            K12 = self.K[i, j]
            K21 = self.K[j, i]
            K22 = self.K[j, j]
            eta = K11 + K22 - 2 * K12
            if eta <= 0:
                # print("eta <= 0")
                return 0

            # update alpha2 and its error in page of 127 Eq.(7.106) and Eq.(7.108)
            self.alphas[j] = old_alpha2 + self.train_label[j]*(Ei - Ej)/eta
            self.alphas[j] = self.updateAlpha2(self.alphas[j], L, H)
            new_alphas2 = self.alphas[j]
            self.upadateError(j)


            # update the alpha1 and its error in page of 127 Eq.(7.109)
            # new_alpha1 = old_alpha1 + y1y2(old_alpha2 - new_alpha2)
            new_alphas1 = old_alpha1 + self.train_label[i] * 
                          self.train_label[j] * (old_alpha2 - new_alphas2)
            self.alphas[i] = new_alphas1
            self.upadateError(i)

            # determine b in page of 130 Eq.(7.115) and Eq.(7.116)
            # new_b1 = -E1 - y1K11(new_alpha1 - old_alpha1) - 
                             y2K21(new_alpha2 - old_alpha2) + old_b
            # new_b2 = -E2 - y1K12(new_alpha1 - old_alpha1) - 
                             y2K22(new_alpha2 - old_alpha2) + old_b
            b1 = - Ei - self.train_label[i] * K11 * (old_alpha1 - self.alphas[i]) - 
                        self.train_label[j] * K21 * (old_alpha2 - self.alphas[j]) + self.b
            b2 = - Ej - self.train_label[i] * K12 * (old_alpha1 - self.alphas[i]) - 
                        self.train_label[j] * K22 * (old_alpha2 - self.alphas[j]) + self.b
            if (self.alphas[i] > 0) and (self.alphas[i] < self.C):
                self.b = b1
            elif (self.alphas[j] > 0) and (self.alphas[j] < self.C):
                self.b = b2
            else:
                self.b = (b1 + b2)/2.0

            return 1
        else:
            return 0

分类

我们可以使用优化后的参数进行预测，计算公式为：

$f\left(x\right)=sign\left(\sum_{i=1}^{N}\alpha_{i}^{*}y_{i}K\left(x,x_{i}\right)+b^{*}\right)$

for i in range(test_num):
    kernel_data = self.kernelTransformation(support_vectors, test_data[i, :], self.kernel)
    probability[i] = np.dot(kernel_data.T, np.multiply
                     (support_vectors_label, support_vectors_alphas)) + self.b
    if probability[i] > 0:
        prediction[i] = 1
    else:
        prediction[i] = -1

结论与分析

SVM比之前的算法更复杂。在本文中，我们简化了搜索过程，使其更容易理解。最后，我们将我们的SVM与Sklearn中的SVM进行比较，并显示检测性能如下：

检测性能比sklearn的略差，这是因为我们的SVM中的SMO比sklearn的更简单。

本文相关的代码和数据集可以在 MachineLearning 中找到。