Sigmoid函數的輸入記為z，由下面公式得出：
如果採用向量的寫法，上述公式可以寫成z = wTx。它表示將這兩個數值向量的對應元素相乘然後全部加起來即得到z值。其中的向量x是分類器的輸入數據，向量w也就是我們要找到的最佳參數（係數），從而使得分類盡可能地精確。為了尋找該最佳參數，需要用到最優化理論的一些知識。
下面首先介紹梯度上升的最優化方法，我們將學習到如何使用該方法求得數據集的最佳參數。接下來，展示如何繪製梯度上升法產生的決策邊界圖，該圖能將梯度上升法的分類效果可視化地呈現出來。最後我們將學習隨機梯度上升算法，以及如何對其進行修改以獲得更好的結果。
5.2.1 梯度上升法
我們介紹的第一個最優化算法叫做梯度上升法。梯度上升法基於的思想是：要找到某函數的最大值，最好的方法是沿著該函數的梯度方向探尋。如果梯度記為∇，則函數f(x,y)的梯度由下式表示：
這是機器學習中最易造成混淆的一個地方，但在數學上並不難，需要做的只是牢記這些符號的意義。這個梯度意味著要沿x的方向移動
，沿y的方向移動
。其中，函數f(x,y)必須要在待計算的點上有定義並且可微。一個具體的函數例子見圖5-2。
圖5-2 梯度上升算法到達每個點後都會重新估計移動的方向。從P0開始，計算完該點的梯度，函數就根據梯度移動到下一點P1。在P1點，梯度再次被重新計算，並沿新的梯度方向移動到P2。如此循環迭代，直到滿足停止條件。迭代的過程中，梯度算子總是保證我們能選取到最佳的移動方向
圖5-2中的梯度上升算法沿梯度方向移動了一步。可以看到，梯度算子總是指向函數值增長最快的方向。這裡所說的是移動方向，而未提到移動量的大小。該量值稱為步長，記做
。用向量來表示的話，梯度上升算法的迭代公式如下：
該公式將一直被迭代執行，直至達到某個停止條件為止，比如迭代次數達到某個指定值或算法達到某個可以允許的誤差範圍。
梯度下降算法
你最經常聽到的應該是梯度下降算法，它與這裡的梯度上升算法是一樣的，只是公式中的加法需要變成減法。因此，對應的公式可以寫成
梯度上升算法用來求函數的最大值，而梯度下降算法用來求函數的最小值。
基於上面的內容，我們來看一個Logistic回歸分類器的應用例子，從圖5-3可以看到我們採用的數據集。
圖5-3 一個簡單數據集，下面將採用梯度上升法找到Logistic回歸分類器在此數據集上的最佳回歸係數
5.2.2 訓練算法：使用梯度上升找到最佳參數
圖5-3中有100個樣本點，每個點包含兩個數值型特徵：X1和X2。在此數據集上，我們將通過使用梯度上升法找到最佳回歸係數，也就是擬合出Logistic回歸模型的最佳參數。
梯度上升法的偽代碼如下：
每個回歸係數初始化為1
重複R次：
計算整個數據集的梯度
使用`alpha × gradient`更新回歸係數的向量
返回回歸係數
下面的代碼是梯度上升算法的具體實現。為瞭解實際效果，打開文本編輯器並創建一個名為logRegres.py的文件，輸入下列代碼：
程序清單5-1 logistic回歸梯度上升優化算法
def loadDataSet:
dataMat = ; labelMat =
fr = open(\'testSet.txt\')
for line in fr.readlines:
lineArr = line.strip.split
dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])
labelMat.append(int(lineArr[2]))
return dataMat,labelMat
def sigmoid(inX):
return 1.0/(1+exp(-inX))
def gradAscent(dataMatIn, classLabels):
# ❶（以下兩行）轉換為NumPy矩陣數據類型
dataMatrix = mat(dataMatIn)
labelMat = mat(classLabels).transpose
m,n = shape(dataMatrix)
alpha = 0.001
maxCycles = 500
weights = ones((n,1))
for k in range(maxCycles):
#❷（以下三行）矩陣相乘
h = sigmoid(dataMatrix*weights)
error = (labelMat - h)
weights = weights + alpha * dataMatrix.transpose* error
return weights
程序清單5-1的代碼在開頭提供了一個便利函數loadDataSet，它的主要功能是打開文本文件testSet.txt並逐行讀取。每行前兩個值分別是X1和X2，第三個值是數據對應的類別標籤。此外，為了方便計算，該函數還將X0的值設為1.0。接下來的函數是5.2節提到的函數sigmoid。
梯度上升算法的實際工作是在函數gradAscent裡完成的，該函數有兩個參數。第一個參數是dataMatIn，它是一個2維NumPy數組，每列分別代表每個不同的特徵，每行則代表每個訓練樣本。我們現在採用的是100個樣本的簡單數據集，它包含了兩個特徵X1和X2，再加上第0維特徵X0，所以dataMath裡存放的將是100×3的矩陣。在❶處，我們獲得輸入數據並將它們轉換成NumPy矩陣。這是本書首次使用NumPy矩陣，如果你對矩陣數學不太熟悉，那麼一些運算可能就會不易理解。比如，NumPy對2維數組和矩陣都提供一些操作支持，如果混淆了數據類型和對應的操作，執行結果將與預期截然不同。對此，本書附錄A給出了對NumPy矩陣的介紹。第二個參數是類別標籤，它是一個1×100的行向量。為了便於矩陣運算，需要將該行向量轉換為列向量，做法是將原向量轉置，再將它賦值給labelMat。接下來的代碼是得到矩陣大小，再設置一些梯度上升算法所需的參數。
變量alpha是向目標移動的步長，maxCycles是迭代次數。在for循環迭代完成後，將返回訓練好的回歸係數。需要強調的是，在❷處的運算是矩陣運算。變量h不是一個數而是一個列向量，列向量的元素個數等於樣本個數，這裡是100。對應地，運算dataMatrix * weights代表的不是一次乘積計算，事實上該運算包含了300次的乘積。
最後還需說明一點，你可能對❷中公式的前兩行覺得陌生。此處本書略去了一個簡單的數學推導，我把它留給有興趣的讀者。定性地說，這裡是在計算真實類別與預測類別的差值，接下來就是按照該差值的方向調整回歸係數。
接下來看看實際效果，打開文本編輯器，添加程序清單5-1的代碼。
在Python提示符下，敲入下面的代碼：
>>> import logRegres
>>> dataArr,labelMat=logRegres.loadDataSet
>>> logRegres.gradAscent(dataArr,labelMat)
matrix([[ 4.12414349],
[ 0.48007329],
[-0.6168482 ]])
5.2.3 分析數據：畫出決策邊界
上面已經解出了一組回歸係數，它確定了不同類別數據之間的分隔線。那麼怎樣能畫出該分隔線，從而使得優化的過程便於理解呢？下面將解決這個問題，打開logRegres.py並添加如下代碼。
程序清單5-2 畫出數據集和logistic回歸最佳擬合直線的函數
def plotBestFit(wei):
import matplotlib.pyplot as plt
weights = wei.getA
dataMat,labelMat=loadDataSet
dataArr = array(dataMat)
n = shape(dataArr)[0]
xcord1 = ; ycord1 =
xcord2 = ; ycord2 =
for i in range(n):
if int(labelMat[i])== 1:
xcord1.append(dataArr[i,1]); ycord1.append(dataArr[i,2])
else:
xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])
fig = plt.figure
ax = fig.add_subplot(111)
ax.scatter(xcord1, ycord1, s=30, c=\'red\', marker=\'s\')
ax.scatter(xcord2, ycord2, s=30, c=\'green\')
x = arange(-3.0, 3.0, 0.1)
# ❶ 最佳擬合直線
y = (-weights[0]-weights[1]*x)/weights[2]
ax.plot(x, y)
plt.xlabel(\'X1\'); plt.ylabel(\'X2\');
plt.show
程序清單5-2中的代碼是直接用Matplotlib畫出來的。唯一要指出的是，❶處設置了sigmoid函數為0。回憶5.2節，0是兩個類別（類別1和類別0）的分界處。因此，我們設定 0 = w0x0 + w1x1 + w2x2，然後解出X2和X1的關係式（即分隔線的方程，注意X0=1）。
運行程序清單5-2的代碼，在Python提示符下輸入：
>>> reload(logRegres)
<module \'logRegres\' from \'logRegres.py\'>
>>> logRegres.plotBestFit(weights.getA)
輸出的結果如圖5-4所示。
圖5-4 梯度上升算法在500次迭代後得到的Logistic回歸最佳擬合直線
這個分類結果相當不錯，從圖上看只錯分了兩到四個點。但是，儘管例子簡單且數據集很小，這個方法卻需要大量的計算（300次乘法）。因此下一節將對該算法稍作改進，從而使它可以用在真實數據集上。
5.2.4 訓練算法：隨機梯度上升
梯度上升算法在每次更新回歸係數時都需要遍歷整個數據集，該方法在處理100個左右的數據集時尚可，但如果有數十億樣本和成千上萬的特徵，那麼該方法的計算複雜度就太高了。一種改進方法是一次僅用一個樣本點來更新回歸係數，該方法稱為隨機梯度上升算法。由於可以在新樣本到來時對分類器進行增量式更新，因而隨機梯度上升算法是一個在線學習算法。與「在線學習」相對應，一次處理所有數據被稱作是「批處理」。
隨機梯度上升算法可以寫成如下的偽代碼：
所有回歸係數初始化為1
對數據集中每個樣本
計算該樣本的梯度
使用alpha × gradient更新回歸係數值
返回回歸係數值
偽代碼最後一行應頂格，原文錯誤，已修改
以下是隨機梯度上升算法的實現代碼。
程序清單5-3 隨機梯度上升算法
def stocGradAscent0(dataMatrix, classLabels):
m,n = shape(dataMatrix)
alpha = 0.01
weights = ones(n)
for i in range(m):
h = sigmoid(sum(dataMatrix[i]*weights))
error = classLabels[i] - h
weights = weights + alpha * error * dataMatrix[i]
return weights
可以看到，隨機梯度上升算法與梯度上升算法在代碼上很相似，但也有一些區別：第一，後者的變量h和誤差error都是向量，而前者則全是數值；第二，前者沒有矩陣的轉換過程，所有變量的數據類型都是NumPy數組。
為了驗證該方法的結果，我們將程序清單5-3的代碼添加到logRegres.py中，並在Python提示符下輸入如下命令：
>>> from numpy import *
>>> reload(logRegres)
<module \'logRegres\' from \'logRegres.py\'>
>>> dataArr,labelMat=logRegres.loadDataSet
>>> weights=logRegres.stocGradAscent0(array(dataArr),labelMat)
>>> logRegres.plotBestFit(weights)
執行完畢後將得到圖5-5所示的最佳擬合直線圖，該圖與圖5-4有一些相似之處。可以看到，擬合出來的直線效果還不錯，但並不像圖5-4那樣完美。這裡的分類器錯分了三分之一的樣本。
直接比較程序清單5-3和程序清單5-1的代碼結果是不公平的，後者的結果是在整個數據集上迭代了500次才得到的。一個判斷優化算法優劣的可靠方法是看它是否收斂，也就是說參數是否達到了穩定值，是否還會不斷地變化？對此，我們在程序清單5-3中隨機梯度上升算法上做了些修改，使其在整個數據集上運行200次。最終繪製的三個回歸係數的變化情況如圖5-6所示。
圖5-5 隨機梯度上升算法在上述數據集上的執行結果，最佳擬合直線並非最佳分類線
圖5-6 運行隨機梯度上升算法，在數據集的一次遍歷中回歸係數與迭代次數的關係圖。回歸係數經過大量迭代才能達到穩定值，並且仍然有局部的波動現象
圖5-6展示了隨機梯度上升算法在200次迭代過程中回歸係數的變化情況。其中的係數2，也就是圖5-5中的X2只經過了50次迭代就達到了穩定值，但係數1和0則需要更多次的迭代。另外值得注意的是，在大的波動停止後，還有一些小的週期性波動。不難理解，產生這種現象的原因是存在一些不能正確分類的樣本點（數據集並非線性可分），在每次迭代時會引發係數的劇烈改變。我們期望算法能避免來回波動，從而收斂到某個值。另外，收斂速度也需要加快。
對於圖5-6存在的問題，可以通過修改程序清單5-3的隨機梯度上升算法來解決，具體代碼如下：
程序清單5-4 改進的隨機梯度上升算法
def stocGradAscent1(dataMatrix, classLabels, numIter=150):
m,n = shape(dataMatrix)
weights = ones(n)
for j in range(numIter): dataIndex = range(m)
for i in range(m):
#❶ alpha每次迭代時需要調整
alpha = 4/(1.0+j+i)+0.01
#❷ 隨機選取更新
randIndex = int(random.uniform(0,len(dataIndex)))
h = sigmoid(sum(dataMatrix[randIndex]*weights))
error = classLabels[randIndex] - h
weights = weights + alpha * error * dataMatrix[randIndex]
del(dataIndex[randIndex])
return weights
程序清單5-4與5-3類似，但增加了兩處代碼來進行改進。
第一處改進在❶處。一方面，alpha在每次迭代的時候都會調整，這會緩解圖5-6上的數據波動或者高頻波動。另外，雖然alpha會隨著迭代次數不斷減小，但永遠不會減小到0，這是因為❶中還存在一個常數項。必須這樣做的原因是為了保證在多次迭代之後新數據仍然具有一定的影響。如果要處理的問題是動態變化的，那麼可以適當加大上述常數項，來確保新的值獲得更大的回歸係數。另一點值得注意的是，在降低alpha的函數中，alpha每次減少 1/(j+i) ，其中j是迭代次數，i是樣本點的下標1 。這樣當j<<max(i)時，alpha就不是嚴格下降的。避免參數的嚴格下降也常見於模擬退火算法等其他優化算法中。
1. 要注意區分這裡的下標與樣本編號，編號表示了樣本在矩陣中的位置（代碼中為randIndex），而這裡的下標i表示本次迭代中第i個選出來的樣本。——譯者注
程序清單5-4第二個改進的地方在❷處，這裡通過隨機選取樣本來更新回歸係數。這種方法將減少週期性的波動（如圖5-6中的波動）。具體實現方法與第3章類似，這種方法每次隨機從列表中選出一個值，然後從列表中刪掉該值（再進行下一次迭代）。
此外，改進算法還增加了一個迭代次數作為第3個參數。如果該參數沒有給定的話，算法將默認迭代150次。如果給定，那麼算法將按照新的參數值進行迭代。
與stocGradAscent1類似，圖5-7顯示了每次迭代時各個回歸係數的變化情況。
圖5-7 使用樣本隨機選擇和alpha動態減少機制的隨機梯度上升算法stocGradAscent1所生成的係數收斂示意圖。該方法比採用固定alpha的方法收斂速度更快
比較圖5-7和圖5-6可以看到兩點不同。第一點是，圖5-7中的係數沒有像圖5-6里那樣出現週期性的波動，這歸功於stocGradAscent1裡的樣本隨機選擇機制；第二點是，圖5-7的水平軸比圖5-6短了很多，這是由於stocGradAscent1可以收斂得更快。這次我們僅僅對數據集做了20次遍歷，而之前的方法是500次。
下面看看在同一個數據集上的分類效果。將程序清單5-4的代碼添加到logRegres.py文件中，並在Python提示符下輸入：
>>> reload(logRegres)
<module \'logRegres\' from \'logRegres.py\'>
>>> dataArr,labelMat=logRegres.loadDataSet
>>> weights=logRegres.stocGradAscent1(array(dataArr),labelMat)
>>> logRegres.plotBestFit(weights)
程序運行結束之後應該可以看到與圖5-8類似的結果圖。該分隔線達到了與GradientAscent差不多的效果，但是所使用的計算量更少。
圖5-8 使用改進的隨機梯度上升算法得到的係數
默認迭代次數是150，可以通過stocGradAscent的第3個參數來對此進行修改，例如：
>>> weights=logRegres.stocGradAscent1(array(dataArr),labelMat, 500)
目前，我們已經學習了幾個優化算法，但還有很多優化算法值得探討，所幸這方面已有大量的文獻可供參考。另外再說明一下，針對給定的數據集，讀者完全可以對算法的各種參數進行調整，從而達到更好的效果。
迄今為止我們分析了回歸係數的變化情況，但還沒有達到本章的最終目標，即完成具體的分類任務。下一節將使用隨機梯度上升算法來解決病馬的生死預測問題。