感知機原理及代碼實現
上篇講完梯度下降,這篇博客我們就來好好整理一下一個非常重要的二分類算法——感知機,這是一種二分類模型,當輸入一系列的數據后,輸出的是一個二分類變量,如0或1
1. 算法原理
1.1 知識引入
說起分類算法,博主想到的另一個算法是邏輯回歸,而感知機從原理上來說和回歸類算法最大的區別就是引入了幾何的思想,將向量放到了高維空間上去想象。首先介紹一個概念——超平面(hyper plane),超平面就是比當前維度還少一個維度的數的集合。打個比方,如果是二維平面,那么超平面就是一維平面,也就是一條直線,形象地來想象,就好比是一條線切割了一張紙,並將一張紙上標記為1和標記為-1的點給分開來了;類比到三維空間,其超平面就是一個二維平面去切割一個立方體,好,下面開始進行原理的數學推導:
高中學數學的時候我們就知道,二維平面上的一條直線內的點可以用一般方程來表示,那么該二維平面中,滿足
的(x, y)表示在此一般方程上方的點,而滿足
的點則表示在此一般方程下方的點,將此情形推廣到高維空間如n維空間(過高維度我們將無法直觀想象),就可得到一超平面方程
,那么推導后可知,滿足等式左邊>0的點自然就在該超平面的上方,<0的點就在該超平面的下方,自然,用這種方式將所有點分成了兩類,因此,原問題就轉換成了找到這么一個超平面將被點分成兩類的問題了,將該方程寫成兩個向量相乘再加上一個常數項,或稱截距項的形式,就可得到下式:
令,
,則原式可轉換為:
,現令
,並令
,則F(x)函數可以用來做線性分類器。
1.2 推導過程
上面的式子得到之后,我們就要想,該如何求得w和b的值呢?模仿博主上兩篇博客里的套路,那自然是分三步走咯,先隨機設置w和b的初始值,然后寫出一個損失函數(又稱成本函數),最后使用梯度下降法通過不斷迭代求得最優解使得損失函數最小。依照這個想法,第一個容易想到的思路就是既然我們的目標是要找到一個超平面方程正確無誤地將兩類點徹底分開,那么最優的解自然就是使得錯誤分類的點降到0個就行了唄,但是該方法有一個明顯的缺點,函數值若是一個兩個的話,就是分類變量了,其函數圖像並不連續,因此第三步梯度下降我們沒法用!於是第二個思路應運而生,既然我們要讓函數連續,我們將求點的個數轉換成誤分類點到超平面方程的距離之和不就行了嗎?那下面,我們首先推導點到直線的距離公式:
在二維直角坐標平面XoY內,有一直線方程,現需要求的是向量
的模
即d的長度:
首先推導d的公式:
設R點的坐標為,既然R點在直線上,那么必滿足:
,可推出
通過原方程可輕松計算出方程的斜率為(-b, a),那么與該方向垂直的向量通過點積為0的原理可算出是(a, b),如此一來,公式中的兩個關鍵向量都已經計算出來,現在開始代入公式
,然后帶入等式
並整理,可得
,由於距離比為正數,加上絕對值符號,最終等式為:
至此,二維平面上的推導過程告一段落,現在將其擴展至高維空間,也就是計算誤分類點到超平面的距離,類比上式,可得:
設誤分類點的集合為M,則損失函數L的公式為:
下一步為對上面的損失函數求導,分兩步,求函數對w的偏導數,以及函數對b的偏導數,求出后運用梯度下降迭代即可,然而,對含有絕對值的函數求導絕非易事,因此需要通過某種方式想方設法去除絕對值符號,
注意到,對於正確分類點來說,正類標記為1,負類標記為-1,而正類位於超平面上方,負類位於超平面下方,因此對於正確分類點,必有不等式,那既然所有的點只能區分為正確分類點以及誤分類點,那么對於誤分類點就必有不等式
,就是一個負負得正的道理,這樣一來,絕對值符號可以去除,原式就變為了:
,我們的目標就是求這個函數的最小值即可。至此,已經可以給出兩個偏導數的計算公式了:
下面開始簡單寫一下梯度下降的偽代碼:
1. 選定一個系數向量和截距的初值
2.在訓練集中進行遍歷,逐個選取
3.遍歷過程中,如果出現的點,那就說明該點是誤分類點,這時就需要更新初值了,設步長為λ,則更新函數為:
4.重復2、3兩步,直到整個數據集沒有誤分類點,即可說明線性分類器已經找到
2.代碼編寫
下面給出上述算法的python代碼實現,並將代碼封裝成了一個類,這個類中定義了兩個不同的梯度下降代碼實現,分別為批量梯度下降以及隨機梯度下降:
import numpy as np import pandas as pd import random from sklearn.datasets import make_blobs class perceptron: def __init__(self, sample_size = 200, centers = [[1,1],[3,4]], cluster_std = 0.6): X,Y = make_blobs(n_samples = sample_size, centers = centers, cluster_std = cluster_std, random_state = 11) Y[Y == 0] = -1 self.X = X self.Y = Y def BGDfit(self, w, b, lam): #This function aims to use Batch Gradient Descent algorithm to find the most suitable function while True: x_mat = np.mat(self.X) tmp = np.ravel(x_mat * w + b) * self.Y if sum(tmp < 0) == 0: break tmp_X = np.array(self.X[tmp < 0]) tmp_Y = np.array(self.Y[tmp < 0]) tmp_Y_final = np.column_stack((tmp_Y, tmp_Y)) w = w + lam * np.mat((tmp_X * tmp_Y_final).sum(axis = 0)).T b = b + lam * sum(tmp_Y) return w, b def SGDfit(self, w, b, lam): #This function aims to use Stochastic Gradient Descent algorithm to find the most suitable function missing_index = 0 while missing_index != -1: missing_index = -1 for i in range(self.X.shape[0]): if self.Y[i] * (np.dot(w, self.X[i]) + b) < 0: missing_index = i break if missing_index != -1: w = w + lam * self.Y[missing_index] * self.X[missing_index] b = b + lam * self.Y[missing_index] return w, b
3.感知機的對偶形式
3.1 原理
使用原始的感知機形式時,在遇到處理海量數據的場景時,運算速度會十分緩慢,那么,存不存在一種方式,我們可以預先將某些矩陣計算好,將其緩存起來而不必每次都進行重復的運算呢?答案是肯定的,這就是感知機的對偶形式,下面開始正式介紹:在1.2章節的時候給出過w和b的更新函數,每遇到一個誤分類點,就將w進行更新:,將b進行更新:
,我們可以從中總結出這樣一個規律,那就是,每遇到一個誤分類點一次,就在w基礎上多累加一次
,在b基礎上多累加一次
,在式子中,λ是確定的,特征值數據集是確定的,特征值對應的分類也是確定的,唯一不確定的就是我們不知道何時會遇到誤分類的點,因此,不確定的變量是誤分類次數,那么,感知機的對偶形式呼之欲出,我們只要維護一個誤分類次數的數組,每當遇到誤分類點的時候,在該數組的對應位置上進行累加即可!
設數據集中一共有n個點,並設誤分類點次數向量為,根據上面的思路,w和b就可以寫成,
,
,我們將w的初始值設置成零向量,那么w就變成了:
,將這個式子代入至判斷某個點是否是誤分類點的判別式中,如果是誤分類點,那就更新a和b,具體的偽代碼如下:
1.選定a和b的初始值:,b = 0
2.在訓練集中進行遍歷,逐個選取
3.如果發現某個點,那就更新a和b,
,
4.重復2、3步直到沒有誤分類點,算法停止
3.2 代碼實現
下面給出感知機對偶形式的python代碼實現,使用的是隨機梯度下降法:
import numpy as np import pandas as pd import matplotlib.pyplot as plt import random from sklearn.datasets import make_blobs sample_size = 200 centers = [[1,1], [3,4]] X, Y = make_blobs(n_samples = sample_size, centers = centers, cluster_std = 0.6, random_state = 11) Y[Y == 0] = -1 gram = np.array(np.mat(X) * np.mat(X).T) a = np.zeros(X.shape[0]) b = 0 lam = 0.1 count = 0 while True: count += 1 missing_index = -1 for i in range(X.shape[0]): checking = Y[i] * ((a * Y * gram[i, :]).sum() + b) if checking <=0 : missing_index = i break if missing_index == -1: break a[missing_index] += lam b += lam * Y[missing_index] theta = np.ravel(np.mat(a * Y) * np.mat(X)) theta, b