**前言：** 生活中，人們經常會遇到各種最優化問題，比如如何在最短時間從一個地點到另外一個地點？如何在投入最少的資金而卻能得到最高的受益？如何設計一款芯片使其功耗最低而性能最好？這一節就要學習一種最優化算法——Logistic回歸，設計最優化算法的目的依然是用于分類。在這里，**Logistic回歸的主要思想是根據現有的數據對分類邊界線建立回歸公式，達到分類的目的。**假設我們有一堆數據，需要劃一條線(最佳直線)對其分類，這就是Logistic回歸的目的。 而“Logistic回歸”中的“回歸”又是代表什么？數學認為，回歸是一個擬合過程，回歸分析本質上就是一個函數估計的問題，就是找出因變量和自變量之間的因果關系。具體到例子，假設我們有一些數據點，現在使用一條直線對這些點進行擬合，使得這條線盡可能地表示數據點的分布，這個擬合過程就稱作“回歸”。 在機器學習任務中，訓練分類器時就是尋找最佳擬合曲線的過程，因此接下來將使用最優化算法。在實現算法之前，先總結Logistic回歸的一些屬性： * 優點：計算代價不高，易于理解和實現 * 缺點：容易欠擬合，分類精度可能不高 * 適用的數據類型：數值型和標稱型數據 **目錄：** 一、基于Sigmoid函數的Logistic回歸 1. Sigmoid函數 2. 基于最優化方法的最佳回歸系數確定 3. 梯度上升法的實現 4. 隨機梯度上升算法的實現 二、一個實例：從疝氣病癥預測病馬的死亡率 1. 準備數據 2. 測試算法，使用Logistic回歸進行分類 三、總結 =============================================== **一、基于Sigmoid函數的Logistic回歸** **1.Sigmoid函數** Logistic回歸想要得到的函數是，能接受所有的輸入然后返回預測的類別，比如，在兩類情況下函數應輸出類別0或1。sigmoid函數可是勝任這一工作，它像是一個階躍函數。其公式如下：  其中：  向量***w***稱為回歸系數，也就是我們要找到的最佳參數，?***x***是`n`維的特征向量，是分類器的輸入數據。 下面附上書本上的圖，該圖是函數在不同坐標尺度下的曲線圖：  為了實現Logistic回歸分類器，我們可以在每個特征上乘以一個回歸系數，然后把所有的結果值相加，將這個總和結果代入sigmoid函數中，進而得到一個范圍在`0~1`之間的數值。任何大于`0.5`的數據被分入`1`類，小于`0.5`的數據被歸入`0`類。所以，Logistic回歸也可以被看成是一種概率估計。 **2.基于最優化方法的最佳回歸系數確定** 上面提到Sigmoid函數里的一部分：  其中向量***w***稱為回歸系數，也就是我們要找到的最佳參數，?***x***是`n`維的特征向量，是分類器的輸入數據。接下來將介紹幾種需找最佳參數的方法： * 梯度上升法：基于的思想是要找到某函數的最大值，最好的方法是沿著該函數的梯度方向探尋。 * 梯度下降法：與梯度上升的思想相似，但方向相反，即如果要找到某函數的最小值，最好的方法是沿著該函數的梯度方向的反方向尋找。 這里使用梯度上升法，對于一個函數`f(x,y)`，其梯度表示方法如下：  該梯度意味著要沿`x`和`y`的方向分別移動一定的距離，這其實是**確立了算法到達每個點后下一步移動的方向**。其中，函數`f(x,y)`?必須要在待計算的點上有定義并且可微。 **移動方向**確定了，這里我們定義移動的大小為步長，用`α`表示，用向量來表示的話，梯度上升算法的迭代公式如下：  該公式表明，最佳參數***w***的結果收到梯度變化和步長`α`的影響，這個公式將一直被迭代執行，直到滿足某個停止條件。 **3.梯度上升法的實現** 以下使用Python實現了梯度上升法，并找到了最佳參數，其中使用的樣本數據有`100`個，每個數據包含兩個數值型的特征：`X1`和`X2`?，這些數據點可以被繪制成二維散點圖。在找出最佳參數***w***后，再編寫函數畫出不同類別數據之間的分割線，觀察最優化算法的效果。 ~~~ # -*- coding: utf-8 -*- """ Created on Sat Sep 12 22:53:07 2015 @author: Herbert """ from numpy import * sys.setrecursionlimit(1500) def loadData(): dataMat = []; labelMat = [] fr = open('testSet.txt') for line in fr.readlines(): lineSplit = line.strip().split() dataMat.append([1.0, float(lineSplit[0]), float(lineSplit[1])]) labelMat.append(int(lineSplit[2])) return dataMat, labelMat def sigmoid(x): return 1.0 / (1 + exp(-x)) def gradAscent(data, label): dataMat = mat(data) labelMat = mat(label).transpose() m, n = shape(dataMat) alpha = 0.001 maxCycles = 500 w = ones((n, 1)) for k in range(maxCycles): p = sigmoid(dataMat * w) error = labelMat - p w = w + alpha * dataMat.transpose() * error return dataMat, labelMat, w data, label = loadData() dataMat, labelMat, w = gradAscent(data, label) print labelMat def plotBestSplit(w): import matplotlib.pyplot as plt dataMat, labelMat = loadData() 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]) # array only 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 = 'blue') x = arange(-3.0, 3.0, 0.1) y = (-w[0] - w[1] * x) / w[2] ax.plot(x, y) plt.xlabel('x1'); plt.ylabel('x2'); plt.show() plotBestSplit(w.getA()) ~~~ 結果如下，可看到使用了Logistic回歸的分類結果還是挺不錯的，盡管有三四個樣本點被錯分：  **4.隨機梯度上升算法的實現** 梯度上升算法在處理`100`個左右的數據集時尚可，但如果有數十億樣本和成千上萬的特征，那么該方法的計算復雜度將變得很恐怖。改進方法為隨機梯度上升算法，該方法一次僅用一個樣本點來更新回歸系數。它占用更少的計算資源，是一種在線算法，可以在數據到來時就完成參數的更新，而不需要重新讀取整個數據集來進行批處理運算。一次處理所有的數據被稱為批處理。以下代碼實現了隨機梯度上升算法： ~~~ # 隨機梯度上升算法 def stocGradAscent(data, label): m,n = shape(data) alpha = 0.01 w = ones(n) for i in range(m): h = sigmoid(sum(data[i]) * w) error = label[i] - h w = w + alpha * error *data[i] return w wNew = stocGradAscent(data, label) plotBestSplit(wNew) ~~~ 結果：  由上圖可以看出，隨機梯度上升算法分錯了很多樣本點，其分類效果并沒有原先實現的普通梯度上升算法好。 書中給出的解釋是，普通梯度上升算法是在整個數據集上迭代`500`次得到的結果，而隨機梯度上升算法只是迭代了`100`次。一個判斷算法優劣的可靠方法是看它是否收斂，也就是說求解的參數是否達到了穩定值，是否還會不斷變化。 書中給出了一個例子，讓隨機梯度上升算法在整個數據集上運行`200`次，迭代過程中?**w向量**?的3個參數的變化如下圖：  圖中，**w**的第二個參數最先達到穩定，而其他兩個參數則還需要更多的迭代次數來達到穩定。而我們也發現，在整個迭代過程中，**w**的三個參數都有不同程度的波動。產生這種現象的原因是存在一些被錯分的樣本點，這些樣本在參與每次迭代的過程中都會引起參數的劇烈變化。為了避免這種劇烈波動，并改善算法的性能，采用以下策略對隨機梯度上升算法做了改進： * 在每次迭代時更新步長alpha的值，使得alpha的值不斷減小，但是不能減少為`0`，這樣做的原因是為了保證在多次迭代后新數據仍然具有一定的影響。 * 通過隨機采取樣本來更新回歸參數，這種方法將減少周期性的波動，**每次隨機從列表中取一個值，然后該值會從列表中刪除。** 實現代碼： ~~~ def stocGradAscentAdvanced(data, label, numIter = 150): m,n = shape(data) w = ones(n) for i in range(numIter): dataIndex = range(m) for j in range(m): alpha = 4 / (1.0 + i + j) + 0.01 randIndex = int(random.uniform(0, len(dataIndex))) h = sigmoid(sum(data[randIndex] * w)) error = label[randIndex] - h w = w + alpha * error* array(data[randIndex]) del(dataIndex[randIndex]) return w wNewAdv = stocGradAscentAdvanced(data, label, numIter = 150) plotBestSplit(wNewAdv) ~~~ 在使用了優化策略后，各個回歸系數參數的變化情況大幅改善，參數的收斂得更快了：  使用改進后的隨機梯度上升算法對樣本點進行劃分，效果與普通梯度上升算法相當，但其中所用的計算量更少了：  **二、一個實例：從疝氣病癥預測病馬的死亡率** **1.準備數據** 該實例使用Logistic回歸來預測患有疝病的馬的存活問題。這里的數據來自2010年1月11日的UCI機器學習數據庫，其中包含`368`個樣本和**28**個特征。**這里的數據集是有`30%`的數據缺失的**。-.- 在實現算法之前，先了解一些關于處理數據中缺失值的方法： * 使用可用特征的均值來填補缺失值； * 使用特殊值來填補缺失值，如-1； * 忽略有缺失值的樣本； * 使用相似樣本的均值來填補缺失值； * 使用另外的機器學習算法預測缺失值； * 對于類別標簽丟失的數據，只能將該數據丟棄。 這里我們使用實數`0`來替換所有缺失的值，這種方法恰好適用于Logistic回歸；而對于類別標簽丟失的數據，則須將該數據丟棄。 **2.測試算法，使用Logistic回歸進行分類** 基于以上工作，以下代碼把測試集上的每個特征向量乘以最優化算法得來的回歸系數w，再將該乘積結果求和，最后輸入到Sigmoid函數，如果對應的Sigmoid函數值大于`0.5`，則將該樣本的類別判定為`1`，否則判定為`0`；最后，統計判定結果與實際結果的誤差，由于誤差具有不確定性，程序最后使用了`10`次分類結果求平均的方法，得出算法的平均分類錯誤率。 ~~~ def classifyVector(x, w): prob = sigmoid(sum(x * w)) if prob > 0.5: return 1.0 else: return 0.0 def horseColicTest(): frTrain = open('horseColicTraining.txt') frTest = open('horseColicTest.txt') trainingData = []; trainingLabel = [] for line in frTrain.readlines(): currLine = line.strip().split('\t') lineArr = [] for i in range(21): lineArr.append(float(currLine[i])) trainingData.append(lineArr) trainingLabel.append(float(currLine[21])) w = stocGradAscentAdvanced(array(trainingData), trainingLabel, numIter = 500) errorCount = 0.0; numOfTest = 0.0 for line in frTest.readlines(): numOfTest += 1.0 currLine = line.strip().split('\t') lineArr = [] for i in range(21): lineArr.append(float(currLine[i])) if int(classifyVector(array(lineArr), w)) != int(currLine[21]): errorCount += 1 errorRate = float(errorCount) / numOfTest print "The error rate of this test is: %f" % errorRate return errorRate def finalTest(): numOfTest = 10; errorSum = 0.0 for i in range(numOfTest): errorSum += horseColicTest() print "After %d iterations the average error rate is: %f" \ % (numOfTest, errorSum / float(numOfTest)) finalTest() ~~~  分類錯誤率是`39.4%`，乍一看是挺高的。但如果你知道所使用的數據集是有`30%`的數據缺失的，你就不會這么想了…-.- **三、總結** Logistic回歸的目的是尋找一個非線性函數sigmoid的最佳擬合參數，求解過程可以由最優化算法來完成。在最優化算法中，最常用的就是梯度上升算法，而梯度上升算法又可以簡化為隨機梯度上升算法。 隨機梯度上升算法和梯度上升算法的效果相當，但占用更少的計算資源。隨機梯度上升算法是一種在線算法，可以在數據到來時就完成參數的更新，而不需要重新讀取整個數據集來進行批處理運算。