邏輯回歸原理（python代碼實現）

Logistic Regression Classifier邏輯回歸主要思想就是用最大似然概率方法構建出方程，為最大化方程，利用牛頓梯度上升求解方程參數。

• 優點：計算代價不高，易於理解和實現。
• 缺點：容易欠擬合，分類精度可能不高。
• 使用數據類型：數值型和標稱型數據。

介紹邏輯回歸之前，我們先看一問題，有個黑箱，里面有白球和黑球，如何判斷它們的比例。

我們從里面抓3個球，2個黑球，1個白球。這時候，有人就直接得出了黑球67%，白球占比33%。這個時候，其實這個人使用了最大似然概率的思想，通俗來講，當黑球是67%的占比的時候，我們抓3個球，出現2黑1白的概率最大。我們直接用公式來說明。

假設黑球占比為P，白球為1-P。於是我們要求解MAX(P*P*(1-P))，顯而易見P=67%（求解方法：對方程求導，使導數為0的P值即為最優解）

我們看邏輯回歸，解決的是二分類問題，是不是和上面黑球白球問題很像，是的，邏輯回歸也是最大似然概率來求解。

假設我們有n個獨立的訓練樣本{(x1, y1) ,(x2, y2),…, (xn, yn)}，y={0, 1}。那每一個觀察到的樣本(xi, yi)出現的概率是： 
 
上面為什么是這樣呢？當y=1的時候，后面那一項是不是沒有了，那就只剩下x屬於1類的概率，當y=0的時候，第一項是不是沒有了，那就只剩下后面那個x屬於0的概率（1減去x屬於1的概率）。所以不管y是0還是1，上面得到的數，都是(x, y)出現的概率。那我們的整個樣本集，也就是n個獨立的樣本出現的似然函數為（因為每個樣本都是獨立的，所以n個樣本出現的概率就是他們各自出現的概率相乘）： 

這里我們稍微變換下L(θ)：取自然對數，然后化簡（不要看到一堆公式就害怕哦，很簡單的哦，只需要耐心一點點，自己動手推推就知道了。注：有xi的時候，表示它是第i個樣本，下面沒有做區分了，相信你的眼睛是雪亮的），得到：

 
其中第三步到第四步使用了下面替換。 
 
這時候為求最大值，對L(θ)對θ求導，得到： 
 
然后我們令該導數為0，即可求出最優解。但是這個方程是無法解析求解（這里就不證明了）。 
最后問題變成了，求解參數使方程L最大化，求解參數的方法梯度上升法（原理這里不解釋了，看詳細的代碼的計算方式應該更容易理解些）。 
根據這個轉換公式 
 
我們代入參數和特征，求P，也就是發生1的概率。 
 
上面這個也就是常提及的sigmoid函數，俗稱激活函數，最后用於分類（若P(y=1|x;ΘΘ )大於0.5，則判定為1）。

下面是詳細的邏輯回歸代碼，代碼比較簡單，主要是要理解上面的算法思想。個人建議，可以結合代碼看一步一步怎么算的，然后對比上面推導公式，可以讓人更加容易理解，並加深印象。

from numpy import *
filename='...\\testSet.txt' #文件目錄
def loadDataSet():   #讀取數據（這里只有兩個特征）
    dataMat = []
    labelMat = []
    fr = open(filename)
    for line in fr.readlines():
        lineArr = line.strip().split()
        dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])   #前面的1，表示方程的常量。比如兩個特征X1,X2，共需要三個參數，W1+W2*X1+W3*X2
        labelMat.append(int(lineArr[2]))
    return dataMat,labelMat

def sigmoid(inX):  #sigmoid函數
    return 1.0/(1+exp(-inX))

def gradAscent(dataMat, labelMat): #梯度上升求最優參數
    dataMatrix=mat(dataMat) #將讀取的數據轉換為矩陣
    classLabels=mat(labelMat).transpose() #將讀取的數據轉換為矩陣
    m,n = shape(dataMatrix)
    alpha = 0.001  #設置梯度的閥值，該值越大梯度上升幅度越大
    maxCycles = 500 #設置迭代的次數，一般看實際數據進行設定，有些可能200次就夠了
    weights = ones((n,1)) #設置初始的參數，並都賦默認值為1。注意這里權重以矩陣形式表示三個參數。
    for k in range(maxCycles):
        h = sigmoid(dataMatrix*weights)
        error = (classLabels - h)     #求導后差值
        weights = weights + alpha * dataMatrix.transpose()* error #迭代更新權重
    return weights

def stocGradAscent0(dataMat, labelMat):  #隨機梯度上升，當數據量比較大時，每次迭代都選擇全量數據進行計算，計算量會非常大。所以采用每次迭代中一次只選擇其中的一行數據進行更新權重。
    dataMatrix=mat(dataMat)
    classLabels=labelMat
    m,n=shape(dataMatrix)
    alpha=0.01
    maxCycles = 500
    weights=ones((n,1))
    for k in range(maxCycles):
        for i in range(m): #遍歷計算每一行
            h = sigmoid(sum(dataMatrix[i] * weights))
            error = classLabels[i] - h
            weights = weights + alpha * error * dataMatrix[i].transpose()
    return weights

def stocGradAscent1(dataMat, labelMat): #改進版隨機梯度上升，在每次迭代中隨機選擇樣本來更新權重，並且隨迭代次數增加，權重變化越小。
    dataMatrix=mat(dataMat)
    classLabels=labelMat
    m,n=shape(dataMatrix)
    weights=ones((n,1))
    maxCycles=500
    for j in range(maxCycles): #迭代
        dataIndex=[i for i in range(m)]
        for i in range(m): #隨機遍歷每一行
            alpha=4/(1+j+i)+0.0001  #隨迭代次數增加，權重變化越小。
            randIndex=int(random.uniform(0,len(dataIndex)))  #隨機抽樣
            h=sigmoid(sum(dataMatrix[randIndex]*weights))
            error=classLabels[randIndex]-h
            weights=weights+alpha*error*dataMatrix[randIndex].transpose()
            del(dataIndex[randIndex]) #去除已經抽取的樣本
    return weights

def plotBestFit(weights):  #畫出最終分類的圖
    import matplotlib.pyplot as plt
    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()

def main():
    dataMat, labelMat = loadDataSet()
    weights=gradAscent(dataMat, labelMat).getA()
    plotBestFit(weights)

if __name__=='__main__':
    main()

 

跑完代碼結果： 
 
當然，還可以換隨機梯度上升和改進的隨機梯度上升算法試試，效果都還不錯。 
下面是代碼使用的數據，可以直接復制本地text里面，跑上面代碼。

-0.017612   14.053064   0
-1.395634   4.662541    1
-0.752157   6.538620    0
-1.322371   7.152853    0
0.423363    11.054677   0
0.406704    7.067335    1
0.667394    12.741452   0
-2.460150   6.866805    1
0.569411    9.548755    0
-0.026632   10.427743   0
0.850433    6.920334    1
1.347183    13.175500   0
1.176813    3.167020    1
-1.781871   9.097953    0
-0.566606   5.749003    1
0.931635    1.589505    1
-0.024205   6.151823    1
-0.036453   2.690988    1
-0.196949   0.444165    1
1.014459    5.754399    1
1.985298    3.230619    1
-1.693453   -0.557540   1
-0.576525   11.778922   0
-0.346811   -1.678730   1
-2.124484   2.672471    1
1.217916    9.597015    0
-0.733928   9.098687    0
-3.642001   -1.618087   1
0.315985    3.523953    1
1.416614    9.619232    0
-0.386323   3.989286    1
0.556921    8.294984    1
1.224863    11.587360   0
-1.347803   -2.406051   1
1.196604    4.951851    1
0.275221    9.543647    0
0.470575    9.332488    0
-1.889567   9.542662    0
-1.527893   12.150579   0
-1.185247   11.309318   0
-0.445678   3.297303    1
1.042222    6.105155    1
-0.618787   10.320986   0
1.152083    0.548467    1
0.828534    2.676045    1
-1.237728   10.549033   0
-0.683565   -2.166125   1
0.229456    5.921938    1
-0.959885   11.555336   0
0.492911    10.993324   0
0.184992    8.721488    0
-0.355715   10.325976   0
-0.397822   8.058397    0
0.824839    13.730343   0
1.507278    5.027866    1
0.099671    6.835839    1
-0.344008   10.717485   0
1.785928    7.718645    1
-0.918801   11.560217   0
-0.364009   4.747300    1
-0.841722   4.119083    1
0.490426    1.960539    1
-0.007194   9.075792    0
0.356107    12.447863   0
0.342578    12.281162   0
-0.810823   -1.466018   1
2.530777    6.476801    1
1.296683    11.607559   0
0.475487    12.040035   0
-0.783277   11.009725   0
0.074798    11.023650   0
-1.337472   0.468339    1
-0.102781   13.763651   0
-0.147324   2.874846    1
0.518389    9.887035    0
1.015399    7.571882    0
-1.658086   -0.027255   1
1.319944    2.171228    1
2.056216    5.019981    1
-0.851633   4.375691    1
-1.510047   6.061992    0
-1.076637   -3.181888   1
1.821096    10.283990   0
3.010150    8.401766    1
-1.099458   1.688274    1
-0.834872   -1.733869   1
-0.846637   3.849075    1
1.400102    12.628781   0
1.752842    5.468166    1
0.078557    0.059736    1
0.089392    -0.715300   1
1.825662    12.693808   0
0.197445    9.744638    0
0.126117    0.922311    1
-0.679797   1.220530    1
0.677983    2.556666    1
0.761349    10.693862   0
-2.168791   0.143632    1
1.388610    9.341997    0
0.317029    14.739025   0

 

參考： 
- http://blog.csdn.net/zouxy09/article/details/20319673 
- Machine Learning in Action 
- 統計學習方法

 

轉載：http://blog.csdn.net/csqazwsxedc/article/details/69690655