本文始發於個人公眾號：TechFlow，原創不易，求個關注

在之前的文章當中，我們推導了線性迴歸的公式，線性迴歸本質是線性函式，模型的原理不難，核心是求解模型引數的過程。通過對線性迴歸的推導和學習，我們基本上了解了機器學習模型學習的過程，這是機器學習的精髓，要比單個模型的原理重要得多。

新關注和有所遺忘的同學可以點選下方的連結回顧一下之前的線性迴歸和梯度下降的內容。

一文講透梯度下降法

詳細推導線性迴歸模型

迴歸與分類

在機器學習當中，模型根據預測結果的不同分為兩類，如果我們希望模型預測一個或者多個連續值，這類問題被稱為是迴歸問題。像是常見的未來股票價格的估計、未來溫度估計等等都算是迴歸問題。還有一類呢是分類問題，模型的預測結果是一個離散值，也就是說只有固定的種類的結果。常見的有垃圾郵件識別、網頁圖片鑑黃等等。

我們之前介紹的邏輯迴歸顧名思義是一個迴歸問題，今天的文章講的呢是如何將這個迴歸模型轉化成分類模型，這個由線性迴歸推導得到的分類模型稱為邏輯迴歸。

邏輯迴歸

邏輯迴歸這個模型很神奇，雖然它的本質也是迴歸，但是它是一個分類模型，並且它的名字當中又包含”迴歸“兩個字，未免讓人覺得莫名其妙。

如果是初學者，覺得頭暈是正常的，沒關係，讓我們一點點捋清楚。

讓我們先回到線性迴歸，我們都知道，線性迴歸當中\(y=WX+b\)。我們通過W和b可以求出X對應的y，這裡的y是一個連續值，是迴歸模型對吧。但如果我們希望這個模型來做分類呢，應該怎麼辦？很容易想到，我們可以人為地設定閾值對吧，比如我們規定y > 0最後的分類是1，y < 0最後的分類是0。從表面上來看，這當然是可以的，但實際上這樣操作會有很多問題。

最大的問題在於如果我們簡單地設計一個閾值來做判斷，那麼會導致最後的y是一個分段函式，而分段函式不連續，使得我們沒有辦法對它求梯度，為了解決這個問題，我們得找到一個平滑的函式使得既可以用來做分類，又可以解決梯度的問題。

很快，資訊學家們找到了這樣一個函式，它就是Sigmoid函式，它的表示式是：

\[\sigma(x)= \frac{1}{1+e^{-x}}\]

它的函式影象如下：

可以看到，sigmoid函式在x=0處取值0.5，在正無窮處極限是1，在負無窮處極限是0，並且函式連續，處處可導。sigmoid的函式值的取值範圍是0-1，非常適合用來反映一個事物發生的概率。我們認為\(\sigma(x)\)表示x發生的概率，那麼x不發生的概率就是\(1-\sigma(x)\)。我們把發生和不發生看成是兩個類別，那麼sigmoid函式就轉化成了分類函式，如果\(\sigma(x)>0.5\)表示類別1，否則表示類別0.

到這裡就很簡單了，通過線性迴歸我們可以得到\(y = WX + b, P(1)=\sigma(y), P(0)=1-\sigma(y)\)。也就是說我們線上性迴歸模型的外面套了一層sigmoid函式，我們通過計算出不同的y，從而獲得不同的概率，最後得到不同的分類結果。

損失函式

下面的推導全程高能，我相信你們看完會三連的（點贊、轉發、關注）。

讓我們開始吧，我們先來確定一下符號，為了區分，我們把訓練樣本當中的真實分類命名為\(y\)，\(y\)的矩陣寫成\(Y\)。同樣，單條樣本寫成\(x\),\(x\)的矩陣寫成\(X\)。單條預測的結果寫成\(\hat{y}\)，所有的預測結果寫成\(\hat{Y}\)。

對於單條樣本來說，y有兩個取值，可能是1，也可能是0，1和0代表兩個不同的分類。我們希望\(y = 1\)的時候，\(\hat{y}\)儘量大，\(y = 0\)時，\(1 - \hat{y}\)儘量大，也就是\(\hat{y}\)儘量小，因為它取值在0-1之間。我們用一個式子來統一這兩種情況：

\[P(y|x)=\hat{y}^y(1-\hat{y})^{1-y}\]

我們代入一下，\(y=0\)時前項為1，表示式就只剩下後項，同理，\(y=1\)時，後項為1，只剩下前項。所以這個式子就可以表示預測準確的概率，我們希望這個概率儘量大。顯然，\(P(y|x)>0\)，所以我們可以對它求對數，因為log函式是單調的。所以\(P(y|x)\)取最值時的取值，就是\(logP(y|x)\)取最值的取值。

\[\log P(y|x)=y\log \hat{y} + (1-y)\log (1-\hat{y})\]

我們期望這個值最大，也就是期望它的相反數最小，我們令\(J=-\log P(y|x)\)，這樣就得到了它的損失函式：

\[J(\theta)=-\frac{1}{m}Y\log \hat{Y} + (1-Y)\log (1-\hat{Y})\]

如果知道交叉熵這個概念的同學，會發現這個損失函式的表示式其實就是交叉熵。交叉熵是用來衡量兩個概率分佈之間的”距離“，交叉熵越小說明兩個概率分佈越接近，所以經常被用來當做分類模型的損失函式。關於交叉熵的概念我們這裡不多贅述，會在之後文章當中詳細介紹。我們隨手推導的損失函式剛好就是交叉熵，這並不是巧合，其實底層是有一套資訊理論的數學邏輯支撐的，我們不多做延伸，感興趣的同學可以瞭解一下。

硬核推導

損失函式有了，接下來就是求梯度來實現梯度下降了。

這個函式看起來非常複雜，要對它直接求偏導過於硬核（危），如果是許久不碰高數的同學直接肝不亞於硬抗葦名一心。

為了簡化難度，我們先來做一些準備工作。首先，我們先來看下\(\sigma\)函式，它本身的形式很複雜，我們先把它的導數搞定。

\[ \begin{aligned} (\frac{1}{1+e^{-x}})'&=\frac{e^{-x}}{(1+e^{-x})^2} \\ &=\frac{1+e^{-x}-1}{(1+e^{-x})^2} \\ &=\frac{1}{1+e^{-x}}\cdot (1-\frac{1}{1+e^{-x}}) \\ &=\sigma(x) (1-\sigma(x)) \end{aligned} \]

因為\(\hat{y}=\sigma(\theta X)\)，我們將它帶入損失函式，可以得到，其中\(\sigma(\theta)=\sigma(\theta X)\)：

\[J(\theta)=-\frac{1}{m}Y\log \sigma(\theta) + (1 - Y)\log( 1- \sigma(\theta))\]

接著我們求\(J(\theta)\)對\(\theta\)的偏導，這裡要代入上面對\(\sigma(x)\)求導的結論：

\[ \begin{aligned} \frac{\partial}{\partial \theta}J(\theta)&=-\frac{1}{m}(Y\cdot \frac{1}{\sigma(\theta)}\cdot \sigma'(\theta)+(1-Y)\cdot \frac{1}{1-\sigma(\theta)}\cdot \sigma'(\theta)\cdot(-1)) \\ &=-\frac{1}{m}(Y\cdot \frac{1}{\sigma(\theta)}\cdot \sigma(\theta)(1-\sigma(\theta))\cdot X-(1-Y)\cdot \frac{1}{1-\sigma(\theta)}\cdot \sigma(\theta)(1-\sigma(\theta))\cdot X) \\ &=-\frac{1}{m}(Y(1-\sigma(\theta)) - \sigma(\theta)+Y\cdot \sigma(\theta))\cdot X\\ &=-\frac{1}{m}(\sigma(\theta) - Y)\cdot X\\ &=-\frac{1}{m}(\hat{Y} - Y)\cdot X \end{aligned} \]

程式碼實戰

梯度的公式都推出來了，離寫程式碼實現還遠嗎？

不過巧婦難為無米之炊，在我們擼模型之前，我們先試著造一批資料。

我們選擇生活中一個很簡單的場景——考試。假設每個學生需要參加兩門考試，兩門考試的成績相加得到最終成績，我們有一批學生是否合格的資料。希望設計一個邏輯迴歸模型，幫助我們直接計算學生是否合格。

為了防止sigmoid函式產生偏差，我們把每門課的成績縮放到(0, 1)的區間內。兩門課成績相加超過140分就認為總體及格。

import numpy as np
def load_data():
    X1 = np.random.rand(50, 1)/2 + 0.5
    X2 = np.random.rand(50, 1)/2 + 0.5
    y = X1 + X2
    Y = y > 1.4
    Y = Y + 0
    x = np.c_[np.ones(len(Y)), X1, X2]
    return x, Y

這樣得到的訓練資料有兩個特徵，分別是學生兩門課的成績，還有一個偏移量1，用來記錄常數的偏移量。

接著，根據上文當中的公式，我們不難（真的不難）實現sigmoid以及梯度下降的函式。

def sigmoid(x):
    return 1.0/(1+np.exp(-x))
    
def gradAscent(data, classLabels):
    X = np.mat(data)
    Y = np.mat(classLabels)
    m, n = np.shape(dataMat)
    # 學習率
    alpha = 0.01
    iterations = 10000
    wei = np.ones((n, 1))
    for i in range(iterations):
        y_hat = sigmoid(X * wei)
        error = Y - y_hat
        wei = wei + alpha * X.transpose()*error
    return wei

這段函式實現的是批量梯度下降，對Numpy熟悉的同學可以看得出來，這就是在直接套公式。

最後，我們把資料集以及邏輯迴歸的分割線繪製出來。

import matplotlib.pyplot as plt
def plotPic():
    dataMat, labelMat = load_data()
    #獲取梯度
    weis = gradAscent(dataMat, labelMat).getA()
    n = np.shape(dataMat)[0]
    xc1 = []
    yc1 = []
    xc2 = []
    yc2 = []
    #資料按照正負類別劃分
    for i in range(n):
        if int(labelMat[i]) == 1:
            xc1.append(dataMat[i][1])
            yc1.append(dataMat[i][2])
        else:
            xc2.append(dataMat[i][1])
            yc2.append(dataMat[i][2])
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(xc1, yc1, s=30, c='red', marker='s')
    ax.scatter(xc2, yc2, s=30, c='green')
    #繪製分割線
    x = np.arange(0.5, 1, 0.1)
    y = (-weis[0]- weis[1]*x)/weis[2]
    ax.plot(x, y)
    plt.xlabel('X1')
    plt.ylabel('X2')
    plt.show()

最後得到的結果如下：

隨機梯度下降版本

可以發現，經過了1萬次的迭代，我們得到的模型已經可以正確識別所有的樣本了。

我們剛剛實現的是全量梯度下降演算法，我們還可以利用隨機梯度下降來進行優化。優化也非常簡單，我們計算梯度的時候不再是針對全量的資料，而是從資料集中選擇一條進行梯度計算。

基本上可以複用梯度下降的程式碼，只需要對樣本選取的部分加入優化。

def stocGradAscent(data, classLab, iter_cnt):
    X = np.mat(data)
    Y = np.mat(classLab)
    m, n = np.shape(dataMat)
    weis = np.ones((n, 1))
    for j in range(iter_cnt):
        for i in range(m):
            # 使用反比例函式優化學習率
            alpha = 4 / (1.0 + j + i) + 0.01
            randIdx = int(np.random.uniform(0, m))
            y_hat = sigmoid(np.sum(X[randIdx] * weis))
            error = Y[randIdx] - y_hat
            weis = weis + alpha * X[randIdx].transpose()*error
    return weis

我們設定迭代次數為2000，最後得到的分隔影象結果如下：

當然上面的程式碼並不完美，只是一個簡單的demo，還有很多改進和優化的空間。只是作為一個例子，讓大家直觀感受一下：其實自己親手寫模型並不難，公式的推導也很有意思。這也是為什麼我會設定高數專題的原因。CS的很多知識也是想通的，在學習的過程當中靈感迸發旁徵博引真的是非常有樂趣的事情，希望大家也都能找到自己的樂趣。

今天的文章就是這些，如果覺得有所收穫，請順手掃碼點個關注吧，你們的舉手之勞對我來說很重要。

相關推薦