機器學習系列 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
相關文章
相關標籤/搜索