机器学习系列 10:支持向量机 02 - SMO(序列最小化)


针对 “支持向量机” 将分为如下几个部分进行介绍:html

  1. 支持向量机 01 - 线性可分支持向量机和线性支持向量机
  2. 支持向量机 02 - SMO(序列最小化)(本篇)
  3. 支持向量机 03 - 非线性支持向量机
  4. 支持向量机 04 - SVR(支持向量回归)


  本内容将介绍 SMO(序列最小化)算法,包含详细公式推导以及 Python 代码实现。

python

  在学习本内容前,须要对支持向量机基本理论知识有必定了解,若是您还不了解,能够参阅 支持向量机 01 - 线性可分支持向量机和线性支持向量机

web

  咱们知道,支持向量机的学习问题能够形式化为求解凸二次规划问题。这样的凸二次规划问题具备全局最优解,而且有许多最优化算法能够用于这一问题的求解。可是当训练样本容量很大时,这些算法每每变得很是低效,以至没法使用。算法

  SMO 算法是一种启发式算法,其基本思路是:若是全部变量的解都知足此最优化问题的 KKT 条件,那么这个最优化问题的解就获得了。由于 KKT 条件是该最优化问题的充要条件。不然,选择两个变量,固定其余变量,针对这两个变量构建一个二次规划问题。这个二次规划问题关于这两个变量的解应该更接近原始二次规划问题的解,由于这会使得原始二次规划问题的目标函数值变得更小。重要的是,这时子问题能够经过解析方法求解,这样就能够大大提升整个算法的计算速度。子问题有两个变量,一个是违反 KKT 条件最严重的那一个,另外一个由约束条件自动肯定。如此,SMO 算法将原问题不断分解为子问题并对子问题求解,进而达到求解原问题的目的。

缓存

1、SMO 算法描述

  整个 SMO 算法包括两个部分:求解两个变量二次规划的解析方法和选择变量的启发式方法。数据结构

1.1 两个变量二次规划的求解方法

1.1.1 问题转化

  在 支持向量机 01 - 线性可分支持向量机和线性支持向量机 中讲解了 SVM 的基本原理,了解到 SMO 算法要解以下凸二次规划的对偶问题app

(1) max α ( i = 1 m α i 1 2 i = 1 m j = 1 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) ) s . t . i = 1 m α i y ( i ) = 0 α i 0 , i = 1 , 2 ,   , m \begin{aligned} & \max_\alpha \left( \sum_{i=1}^{m}\alpha_i - \frac{1}{2} \sum_{i=1}^{m} \sum_{j=1}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} \right) \\\\ & \begin{aligned} s.t. \quad &\sum_{i=1}^{m}\alpha_i y^{(i)} = 0 \\\\ & \alpha_i \geq 0,\quad i=1,2,\cdots,m \end{aligned} \end{aligned} \tag{1} dom

  添加一个负号,将最大值问题转换成最小值问题机器学习

(2) min α ( 1 2 i = 1 m j = 1 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) i = 1 m α i ) s . t . i = 1 m α i y ( i ) = 0 α i 0 , i = 1 , 2 ,   , m \begin{aligned} & \min_\alpha \left( \frac{1}{2} \sum_{i=1}^{m} \sum_{j=1}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} - \sum_{i=1}^{m}\alpha_i \right) \\\\ & \begin{aligned} s.t. \quad & \sum_{i=1}^{m}\alpha_i y^{(i)} = 0 \\\\ & \alpha_i \geq 0,\quad i=1,2,\cdots,m \end{aligned} \end{aligned} \tag{2} ide

1.1.2 转换为二元函数

  假设选择的两个变量为 α 1 \alpha_1 α 2 \alpha_2 ,其余变量 α i ( i = 3 , 4 ,   , m ) \alpha_i(i=3,4, \cdots ,m) 是固定的,式(2)可写为

(3) W ( α 1 , α 2 ) = 1 2 i = 1 m j = 1 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) i = 1 m α i = 1 2 i = 1 m ( j = 1 2 α i α j y ( i ) y ( j ) x ( i ) x ( j ) + j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) ) α i α 2 i = 3 m α i = 1 2 i = 1 2 ( j = 1 2 α i α j y ( i ) y ( j ) x ( i ) x ( j ) + j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) ) + 1 2 i = 3 m ( j = 1 2 α i α j y ( i ) y ( j ) x ( i ) x ( j ) + j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) ) α i α 2 i = 3 m α i = 1 2 i = 1 2 j = 1 2 α i α j y ( i ) y ( j ) x ( i ) x ( j ) + i = 1 2 j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) + 1 2 i = 3 m j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) α 1 α 2 i = 3 m α i = 1 2 α 1 2 x ( 1 ) x ( 1 ) + 1 2 α 2 2 x ( 2 ) x ( 2 ) + y ( 1 ) y ( 2 ) α 1 α 2 x ( 1 ) x ( 2 ) + y ( 1 ) α 1 j = 3 m α j y ( j ) x ( 1 ) x ( j ) + y ( 2 ) α 2 j = 3 m α j y ( j ) x ( 2 ) x ( j ) α 1 α 2 + 1 2 i = 3 m j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) i = 3 m α i \begin{aligned} W(\alpha_1, \alpha_2) & = \frac{1}{2} \sum_{i=1}^{m} \sum_{j=1}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} - \sum_{i=1}^{m}\alpha_i \\\\ & = \frac{1}{2} \sum_{i=1}^{m} \left( \sum_{j=1}^{2} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} + \sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} \right) - \alpha_i - \alpha_2 - \sum_{i=3}^{m} \alpha_i \\\\ & = \frac{1}{2} \sum_{i=1}^{2} \left( \sum_{j=1}^{2} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} + \sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} \right) \\ & \quad\quad + \frac{1}{2} \sum_{i=3}^{m} \left( \sum_{j=1}^{2} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} + \sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} \right) - \alpha_i - \alpha_2 - \sum_{i=3}^{m} \alpha_i \\\\ & = \frac{1}{2}\sum_{i=1}^{2}\sum_{j=1}^{2} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} + \sum_{i=1}^{2}\sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} \\ & \quad\quad +\frac{1}{2}\sum_{i=3}^{m}\sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} - \alpha_1 - \alpha_2 - \sum_{i=3}^m \alpha_i \\\\ & = \frac{1}{2} \alpha_1^2 \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} + \frac{1}{2} \alpha_2^2 \mathbf{x}^{(2)} \mathbf{x}^{(2)} + y^{(1)} y^{(2)} \alpha_1 \alpha_2 \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} + y^{(1)} \alpha_1 \sum_{j=3}^m \alpha_j y^{(j)} \mathbf{x}^{(1)} \cdot \mathbf{x}^{(j)} \\ & \quad\quad + y^{(2)} \alpha_2 \sum_{j=3}^m \alpha_j y^{(j)} \mathbf{x}^{(2)} \cdot \mathbf{x}^{(j)} - \alpha_1 - \alpha_2 + \frac{1}{2}\sum_{i=3}^{m}\sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} - \sum_{i=3}^m \alpha_i \end{aligned} \tag{3}

  为了描述方便,咱们定义以下符号:

(4) f ( x ( i ) ) = j = 1 m α j y ( j ) x ( j ) x ( i ) + b = j = 1 m α j y ( j ) x ( i ) x ( j ) + b f(\mathbf{x}^{(i)}) = \sum_{j=1}^{m} \alpha_j y^{(j)} \mathbf{x}^{(j)} \cdot \mathbf{x}^{(i)} + b = \sum_{j=1}^{m} \alpha_j y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} + b \tag{4}

(5) v i = j = 3 m α j y ( j ) x ( i ) x ( j ) = f ( x ( i ) ) j = 1 2 α j y ( j ) x ( i ) x ( j ) b v_i = \sum_{j=3}^m \alpha_j y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} = f(\mathbf{x}^{(i)}) - \sum_{j=1}^2 \alpha_j y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)}- b \tag{5}

  根据式(4)和(5),目标函数式(3)可写为

(6) W ( α 1 , α 2 ) = 1 2 α 1 2 x ( 1 ) x ( 1 ) + 1 2 α 2 2 x ( 2 ) x ( 2 ) + y ( 1 ) y ( 2 ) α 1 α 2 x ( 1 ) x ( 2 ) + y ( 1 ) α 1 v 1 + y ( 2 ) α 2 v 2 α 1 α 2 + c o n s t a n t \begin{aligned} W(\alpha_1,\alpha_2) & = \frac{1}{2} \alpha_1^2 \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} + \frac{1}{2} \alpha_2^2 \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)} + y^{(1)} y^{(2)} \alpha_1 \alpha_2 \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} \\ & + y^{(1)} \alpha_1 v_1 + y^{(2)} \alpha_2 v_2 - \alpha_1 - \alpha_2 + constant \end{aligned} \tag{6}

其中

(7) c o n s t a n t = 1 2 i = 3 m j = 3 m α i α j y ( i ) y ( j ) x ( i ) x ( j ) i = 3 m α i constant = \frac{1}{2}\sum_{i=3}^{m}\sum_{j=3}^{m} \alpha_i \alpha_j y^{(i)} y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} - \sum_{i=3}^m \alpha_i \tag{7}

  由于 c o n s t a n t constant 部分对 α 1 \alpha_1 α 2 \alpha_2​ 来讲,属于常数项,在求导的时候,直接变为 0,因此咱们不须要关心这部分的内容。

1.1.3 转换为一元函数

  由于有 i = 1 m α i y ( i ) = 0 \sum_{i=1}^{m}\alpha_i y^{(i)} = 0 ,可得

(8) α 1 y ( 1 ) + α 2 y ( 2 ) = i = 3 m α i y ( i ) \alpha_1 y^{(1)} + \alpha_2 y^{(2)} = -\sum_{i=3}^m \alpha_i y^{(i)} \tag{8}

  对变量 α 1 \alpha_1 α 2 \alpha_2 来讲, α i \alpha_i y ( i ) y^{(i)} i = 3 , 4 ,   , m i = 3,4,\cdots,m )可看做常数项,所以有

(9) α 1 y ( 1 ) + α 2 y ( 2 ) = B \alpha_1 y^{(1)} + \alpha_2 y^{(2)} = B \tag{9}

其中, B B 为一个常数。将等式(9)两边同时乘以 y ( 1 ) y^{(1)} ,得

(10) α 1 = B y ( 1 ) α 2 y ( 1 ) y ( 2 ) = ( B α 2 y ( 2 ) ) y ( 1 ) \alpha_1 = By^{(1)} - \alpha_2 y^{(1)} y^{(2)} = \left(B - \alpha_2 y^{(2)}\right)y^{(1)} \tag{10}

将其带入式(6),获得只含变量 α 2 \alpha_2 的目标函数

(11) W ( α 2 ) = 1 2 ( B α 2 y ( 2 ) ) 2 x ( 1 ) x ( 1 ) + 1 2 α 2 2 x ( 2 ) x ( 2 ) + y ( 2 ) ( B α 2 y ( 2 ) ) α 2 x ( 1 ) x ( 2 ) + ( B α 2 y ( 2 ) ) v 1 + y ( 2 ) α 2 v 2 ( B α 2 y ( 2 ) ) y ( 1 ) α 2 + c o n s t a n t \begin{aligned} W(\alpha_2) & = \frac{1}{2} (B - \alpha_2 y^{(2)})^{2} \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} + \frac{1}{2} \alpha_2^{2} \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)} + y^{(2)} (B - \alpha_2 y^{(2)}) \alpha_2 \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} \\ & \quad\quad + (B - \alpha_2 y^{(2)}) v_1 + y^{(2)} \alpha_2 v_2 - (B - \alpha_2 y^{(2)})y^{(1)} - \alpha_2 + constant \end{aligned} \tag{11}

1.1.4 求解 α 2 n e w , u n c l i p p e d \alpha_2^{new,unclipped}

  将式(11)对 α 2 \alpha_2 求导,得

(12) W ( α 2 ) α 2 = x ( 1 ) x ( 1 ) α 2 + x ( 2 ) x ( 2 ) α 2 2 x ( 1 ) x ( 2 ) α 2 B x ( 1 ) x ( 1 ) y ( 2 ) + B x ( 1 ) x ( 2 ) y ( 2 ) + y ( 1 ) y ( 2 ) v 1 y ( 2 ) + v 2 y ( 2 ) 1 \begin{aligned} \frac{\partial W(\alpha_2)}{\partial \alpha_2} &= \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} \alpha_2 + \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)} \alpha_2 - 2 \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} \alpha_2 \\ &\quad\quad - B \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} y^{(2)} + B \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} y^{(2)} + y^{(1)} y^{(2)} - v_1 y^{(2)} + v_2 y^{(2)} - 1 \end{aligned} \tag{12}

令其为 0 0​ ,得

(13) α 2 n e w = y ( 2 ) ( y ( 2 ) y ( 1 ) + B ( x ( 1 ) x ( 1 ) x ( 1 ) x ( 2 ) ) + v 1 v 2 ) x ( 1 ) x ( 1 ) + x ( 2 ) x ( 2 ) 2 x ( 1 ) x ( 2 ) \alpha_2^{new} = \frac {y^{(2)}\left(y^{(2)} - y^{(1)} + B(\mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} - \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)}) + v_1 - v_2\right)} {\mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} + \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)} - 2\mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)}} \tag{13}

  定义以下符号

(14) E i = f ( x ( i ) ) y ( i ) E_i = f(\mathbf{x}^{(i)}) - y^{(i)} \tag{14}

(15) η = x ( 1 ) x ( 1 ) + x ( 2 ) x ( 2 ) 2 x ( 1 ) x ( 2 ) \eta = \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} + \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)} - 2\mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} \tag{15}

其中 E i E_i 为偏差项, η \eta 为学习速率。

  因为已知

(16) B = α 1 o l d y ( 1 ) + α 2 o l d y ( 2 ) B = \alpha_1^{old} y^{(1)} + \alpha_2^{old} y^{(2)} \tag{16}

(17) v i = f ( x ( i ) ) j = 1 2 α j y ( j ) x ( i ) x ( j ) b v_i = f\left(\mathbf{x}^{(i)}\right) - \sum_{j=1}^2 \alpha_j y^{(j)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(j)} - b \tag{17}

将式(16)和(17)带入式(13), 将 α 2 n e w , u n c l i p p e d \alpha_2^{new,unclipped}​ 化简后得

(18) α 2 n e w , u n c l i p p e d = α 2 o l d + y ( 2 ) ( E 1 E 2 ) η \alpha_2^{new,unclipped} = \alpha_2^{old} + \frac{y^{(2)}(E_1-E_2)}{\eta} \tag{18}

1.1.5 求解 α 2 n e w \alpha_2^{new} (对原始 α 2 n e w , u n c l i p p e d \alpha_2^{new,unclipped} 解进行修剪)

  上面求出的 α 2 n e w , u n c l i p p e d \alpha_2^{new, unclipped} 解是没有通过修剪的。咱们知道 α i \alpha_i 存在如下约束条件

(19) { 0 < α i < C α 1 y ( 1 ) + α 2 y ( 2 ) = B \left \{ \begin{array}{cc} \begin{aligned} & 0 < \alpha_i < C \\\\ & \alpha_1 y^{(1)} + \alpha_2 y^{(2)} = B \end{aligned} \end{array} \right. \tag{19}

可知 α 2 n e w \alpha_2^{new}​ 的取值须要知足必定条件,假设为

(20) L α 2 n e w H L \leq \alpha_2^{new} \leq H \tag{20}

  图-1 在二维平面上直观地表达了式(19)中的约束条件,从而可知

  • y ( 1 ) y ( 2 ) y^{(1)} \neq y^{(2)} 时,存在 α 2 o l d α 1 o l d = k \alpha_2^{old} - \alpha_1^{old}=k ,得

(21) L = m a x ( 0 , α 2 o l d α 1 o l d ) , H = m i n ( C , C + α 2 o l d α 1 o l d ) L = max(0, \alpha_2^{old} - \alpha_1^{old}),\quad H = min(C, C + \alpha_2^{old} - \alpha_1^{old}) \tag{21}

  • y ( 1 ) = y ( 2 ) y^{(1)} = y^{(2)} 时,存在 α 2 o l d + α 1 o l d = k \alpha_2^{old} + \alpha_1^{old}=k​ ,得

(22) L = m a x ( 0 , α 2 o l d + α 1 o l d C ) , H = m i n ( C , α 2 o l d + α 1 o l d ) L = max(0, \alpha_2^{old} + \alpha_1^{old} - C), \quad H = min(C, \alpha_2^{old} + \alpha_1^{old}) \tag{22}

图-1

  则通过修剪后, α 2 n e w \alpha_2^{new} 的解为

(23) α 2 n e w = { H , α 2 n e w , u n c l i p p e d > H α 2 n e w , u n c l i p p e d , L α 2 n e w , u n c l i p p e d H L , α 2 n e w , u n c l i p p e d < L \alpha_2^{new} = \left \{ \begin{array}{cc} \begin{aligned} & H, & \alpha_2^{new,unclipped} > H \\\\ & \alpha_2^{new,unclipped}, & L \leq \alpha_2^{new,unclipped} \leq H \\\\ & L, & \alpha_2^{new,unclipped} < L \end{aligned} \end{array} \right. \tag{23}

1.1.6 求解 α 1 n e w \alpha_1^{new}

  由于存在如下等式条件

(24) α 1 o l d y ( 1 ) + α 2 o l d y ( 2 ) = α 1 n e w y ( 1 ) + α 2 n e w y ( 2 ) \alpha_1^{old} y^{(1)} + \alpha_2^{old} y^{(2)} = \alpha_1^{new} y^{(1)} + \alpha_2^{new} y^{(2)} \tag{24}

从而得出

(25) α 1 n e w = α 1 o l d + y ( 1 ) y ( 2 ) ( α 2 o l d α 2 n e w ) \alpha_1^{new} = \alpha_1^{old} + y^{(1)} y^{(2)} (\alpha_2^{old} - \alpha_2^{new}) \tag{25}

1.1.7 求解阈值 b n e w b^{new}

  在求得 α 1 n e w \alpha_1^{new} α 2 n e w \alpha_2^{new} 后,须要根据它们更新 b b

  1. 若是 0 < α 1 n e w < C 0 < \alpha_1^{new} < C
      由 KKT 条件可知(在 支持向量机 01 - 线性可分支持向量机和线性支持向量机 中 2.2 节有进行介绍),当 0 < α i < C 0 < \alpha_i < C 时,这个样本点是支持向量。所以,知足
    (26) y ( 1 ) ( w T x ( 1 ) + b ) = 1 y^{(1)}(\mathbf{w}^T \mathbf{x}^{(1)} + b) = 1 \tag{26}

    上面公式两边同时乘以 y ( 1 ) y^{(1)}​
    (27) w T x ( 1 ) + b = y ( 1 ) \mathbf{w}^T \mathbf{x}^{(1)} + b = y^{(1)} \tag{27}

    w = i = 1 m α i y ( i ) x ( i ) \mathbf{w} = \sum_{i = 1}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)}​ 代入,得

    (28) i = 1 m α i y ( i ) x ( i ) x ( 1 ) + b = y ( 1 ) \sum_{i=1}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(1)} + b = y^{(1)} \tag{28}

      由于咱们是根据 α 1 n e w \alpha_1^{new}​ α 2 n e w \alpha_2^{new}​ 的值更新 b b​ ,整理可得:

    (29) b 1 n e w = y ( 1 ) i = 3 m α i y ( i ) x ( i ) x ( 1 ) α 1 n e w y ( 1 ) x ( 1 ) x ( 1 ) α 2 n e w y ( 2 ) x ( 2 ) x ( 1 ) b_1^{new} =y^{(1)} -\sum_{i=3}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(1)} -\alpha_1^{new} y^{(1)} \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} -\alpha_2^{new} y^{(2)} \mathbf{x}^{(2)} \cdot \mathbf{x}^{(1)} \tag{29}

    其中有

    (30) y ( 1 ) i = 3 m α i y ( i ) x ( i ) x ( 1 ) = y ( 1 ) i = 1 m α i y ( i ) x ( i ) x ( 1 ) b o l d + α 1 o l d y ( 1 ) x ( 1 ) x ( 1 ) + α 2 o l d y ( 2 ) x ( 2 ) x ( 1 ) + b o l d = y ( 1 ) f ( x ( 1 ) ) + α 1 o l d y ( 1 ) x ( 1 ) x ( 1 ) + α 2 o l d y ( 2 ) x ( 2 ) x ( 1 ) + b o l d = E 1 + α 1 o l d y ( 1 ) x ( 1 ) x ( 1 ) + α 2 o l d y ( 2 ) x ( 2 ) x ( 1 ) + b o l d \begin{aligned} y^{(1)} - \sum_{i=3}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(1)} &= y^{(1)} - \sum_{i=1}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(1)} -b^{old} \\ & \quad\quad +\alpha_1^{old} y^{(1)} \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} +\alpha_2^{old} y^{(2)} \mathbf{x}^{(2)} \cdot \mathbf{x}^{(1)} +b^{old} \\\\ &= y^{(1)} - f\left(\mathbf{x}^{(1)}\right) +\alpha_1^{old} y^{(1)} \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} +\alpha_2^{old} y^{(2)} \mathbf{x}^{(2)} \cdot \mathbf{x}^{(1)} +b^{old} \\\\ &= - E_1 +\alpha_1^{old} y^{(1)} \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} +\alpha_2^{old} y^{(2)} \mathbf{x}^{(2)} \cdot \mathbf{x}^{(1)} +b^{old} \end{aligned} \tag{30}

    将式(30)代入式(29)得

    (31) b 1 n e w = b o l d E 1 y ( 1 ) ( α 1 n e w α 1 o l d ) x ( 1 ) x ( 1 ) y ( 2 ) ( α 2 n e w α 2 o l d ) x ( 2 ) x ( 1 ) b_1^{new} = b^{old} - E_1 -y^{(1)}(\alpha_1^{new} - \alpha_1^{old}) \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} -y^{(2)}(\alpha_2^{new} - \alpha_2^{old}) \mathbf{x}^{(2)} \mathbf{x}^{(1)} \tag{31}

  2. 若是 0 < α 2 n e w < C 0 < \alpha_2^{new} < C
      按照上面的步骤,一样可求得 b 2 n e w b_2^{new}​
    (32) b 2 n e w = b o l d E 2 y ( 1 ) ( α 1 n e w α 1 o l d ) x ( 1 ) x ( 2 ) y ( 2 ) ( α 2 n e w α 2 o l d ) x ( 2 ) x ( 2 ) b_2^{new} = b^{old} - E_2 -y^{(1)}(\alpha_1^{new} - \alpha_1^{old}) \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} -y^{(2)}(\alpha_2^{new} - \alpha_2^{old}) \mathbf{x}^{(2)} \mathbf{x}^{(2)} \tag{32}

      若是 α 1 n e w \alpha_1^{new} α 2 n e w \alpha_2^{new} 同时知足条件 0 < α i n e w < C 0 < \alpha_i^{new} < C ,则 b n e w = b 1 n e w = b 2 n e w b^{new} = b_1^{new} = b_2^{new} ;若是 α 1 n e w \alpha_1^{new} α 2 n e w \alpha_2^{new} 0 0 或者 C C ,那么 b 1 n e w b_1^{new} b 2 n e w b_2^{new} 以及它们之间的数都是符合 KKT 条件的阈值,这是选择它们的重点做为 b n e w b^{new} 。则 b n e w b^{new}

(33) b n e w = { b 1 n e w , 0 < α 1 n e w < C b 2 n e w , 0 < α 2 n e w < C ( b 1 n e w + b 2 n e w ) 2 , o t h e r w i s e b^{new} = \left \{ \begin{array}{cc} & b_1^{new}, & 0 < \alpha_1^{new} < C \\\\ & b_2^{new}, & 0 < \alpha_2^{new} < C \\\\ & \frac{(b_1^{new} + b_2^{new})}{2}, & otherwise \end{array} \right. \tag{33}

1.1.8 变量更新算法的步骤

  根据上面讲解的内容,咱们来整理一下变量更新算法的步骤:

  • 步骤 1:计算偏差

E 1 = f ( x ( 1 ) ) y ( 1 ) = i = 1 m α i y ( i ) x ( i ) x ( 1 ) + b y ( 1 ) E_1 = f(\mathbf{x}^{(1)}) - y^{(1)} = \sum_{i=1}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(1)} + b - y^{(1)}

E 2 = f ( x ( 2 ) ) y ( 2 ) = i = 1 m α i y ( i ) x ( i ) x ( 2 ) + b y ( 2 ) E_2 = f(\mathbf{x}^{(2)}) - y^{(2)} = \sum_{i=1}^{m} \alpha_i y^{(i)} \mathbf{x}^{(i)} \cdot \mathbf{x}^{(2)} + b - y^{(2)}

  • 步骤 2:计算 α 2 n e w \alpha_2^{new} 取值范围

{ i f ( y ( 1 ) y ( 2 ) ) L = m a x ( 0 , α 2 o l d α 1 o l d ) , H = m i n ( C , C + α 2 o l d α 1 o l d ) i f ( y ( 1 ) = y ( 2 ) ) L = m a x ( 0 , α 2 o l d + α 1 o l d C ) , H = m i n ( C , α 2 o l d + α 1 o l d ) \left \{ \begin{array}{cc} \begin{aligned} & if(y^{(1)} \neq y^{(2)}) \quad L = max(0, \alpha_2^{old} - \alpha_1^{old}),\quad H = min(C, C + \alpha_2^{old} - \alpha_1^{old}) \\\\ & if(y^{(1)} = y^{(2)}) \quad L = max(0, \alpha_2^{old} + \alpha_1^{old} - C), \quad H = min(C, \alpha_2^{old} + \alpha_1^{old}) \end{aligned} \end{array} \right.

  • 步骤 3:计算 η \eta​

η = x ( 1 ) x ( 1 ) + x ( 2 ) x ( 2 ) 2 x ( 1 ) x ( 2 ) \eta = \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} + \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)} - 2\mathbf{x}^{(1)} \mathbf{x}^{(2)}

  • 步骤 4:求解 α 2 n e w , u n c l i p p e d \alpha_2^{new,unclipped}​

α 2 n e w , u n c l i p p e d = α 2 o l d + y ( 2 ) ( E 1 E 2 ) η \alpha_2^{new,unclipped} = \alpha_2^{old} + \frac{y^{(2)}(E_1-E_2)}{\eta}

  • 步骤 5:求解 α 2 n e w \alpha_2^{new}​

α 2 n e w = { H , α 2 n e w , u n c l i p p e d > H α 2 n e w , u n c l i p p e d , L α 2 n e w , u n c l i p p e d H L , α 2 n e w , u n c l i p p e d < L \alpha_2^{new} = \left \{ \begin{array}{cc} \begin{aligned} & H, & \alpha_2^{new,unclipped} > H \\\\ & \alpha_2^{new,unclipped}, & L \leq \alpha_2^{new,unclipped} \leq H \\\\ & L, & \alpha_2^{new,unclipped} < L \end{aligned} \end{array} \right.

  • 步骤 6:求解 α 1 n e w \alpha_1^{new}

α 1 n e w = α 1 o l d + y ( 1 ) y ( 2 ) ( α 2 o l d α 2 n e w ) \alpha_1^{new} = \alpha_1^{old} + y^{(1)} y^{(2)} (\alpha_2^{old} - \alpha_2^{new})

  • 步骤 7:求解 b 1 n e w b_1^{new} b 2 n e w b_2^{new}

b 1 n e w = b o l d E 1 y ( 1 ) ( α 1 n e w α 1 o l d ) x ( 1 ) x ( 1 ) y ( 2 ) ( α 2 n e w α 2 o l d ) x ( 2 ) x ( 1 ) b_1^{new} = b^{old} - E_1 - y^{(1)}(\alpha_1^{new} - \alpha_1^{old}) \mathbf{x}^{(1)} \cdot \mathbf{x}^{(1)} - y^{(2)}(\alpha_2^{new} - \alpha_2^{old}) \mathbf{x}^{(2)} \cdot \mathbf{x}^{(1)}

b 2 n e w = b o l d E 2 y ( 1 ) ( α 1 n e w α 1 o l d ) x ( 1 ) x ( 2 ) y ( 2 ) ( α 2 n e w α 2 o l d ) x ( 2 ) x ( 2 ) b_2^{new} = b^{old} - E_2 - y^{(1)}(\alpha_1^{new} - \alpha_1^{old}) \mathbf{x}^{(1)} \cdot \mathbf{x}^{(2)} - y^{(2)}(\alpha_2^{new} - \alpha_2^{old}) \mathbf{x}^{(2)} \cdot \mathbf{x}^{(2)}

  • 步骤 8:求解 b n e w b^{new}​

b n e w = { b 1 n e w , 0 < α 1 n e w < C b 2 n e w , 0 < α 2 n e w < C ( b 1 n e w + b 2 n e w ) 2 , o t h e r w i s e b^{new} = \left \{ \begin{array}{cc} & b_1^{new}, & 0 < \alpha_1^{new} < C \\\\ & b_2^{new}, & 0 < \alpha_2^{new} < C \\\\ & \frac{(b_1^{new} + b_2^{new})}{2}, & otherwise \end{array} \right.


  前面进行介绍时,咱们有假设选择两个变量 α 1 \alpha_1​ α 2 \alpha_2​ ,可是在实际实现算法时,咱们应该如何选择这两个变量呢?下面就来介绍一下。

1.2 变量 α 1 \alpha_1 α 2 \alpha_2 的选择方法

1.2.1 第 1 个变量的选择

  SMO 称选择第 1 个变量的过程为外层循环。外层循环在训练样本中选取违反 KKT 条件最严重的样本点,并将其对应的变量做为第 1 个变量。具体地,检验训练样本点是否知足 KKT 条件,即(在 支持向量机 01 - 线性可分支持向量机和线性支持向量机 中 2.2 节有进行介绍)

(34) α i = 0 y ( i ) f ( x ( i ) ) 1 0 < α i < C y ( i ) f ( x ( i ) ) = 1 α i = C y ( i ) f ( x ( i ) ) 1 \begin{aligned} \alpha_i = 0 \quad \Leftrightarrow \quad y^{(i)} f(\mathbf{x}^{(i)}) \geq 1 \\\\ 0 < \alpha_i < C \quad \Leftrightarrow \quad y^{(i)} f(\mathbf{x}^{(i)}) = 1 \\\\ \alpha_i = C \quad \Leftrightarrow \quad y^{(i)} f(\mathbf{x}^{(i)}) \leq 1 \end{aligned} \tag{34}

  从而可得知,如下几种状况将会不知足 KKT 条件

(35) y ( i ) f ( x ( i ) ) 1 a n d α i > 0 y ( i ) f ( x ( i ) ) = 1 a n d ( α i = 0 o r α i = C ) y ( i ) f ( x ( i ) ) 1 a n d α i < C \begin{aligned} & y^{(i)} f\left(\mathbf{x}^{(i)}\right) \geq 1 \quad and \quad \alpha_i > 0 \\\\ & y^{(i)} f\left(\mathbf{x}^{(i)}\right) = 1 \quad and \quad (\alpha_i = 0 \quad or \quad \alpha_i = C) \\\\ & y^{(i)} f\left(\mathbf{x}^{(i)}\right) \leq 1 \quad and \quad \alpha_i < C \end{aligned} \tag{35}

  在检验过程当中,外层循环首先遍历全部知足条件 0 < α i < C 0 < \alpha_i < C​ 的样本点,即在间隔边界上的支持向量点,检验它们是否知足 KKT 条件。若是这些样本点都知足 KKT 条件,那么遍历整个训练集,检验它们是否知足 KKT 条件。

1.2.2 第 2 个变量的选择

  SMO 称选择第 2 个变量的过程为内循环。假设在外层循环中已经找到第 1 个变量 α 1 \alpha_1 ,如今要在内层循环中找第 2 个变量 α 2 \alpha_2 。第 2 个变量选择的标准是但愿能使 α 2 \alpha_2 有足够大的变化。

  由式(18)可知, α 2 n e w \alpha_2^{new} 是依赖于 E 1 E 2 |E_1 - E_2| 的,为了加快计算速度,一种简单的作法是选择 α 2 \alpha_2 ,使其对应的 E 1 E 2 |E_1 - E_2| 最大。为了节省计算时间,能够将全部 E i E_i 值保存在一个列表中。

  在特殊状况下,若是内层循环经过以上方法选择的 α 2 \alpha_2 不能使目标函数有足够的降低(等价于 α 2 n e w α 2 |\alpha_2^{new} - \alpha_2| 很小),那么能够采用如下启发式规则继续选择 α 2 \alpha_2 。遍历在间隔边界上的支持向量点,依次将其对应的变量做为 α 2 \alpha_2 试用,知道目标函数有足够的降低。若仍然找不到合适的 α 2 \alpha_2 ,那么遍历训练数据集;若仍是找不到合适的 α 2 \alpha_2 ,则放弃第 1 个 α 2 \alpha_2 ,再经过外层循环寻求另外的 α 1 \alpha_1​

  由 Osuna 定理可知,只需选取的 α 1 \alpha_1 α 2 \alpha_2​ 中有一个不知足 KKT 条件,目标函数就会在迭代后减少。

1.3 SMO 算法的步骤

咱们总结一下 SMO 算法的步骤:

  • 步骤 1:初始化 α \mathbf{\alpha} b b ,好比初始化 α \mathbf{\alpha} 为零向量, b b 为 0。(注意这里的 α \mathbf{\alpha} 是一个列向量)
  • 步骤 2:选取优化变量 α 1 \alpha_1 α 2 \alpha_2 ,而后更新 α 1 \alpha_1 α 2 \alpha_2 b b
  • 步骤 3:若是知足中止条件(即前面提到的“全部变量的解都知足此最优化问题的 KKT 条件”)则结束;不然,跳到步骤 2。

2、Python 代码实现

  如下是 Python 3 的代码实现:

import numpy as np
import random
import matplotlib.pyplot as plt


class OptStruct:
    """ 数据结构,存储须要操做的数据 """
    def __init__(self, input_mat, label_mat, c, toler):
        """ :param input_mat: 样本特征值 :param label_mat: 样本标签值 :param c: 惩罚因子 :param toler: 容错率 """
        self.x = input_mat
        self.label_mat = label_mat
        self.c = c
        self.toler = toler
        self.m = np.shape(input_mat)[0]
        # 初始化 alphas 为零向量,初始化 b 为 0
        self.alphas = np.mat(np.zeros((self.m, 1)))
        self.b = 0
        # e_cache 用于缓存偏差值
        # 第一列给出 e_cache 是否有效的标志位,第二列给出实际的偏差值
        self.e_cache = np.mat(np.zeros((self.m, 2)))


class SVM:
    def __init__(self):
        self.alphas = None
        self.b = None
        self.w = None
        pass

    @staticmethod
    def calc_ek(os, k):
        """ 计算偏差 E 值 :param os: :param k: :return: """
        fxk = float(np.multiply(os.alphas, os.label_mat).T *\
                    (os.x * os.x[k, :].T)) + os.b
        ek = fxk - float(os.label_mat[k])
        return ek

    @staticmethod
    def select_j_rand(i, m):
        """ 随机选择第二个 alpha :param i: 第一个 alpha 对应的索引值 :param m: alpha 参数的个数(即全部训练样本的总数) :return: 第二个 alpha 对应的索引值 """
        j = i
        while j == i:
            j = int(random.uniform(0, m))
        return j

    @staticmethod
    def select_j(i, os, ei):
        """ 选择第二个 alpha :param i: 第一个 alpha 的索引值 :param os: 数据结构 :param ei: 第一个 alpha 对应样本点的偏差 :return: 返回第二个 alpha 的索引值 和 对应的偏差值 """
        max_k = -1
        # 保存最大的 |ei - ej|
        max_delta_e = 0
        ej = 0
        # 将 ei 更新到缓存 e_cache 中
        os.e_cache[i] = [1, ei]
        # 返回 os.e_cache 中 e 为非零的索引值列表
        valid_e_cache_list = np.nonzero(os.e_cache[:, 0].A)[0]
        # 第一次执行这个函数时,会执行 else 语句
        if (len(valid_e_cache_list)) > 1:
            # 经过遍历找到使 |ei - ej| 最大的第二个 alpha
            for k in valid_e_cache_list:
                # 若是是其自己,不进行比较
                if k == i:
                    continue
                ek = SVM.calc_ek(os, k)
                delta_e = abs(ei - ek)
                if delta_e > max_delta_e:
                    max_k = k
                    max_delta_e = delta_e
                    ej = ek
            return max_k, ej
        else:
            # 随机选择第二个 alpha
            j = SVM.select_j_rand(i, os.m)
            ej = SVM.calc_ek(os, j)
        return j, ej

    @staticmethod
    def update_ek(os, k):
        """ 更新偏差缓存 :param os: 数据结构 :param k: 须要更新偏差项对应的索引值 """
        ek = SVM.calc_ek(os, k)
        os.e_cache[k] = [1, ek]

    @staticmethod
    def clip_alpha(alpha_j, high, low):
        """ 修剪 alpha_j :param alpha_j: alpha_j 的未修剪的值 :param high: alpha_j 的上限 :param low: alpha_j 的下限 :return: alpha_j 修剪后的值 """
        if alpha_j > high:
            alpha_j = high
        elif alpha_j < low:
            alpha_j = low
        return alpha_j

    def inner_l(self, i, os):
        """ 选择第二个 alpha,并更新 alpha 和 b :param i: 第一个 alpha 对应样本点的索引值 :param os: 数据结构 :return: 是否有更新 alpha 对 """
        # 步骤 1:计算第一个 alpha 对应的样本点的偏差
        ei = SVM.calc_ek(os, i)
        # 若是 (yi * f(xi) < 1 && alpha_i < c) 或者 (yi * f(xi) > 1 && alpha_i > 0)
        # 说明不知足 KKT 条件,能够将其做为第一个 alpha
        if ((os.label_mat[i] * ei < -os.toler) and (os.alphas[i] < os.c)) or \
                ((os.label_mat[i] * ei > os.toler) and (os.alphas[i] > 0)):
            # 步骤 1:选择第二个 alpha,并计算偏差
            j, ej = SVM.select_j(i, os, ei)
            alpha_iold = os.alphas[i].copy()
            alpha_jold = os.alphas[j].copy()
            # 步骤 2:求取 alphas_j 取值的上下限
            if os.label_mat[i] != os.label_mat[j]:
                low = max(0, os.alphas[j] - os.alphas[i])
                high = min(os.c, os.c + os.alphas[j] - os.alphas[i])
            else:
                low = max(0, os.alphas[j] + os.alphas[i] - os.c)
                high = min(os.c, os.alphas[j] + os.alphas[i])
            if low == high:
                print('low == high')
                return 0
            # 步骤 3:计算 2*x1*x2 - x1*x2 - x2*x2
            eta = 2.0 * os.x[i, :] * os.x[j, :].T - os.x[i, :] * os.x[i, :].T -\
                os.x[j, :] * os.x[j, :].T
            if eta >= 0:
                print('eta >= 0')
                return 0
            # 步骤 4:求解未修剪的 alpha_2
            os.alphas[j] -= os.label_mat[j] * (ei - ej) / eta
            # 步骤 5:求解修剪后的 alpha_2
            os.alphas[j] = SVM.clip_alpha(os.alphas[j], high, low)
            SVM.update_ek(os, j)
            if abs(os.alphas[j] - alpha_jold) < 0.00001:
                print('j not moving enough')
                return 0
            # 步骤 6:求解 alpha_1
            os.alphas[i] += os.label_mat[j] * os.label_mat[i] * \
                            (alpha_jold - os.alphas[j])
            SVM.update_ek(os, i)
            # 步骤 7:求解 b1 和 b2
            b1 = os.b - ei - os.label_mat[i] * (os.alphas[i] - alpha_iold) * \
                os.x[i, :] * os.x[i, :].T
相关文章
相关标签/搜索