关于算法:机器学习算法系列十六非线性支持向量机算法NonLinear-Support-Vector-Machine

41次阅读

共计 9311 个字符，预计需要花费 24 分钟才能阅读完成。

浏览本文须要的背景知识点：线性反对向量机、一丢丢编程常识

一、引言

  后面咱们用两节介绍了两种反对向量机模型——硬距离反对向量机、软距离反对向量机，这两种模型能够统称为线性反对向量机，上面来介绍另一种反对向量机模型——非线性反对向量机 1（Non-Linear Support Vector Machine）。

二、模型介绍

  在介绍非线性反对向量机之前，让咱们来看如下的线性不可分的数据集：


<center> 图 2 -1 </center>

  图 2 - 1 中别离用圆和叉来代表不同的标签分类，很显著该数据集无奈通过一条直线正确分类，然而如图 2 - 2 所示，该数据集能够被一条椭圆曲线正确的分类。


<center> 图 2 -2</center>

  然而想要求解如图 2 - 2 中的椭圆曲线这种非线性分类的问题，相对来说难度较大。既然线性分类问题绝对容易求解，那么该数据是否通过肯定的非线性变动转化为一个线性分类的数据呢，答案是能够的。


<center> 图 2 -3</center>

  如图 2 - 3 所示，将数据集进行非线性转换（Z = X * X），这时能够看到转换后的数据能够通过一条直线正确分类，原来图 2 - 2 中的椭圆曲线变换成了图 2 - 3 中的直线，原来非线性分类问题变换成了线性分类问题。

原始模型

  先来回顾一下后面软距离反对向量机的原始模型：

$$
\begin{array}{c}
\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 -1</center>

  假如存在一个 φ 函数能对 x 进行变换，例如下面例子中的将 x 映射成 z，带入下面的软距离反对向量机的原始模型中得：

$$
\begin{array}{c}
\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} \phi(x_{i})+b\right) \quad i=1,2, \cdots, N
\end{array}
$$

<center> 式 2 -2</center>

  上式 2 - 2 就是非线性反对向量机的原始模型。

对偶模型

  同后面两节中所求解的对偶模型，能够失去非线性反对向量机的对偶模型：

$$
\begin{array}{c}
\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} \phi(x_{i})^{T} \phi(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 -3</center>

核技巧

  察看式 2 - 3 中对软距离反对向量机的对偶模型的变动局部，特征向量通过 φ 函数变换后，其维度数可能会很高，再求内积通常比拟艰难。为了绕过这个内积计算，咱们能够找到一个如式 2 - 4 所示的函数，使得通过该函数计算后的值等于特色转换后的内积，这样的函数被称为核函数（Kernel Function），这种替换的办法被称为核技巧 2（Kernel Trick）

$$
K(x_i, x_j) = \phi(x_i)^T\phi(x_j)
$$

<center> 式 2 -4</center>

  应用核技巧带入式 2 - 3 后失去最终的非线性反对向量机的对偶模型：

$$
\begin{array}{c}
\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} K(x_i, 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 -5</center>

  核技巧下的决策面：

（1）原始的决策面函数

（2）带入权重系数 w

（3）应用核技巧替换后失去最初的决策面函数

$$
\begin{aligned}
f(x) &=w^{T} \phi(x)+b & (1)\\
&=\sum_{i=1}^{N} \lambda_{i} y_{i} \phi\left(x_{i}\right) \phi(x)+b & (2)\\
&=\sum_{i=1}^{N} \lambda_{i} y_{i} K\left(x_{i}, x\right)+b & (3)
\end{aligned}
$$

<center> 式 2 -6</center>

  核技巧不仅仅实用于反对向量机，在其余机器学习算法中也有宽泛的利用。

核函数

  罕用核函数：

（1）线性核函数（Linear Kernel）：应用线性核函数的反对向量机分类后果与软距离反对向量机的后果雷同

（2）多项式核函数（Polynomial Kernel）：当 γ = 1，ε = 0，d = 1 时，多项式核函数进化为线性核函数

（3）径向基核函数（Radial basis function/RBF Kernel）：常见的高斯核函数（Gaussian Kernel）是径向基核函数的一种

（4）Sigmoid 核函数（Sigmoid Kernel）：tanh 为双曲正切函数

$$
\begin{aligned}
K\left(x_{i}, x_{j}\right)&=x_{i}^{T} x_{j} & (1) \\
K\left(x_{i}, x_{j}\right)&=\left(\gamma x_{i}^{T} x_{j}+\varepsilon\right)^{d} & (2)\\
K\left(x_{i}, x_{j}\right)&=e^{-\gamma\left\|x_{i}-x_{j}\right\|^{2}} & (3)\\
K\left(x_{i}, x_{j}\right)&=\tanh \left(\gamma x_{i}^{T} x_{j}+\varepsilon\right) & (4)
\end{aligned}
$$

<center> 式 2 -7</center>

三、算法步骤

  比照软距离反对向量机的对偶模型，非线性反对向量机只有特征向量的内积应用核函数替换，其余没有变动，所以对应的序列最小优化算法（SMO）变动也不大。

  算法步骤请参考后面两节中的内容和上面的代码实现，同时能够参看 SMO 算法对应的论文 3 实现。

四、代码实现

应用 Python 实现非线性反对向量机（SMO 算法）

import numpy as np

class SMO:
    """
    反对向量机
    序列最小优化算法实现（Sequential minimal optimization/SMO）"""

    def __init__(self, X, y, kernel, degree = 3, coef0 = 0.0, gamma = 1.0):
        # 训练样本特色矩阵（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
        # 代价值
        self.cost = -np.inf
        # 核函数
        self.kernel = kernel
        # 核函数相干参数
        self.degree = degree
        self.coef0 = coef0
        self.gamma = gamma

    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]
        
        s = 0
        for idx in range(X_gt.shape[0]):
            ss = 0
            for jdx in range(X_gt.shape[0]):
                ss += alpha_gt[jdx] * y_gt[jdx] * self.kernel(X_gt[jdx], X_gt[idx], self.degree, self.coef0, self.gamma)
            s += y_gt[idx] - ss
        return 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] * self.kernel(X[idx], X[jdx], self.degree, self.coef0, self.gamma) * 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 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]
        U = self.calcU(idx, jdx, C)
        V = self.calcV(idx, jdx, C)
        if U == V:
            return False
        Kii = self.kernel(X[idx], X[idx], self.degree, self.coef0, self.gamma)
        Kjj = self.kernel(X[jdx], X[jdx], self.degree, self.coef0, self.gamma)
        Kij = self.kernel(X[idx], X[jdx], self.degree, self.coef0, self.gamma)
        # 计算 K
        K = Kii + Kjj - 2 * Kij
        oldAlphaIdx = alpha[idx]
        oldAlphaJdx = alpha[jdx]
        oldB = self.b
        s = y[idx] * y[jdx]
        if K > 0:
            # 计算 jdx 的新拉格朗日乘子
            newAlphaJdx = oldAlphaJdx + y[jdx] * (Ei - Ej) / K
            if newAlphaJdx < U:
                # 当新值超过上界时，批改其为上界
                newAlphaJdx = U
            if newAlphaJdx > V:
                # 当新值低于下界时，批改其为下界
                newAlphaJdx = V
        else:
            fi = y[idx] * (Ei + oldB) - oldAlphaIdx * Kii - s * oldAlphaJdx * Kij
            fj = y[jdx] * (Ej + oldB) - s * oldAlphaIdx * Kij - oldAlphaJdx * Kjj
            Vv = oldAlphaIdx + s * (oldAlphaJdx - V)
            Uu = oldAlphaIdx + s * (oldAlphaJdx - U)
            Vv = Vv * fi + V * fj + 0.5 * (Vv ** 2) * Kii + 0.5 * (V ** 2) * Kjj + s * V * Vv * Kij
            Uu = Uu * fi + U * fj + 0.5 * (Uu ** 2) * Kii + 0.5 * (U ** 2) * Kjj + s * U * Uu * Kij
            if Vv < Uu:
                newAlphaJdx = V
            elif Vv > Uu:
                newAlphaJdx = U
            else:
                newAlphaJdx = oldAlphaJdx
        if oldAlphaJdx == newAlphaJdx:
            # 当新值与旧值相等时，判断为未更新，间接返回
            return False
        # 计算 idx 的新拉格朗日乘子
        newAlphaIdx = oldAlphaIdx + s * (oldAlphaJdx - newAlphaJdx)
        # 更新拉格朗日乘子向量
        alpha[idx] = newAlphaIdx
        alpha[jdx] = newAlphaJdx
        oldB = self.b
        # 从新计算偏移量
        self.b = self.calcB()
        # 从新计算误差向量
        newErrors = []
        for i in range(X.shape[0]):
            newError = errors[i] + y[idx] * (newAlphaIdx - oldAlphaIdx) * self.kernel(X[idx], X[i], self.degree, self.coef0, self.gamma) + y[jdx] * (newAlphaJdx - oldAlphaJdx) * self.kernel(X[jdx], X[i]) - oldB + self.b
            newErrors.append(newError)
        self.errors = newErrors
        return True
    
    def predict(self, X):
        fxs = []
        for idx in range(len(X)):
            fx = 0
            for jdx in range(self.X.shape[0]):
                fx += self.y[jdx] * self.alpha[jdx] * self.kernel(self.X[jdx], X[idx], self.degree, self.coef0, self.gamma)
            fxs.append(fx + self.b)
        return np.sign(fxs)

线性核函数反对向量机 

# 线性核函数反对向量机
def linearKernel(x, z, degree = 3, coef0 = 0.0, gamma = 1.0):
    return x.dot(z)

linearSmo = SMO(X, y, kernel = linearKernel)
linearSmo.fit()

多项式核函数反对向量机 

# 多项式核函数反对向量机
def polynomialKernel(x, z, degree = 3, coef0 = 0.0, gamma = 1.0):
    return np.power(gamma * x.dot(z) + coef0, degree)

polynomialSmo = SMO(X, y, kernel = polynomialKernel)
polynomialSmo.fit()

径向基核函数反对向量机 

# 径向基核函数反对向量机
def rbfKernel(x, z, degree = 3, coef0 = 0.0, gamma = 1.0):
    return np.exp(-gamma * np.power(np.linalg.norm(x - z), 2))

rbfSmo = SMO(X, y, kernel = rbfKernel)
rbfSmo.fit()

五、第三方库实现

scikit-learn4 实现 

from sklearn.svm import SVC

# 线性核函数反对向量机
svc = SVC(kernel = "linear")
# 多项式核函数反对向量机
svc = SVC(kernel = "poly", gamma = 1)
# 径向基核函数反对向量机
svc = SVC(kernel = "rbf", gamma = 1)

# 拟合数据
svc.fit(X, y)

scikit-learn 外部应用的是 LIBSVM5 库，对于该库实现的细节能够参看这篇文章 6。

六、动画演示

  上面三个图别离展现了线性不可分的数据集针对不同的核函数的分类后果，红色示意标签值为 - 1 的样本点，蓝色代表标签值为 1 的样本点。浅红色的区域为预测值为 - 1 的局部，浅蓝色的区域则为预测值为 1 的局部：


<center> 图 6 -1</center>


<center> 图 6 -2</center>


<center> 图 6 -3</center>

  能够看到不同的核函数对于分类后果的影响很大，当不晓得特色转换的的模式时，也就无奈抉择出适合的核函数，但通常状况下，非线性分类问题应用径向基核函数通过穿插验证的形式找出适合的参数即可。

七、思维导图


<center> 图 7 -1</center>

八、参考文献

  1. https://en.wikipedia.org/wiki…
  2. https://en.wikipedia.org/wiki…
  3. https://www.microsoft.com/en-…
  4. https://scikit-learn.org/stab…
  5. https://www.csie.ntu.edu.tw/~…
  6. https://github.com/Kaslanaria…

残缺演示请点击这里

注：本文力求精确并通俗易懂，但因为笔者也是初学者，程度无限，如文中存在谬误或脱漏之处，恳请读者通过留言的形式批评指正

本文首发于——AI 导图 ，欢送关注

正文完
算法
发表至： 算法
2022-02-07
 0
关于算法:AI在出行场景的应用实践路线规划ETA动态事件挖掘…
关于算法:力扣之两数之和
关于算法:奇妙的安全旅行之SHA算法
关于算法:Conways-Problem