寫程式學ML：Logistic迴歸演算法原理及實現（二）
2、Logistic迴歸演算法的實現
2.1 Logistic演算法的實現
首先，我們實現梯度上升演算法。
Sigmoid函式的定義如下：
#sigmoid函式的實現
def sigmoid(inX):
return 1.0 / (1 + exp(-inX))
定義函式gradAscent(dataMatIn, classLabels)實現梯度上升演算法。它有兩個形參，形參dataMatIn是一個二維陣列，包含訓練樣本，每個樣本有三個特徵值。形參classLabels儲存每個樣本的分類情況。
此函式中呼叫函式mat()將兩個形參轉換成了NumPy矩陣，此處使用到了矩陣操作。為了獲得最佳迴歸係數，使用了500次迭代。每次迭代時，將每個樣本3個特徵與迴歸係數的乘積累加和作為sigmoid函式的輸入引數，求得每個樣本的分類結果。然後與每個樣本真實分類情況作差。最後利用梯度上升演算法的迭代公式更新每個樣本的各個特徵值，更新步長為0.001。如果之前分類情況作差結果為1，則特徵值增加原來的0.001倍；如果是-1，則特徵值減少原來的0.001倍；如果是0，否則不做調整。
該函式的具體實現如下：
#dataMatIn：訓練樣本矩陣，每個樣本包含3個特徵值 #classLabels：訓練樣本對應的分類矩陣 #該函式利用梯度上升演算法，產生經過多次調整後的迴歸係數 def gradAscent(dataMatIn, classLabels): dataMatrix = mat(dataMatIn) #呼叫mat()函式將陣列轉化為矩陣 labelMat = mat(classLabels).transpose() #呼叫transpose()函式對矩陣進行轉置 m, n = shape(dataMatrix) #獲取矩陣或者陣列的維數 alpha = 0.001 #向目標移動的步長 maxCycles = 500 #迭代次數 weights = ones((n, 1)) #產生n*1的矩陣，初始化為全1；該變數儲存迴歸係數 #迴圈500次，每次按照當前的迴歸係數，計算每個樣本的預測類別；計算出其與真實類別的差值，進而調整迴歸係數，然後再迴圈處理。 for k in range(maxCycles): #dataMatrix * weights為矩陣乘法，dataMatrix為mxn的矩陣，weights為nx1的矩陣，它們乘積為mx1的矩陣 h = sigmoid(dataMatrix * weights) #針對dataMatrix * weights的乘積結果，獲取各個元素對應的sigmoid值，h為mx1的矩陣 error = (labelMat - h) #獲取各個樣本預測類別與真實類別資訊資料的差值 #根據差值error的方向調整迴歸係數，如果error為負值 ，則減小系數；反之，增大系數 weights = weights + alpha * dataMatrix.transpose() * error return weights
從該函式的實現可以看到，每次迭代時，對於有N個特徵的M個樣本，需要進行M*N次乘法，這個計算量是很大的，特別是樣本和特徵比較多的時候，計算量會飛速增長，函式效能會非常差。
為了提升效能，可以使用隨機梯度上升演算法，即每次迭代時僅用一個樣本點來更新迴歸係數。由於可以在新樣本到來時對分類器進行增量式更新，因而隨機梯度上升演算法是一個線上學習演算法。與“線上學習”相對應，一次迭代中處理所有資料被稱為“批處理”。
隨機梯度上升演算法的具體實現如下：
#dataMatIn：訓練樣本矩陣，每個樣本包含3個特徵值 #classLabels：訓練樣本對應的分類矩陣 #該函式利用隨機梯度上升演算法，產生按照樣本數進行多次調整後的迴歸係數 def stocGradAscent0(dataMatrix, claaLabels): m, n = shape(dataMatrix) alpha = 0.01 #向目標移動的步長 weights = ones(n) #建立有n個元素的迴歸係數陣列，並初始化為1 for i in range(m): #按照樣本數目進行迴圈迭代 h = sigmoid(sum(dataMatrix[i] * weights)) #求得每個樣本的特徵值與迴歸係數乘積之和，然後呼叫函式sigmoid獲得分類結果 error = classLabels[i] - h #計算計算的分類結果與真實結果差值 weights = weights + alpha * error * dataMatrix[i] #對迴歸係數進行調整 return weights
從程式實現可以看出，隨機梯度上升演算法與梯度上升演算法在程式碼上很相似，但也有區別：第一，後者的變數h和誤差error都是向量，而前者則全是數值；第二，前者沒有矩陣的轉換過程，所有變數的資料型別都是NumPy陣列。
隨機梯度上升演算法中迭代次數與樣本個數相等，本章中使用樣本個數為100，即迭代100次，而梯度上升演算法中迭代了500次。雖然前者程式執行速度快了很多，但演算法結果卻差強人意。為了提高隨機梯度上升演算法的質量，需要對該演算法做一些修改。
圖5-6展示了隨機梯度上升演算法在200次迭代過程中迴歸係數的變化情況。其中係數2只經過50次迭代就達到了穩定值，但係數1和0則需要更多次的迭代。另外值得注意的是，在大的波動停止後，還有一些小的週期性波動。產生這種現象的原因是存在一些不能正確分類的樣本點（資料集並非線性可分），在每次迭代時會引發係數的劇烈改變。我們期望演算法能避免來回波動，從而收斂到某個值。另外，收斂速度也需要加快。
對於隨機梯度上升演算法中存在的問題，我們使用改進的演算法來解決。改進的梯度上升演算法具體實現如下：
#dataMatIn：訓練樣本矩陣，每個樣本包含3個特徵值
#classLabels：訓練樣本對應的分類矩陣
#numIter：迭代次數，預設值為150
#該函式利用隨機梯度上升演算法，多次迭代，產生多次調整後的迴歸係數
def stocGradAscent1(dataMatrix, classLabels, numIter = 150):
m, n = shape(dataMatrix)
weights = ones(n) #建立有n個元素的迴歸係數陣列，並初始化為1
for j in range(numIter): #迭代處理，預設迭代150次
#range函式用來建立序列，rang(5)結果為[0,1,2,3,4]
dataIndex = range(m)
for i in range(m): #每次迭代，需要更新m次係數，m等於樣本次數
alpha = 4 / (1.0 + j + i) + 0.01 #每次迭代時係數調整的步長
#numpy.random.uniform(low,high,size)：
#從一個均勻分佈[low,high)中隨機取樣，注意定義域是左閉右開，即包含low，不包含high.
#從樣本中隨機選擇一個樣本對應的index，作為此次迭代的樣本id
randIndex = int(random.uniform(0, len(dataIndex)))
#求得第randIndex個樣本的特徵值與迴歸係數乘積之和，然後呼叫函式sigmoid獲得分類結果
h = sigmoid(sum(dataMatrix[randIndex] * weights))
error = classLabels[randIndex] - h #計算計算的分類結果與真實結果差值
weights = weights + alpha * error * dataMatrix[randIndex] #對迴歸係數進行調整
del(dataIndex[randIndex]) #從樣本ID集中刪除此次迭代所用樣本ID
return weights程式中alpha在每次迭代的時候都會調整，這會緩解圖5-6上的資料波動或者高頻波動。雖然alpha會隨著迭代次數不斷減小，但永遠不會減小到0，這是因為存在一個常數項。必須這樣做的原因是為了保證在多次迭代之後新資料仍然具有一定的影響。如果要處理的問題是動態變化的，那麼可以適當加大上述常數項，來確保新的值獲得更大的迴歸係數。
另一點值得注意的是，在降低alpha的函式中，alpha每次減少1/(j+1)，其中j是迭代次數，i是樣本點的下標。這樣當j<<max(i)時，alpha就不是嚴格下降的。避免參數的嚴格下降也常見於模擬退火演算法等其他優化演算法中。
改進後的演算法中通過隨機選取樣本來更新迴歸係數。這種方法將減少週期性的波動。這種方法每次隨機從列表中選出一個值，然後從列表中刪掉該值（再進行下一次迭代）。
比較圖5-7和圖5-6可以看到兩點不同。第一點是，圖5-7中的係數沒有像圖5-6裡那樣出現週期性的波動，這歸功於stocGrdAscent1()裡的樣本隨機選擇機制；第二點是，圖5-7的水平軸比圖5-6短了很多，這是由於stocGrdAscent1()可以收斂得更快。
對於訓練出來的迴歸係數，如果要對測試樣本進行測試，可以通過下面的函式來實現：#inX：測試樣本
#weights：訓練好的係數
#獲得測試樣本的分類情況
def classifyVector(inX, weights):
#求得測試樣本的特徵值與迴歸係數乘積之和，然後呼叫函式sigmoid獲得分類結果
prob = sigmoid(sum(inX * weights))
if prob > 0.5: #如果分類結果大於0.5，則認為類別1；否則分為類別0
return 1.0
else:
return 0.0