大叔學ML第五:邏輯回歸


基本形式

邏輯回歸是最常用的分類模型,在線性回歸基礎之上擴展而來,是一種廣義線性回歸。下面舉例說明什么是邏輯回歸:假設我們有樣本如下(是我編程生成的數據):

image.png-9.7kB

我們要做的是找到一個決策邊界,把兩類樣本給分開,當有新數據進來時,就判斷它在決策邊界的哪一邊。設邊界線為線性函數

\[h_\theta(\vec x) = \theta_0 + \theta_1x_1 + \theta_2x_2 \tag {1}$$取0時的直線,如下圖: ![image.png-12.4kB][2] **我們的目標就是根據已知的樣本來確定$\vec\theta$取值**。 上圖中,處在邊界線左邊的為負例(帶入(1)式結果小於0),邊界線右邊的為正例(帶入(1)式結果大於0)。**從概率的角度考慮**:把一個點代入(1)式,如果為正,且越大越離邊界線遠,它是正例的概率就越大,直到接近1;相反,把一個點代入(1)式,如果為負,且越小離邊界線越遠,它是正例的概率就越小,直到接近0;此外,把一個點代入(1)式,如果結果正好等於0,那么它在邊界線上,為正例和為負例的概率都是0.5。 為了用數學的方式精確表達上面的概率論述,前人找到一個好用的函數: $$s(z) = \frac{1}{1+e^{-z}} \tag{2}\]

這個函數叫做sigmoid(sigmoid:S形狀的)函數(下文用S(z)或s(z)表示sigmoid函數),樣子如下:

image.png-11kB

\(z=0\)時,函數值為0.5,當\(z > 0\)時,函數取值\((0.5, 1)\);當\(z < 0\)時,函數取值\((0, 0.5)\),如果我們把欲判斷點代入決策邊界(1)式后得到的結果作為sigmoid函數的輸入,那么輸出就可以表示該點是正例的概率,簡直完美。其實S型的函數應該還有別的,為什么前人獨愛這個呢?那是因為,這個函求導比較簡單,用鏈式法則可以非常容易算出其導數式為:

\[\frac{d}{z}s(z)=s(z)(1 - s(z)) \tag{3} \]

simgiod函數求導過程:

\[\begin{align}\frac{d}{z}s(z)&=-\frac{1}{(1+e^{-z})^2}\cdot e^{-z} \cdot -1 \\ &=\frac{e^{-z}}{(1+e^{-z})^2}\\ &=\frac{1+e^{-z}-1}{(1+e^{-z})^2}\\ &=\frac{1+e^{-z}}{(1+e^{-z})^2} - \frac{1}{(1+e^{-z})^2}\\ &=s(z) - s(z)^2\\ &=s(z)(1-s(z)) \end{align}\]

第一步用了鏈式法則。

代價函數

邏輯回歸的代價函數可由極大似然估計法得出。我們暫且不管極大似然估計法的證明,直觀地理解非常容易:如果你是一個班級的新老師,發現有個孩子考了95.5分,你肯定會認為這個孩子很可能是學霸,雖然學霸有時也會考低分,學渣有時也考高分,但是發生概率很小。更一般的敘述是:有若干事件A、B和C,已知其發生概率為分別為\(P(A|\theta)\)\(P(B|\theta)\)\(P(C|\theta)\),如果我們觀察到A、B和C已經發生了,那么我們就認為\(P(ABC|\theta)\)是個盡可能大的值,如果我們要求\(\theta\),那么\(\theta\)應該是使得\(P(ABC|\theta)\)最大的那個值。如果A、B和C互相獨立,我們所求的就是使得\(P(A|\theta)P(B|\theta)P(C|\theta)\)最大化的\(\theta\)

\(t_\theta(\vec x)=s(h_\theta(\vec x))\)已知:$$P(Y=y) = \begin{cases}
t_\theta(\vec x), &\text{如果 y = 1} \
1 - t_\theta(\vec x), &\text{如果 y = 0}
\end{cases}$$,寫到一起:

\[P(Y=y)=t_\theta(\vec x)^y(1-t_\theta(\vec x))^{1-y} \tag{4}$$,根據(4)式寫出極大似然函數: $$l(\vec\theta)=\prod_{i=1}^mt_\theta(\vec x^{(i)})^{y^{(i)}}(1-t_\theta(\vec x^{(i)}))^{1-y^{(i)}} $$,設代價函數$$j(\vec\theta)=-\frac{1}{m}\ln l(\vec\theta)$$ ,最大化$l(\vec\theta)$即最小化$j(\vec\theta)$。求對數是為了方便計算,把乘法轉換為加法,把冪運算轉換為乘法,單調性不變;前面的負號是為了把求最大值問題轉換為求最小值問題,習慣上,應用梯度下降法時都是求最小值,不然叫“梯度上升了法”了,手動滑稽。當然,使用梯度下降法的前提條件是函數是凸的,大叔懶得證明了,這個函數就是個凸函數,不管你們信不信,我反正是信了。進一步對上式化解得: $$j(\vec\theta)=-\frac{1}{m}\sum_{i=1}^m\left[y^{(i)}\ln t_\theta(\vec x^{(i)}) + (1-y^{(i)})\ln(1-t_\theta(\vec x^{(i)}))\right] \tag{5}\]

用梯度下降法求\(\vec\theta\)

將(5)式對\(\vec\theta\)求偏導:

  • \(\frac{\partial}{\partial\theta_0}j(\vec\theta)=\frac{1}{m}\sum_{i=1}^m(t_\theta(\vec x^{(i)})-y^{(i)})x_0^{(i)}\)
  • \(\frac{\partial}{\partial\theta_1}j(\vec\theta)=\frac{1}{m}\sum_{i=1}^m(t_\theta(\vec x^{(i)})-y^{(i)})x_1^{(i)}\)
  • \(\cdots\)
  • \(\frac{\partial}{\partial\theta_n}j(\vec\theta)=\frac{1}{m}\sum_{i=1}^m(t_\theta(\vec x^{(i)})-y^{(i)})x_n^{(i)}\)

如果讀者對矩陣計算非常熟悉的話,應該可以看出,上式可以寫成矩陣形式,這樣計算更方便:

\[\frac{d}{d\vec\theta}=\frac{1}{m}X^T(T-Y),其中,T=\begin{pmatrix}t_\theta(\vec x^{(1)})\\ t_\theta(\vec x^{(2)})\\ \vdots \\ t_\theta(\vec x^{(m)}) \end{pmatrix}\]

\(j(\vec\theta)\)求導的過程如下,多次應用鏈式法則即可:

\[\begin{align} \frac{\partial}{\partial\theta_n}j(\vec\theta)&=-\frac{1}{m}\sum_{i=1}^m\left[y^{(i)}\frac{1}{t_\theta(\vec x^{(i)})}t_\theta(\vec x)^\prime+(1-y^{(i)})\frac{1}{1-t_\theta(\vec x^{(i)})}t_\theta(\vec x)^\prime\right]\\ &=-\frac{1}{m}\sum_{i=1}^m\left[y^{(i)}\frac{1}{t_\theta(\vec x^{(i)})}t_\theta(\vec x^{(i)})(1-t_\theta(\vec x^{(i)}))h_\theta(\vec x)^\prime+(1-y^{(i)})\frac{1}{1-t_\theta(\vec x^{(i)})}t_\theta(\vec x^{(i)})(t_\theta(\vec x^{(i)}) - 1)h_\theta(\vec x)^\prime\right]\\ &=-\frac{1}{m}\sum_{i=1}^m\left[y^{(i)}(1-t_\theta(\vec x^{(i)}))h_\theta(\vec x)^\prime + (y^{(i)})t_\theta(\vec x^{(i)} - 1)h_\theta(\vec x)^\prime\right]\\ &=-\frac{1}{m}\sum_{i=1}^m\left[y^{(i)}(1-t_\theta(\vec x^{(i)}))x_n^{(i)} + (y^{(i)})t_\theta(\vec x^{(i)} - 1)x_n^{(i)}\right]\\ &=-\frac{1}{m}\sum_{i=1}^m(y^{(i)} - t_\theta(\vec x^{(i)}))x_n^{(i)}\\ &=\frac{1}{m}\sum_{i=1}^m(t_\theta(\vec x^{(i)}) - y^{(i)})x_n^{(i)} \end{align}\]

  • 步驟1用到到了公式:\(\ln x^\prime = \frac{1}{x}\)
  • 步驟2用到式(3)

有了偏導式后就可以編程了:

import matplotlib.pyplot as plt
import numpy as np

def sigmoid(z):
    return 1 / (1 + np.e**(-z))

def draw_samples(X, Y, sample_count):
    ''' 繪制正負例. '''
    positiveX1 = []
    positiveX2 = []

    negativeX1 = []
    negativeX2 = []

    for i in range(sample_count):
        if Y[i, 0] == 1:
            positiveX1.append(X[i, 0])
            positiveX2.append(X[i, 1])
        else:
            negativeX1.append(X[i, 0])
            negativeX2.append(X[i, 1])
    
    plt.scatter(positiveX1, positiveX2, marker='+')
    plt.scatter(negativeX1, negativeX2, marker='.')

def draw_border(theta):
    '''繪制邊界線'''
    X = []
    Y = []
    for x in range(-2, 12):
        X.append(x)
        Y.append(-theta[0] / theta[2] - theta[1] / theta[2] * x )
    plt.plot(X, Y)

def create_samples(samples_count):
    ''' 生成樣本數據'''
    X = np.empty((samples_count, 2))
    Y = np.empty((samples_count, 1))

    for i in range(samples_count):
        x1 = np.random.randint(0, 10)
        x2 = np.random.randint(0, 10)
        y = 0
        if x1 + x2 - 10 > 0:
            y =1
        X[i, 0] = x1
        X[i, 1] = x2
        Y[i, 0] = y

        noise = np.random.normal(0, 0.1, (samples_count, 2))
        X += noise
    return X, Y

def grad(X, Y, samples_count, theta):
    ''' 求代價函數在theta處的梯度. '''
    T = sigmoid(np.dot(X, theta))
    g =  np.dot(X.T, (T - Y)) / samples_count
    return g

def descend(X, Y, samples_count, theta = np.array([[1],[1],[1]]), step = 0.01, threshold = 0.05):
    ''' 梯度下降.
        Args: 
            X: 樣本
            Y:樣本標記
            theta:初始值
            step:步長
            threshold:閾值
        Returns:
            theta:求出來的最優theta
    '''
    g = grad(X, Y, samples_count, theta)
    norm = np.linalg.norm(g)

    while norm > threshold:
        g = grad(X, Y, samples_count, theta)
        norm = np.linalg.norm(g)
        theta = theta - g * step
        print(norm)

    return theta

samples_count = 100
X, Y = create_samples(samples_count)
MatrixOnes = np.ones((100, 1))
X_with_noes = np.hstack((MatrixOnes, X)) # 添加等於1的x0
theta = descend(X_with_noes, Y, samples_count)

plt.xlabel('$x_1$')
plt.ylabel('$x_2$')

draw_samples(X, Y, samples_count)
draw_border(theta.flatten())
plt.show()

運行結果:
image.png-14.9kB

上面的代碼中,我在樣本矩陣中增加了一列為1的元素,這樣做是為了計算方便,使得樣本矩陣的列數等於欲求\(\vec\theta\)的行數,滿足矩陣乘法的要求,加上這一列對結果沒有影響,請參考大叔學ML第二:線性回歸

擴展

  • 和線性回歸一樣,邏輯回歸也存在過擬合的情況,可以通過添加一范數、二范數正則化來解決
  • 和線性回歸一樣,邏輯回歸的決策邊界可以是曲折的多項式,可參考大叔學ML第三:多項式回歸,依葫蘆畫瓢來搞定決曲折的決策邊界
  • 和線性回歸一樣,邏輯回歸也可以擴展到多余3維的情況,只是不可以做可視化了,代數原理是一致的
  • 和線性回歸不一樣,邏輯回歸應該沒有“正規方程”,大叔昨天按照線性回歸的正規方程的推導思路推導並未得出結果,查資料也未果
  • 除了梯度下降法,牛頓法和擬牛頓法也是迭代下降算法,而且下降效率更好,大叔稍后將會寫一篇博文介紹之
  • 多分類,多分類非常簡單,假如樣本數據共有3類,記作A、B和C,先把A分為正例,B和C分為負例,會得到一條決策邊界,然后再把B分為正例,C分為正例做同樣的操作,最后得到3條決策邊界,也就是得到3個不同的概率計算函數,如果有新數據進來,分別帶入這三個函數,取概率值最大的那個函數所代表的分類即可。

祝元旦快樂!


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM