1. 程式人生 > >Logistic迴歸基礎篇之梯度上升演算法

Logistic迴歸基礎篇之梯度上升演算法

作者:崔家華

編輯:趙一帆

一、前言

        本文從Logistic迴歸的原理開始講起,補充了書上省略的數學推導。本文可能會略顯枯燥,理論居多,Sklearn實戰內容會放在下一篇文章。自己慢慢推導完公式,還是蠻開心的一件事。

二、Logistic迴歸與梯度上升演算法

        Logistic迴歸是眾多回歸演算法中的一員。迴歸演算法有很多,比如:線性迴歸、Logistic迴歸、多項式迴歸、逐步迴歸、令迴歸、Lasso迴歸等。我們常用Logistic迴歸模型做預測。通常,Logistic迴歸用於二分類問題,例如預測明天是否會下雨。當然它也可以用於多分類問題,不過為了簡單起見,本文暫先討論二分類問題。首先,讓我們來了解一下,什麼是Logistic迴歸。

1、Logistic迴歸

       假設現在有一些資料點,我們利用一條直線對這些點進行擬合(該線稱為最佳擬合直線),這個擬合過程就稱作為迴歸,如下圖所示:

640?wx_fmt=png&wxfrom=5&wx_lazy=1

        Logistic迴歸是迴歸的一種方法,它利用的是Sigmoid函式閾值在[0,1]這個特性。Logistic迴歸進行分類的主要思想是:根據現有資料對分類邊界線建立迴歸公式,以此進行分類。其實,Logistic本質上是一個基於條件概率的判別模型(Discriminative Model)。

        所以要想了解Logistic迴歸,我們必須先看一看Sigmoid函式 ,我們也可以稱它為Logistic函式。它的公式如下:

0?wx_fmt=png0?wx_fmt=png0?wx_fmt=png

        整合成一個公式,就變成了如下公式:

0?wx_fmt=jpeg

        下面這張圖片,為我們展示了Sigmoid函式的樣子。

0?wx_fmt=png

        z是一個矩陣,θ是引數列向量(要求解的),x是樣本列向量(給定的資料集)。θ^T表示θ的轉置。g(z)函式實現了任意實數到[0,1]的對映,這樣我們的資料集([x0,x1,…,xn]),不管是大於1或者小於0,都可以對映到[0,1]區間進行分類。hθ(x)給出了輸出為1的概率。比如當hθ(x)=0.7,那麼說明有70%的概率輸出為1。輸出為0的概率是輸出為1的補集,也就是30%。

        如果我們有合適的引數列向量θ([θ0,θ1,…θn]^T),以及樣本列向量x([x0,x1,…,xn]),那麼我們對樣本x分類就可以通過上述公式計算出一個概率,如果這個概率大於0.5,我們就可以說樣本是正樣本,否則樣本是負樣本。

        舉個例子,對於”垃圾郵件判別問題”,對於給定的郵件(樣本),我們定義非垃圾郵件為正類,垃圾郵件為負類。我們通過計算出的概率值即可判定郵件是否是垃圾郵件。

        那麼問題來了!如何得到合適的引數向量θ?

        根據sigmoid函式的特性,我們可以做出如下的假設:

0?wx_fmt=jpeg

        上式即為在已知樣本x和引數θ的情況下,樣本x屬性正樣本(y=1)和負樣本(y=0)的條件概率。理想狀態下,根據上述公式,求出各個點的概率均為1,也就是完全分類都正確。但是考慮到實際情況,樣本點的概率越接近於1,其分類效果越好。比如一個樣本屬於正樣本的概率為0.51,那麼我們就可以說明這個樣本屬於正樣本。另一個樣本屬於正樣本的概率為0.99,那麼我們也可以說明這個樣本屬於正樣本。但是顯然,第二個樣本概率更高,更具說服力。我們可以把上述兩個概率公式合二為一:

0?wx_fmt=png

        合併出來的Cost,我們稱之為代價函式(Cost Function)。當y等於1時,(1-y)項(第二項)為0;當y等於0時,y項(第一項)為0。為了簡化問題,我們對整個表示式求對數,(將指數問題對數化是處理數學問題常見的方法):

0?wx_fmt=png

        這個代價函式,是對於一個樣本而言的。給定一個樣本,我們就可以通過這個代價函式求出,樣本所屬類別的概率,而這個概率越大越好,所以也就是求解這個代價函式的最大值。既然概率出來了,那麼最大似然估計也該出場了。假定樣本與樣本之間相互獨立,那麼整個樣本集生成的概率即為所有樣本生成概率的乘積,再將公式對數化,便可得到如下公式:

0?wx_fmt=png

        其中,m為樣本的總數,y(i)表示第i個樣本的類別,x(i)表示第i個樣本,需要注意的是θ是多維向量,x(i)也是多維向量。

        綜上所述,滿足J(θ)的最大的θ值即是我們需要求解的模型。

        怎麼求解使J(θ)最大的θ值呢?因為是求最大值,所以我們需要使用梯度上升演算法。如果面對的問題是求解使J(θ)最小的θ值,那麼我們就需要使用梯度下降演算法。面對我們這個問題,如果使J(θ) := -J(θ),那麼問題就從求極大值轉換成求極小值了,使用的演算法就從梯度上升演算法變成了梯度下降演算法,它們的思想都是相同的,學會其一,就也會了另一個。本文使用梯度上升演算法進行求解。

2、梯度上升演算法

        說了半天,梯度上升演算法又是啥?J(θ)太複雜,我們先看個簡單的求極大值的例子。一個看了就會想到高中生活的函式:

0?wx_fmt=png

        來吧,做高中題。這個函式的極值怎麼求?顯然這個函式開口向下,存在極大值,它的函式影象為:

0?wx_fmt=png

        求極值,先求函式的導數:

0?wx_fmt=png

        令導數為0,可求出x=2即取得函式f(x)的極大值。極大值等於f(2)=4

        但是真實環境中的函式不會像上面這麼簡單,就算求出了函式的導數,也很難精確計算出函式的極值。此時我們就可以用迭代的方法來做。就像爬坡一樣,一點一點逼近極值。這種尋找最佳擬合引數的方法,就是最優化演算法。爬坡這個動作用數學公式表達即為:

0?wx_fmt=png

        其中,α為步長,也就是學習速率,控制更新的幅度。效果如下圖所示:

0?wx_fmt=png

        比如從(0,0)開始,迭代路徑就是1->2->3->4->…->n,直到求出的x為函式極大值的近似值,停止迭代。我們可以編寫Python3程式碼,來實現這一過程

# -*- coding:UTF-8 -*-

""" 

函式說明:梯度上升演算法測試函式

求函式f(x) = -x^2 + 4x的極大值

Parameters:    無

Returns:    無

Author:    Jack Cui

Blog:    http://blog.csdn.net/c406495762

Zhihu:    https://www.zhihu.com/people/Jack--Cui/

Modify:    2017-08-28

"""

def Gradient_Ascent_test():

    def f_prime(x_old):           #f(x)的導數

        return -2 * x_old + 4

    x_old = -1         #初始值,給一個小於x_new的值

    x_new = 0        #梯度上升演算法初始值,即從(0,0)開始

    alpha = 0.01    #步長,也就是學習速率,控制更新的幅度

    presision = 0.00000001           #精度,也就是更新閾值 

    while abs(x_new - x_old) > presision:

        x_old = x_new

        x_new = x_old + alpha * f_prime(x_old)

                                             #上面提到的公式

        print(x_new)                 #列印最終求解的極值近似值

if __name__ == '__main__':

     Gradient_Ascent_test()

        程式碼執行結果如下:

0?wx_fmt=png

        結果很顯然,已經非常接近我們的真實極值了。這一過程,就是梯度上升演算法。那麼同理,J(θ)這個函式的極值,也可以這麼求解。公式可以這麼寫:

0?wx_fmt=png

        由上小節可知J(θ)為:

0?wx_fmt=png

        sigmoid函式為:

0?wx_fmt=jpeg

那麼,現在我只要求出J(θ)的偏導,就可以利用梯度上升演算法,求解J(θ)的極大值了。

那麼現在開始求解J(θ)對θ的偏導,求解如下(數學推導):

0?wx_fmt=jpeg0?wx_fmt=jpeg

        知道了,梯度上升迭代公式,我們就可以自己編寫程式碼,計算最佳擬合引數了。

三、Python3實戰

1、資料準備

        資料集已經為大家準備好,下載地址:

https://github.com/Jack-Cherish/Machine-Learning/blob/master/Logistic/testSet.txt

        這就是一個簡單的資料集,沒什麼實際意義。讓我們先從這個簡單的資料集開始學習。先看下資料集有哪些資料:

0?wx_fmt=png

         這個資料有兩維特徵,因此可以將資料在一個二維平面上展示出來。我們可以將第一列資料(X1)看作x軸上的值,第二列資料(X2)看作y軸上的值。而最後一列資料即為分類標籤。根據標籤的不同,對這些點進行分類。

        那麼,先讓我們編寫程式碼,看下資料集的分佈情況:

# -*- coding:UTF-8 -*-

import matplotlib.pyplot as plt

import numpy as np

"""

函式說明:載入資料

Parameters:    無

Returns:    dataMat - 資料列表

                  labelMat - 標籤列表

Author:    Jack Cui 

Blog:    http://blog.csdn.net/c406495762

Zhihu:    https://www.zhihu.com/people/Jack--Cui/ 

Modify:    2017-08-28

"""

def loadDataSet():

    dataMat = []                                 #建立資料列表

    labelMat = []                                 #建立標籤列表

    fr = open('testSet.txt')                  #開啟檔案

    for line in fr.readlines():               #逐行讀取

    lineArr = line.strip().split()            #去回車,放入列表

    dataMat.append([1.0, float(lineArr[0]),

                                 float(lineArr[1])])       #新增資料

    labelMat.append(int(lineArr[2]))            #新增標籤

    fr.close()                                                #關閉檔案

    return dataMat, labelMat         

"""

函式說明:繪製資料集

Parameters:    無

Returns:    無

Author:    Jack Cui

Blog:    http://blog.csdn.net/c406495762

Zhihu:    https://www.zhihu.com/people/Jack--Cui/

Modify:    2017-08-28

"""

def plotDataSet():

    dataMat, labelMat = loadDataSet()       #載入資料集

    dataArr = np.array(dataMat)                                                                       #轉換成numpy的array陣列

    n = np.shape(dataMat)[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])    #1為正樣本

        else:

            xcord2.append(dataArr[i,1])

            ycord2.append(dataArr[i,2])    #0為負樣本

    fig = plt.figure()

    ax = fig.add_subplot(111)                #新增subplot

    ax.scatter(xcord1, ycord1, s = 20, c = 'red'

                        marker = 's',alpha=.5)     #繪製正樣本

    ax.scatter(xcord2, ycord2, s = 20,

                     c = 'green',alpha=.5)          #繪製負樣本

    plt.title('DataSet')                                #繪製title

    plt.xlabel('x'); plt.ylabel('y')                  #繪製label

    plt.show()                                            #顯示

if __name__ == '__main__'

    plotDataSet()

        執行結果如下:

0?wx_fmt=png

        從上圖可以看出資料的分佈情況。假設Sigmoid函式的輸入記為z,那麼z=w0x0 + w1x1 + w2x2,即可將資料分割開。其中,x0為全是1的向量,x1為資料集的第一列資料,x2為資料集的第二列資料。另z=0,則0=w0 + w1x1 + w2x2。橫座標為x1,縱座標為x2。這個方程未知的引數為w0,w1,w2,也就是我們需要求的迴歸係數(最優引數)。

2、訓練演算法

        在編寫程式碼之前,讓我們回顧下梯度上升迭代公式:

0?wx_fmt=png

        將上述公式向量化:

0?wx_fmt=png

        根據向量化的公式,編寫程式碼如下:

0?wx_fmt=png

        執行結果如圖所示:

0?wx_fmt=png

         可以看出,我們已經求解出迴歸係數[w0,w1,w2]。

         通過求解出的引數,我們就可以確定不同類別資料之間的分隔線,畫出決策邊界。

3、繪製決策邊界

        我們已經解出了一組迴歸係數,它確定了不同類別資料之間的分隔線。現在開始繪製這個分隔線,編寫程式碼如下:

0?wx_fmt=png0?wx_fmt=png0?wx_fmt=png0?wx_fmt=png0?wx_fmt=png

        執行結果如下:

0?wx_fmt=png

        這個分類結果相當不錯,從上圖可以看出,只分錯了幾個點而已。但是,儘管例子簡單切資料集很小,但是這個方法卻需要大量的計算(300次乘法)。因此下篇文章將對改演算法稍作改進,從而減少計算量,使其可以應用於大資料集上。

四、總結

Logistic迴歸的一般過程:

     收集資料:採用任意方法收集資料。

    準備資料:由於需要進行距離計算,因此要求資料型別為數值型。另外,結構化資料格式則最佳。

     分析資料:採用任意方法對資料進行分析。

    訓練演算法:大部分時間將用於訓練,訓練的目的是為了找到最佳的分類迴歸係數。

     測試演算法:一旦訓練步驟完成,分類將會很快。

    使用演算法:首先,我們需要輸入一些資料,並將其轉換成對應的結構化數值;接著,基於訓練好的迴歸係數,就可以對這些數值進行簡單的迴歸計算,判定它們屬於哪個類別;在這之後,我們就可以在輸出的類別上做一些其他分析工作。

其他:

    Logistic迴歸的目的是尋找一個非線性函式Sigmoid的最佳擬合引數,求解過程可以由最優化演算法完成。

    本文講述了Logistic迴歸原理以及數學推導過程。

下篇文章將講解Logistic迴歸的改進以及Sklearn實戰內容。

    如有問題,請留言。如有錯誤,還望指正,謝謝!

        PS: 如果覺得本篇本章對您有所幫助,歡迎關注、評論、頂!

         本文出現的所有程式碼和資料集,均可在我的github上下載,歡迎Follow、Star:

https://github.com/Jack-Cherish/Machine-Learning

參考文獻:

  1. 斯坦福大學的吳恩達《機器學習》:https://www.coursera.org/learn/machine-learning

  2. 《機器學習實戰》第五章內容

燕哥開課了,趕緊報名!

還有機會獲取簽名機器學習圖書贈送!!!

0?wx_fmt=jpeg

掃描燕哥微訊號,

拉你進機器學習大牛群。

福利滿滿,名額已不多…

640.jpeg?

群裡目前包括:

清華張長水教授,

清華顧險峰教授,

北大黃鐵軍教授,

西安電子科技大學焦李成教授,

新加坡南洋理工大學黃廣斌教授,

北交李清勇教授

等等……

0?wx_fmt=gif