SVM
支持向量機(jī)（support vector machines, SVM）,是一種二分類模型姨拥， 屬于判別模型。
基本原理是在特征空間中找到一個(gè)間隔最大的分離超平面，將正負(fù)類分開疯暑。
1. 基本模型定義：
（硬間隔）支持向量機(jī)：
假設(shè)給定線性可分訓(xùn)練集
通過間隔最大化或者等價(jià)地求解相應(yīng)的凸二次規(guī)劃問題，從而學(xué)習(xí)到的分隔超平面為：
以及相應(yīng)的分類決策函數(shù)為
以上模型即（硬間隔）支持向量機(jī)
其中
均為模型參數(shù)
表示分離超平面的法向量
為 分離超平面的截距
支持向量機(jī)的目標(biāo)是最大化分超平面與最近的二分類點(diǎn)的距離
樣本空間中任意一點(diǎn)到分離超平面的距離推導(dǎo)：
設(shè)一點(diǎn)
, 求它到超平面
的距離哑舒，
設(shè)上的投影點(diǎn)為
, 則有
妇拯，
同時(shí)向量與 超平面的法向量
平行，故有洗鸵，
其中 d 是我們要求的距離越锈，而
代入得到
即
因此 支持向量機(jī)的目標(biāo)就是 通過特征空間，找到來最大化這個(gè)距離d
2. 函數(shù)間隔和幾何間隔
函數(shù)間隔定義為
對(duì)于給定的訓(xùn)練數(shù)據(jù)集
和超平面
, 定義樣本點(diǎn)
到超平面的函數(shù)間隔為：
PS：定義訓(xùn)練數(shù)據(jù)集
的函數(shù)間隔為所有樣本點(diǎn)
的函數(shù)間隔中的最小值 即
幾何間隔定義為：
對(duì)于給定的訓(xùn)練數(shù)據(jù)集
定義訓(xùn)練數(shù)據(jù)集
所以 當(dāng)||w|| = 1 時(shí) 函數(shù)間隔與幾何間隔等價(jià)瞪浸。
3. 幾何間隔最大化
支持向量機(jī)的基本目標(biāo)是找到能夠正確劃分?jǐn)?shù)據(jù)集并且最大化幾何間隔的超平面。
和感知機(jī)原理不同吏祸，感知機(jī)中能夠正確劃分?jǐn)?shù)據(jù)集的分離超平面可能會(huì)有很多对蒲。
SVM中只有一個(gè)符合目標(biāo)的超平面
數(shù)學(xué)表示為約束最優(yōu)化問題：
考慮到幾何間隔和函數(shù)間隔的關(guān)系, 有：
假設(shè)將和
按比例改變?yōu)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Clambda%20w" alt="\lambda w" mathimg="1"> 和
,此時(shí) 函數(shù)間隔將變?yōu)?
,因此函數(shù)間隔
的取值對(duì)于上述最優(yōu)化問題沒有影響贡翘。
不妨取
意識(shí)到最大化 和最小化
是等價(jià)的
則上述最優(yōu)化問題 變?yōu)?/p>
這是個(gè)凸二次規(guī)劃問題蹈矮。
可以采用拉格朗日乘子法，將其變換成對(duì)偶問題鸣驱。
4. 拉格朗日法求解凸二次規(guī)劃問題
定義拉格朗日函數(shù)為：
其中
為拉格朗日乘子向量
分別對(duì) 并令其 = 0得：
得
代回拉格朗日函數(shù)
即轉(zhuǎn)換后得對(duì)偶問題為
將目標(biāo)函數(shù)由求極大轉(zhuǎn)為求極小泛鸟，得到與之等價(jià)的對(duì)偶問題為
相應(yīng)的KKT條件為：
5. 軟間隔支持向量機(jī)
如果訓(xùn)練樣本含有噪聲，這個(gè)時(shí)候支持向量機(jī)可能就不存在這個(gè) 超平面可以將兩個(gè)樣本分開或者達(dá)到比較好的最大間隔踊东，為了解決這個(gè)問題引入 一個(gè)新概念：軟間隔北滥。
即：
同時(shí)我們的約束變?yōu)?/p>
?
同之前硬間隔一樣，變換為對(duì)偶問題闸翅，最后變成：
所以我們只要求解滿足約束的最優(yōu)解的 出來就計(jì)算
就可以找到超平面了再芋。
6. 核函數(shù)
實(shí)際上到此我們已經(jīng)掌握了SVM，但是如果我們就這樣直接去尋找最優(yōu)解坚冀，速度會(huì)非常慢, 這里的核技巧 和 后面序列最優(yōu)化算法（SMO）都是為了改善SVM的速度而提出的济赎。
?
上式中，如果直接計(jì)算，對(duì)于整個(gè)訓(xùn)練集而言司训，計(jì)算起來會(huì)有巨大的時(shí)間開銷构捡，核函數(shù)可以減少這里時(shí)間開銷。
即向量 的點(diǎn)積 可以用
來表示
比如高斯核：
因此使用的時(shí)候我們把 替換成
就好
7. 序列最優(yōu)化算法
將上一部分的核函數(shù)替換掉點(diǎn)積壳猜，我們得到:
我們需要做的就是計(jì)算
由于 ,直接求解N個(gè)未知量難度很大
序列最優(yōu)化算法就是為了解決這個(gè)問題勾徽。
其核心思想就是，先將視為變量, 將其他的
全部當(dāng)成參數(shù)捂蕴，然后迭代使得
滿足KKT條件。
由第一個(gè)約束有：
所以將原式中的變量顯式寫出闪幽，最終整理得如下：
根據(jù)書中表述啥辨，假設(shè)上述問題得初始可行解為 ，最優(yōu)解為
,則有：
其中和
是
所在對(duì)角線端點(diǎn)的界盯腌。
如果 則
如果,則
假設(shè)在沿著約束方向未剪輯時(shí)的最優(yōu)解為
為了表示 ,記
8.變量的選擇方法
第一個(gè)變量的選擇
從訓(xùn)練樣本中選擇違反KKT條件最嚴(yán)重的點(diǎn) 即檢查KKT條件
第二個(gè)變量的選擇
選擇最大化的點(diǎn)。
9. python代碼實(shí)現(xiàn)
實(shí)現(xiàn)遵循文中表述
import time
import numpy as np
import math
import random
'''
數(shù)據(jù)集: MNIST
訓(xùn)練集大小: 60000（實(shí)際使用1000）
測(cè)試集大小: 10000（實(shí)際使用100）
-----
運(yùn)行結(jié)果：
the accuracy is : 0.92
time span: 62.679516553878784 s
'''
def loadData(fileName):
'''
加載數(shù)據(jù)集
:param fileName: 文件路徑
:return: dataList腕够，labelList分別為特征集X和標(biāo)記Y. 均為list
'''
dataList = [] #
labelList = []
f = open(fileName, 'r')
for line in f.readlines():
curline = line.strip().split(',')
'''
這里考慮到我用的文件是csv格式级乍，所以用split(',')
Mnsit有0-9十個(gè)標(biāo)記
文件每行開頭第一個(gè)數(shù)字為該行的label標(biāo)記
這里為了簡(jiǎn)化，限定二分類任務(wù)帚湘，所以將標(biāo)記0的作為1(正例)玫荣，其余為0(反例)
'''
if (int(curline[0]) == 0):
labelList.append(1)
else:
labelList.append(0)
'''
加入特征集X
由于開頭第一個(gè)數(shù)字為標(biāo)記，故從下標(biāo) 1 開始
這里轉(zhuǎn)為int類型
/255 是為了歸一化大诸，有效減少數(shù)字爆炸捅厂。
'''
dataList.append([int(num) / 255 for num in curline[1:]])
# 讀取完畢關(guān)閉文件
f.close()
# 返回?cái)?shù)據(jù)集和標(biāo)記
return dataList, labelList
class SVM_CLF:
def __init__(self, X_train,y_train, sigma = 10, C = 200, toler = 0.001):
"""
相關(guān)參數(shù)初始化
:param X_train: 訓(xùn)練數(shù)據(jù)集
:param y_train: 訓(xùn)練數(shù)據(jù)集標(biāo)記
:param sigma: 高斯核中的\sigma
:param C: 懲罰項(xiàng)系數(shù)
:param toler: 松弛變量
"""
self.trainData = np.mat(X_train)
self.trainLabel = np.mat(y_train).T
self.m, self.n = np.shape(self.trainData) #訓(xùn)練及大小，循環(huán)中會(huì)用到
self.sigma = sigma #高斯核參數(shù)
self.C = C
self.toler = toler
self.k = self.calcKernel() #初始化時(shí)提前計(jì)算高斯核
self.b = 0 #初始化偏置項(xiàng)為0
self.alpha = [0] * self.m
self.E = [0 * self.trainLabel[i,0] for i in range(self.m)]
self.supportVectorIndex = []
def calcKernel(self):
"""
計(jì)算核函數(shù)资柔，本例使用的是高斯核
:return: 高斯核矩陣
"""
#初始化高斯核矩陣 大小=[m * m] 焙贷，m為訓(xùn)練集樣本數(shù)量
k_matrix = np.zeros((self.m, self.m))
for i in range(self.m):
X = self.trainData[i,:]
for j in range(i, self.m):
Z = self.trainData[j,:]
XZ = (X - Z) * (X - Z).T
result = np.exp(-1 * XZ / (2 * self.sigma **2))
k_matrix[i][j] = result
k_matrix[j][i] = result
return k_matrix
def isSatisfyKKT(self, i):
"""
判斷alpha i 是否符合KKT條件
i為下標(biāo)
:return: True (滿足)或者 False (不滿足)
"""
gxi = self.calc_gxi(i)
yi = self.trainLabel[i]
if (math.fabs(self.alpha[i] < self.toler)) and (yi * gxi >= 1):
return True
elif(math.fabs(self.alpha[i] - self.C) < self.toler) and (yi * gxi )<=1:
return True
elif(math.fabs(self.alpha[i]) > -self.toler) and (self.alpha[i] < self.C)\
and (math.fabs(yi * gxi) - 1) < self.toler:
return True
return False
def calc_gxi(self, i):
gxi = 0
index = [i for i,alpha in enumerate(self.alpha) if alpha != 0]
for j in index:
gxi += self.alpha[j] * self.trainLabel[j] * self.k[j][i]
gxi += self.b
return gxi
def calc_EI(self, i):
'''
計(jì)算Ei
根據(jù)書中章節(jié)7.4.1 兩個(gè)變量二次規(guī)劃的求解方法 中 式7.105
:param i: E的下標(biāo)
:return: E2,和下標(biāo)
'''
gxi = self.calc_gxi(i)
return gxi - self.trainLabel[i]
def get_Alpha2(self, E1, i):
'''
選擇第二個(gè)alpha
:param E1:
:param i:
:return:
'''
E2 = 0
maxE1_E2 = -1
E2_index = -1
noZero_E = [i for i ,Ei in enumerate(self.E) if Ei != 0]
for j in noZero_E:
E2_tmp = self.calc_EI(j)
if (math.fabs(E1 - E2_tmp) > maxE1_E2):
maxE1_E2 = math.fabs(E1 - E2_tmp)
E2 = E2_tmp
E2_index = j
if E2_index == -1:
E2_index = i
while E2_index == i:
E2_index = random.randint(0, self.m)
E2 = self.calc_EI(E2_index)
return E2, E2_index
def train(self, Maxiter = 100):
'''
用訓(xùn)練數(shù)據(jù)集學(xué)習(xí)模型
:param Maxiter: 最大迭代次數(shù)
:return: w， b
'''
iterstep = 0
is_Alpha_changed = 1
while(iterstep < Maxiter) and (is_Alpha_changed > 0):
print('iterstep:%d:%d'%(iterstep, Maxiter))
iterstep += 1
is_Alpha_changed = 0
for i in range(self.m):
if self.isSatisfyKKT(i) == False:
E1 = self.calc_EI(i)
E2, j = self.get_Alpha2(E1, i )
y1 = self.trainData[i]
y2 = self.trainLabel[j]
alpha_old_1 = self.alpha[i]
alpha_old_2 = self.alpha[j]
if y1 != y2:
L = max(0, alpha_old_2 - alpha_old_1)
H = min(self.C, self.C+ alpha_old_2 - alpha_old_1)
else:
L = max(0, alpha_old_2 + alpha_old_1 - self.C)
H = min(self.C, alpha_old_2 + alpha_old_1)
#如果L 和 H相等 贿堰，說明變量無法優(yōu)化辙芍， 直接跳過
if L == H: continue
#根據(jù)書中式7.106
# 計(jì)算alpha的新值
K11 = self.k[i][i]
K22 = self.k[j][j]
K12 = self.k[i][j]
alpha_new_2 = alpha_old_2 + y2 * (E1 - E2) / (k11 + k22 - 2 * k12)
#剪切alpha2
if alpha_new_2 < L : alpha_new_2 = L
elif alpha_new_2 > H: alpha_new_2 = H
#根據(jù) 7.109式 更新alpha1
alpha_new_1 = alpha_old_1 + y1 * y2 * (alpha_old_2 - alpha_new_2)
#根據(jù)書中“7.4.2 變量的選擇方法” 第三步式 7.115 和 7.116
b1_new = -E1 - y1 * K11 * (alpha_new_1 - alpha_old_1) \
- y2* K21 * (alpha_new_2 - alpha_old_2) + self.b
b2_new = -E2 - y1 * K12 * (alpha_new_1 - alpha_old_1) \
- y2* K22 * (alpha_new_2 - alpha_old_2) + self.b
if(alpha_new_1 > 0) and (alpha_new_1 < self.C):
b_new = b1_new
elif(alpha_new_2 > 0) and (alpha_new_2 < self.C):
b_new = b2_new
else:
b_new = (b1_new + b2_new) / 2
#將各參數(shù)更新
self.alpha[i] = alpha_new_1
self.alpha[j] = alpha_new_2
self.b = b_new
self.E[i] = self.calc_EI(i)
self.E[j] = self.calc_EI(j)
if math.fabs(alpha_new_2 - alpha_old_2) > 0.0001:
is_Alpha_changed += 1
#print("iter: %d i: %d, pairs changed %d" %(iterstep, i, is_Alpha_changed))
#全部計(jì)算結(jié)束后遍歷一遍alpha ,找出支持向量。
for i in range(self.m):
if self.alpha[i] > 0:
self.supportVectorIndex.append(i)
def calc_Single_Kernel(self, x1, x2):
'''
單獨(dú)計(jì)算核函數(shù)
在predict的時(shí)候用到
:param x1:
:param x2:
:return:
'''
x1_x2 = (x1 - x2) * (x1 - x2).T
result = np.exp(- x1_x2 / (2 * self.sigma ** 2))
return result
def predict(self, x):
'''
預(yù)測(cè)樣本
:param x: 預(yù)測(cè)的樣本
:return: 預(yù)測(cè)結(jié)果
'''
result = 0
for i in self.supportVectorIndex:
tmp = self.calc_Single_Kernel(np.mat(x), self.trainData[i,:])
result += self.alpha[i] * self.trainLabel[i] * tmp
result += self.b
return np.sign(result)
def test(self, testData, testLabel):
'''
測(cè)試測(cè)試集
:param testData: 測(cè)試數(shù)據(jù)集
:param testLabel: 測(cè)試標(biāo)記集
:return: 正確率
'''
errorCount = 0
for i in range(len(testData)):
print('testing : %d: %d'%(i, len(testData)))
result = self.predict(testData[i])
if result != testLabel[i]:
errorCount += 1
return 1 - errorCount/len(testLabel)
if __name__ == '__main__':
start = time.time()
print('read trainSet')
trainData, trainLabel = loadData('D:/PythonLearn/MLA/mnist_train.csv')
print('read testSet')
testData, testLabel = loadData('D:/PythonLearn/MLA/mnist_test.csv')
print('init SVM ')
svm = SVM_CLF(trainData[:1000], trainLabel[:1000], sigma=10, C = 200,toler = 0.001)
print('training')
svm.train()
print('testing')
accuracy = svm.test(testData[:100], testLabel[:100])
end = time.time()
print('the accuracy is :', accuracy)
print('time span:', end - start, 's')
