目錄
0. 前言
1. 尋找最大間隔
2. 拉格朗日乘子法和KKT條件
3. 松弛變量
4. 帶松弛變量的拉格朗日乘子法和KKT條件
5. 序列最小優化SMO(Sequential Minimal Optimization)
6. 高斯核函數(Gaussian Kernel)
7. 實戰案例
7.1. 簡化版SMO算法案例
7.2. 完整版SMO算法案例
7.3. 高斯核函數案例
學習完機器學習實戰的支援向量機,簡單的做個筆記。文中部分描述屬于個人消化後的了解,僅供參考。
本篇綜合了先前的文章,如有不了解,可參考:
吳恩達機器學習(十)支援向量機
所有代碼和資料可以通路 我的 github
如果這篇文章對你有一點小小的幫助,請給個關注喔~我會非常開心的~
0. 前言
支援向量機(Support Vector Machine)是一類監督學習的二分類算法。
通過建構超平面,将資料分割開來,超平面一邊的資料屬于同一類别。
支援向量是指距離分割超平面最近的那些點,支援向量機目的是使得支援向量到超平面的間隔最大。
- 優點:泛化錯誤率低,結果易于解釋
- 對參數調節和核函數的選擇敏感
- 适用資料類型:數值型和标稱型資料
如下圖所示(圖源:百度百科):
為分割超平面,在超平面上下方各建立一個界面
和
,在
之上的為一類,在
之下的為一類。
SVM的目的是為了最大化間隔,即
到
的距離。
1. 尋找最大間隔
在SVM中,定義正類
,定義反類
,這是為了友善後面計算。
通過
拟合資料的分界線,即超平面,則可定義上界面為
,下界面為
。因可以通過改變系數的倍數表示距離,是以可粗略定義為
。
上界面以上的點表示為
,下界面以下的點表示為
,是以綜合可得
。
需最大化上界面到下界面的距離:
是以,可得目标函數和限制條件:
2. 拉格朗日乘子法和KKT條件
- 拉格朗日乘子法:求在限制條件下,函數的極值問題
- KKT條件:若限制條件為不等式,則需滿足一定條件
在拉格朗日乘子法中,引入了其他的參數,并将原問題轉換為對偶問題(求最小值轉換為求最大值)。
首先修改目标函數和限制條件為:
應用拉格朗日乘子法:
應滿足的KKT條件為:
将KKT條件代回拉格朗日乘子法:
因轉換為對偶問題,是以此時目标函數和限制條件:
3. 松弛變量
有時,資料并非
線性可分,有一部分資料會在上下界面之間,或者處于錯誤的類别。
此時可增加一個松弛變量
,允許有資料點位于錯誤的一側。
定義上界面以上為
,下界面以下為
,可綜合為
是以,可得目标函數和限制條件,其中
是控制權重:
4. 帶松弛變量的拉格朗日乘子法和KKT條件
同理,根據拉格朗日乘子法和KKT條件,可得:
是以,目标函數和限制條件為:
5. 序列最小優化SMO(Sequential Minimal Optimization)
SMO算法的目标是求出一系列的
和
,然後通過一開始的KKT條件,得出權重
:
SMO算法的工作原理是,每次循環都選擇兩個
,同時對兩個
進行優化。
因需滿足
,是以增大其中一個
的同時,減小另一個
。
6. 高斯核函數(Gaussian Kernel)
若資料是非線性可分,就無法使用線性的SVM進行解決。
可通過核函數,将低維空間中的資料映射到高維空間中,可将一個在低維空間中的非線性問題轉換為高維空間下的線性問題求解。
高斯核函數的定義如下,其中
表示到達率(函數值跌落到
的速度參數),
越小,支援向量越多:
其中,
稱為标記點,
,每一個标記點與每一個樣本資料在空間中位于相同位置。是以有:
- 如果 與 相近
- 如果 與 相隔遠
超平面的定義如下:
可使用核函數,替換向量的内積:
高斯核函數的SVM流程可表示為:
- 給定資料集
- SMO算法計算 和
- 高斯核函數計算
- 代入 ,判斷大于 則屬于正類,小于 則屬于反類
7. 實戰案例
以下将展示書中案例的代碼段,所有代碼和資料可以在github中下載下傳:
7.1. 簡化版SMO算法案例
# coding:utf-8
from numpy import *
"""
簡化版SMO算法案例
"""
# 加載資料集
def loadDataSet(fileName):
dataMat = []
labelMat = []
fr = open(fileName)
for line in fr.readlines():
lineArr = line.strip().split('\t')
dataMat.append([float(lineArr[0]), float(lineArr[1])])
labelMat.append(float(lineArr[2]))
return dataMat, labelMat
# 根據第i個alpha,從m中選擇另一個不同的alpha
def selectJrand(i, m):
j = i
while (j == i):
j = int(random.uniform(0, m))
return j
# 設定alpha的上界和下界
def clipAlpha(aj, H, L):
if aj > H:
aj = H
if L > aj:
aj = L
return aj
# 簡化版本SMO算法,并沒有确定優化的最佳alpha對,而是随機選擇
# dataMatIn 資料集
# classLabels 類别标簽
# C 常數C
# toler 容錯率
# maxIter 最大疊代次數
def smoSimple(dataMatIn, classLabels, C, toler, maxIter):
dataMatrix = mat(dataMatIn)
labelMat = mat(classLabels).transpose()
b = 0
m, n = shape(dataMatrix)
# alphas的次元為m*1,m為資料集大小
alphas = mat(zeros((m, 1)))
iter = 0
while (iter < maxIter):
alphaPairsChanged = 0
for i in range(m):
# 預測第i個向量屬于的類别
fXi = float(multiply(alphas, labelMat).T * (dataMatrix * dataMatrix[i, :].T)) + b
# 判斷這個預測的誤差
Ei = fXi - float(labelMat[i])
# 滿足KKT條件
if ((labelMat[i] * Ei < -toler) and (alphas[i] < C)) or ((labelMat[i] * Ei > toler) and (alphas[i] > 0)):
# 選擇另一個向量j
j = selectJrand(i, m)
# 預測第j個向量屬于的類别
fXj = float(multiply(alphas, labelMat).T * (dataMatrix * dataMatrix[j, :].T)) + b
# 判斷這個預測的誤差
Ej = fXj - float(labelMat[j])
# 建立alpha i 和alpha j的副本
alphaIold = alphas[i].copy()
alphaJold = alphas[j].copy()
# 判斷兩個類别标簽是否相等
# 設定最大值為0,最小值為C
if labelMat[i] != labelMat[j]:
L = max(0, alphas[j] - alphas[i])
H = min(C, C + alphas[j] - alphas[i])
else:
L = max(0, alphas[j] + alphas[i] - C)
H = min(C, alphas[j] + alphas[i])
if L == H:
print("L==H")
continue
# 計算最優修改量
eta = 2.0 * dataMatrix[i, :] * dataMatrix[j, :].T - \
dataMatrix[i, :] * dataMatrix[i, :].T - \
dataMatrix[j, :] * dataMatrix[j, :].T
if eta >= 0:
print("eta>=0")
continue
# 對alpha j作調整
alphas[j] -= labelMat[j] * (Ei - Ej) / eta
alphas[j] = clipAlpha(alphas[j], H, L)
# 如果alpha j隻是輕微的調整
if abs(alphas[j] - alphaJold) < 0.00001:
print("j not moving enough")
continue
# alpha a改變大小與alpha j一樣,但是方向相反
alphas[i] += labelMat[j] * labelMat[i] * (alphaJold - alphas[j])
# 調整常數項b
b1 = b - Ei - labelMat[i] * (alphas[i] - alphaIold) * dataMatrix[i, :] * dataMatrix[i, :].T \
- labelMat[j] * (alphas[j] - alphaJold) * dataMatrix[i, :] * dataMatrix[j, :].T
b2 = b - Ej - labelMat[i] * (alphas[i] - alphaIold) * dataMatrix[i, :] * dataMatrix[j, :].T \
- labelMat[j] * (alphas[j] - alphaJold) * dataMatrix[j, :] * dataMatrix[j, :].T
if (0 < alphas[i]) and (C > alphas[i]):
b = b1
elif (0 < alphas[j]) and (C > alphas[j]):
b = b2
else:
b = (b1 + b2) / 2.0
alphaPairsChanged += 1
print("iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
# 如果沒有對任何向量進行修改,則疊代次數+1,否則疊代次數清零
# 隻有在所有疊代中,都沒有對alpha進行修改,才能說明所有alpha已經最優
if alphaPairsChanged == 0:
iter += 1
else:
iter = 0
print("iteration number: %d" % iter)
return b, alphas
if __name__ == '__main__':
dataArr, labelArr = loadDataSet('testSet.txt')
b, alphas = smoSimple(dataArr, labelArr, 0.6, 0.001, 40)
print(b)
print(alphas)
7.2. 完整版SMO算法案例
# coding:utf-8
from numpy import *
"""
完整版SMO算法案例
"""
# 加載資料集
def loadDataSet(fileName):
dataMat = []
labelMat = []
fr = open(fileName)
for line in fr.readlines():
lineArr = line.strip().split('\t')
dataMat.append([float(lineArr[0]), float(lineArr[1])])
labelMat.append(float(lineArr[2]))
return dataMat, labelMat
# 建立存儲各參數的類
class optStruct:
def __init__(self, dataMatIn, classLabels, C, toler):
self.X = dataMatIn
self.labelMat = classLabels
self.C = C
self.tol = toler
self.m = shape(dataMatIn)[0]
self.alphas = mat(zeros((self.m, 1)))
self.b = 0
self.eCache = mat(zeros((self.m, 2)))
# 根據第i個alpha,從m中選擇另一個不同的alpha
def selectJrand(i, m):
j = i
while (j == i):
j = int(random.uniform(0, m))
return j
# 設定alpha的上界和下界
def clipAlpha(aj, H, L):
if aj > H:
aj = H
if L > aj:
aj = L
return aj
# 計算預測誤差函數
def calcEk(oS, k):
fXk = float(multiply(oS.alphas, oS.labelMat).T * (oS.X * oS.X[k, :].T)) + oS.b
Ek = fXk - float(oS.labelMat[k])
return Ek
# 内循環尋找第二個alpha
# 啟發式方法,最大化步長尋找第二個alpha
def selectJ(i, oS, Ei):
maxK = -1
maxDeltaE = 0
Ej = 0
# 誤差全局緩存指派
# 第0個元素表示是否指派過
# 第1個元素表示對應誤差
oS.eCache[i] = [1, Ei]
# 擷取所有指派過的誤差全局緩存
validEcacheList = nonzero(oS.eCache[:, 0].A)[0]
# 因不能選擇目前i,是以需判斷數量是否大于1
# 若是則最大化步長, 若否則随機選擇alpha j
if (len(validEcacheList)) > 1:
for k in validEcacheList:
if k == i:
continue
Ek = calcEk(oS, k)
deltaE = abs(Ei - Ek)
# 選擇最大步長的
if (deltaE > maxDeltaE):
maxK = k
maxDeltaE = deltaE
Ej = Ek
return maxK, Ej
else:
j = selectJrand(i, oS.m)
Ej = calcEk(oS, j)
return j, Ej
# 更新誤差全局緩存
def updateEk(oS, k):
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]
# 完整版SMO,内循環優化
def innerL(i, oS):
# alpha_i的誤差
Ei = calcEk(oS, i)
if ((oS.labelMat[i] * Ei < -oS.tol) and (oS.alphas[i] < oS.C)) \
or ((oS.labelMat[i] * Ei > oS.tol) and (oS.alphas[i] > 0)):
# 選擇alpha_j
j, Ej = selectJ(i, oS, Ei)
# 記錄原始的alpha
alphaIold = oS.alphas[i].copy()
alphaJold = oS.alphas[j].copy()
if oS.labelMat[i] != oS.labelMat[j]:
L = max(0, oS.alphas[j] - oS.alphas[i])
H = min(oS.C, oS.C + oS.alphas[j] - oS.alphas[i])
else:
L = max(0, oS.alphas[j] + oS.alphas[i] - oS.C)
H = min(oS.C, oS.alphas[j] + oS.alphas[i])
if L == H:
print("L==H")
return 0
# 計算最優修改量
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
# 更新alpha_j
oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej) / eta
oS.alphas[j] = clipAlpha(oS.alphas[j], H, L)
# 更新全局緩存
updateEk(oS, j)
if abs(oS.alphas[j] - alphaJold) < 0.00001:
print("j not moving enough")
return 0
# 更新alpha_i
oS.alphas[i] += oS.labelMat[j] * oS.labelMat[i] * (alphaJold - oS.alphas[j])
# 更新全局緩存
updateEk(oS, i)
# 更新常數b
b1 = oS.b - Ei - oS.labelMat[i] * (oS.alphas[i] - alphaIold) * oS.X[i, :] * oS.X[i, :].T - oS.labelMat[j] * (
oS.alphas[j] - alphaJold) * oS.X[i, :] * oS.X[j, :].T
b2 = oS.b - Ej - oS.labelMat[i] * (oS.alphas[i] - alphaIold) * oS.X[i, :] * oS.X[j, :].T - oS.labelMat[j] * (
oS.alphas[j] - alphaJold) * oS.X[j, :] * oS.X[j, :].T
if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]):
oS.b = b1
elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]):
oS.b = b2
else:
oS.b = (b1 + b2) / 2.0
return 1
else:
return 0
# 完整版SMO,外循環
def smoP(dataMatIn, classLabels, C, toler, maxIter, kTup=('lin', 0)):
oS = optStruct(mat(dataMatIn), mat(classLabels).transpose(), C, toler)
iter = 0
entireSet = True
alphaPairsChanged = 0
while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):
alphaPairsChanged = 0
# 完整周遊
if entireSet:
for i in range(oS.m):
alphaPairsChanged += innerL(i, oS)
print("fullSet, iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
iter += 1
# 非邊界值周遊
else:
nonBoundIs = nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
for i in nonBoundIs:
alphaPairsChanged += innerL(i, oS)
print("non-bound, iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
iter += 1
# 切換周遊的地方
if entireSet:
entireSet = False
elif alphaPairsChanged == 0:
entireSet = True
print("iteration number: %d" % iter)
return oS.b, oS.alphas
# 根據alpha計算權重系數w
def calcWs(alphas, dataArr, classLabels):
X = mat(dataArr)
labelMat = mat(classLabels).transpose()
m, n = shape(X)
w = zeros((n, 1))
for i in range(m):
w += multiply(alphas[i] * labelMat[i], X[i, :].T)
return w
if __name__ == '__main__':
dataArr, labelArr = loadDataSet('testSet.txt')
b, alphas = smoP(dataArr, labelArr, 0.6, 0.001, 40)
w = calcWs(alphas, dataArr, labelArr)
print(w)
7.3. 高斯核函數案例
# coding:utf-8
from numpy import *
"""
高斯核函數案例
"""
# 加載資料集
def loadDataSet(fileName):
dataMat = []
labelMat = []
fr = open(fileName)
for line in fr.readlines():
lineArr = line.strip().split('\t')
dataMat.append([float(lineArr[0]), float(lineArr[1])])
labelMat.append(float(lineArr[2]))
return dataMat, labelMat
# 核函數轉換
def kernelTrans(X, A, kTup):
m, n = shape(X)
K = mat(zeros((m, 1)))
# 線性核心
if kTup[0] == 'lin':
K = X * A.T
# 高斯核心
elif kTup[0] == 'rbf':
for j in range(m):
deltaRow = X[j, :] - A
K[j] = deltaRow * deltaRow.T
K = exp(K / (-1 * kTup[1] ** 2))
else:
raise NameError('Houston We Have a Problem -- That Kernel is not recognized')
return K
# 存儲參數資料結構
class optStruct:
def __init__(self, dataMatIn, classLabels, C, toler, kTup):
self.X = dataMatIn
self.labelMat = classLabels
self.C = C
self.tol = toler
self.m = shape(dataMatIn)[0]
self.alphas = mat(zeros((self.m, 1)))
self.b = 0
self.eCache = mat(zeros((self.m, 2)))
self.K = mat(zeros((self.m, self.m)))
for i in range(self.m):
self.K[:, i] = kernelTrans(self.X, self.X[i, :], kTup)
# 計算預測的誤差
def calcEk(oS, k):
fXk = float(multiply(oS.alphas, oS.labelMat).T * oS.K[:, k] + oS.b)
Ek = fXk - float(oS.labelMat[k])
return Ek
# 根據第i個alpha,從m中選擇另一個不同的alpha
def selectJrand(i, m):
j = i
while (j == i):
j = int(random.uniform(0, m))
return j
# 設定alpha的上界和下界
def clipAlpha(aj, H, L):
if aj > H:
aj = H
if L > aj:
aj = L
return aj
# 内循環尋找第二個alpha
# 啟發式方法,最大化步長尋找第二個alpha
def selectJ(i, oS, Ei):
maxK = -1
maxDeltaE = 0
Ej = 0
# 誤差全局緩存指派
# 第0個元素表示是否指派過
# 第1個元素表示對應誤差
oS.eCache[i] = [1, Ei]
# 擷取所有指派過的誤差全局緩存
validEcacheList = nonzero(oS.eCache[:, 0].A)[0]
# 因不能選擇目前i,是以需判斷數量是否大于1
# 若是則最大化步長, 若否則随機選擇alpha j
if (len(validEcacheList)) > 1:
for k in validEcacheList:
if k == i:
continue
Ek = calcEk(oS, k)
deltaE = abs(Ei - Ek)
# 選擇最大步長的
if (deltaE > maxDeltaE):
maxK = k
maxDeltaE = deltaE
Ej = Ek
return maxK, Ej
else:
j = selectJrand(i, oS.m)
Ej = calcEk(oS, j)
return j, Ej
# 更新誤差全局緩存
def updateEk(oS, k):
Ek = calcEk(oS, k)
oS.eCache[k] = [1, Ek]
# 完整版SMO,内循環優化
def innerL(i, oS):
# alpha_i的誤差
Ei = calcEk(oS, i)
if ((oS.labelMat[i] * Ei < -oS.tol) and (oS.alphas[i] < oS.C)) \
or ((oS.labelMat[i] * Ei > oS.tol) and (oS.alphas[i] > 0)):
# 選擇alpha_j
j, Ej = selectJ(i, oS, Ei)
# 記錄原始的alpha
alphaIold = oS.alphas[i].copy()
alphaJold = oS.alphas[j].copy();
if oS.labelMat[i] != oS.labelMat[j]:
L = max(0, oS.alphas[j] - oS.alphas[i])
H = min(oS.C, oS.C + oS.alphas[j] - oS.alphas[i])
else:
L = max(0, oS.alphas[j] + oS.alphas[i] - oS.C)
H = min(oS.C, oS.alphas[j] + oS.alphas[i])
if L == H:
print("L==H")
return 0
# 計算最優修改量
eta = 2.0 * oS.K[i, j] - oS.K[i, i] - oS.K[j, j]
if eta >= 0:
print("eta>=0")
return 0
# 更新alpha_j
oS.alphas[j] -= oS.labelMat[j] * (Ei - Ej) / eta
oS.alphas[j] = clipAlpha(oS.alphas[j], H, L)
# 更新全局緩存
updateEk(oS, j)
if abs(oS.alphas[j] - alphaJold) < 0.00001:
print("j not moving enough")
return 0
# 更新alpha_i
oS.alphas[i] += oS.labelMat[j] * oS.labelMat[i] * (alphaJold - oS.alphas[j])
# 更新全局緩存
updateEk(oS, i)
# 更新常數b
b1 = oS.b - Ei - oS.labelMat[i] * (oS.alphas[i] - alphaIold) * oS.K[i, i] - oS.labelMat[j] * (
oS.alphas[j] - alphaJold) * oS.K[i, j]
b2 = oS.b - Ej - oS.labelMat[i] * (oS.alphas[i] - alphaIold) * oS.K[i, j] - oS.labelMat[j] * (
oS.alphas[j] - alphaJold) * oS.K[j, j]
if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]):
oS.b = b1
elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]):
oS.b = b2
else:
oS.b = (b1 + b2) / 2.0
return 1
else:
return 0
# 完整版SMO,外循環
def smoP(dataMatIn, classLabels, C, toler, maxIter, kTup=('lin', 0)):
oS = optStruct(mat(dataMatIn), mat(classLabels).transpose(), C, toler, kTup)
iter = 0
entireSet = True
alphaPairsChanged = 0
while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):
alphaPairsChanged = 0
# 完整周遊
if entireSet:
for i in range(oS.m):
alphaPairsChanged += innerL(i, oS)
print("fullSet, iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
iter += 1
# 非邊界值周遊
else:
nonBoundIs = nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
for i in nonBoundIs:
alphaPairsChanged += innerL(i, oS)
print("non-bound, iter: %d i:%d, pairs changed %d" % (iter, i, alphaPairsChanged))
iter += 1
# 切換周遊的地方
if entireSet:
entireSet = False
elif alphaPairsChanged == 0:
entireSet = True
print("iteration number: %d" % iter)
return oS.b, oS.alphas
# 文本轉換為向量
def img2vector(filename):
returnVect = zeros((1, 1024))
fr = open(filename)
for i in range(32):
lineStr = fr.readline()
for j in range(32):
returnVect[0, 32 * i + j] = int(lineStr[j])
return returnVect
# 加載資料
def loadImages(dirName):
from os import listdir
hwLabels = []
trainingFileList = listdir(dirName)
m = len(trainingFileList)
trainingMat = zeros((m, 1024))
for i in range(m):
fileNameStr = trainingFileList[i]
fileStr = fileNameStr.split('.')[0]
classNumStr = int(fileStr.split('_')[0])
if classNumStr == 9:
hwLabels.append(-1)
else:
hwLabels.append(1)
trainingMat[i, :] = img2vector('%s/%s' % (dirName, fileNameStr))
return trainingMat, hwLabels
# 測試函數
def testDigits(kTup=('rbf', 10)):
dataArr, labelArr = loadImages('trainingDigits')
b, alphas = smoP(dataArr, labelArr, 200, 0.0001, 10000, kTup)
datMat = mat(dataArr)
labelMat = mat(labelArr).transpose()
# 擷取所有支援向量
svInd = nonzero(alphas.A > 0)[0]
sVs = datMat[svInd]
labelSV = labelMat[svInd]
print("there are %d Support Vectors" % shape(sVs)[0])
m, n = shape(datMat)
errorCount = 0
for i in range(m):
# 建構核函數
kernelEval = kernelTrans(sVs, datMat[i, :], kTup)
# 利用核函數進行預測
predict = kernelEval.T * multiply(labelSV, alphas[svInd]) + b
if sign(predict) != sign(labelArr[i]): errorCount += 1
print("the training error rate is: %f" % (float(errorCount) / m))
dataArr, labelArr = loadImages('testDigits')
errorCount = 0
datMat = mat(dataArr)
labelMat = mat(labelArr).transpose()
m, n = shape(datMat)
for i in range(m):
kernelEval = kernelTrans(sVs, datMat[i, :], kTup)
predict = kernelEval.T * multiply(labelSV, alphas[svInd]) + b
if sign(predict) != sign(labelArr[i]): errorCount += 1
print("the test error rate is: %f" % (float(errorCount) / m))
if __name__ == '__main__':
testDigits()
如果這篇文章對你有一點小小的幫助,請給個關注喔~我會非常開心的~