機器學習之邏輯迴歸(純python實現)

logistic迴歸是一種廣義的線性迴歸，通過構造迴歸函式，利用機器學習來實現分類或者預測。

原理

上一文簡單介紹了線性迴歸，與邏輯迴歸的原理是類似的。

1. 預測函式(h)。該函式就是分類函式，用來預測輸入資料的判斷結果。過程非常關鍵，需要預測函式的“大概形式”， 比如是線性還是非線性的。 本文參考機器學習實戰的相應部分，看一下資料集。

// 兩個特徵
-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
複製程式碼

如上圖，紅綠代表兩種不同的分類。可以預測分類函式大概是一條直線。

2. Cost函式(損失函式)：該函式預測的輸出h和訓練資料類別y之間的偏差，(h-y)或者其他形式。綜合考慮所有訓練資料的cost， 將其求和或者求平均，極為J函式， 表示所有訓練資料預測值和實際值的偏差。

3. 顯然，J函式的值越小，表示預測的函式越準確(即h函式越準確)，因此需要找到J函式的最小值。有時需要用到梯度下降。

具體過程

構造預測函式

邏輯迴歸名為迴歸，實際為分類，用於兩分類問題。 這裡直接給出sigmoid函式。

接下來確定分類的邊界，上面有提到，該資料集需要一個線性的邊界。 不同資料需要不同的邊界。

確定了分類函式，將其輸入記做z ，那麼

向量x是特徵變數, 是輸入資料。此資料有兩個特徵，可以表示為z = w0x0 + w1x1 + w2x2。w0是常數項，需要構造x0等於1(見後面程式碼)。 向量W是迴歸係數特徵，T表示為列向量。 之後就是確定最佳迴歸係數w(w0, w1, w2)。

cost函式

綜合以上，預測函式為：

這裡不做推導，可以參考文章 Logistic迴歸總結

有了上述的cost函式，可以使用梯度上升法求函式J的最小值。推導見上述連結。

綜上：梯度更新公式如下：

接下來是python程式碼實現：

# sigmoid函式和初始化資料
def sigmoid(z):
    return 1 / (1 + np.exp(-z))

def init_data():
    data = np.loadtxt('data.csv')
    dataMatIn = data[:, 0:-1]
    classLabels = data[:, -1]
    dataMatIn = np.insert(dataMatIn, 0, 1, axis=1)  #特徵資料集，新增1是構造常數項x0
    return dataMatIn, classLabels
複製程式碼

//  梯度上升
def grad_descent(dataMatIn, classLabels):
    dataMatrix = np.mat(dataMatIn)  #(m,n)
    labelMat = np.mat(classLabels).transpose()
    m, n = np.shape(dataMatrix)
    weights = np.ones((n, 1))  #初始化迴歸係數（n, 1)
    alpha = 0.001 #步長
    maxCycle = 500  #最大迴圈次數

    for i in range(maxCycle):
        h = sigmoid(dataMatrix * weights)  #sigmoid 函式
        weights = weights + alpha * dataMatrix.transpose() * (labelMat - h)  #梯度
    return weights
複製程式碼

// 計算結果
if __name__ == '__main__':
    dataMatIn, classLabels = init_data()
    r = grad_descent(dataMatIn, classLabels)
    print(r)
複製程式碼

輸入如下：

[[ 4.12414349]
 [ 0.48007329]
 [-0.6168482 ]]
複製程式碼

上述w就是所求的迴歸係數。w0 = 4.12414349, w1 = 0.4800, w2=-0.6168 之前預測的直線方程0 = w0x0 + w1x1 + w2x2， 帶入迴歸係數，可以確定邊界。 x2 = (-w0 - w1*x1) / w2

畫出函式影象：

def plotBestFIt(weights):
    dataMatIn, classLabels = init_data()
    n = np.shape(dataMatIn)[0]
    xcord1 = []
    ycord1 = []
    xcord2 = []
    ycord2 = []
    for i in range(n):
        if classLabels[i] == 1:
            xcord1.append(dataMatIn[i][1])
            ycord1.append(dataMatIn[i][2])
        else:
            xcord2.append(dataMatIn[i][1])
            ycord2.append(dataMatIn[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 = np.arange(-3, 3, 0.1)
    y = (-weights[0, 0] - weights[1, 0] * x) / weights[2, 0]  #matix
    ax.plot(x, y)
    plt.xlabel('X1')
    plt.ylabel('X2')
    plt.show()
複製程式碼

如下：

演算法改進

隨機梯度上升

上述演算法中，每次迴圈矩陣都會進行m * n次乘法計算，時間複雜度是maxCycles* m * n。當資料量很大時， 時間複雜度是很大。 這裡嘗試使用隨機梯度上升法來進行改進。 隨機梯度上升法的思想是，每次只使用一個資料樣本點來更新迴歸係數。這樣就大大減小計算開銷。 演算法如下：

def stoc_grad_ascent(dataMatIn, classLabels):
    m, n = np.shape(dataMatIn)
    alpha = 0.01
    weights = np.ones(n)
    for i in range(m):
        h = sigmoid(sum(dataMatIn[i] * weights))  #數值計算
        error = classLabels[i] - h
        weights = weights + alpha * error * dataMatIn[i]
    return weights
複製程式碼

進行測試：

隨機梯度上升的改進

def stoc_grad_ascent_one(dataMatIn, classLabels, numIter=150):
    m, n = np.shape(dataMatIn)
    weights = np.ones(n)
    for j in range(numIter):
        dataIndex = list(range(m))
        for i in range(m):
            alpha = 4 / (1 + i + j) + 0.01 #保證多次迭代後新資料仍然有影響力
            randIndex = int(np.random.uniform(0, len(dataIndex)))
            h = sigmoid(sum(dataMatIn[i] * weights))  # 數值計算
            error = classLabels[i] - h
            weights = weights + alpha * error * dataMatIn[i]
            del(dataIndex[randIndex])
    return weights
複製程式碼

可以對上述三種情況的迴歸係數做個波動圖。 可以發現第三種方法收斂更快。 評價演算法優劣勢看它是或否收斂，是否達到穩定值，收斂越快，演算法越優。

總結

這裡用到的梯度上升和梯度下降是一樣的，都是求函式的最值， 符號需要變一下。 梯度意味著分別沿著x， y的方向移動一段距離。(cost分別對x， y)的導數。

完整程式碼請檢視： github: logistic regression

參考文章： 機器學習之Logistic迴歸與Python實現
機器學習筆記:Logistic迴歸總結
機器學習基本演算法系列之邏輯迴歸