機器學習算法——支援向量機（SVM）
0 前序
本人是生物類的一名研究所學生，學習Python也有一段時間了，從資料處理到人工智能，機器學習， Python幾乎是無所不能的，用其處理生信資料想必也是很好的工具。這是本人寫的第一篇部落格，某法則曾說：能把問題清楚地寫下來，就已解決了一半。是以打算把學習過程中遇到問題記錄下來，以便日後回顧解決，同時打算把最近幾天的學習心得寫下來，也算是對期間學習的一個總結。
1 介紹
SVM（Support Vecor Machine）是一個二進制分類算法，不僅可以用于資料分類，也可以應用于回歸問題。對于沒有高數基礎的同學（比如我），在了解SVM之前，還需了解一下基本的高數常識，比如點到面的距離、拉格朗日乘子法等，以便後續了解。
1.1 距離公式介紹
點到直線的距離公式推導：
其中向量n=(A,B,C)的模:
1.2 拉格朗日乘子法
關于拉格朗日乘子法在知乎上看到了一個很直覺的解釋：https://www.zhihu.com/question/38586401/answer/105588901 ，可供參考。
目标函數f(x,y)是一座山的高度，限制條件g(x,y)=C是鑲嵌在山上的一條曲線如下圖：
為了找到曲線上的最低點，就從最低的等高線（0那條）開始網上數。數到第三條，等高線終于和曲線有交點了（如上圖所示）。因為比這條等高線低的地方都不在限制範圍内，是以這肯定是這條限制曲線的最低點了。
那麼什麼是等高線呢？下圖還是比較直覺的。
是以當等高線與限制條件相切時，兩者的法向量相平行，是以可以列出二者的關系式：
我們引入任意的常數乘子（取為λ）：, 我們把這個式子的右邊移到左邊，并把常數移進微分算子，就得到 ：
是以根據上式構造輔助函數：
然後令
就可求出極值點（ 意為求偏導）。
2 線性可分支援向量機
對于線性支援向量機我們需要考慮的是一個二分類問題，即資料是線性可分的。如下圖所示：
對于一個資料集，我們可以構造一個超平面（在二維空間裡表示為直線），将資料集分為兩類。如圖中紫色與黃色的點，每個點都為一個資料點，用X(x,y)表示，該資料集有明顯的區分度，圖中三條超平面（直線）都可以将其區分開，甚至還有無數個超平面（隻是未畫出），可将其區分，但哪一個超平面的區分效果更好呢？超平面滿足方程：
其中x為資料點， w T x w^T x wTx為超平面的法向量。令 f ( x ) = w T x + b f(x)=w^T x+b f(x)=wTx+b ，若f(x)=0，則表示資料點落在超平面上，無法被分類；若f(x)>0，其對應的類别為 y i = + 1 y_i=+1 yi=+1；相反若f(x)<0，其對應類别 y i = − 1 y_i=-1 yi=−1（取±1為了友善運算），是以無論f(x)是否大于0，都有 y i f ( x ) > 0 y_i f(x)>0 yif(x)>0。
根據點到平面的距離公式可得，資料點到超平面的間隔為：
間隔如下圖灰色部分，可以看出綠線和藍線分類效果沒有橙線效果好。
由于 ( w T x + b ) ⋅ y i = y i f ( x ) > 0 (w^T x+b)·y_i=y_i f(x)>0 (wTx+b)⋅yi=yif(x)>0。我們可以成倍的縮放 w w w和 b b b使得 y i f ( x ) ≥ 1 y_i f(x)≥1 yif(x)≥1，我們的目标是要找到一個超平面，使得離超平面最近的資料點間的間隔最大。即：
在支援向量機中，距離超平面最近的且滿足一定條件的幾個樣本點稱之為支援向量。
因：
即求
要使得 1 / ‖ w ‖ 1/‖w‖ 1/‖w‖最大，那麼 ‖ w ‖ ‖w‖ ‖w‖就要最小，上式展開轉化為求：
對于這個式子，我們回顧上面的拉格朗日乘子法發現，就和其标準格式一樣
将①轉化為拉格朗日乘子法形式：
因為引入了拉格朗日乘子α，我們優化目标就轉化為要求：
由于我們求α的極大值并不容易，是以我們可以根據對偶問題，将
這樣求 L ( w , b , a ) L(w,b,a) L(w,b,a)的極小值就比較容易了，隻要對 w w w和 b b b求偏導就可以得到結果。
将③④帶入②中，計算過程省略，最後得出結果：
接着對φ(α)求極大值，将1/2的負号提出，負數的極大值即為正數的極小值，等價：
其中 x i T x j x_i^T x_j xiTxj是兩個資料點的内積。以上即為線性支援向量機目标函數建構推導過程了，最後根據已知的資料點 （ x i , y i ） （x_i,y_i） （xi,yi）計算：
的極小值，通常采用SMO算法進行求解，求出 α α α後，帶入即可以得到 w w w和 b b b，以及分隔超平面 w T x + b = 0 w^T x+b=0 wTx+b=0。
SMO算法可見：https://www.jianshu.com/p/eef51f939ace
3 軟間隔線性支援向量機
假如一個資料點絕大部分都是可分的，但有一些異常點落入了超平面的間隔當中，如圖：
如果硬要将其分開，可能會導緻超平面的間隔小于1，無法使分類效果達到最佳。是以為了照顧到絕大多數點都能被正确分類，我們要容忍這些異常的點，出現在不該不出現的地方。是以為了使得異常點的函數間隔 ( w T x + b ) ⋅ y i (w^T x+b)·y_i (wTx+b)⋅yi也能≥1，引入了松弛因子 ξ i ξ_i ξi，即： ( w T x + b ) ⋅ y i + ξ i ≥ 1 (w^T x+b)·y_i+ξ_i≥1 (wTx+b)⋅yi+ξi≥1
原始的目标函數也會因容忍這些異常點而調整，需在原始的目标函數後面加上因為容忍異常點産生的代價，新的目标函數為：
其中C為懲罰系數，C越大，對誤分類的懲罰越大，而函數模型為了避免懲罰，就會盡可能的減少誤分類的樣本點，是以函數的帶寬（間隔）就會越窄；相反，C越小，對誤分類的懲罰就越小，函數模型可以較多地包容那些誤分類的點，函數的間隔也就會越寬。就像父母教育小孩，因為懲罰系數大，小孩盡可能的不犯任何錯誤，但這會導緻“過拟合”的風險；但如果過分溺愛，小孩容易無法無天。
其步驟和前面的線性SVM類似，現将該目标函數用拉格朗日乘子法表示：
其中 μ i ≥ 0 μ_i≥0 μi≥0， a i ≥ 0 a_i≥0 ai≥0，均為拉格朗日乘子系數
接着對 L ( w , b , a , ξ i , μ i ) L(w,b,a,ξ_i,μ_i ) L(w,b,a,ξi,μi)中的參數 w 、 b 、 ξ i w、b、ξ_i w、b、ξi求偏導，并使導函數為0：
将⑥⑦⑧帶入 L ( w , b , a , ξ i , μ i ) L(w,b,a,ξ_i,μ_i ) L(w,b,a,ξi,μi)式中，得到：
我們發現最後結果和未引入松弛因子的線性SVM是一樣的，隻是限制條件發生了改變。因為 C − α i = μ i ≥ 0 C-α_i=μ_i≥0 C−αi=μi≥0，我們可以得到： 0 ≤ α i ≤ C 0≤α_i≤C 0≤αi≤C，整理完限制條件軟間隔線性SVM的最終結果為：
最後我們可以通過SMO算法，求出 w w w和 b b b。
不過在進行資料分類時，為了使得每個樣本點的函數間隔大于等于1，我們引入了松弛因子 ξ i ξ_i ξi。如果所有樣本點的函數間隔都大于等于1，那麼采用線性可分SVM就能正确分類；如果存在樣本點的函數間隔小于1，這時則需要借助松弛因子 ξ i ξ_i ξi，讓其間隔 d + ξ i ≥ 1 d+ξ_i≥1 d+ξi≥1，間隔越小則 ξ i ξ_i ξi越大，對應函數模型的損失也越大。如果 ξ i = 0 ξ_i=0 ξi=0，表示該樣本點的間隔為1，在邊界上為支援向量，如下圖所示X1；如果 0 ＜ ξ i ＜ 1 0＜ξ_i＜1 0＜ξi＜1，表示該樣本點在超平面和自己類别的間隔邊界之間，如X2；如果 ξ i = 1 ξ_i=1 ξi=1，表示樣本點在超平面上，無法被分類；如果 ξ i ＞ 1 ξ_i＞1 ξi＞1，表示樣本點的間隔為負數，即在超平面的另一側，說明這個樣本點無法被正确分類，如X3。
4 非線性可分支援向量機
前兩節所針對的資料均為線性可分的，采用線性超平面就可以很好地将樣本點區分開來。但采用線性超平面分割非線性可分資料，将産生很差的效果，如下圖：
這時，就需要想樣本點資料映射到高維空間中，最終總會在一個高維空間中有一個超平面使其線性可分。是以建立非線性可分SVM需要：（1）将原始資料映射到高維；（2）在高維空間中找到一個超平面，将樣本點分類。
我們将樣本點x映射到高維空間的變化設為 ϕ ( x i ) ϕ(x_i ) ϕ(xi)，這種映射在二維空間中即為内積，就是前文的 x i T x j x_i^T x_j xiTxj。非線性SVM的目标函數可以表示為：
對于 ϕ ( x i ) ⋅ ϕ ( x j ) ϕ(x_i )·ϕ(x_j ) ϕ(xi)⋅ϕ(xj)可以用核函數 K ( x i , x j ) K(x_i,x_j ) K(xi,xj)表示，即 K ( x i , x j ) = ϕ ( x i ) ⋅ ϕ ( x j ) K(x_i,x_j )= ϕ(x_i )·ϕ(x_j ) K(xi,xj)=ϕ(xi)⋅ϕ(xj)。對于這個函數 K ( x i , x j ) K(x_i,x_j ) K(xi,xj)，我們稱之為核函數。核函數低維上進行計算，而将實質上的分類效果（利用了内積）表現在了高維上，這樣避免了直接在高維空間中的複雜計算，真正解決了SVM線性不可分的問題。
前人已經發現了很多核函數，在實際應用中可供我們選擇，具體如下：
（1）線性核函數
線性核函數實際上就是線性可分SVM的函數模型，但其缺點是隻能解決線性可分的問題。
（2） 多項式核函數
（3） 高斯核函數
對于非線性可分的資料常常使用高斯核函數，因為高斯核函數是指數函數，可以将低維資料映射到無窮維。下圖即為對上面的資料集采用高斯核函數的分類結果：
5 一個執行個體
利用sklearn提供的svm方法對于鸢尾花資料集進行分類，也算是一個很經典的實戰案例。
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import os
%matplotlib inline
#os.getcwd()
#os.chdir('C:\\Users\\HP')
df=pd.read_csv('iris.csv')
#print(df)
#添加列标簽
df.columns = ['s_len', 's_wid', 'p_len', 'p_wid', 'class']
df.head()
#資料描述
#df.describe()
輸出的是該資料集的前五條資料，s_len和s_wid分别為鸢尾花的花萼長度和寬度，p_len和p_wid為鸢尾花花瓣的長度和寬度，class為鸢尾花的種類，該資料中有三種類别。
#可視化，畫矩陣圖
import seaborn as sns
sns.set()
sns.pairplot(df,hue='class',size=2.5)
從圖中可以看出，在花瓣長度寬度上，三種花的差別還是比較大的，是以采用p_len和p_wid這兩列資料進行分析。
#小提琴圖
plt.figure(figsize=(20, 14))
for column_index, column in enumerate(df.columns):
if column == 'class':
continue
plt.subplot(2, 2, column_index + 1)
sns.violinplot(x='class', y=column, data=df)
同樣是可視化一種對資料的表現形式
from sklearn.cross_validation import train_test_split
X = df[['s_len', 's_wid','p_len', 'p_wid']].values
y = df['class'].values
X=X[:,2:4]
x=X
y = pd.Categorical(y).codes #将類别資訊轉化為數值資訊
#print(y)
(x_train, x_test, y_train, y_test) = train_test_split(X,y,train_size=0.7, random_state=0)
#train_test_split為簡單的交叉驗證：随機将資料劃分為兩部分，訓練集和測試集。一般 70%的資料為訓練集，30%為測試集。
#random_state填1的話，其他參數一樣的情況下你得到的随機數組是一樣的。但填0或不填，每次都會不一樣。
#print(X_train)
x_train，y_train為那70%的資料集，用來給機器學習。
x_test，y_test為用來測試機器學習準不準确的那30%資料，其中y_test是機器根據之前的學習然後預測的值。
#建構函數模型
from sklearn.svm import SVC
clf=SVC(kernel='rbf',gamma=0.1,decision_function_shape='ovo',C=0.8)
# kernel='rbf'（default）時，為高斯核函數，gamma值越小，分類界面越連續；gamma值越大，分類界面越“散”，分類效果越好，但有可能會過拟合。
# decision_function_shape='ovo'時，為one v one分類問題，即将類别兩兩之間進行劃分，用二分類的方法模拟多分類的結果。
# decision_function_shape='ovr'時，為one v rest分類問題，即一個類别與其他類别進行劃分。
采用高斯核函數
#訓練
from sklearn.metrics import accuracy_score
clf.fit(x_train,y_train)
print ('測試集準确率：',accuracy_score(y_test,clf.predict(x_test)))
輸出結果：測試集準确率： 0.9555555555555556
接下來對預測的資料進行繪圖，首先建立網格
# 建立網格，以繪制圖像
N=500
x1_min, x2_min = x[:,0].min(),x[:,1].min()
x1_max, x2_max = x[:,0].max(),x[:,1].max()
#print(x1_min, x2_min)
#print(x1_max, x2_max)
t1 = np.linspace(x1_min, x1_max, N)
t2 = np.linspace(x2_min, x2_max, N)
x1, x2 = np.meshgrid(t1, t2) # 生成網格采樣點
grid_show = np.dstack((x1.flat, x2.flat))[0] # 測試點
#print(grid_show)
grid_hat = clf.predict(grid_show) # 預測分類值
grid_hat = grid_hat.reshape(x1.shape) # 使之與輸入的形狀相同
#print(grid_hat)
#準備繪圖
import matplotlib as mpl
cm_light = mpl.colors.ListedColormap(['yellow', 'pink', '#A0A0FF'])#背景顔色
cm_dark = mpl.colors.ListedColormap(['g', 'r', 'b'])#點的顔色
plt.figure(facecolor='w')
定義各label值
plt.pcolormesh(x1, x2, grid_hat, cmap=cm_light)#畫背景圖
plt.scatter(x[:,0], x[:,1], c=y, edgecolors='k', s=66, cmap=cm_dark) # 樣本
plt.scatter(x_test[0], x_test[1], s=120, facecolors='none', zorder=10) # 圈中測試集樣本
plt.xlabel('petal_length', fontsize=13)#坐标軸名稱
plt.ylabel('petal_width', fontsize=13)
plt.xlim(x1_min, x1_max)#坐标軸範圍
plt.ylim(x2_min, x2_max)
plt.title('Feature classification of iris by SVM', fontsize=20)
#plt.grid()
plt.show()
最後結果如圖：
在資料分類可視化這一方面搞了半天，還有待加強^ ^