Python機器學習筆記:Logistic Regression


Logistic回歸公式推導和代碼實現

1,引言

  logistic回歸是機器學習中最常用最經典的分類方法之一,有人稱之為邏輯回歸或者邏輯斯蒂回歸。雖然他稱為回歸模型,但是卻處理的是分類問題,這主要是因為它的本質是一個線性模型加上一個映射函數Sigmoid,將線性模型得到的連續結果映射到離散型上。它常用於二分類問題,在多分類問題的推廣叫softmax。

  本文首先闡述Logistic回歸的定義,然后介紹一些最優化算法,其中包括基本的梯度上升法和一個改進的隨機梯度上升法,這些最優化算法將用於分類器的訓練,最好本文將給出一個Logistic回歸的實例,預測一匹病馬是否能被治愈。

  在我們的日常生活中遇到過很多最優化問題,比如如何在最短時間內從A點到達B點?如何投入最少工作量卻獲得最大的效益?如何設計發動機使得油耗最少而功率最大?可見,最優化的作用十分強大,所以此處我們介紹幾個最優化算法,並利用它們訓練出一個非線性函數用於分類。

  現在假設有一些數據點,我們用一條直線對這些點進行擬合(該線稱為最佳擬合直線),這個擬合過程就稱作回歸。利用logistic回歸進行分類的主要思想是:根據現有數據對分類邊界線建立回歸公式,以此進行分類,這里的“回歸”一詞源於最佳擬合,表示要找到最佳擬合參數集。訓練分類器時的做法就是尋找最佳擬合參數,使用的是最優化算法,下面我們首先介紹一下這個二值型輸出分類器的數學原理。

2,Logistic回歸的一般過程

  • (1)收集數據:采用任意方法收集數據
  • (2)准備數據:由於需要進行距離計算,因此要求數據類型為數值型。另外,結構化數據格式則最佳
  • (3)分析數據:采用任意方法對數據進行分析
  • (4)訓練算法:大部分時間將用於訓練,訓練的目的是為了找到最佳的分類回歸系數
  • (5)使用算法:首先,我們需要輸入一些數據,並將其轉換成對應的結構化數值;接着,基於訓練好的回歸系數就可以對這些數值進行簡單的回歸計算,判定他們屬於哪個類別;在這之后,我們就可以在輸出的類別上做一些其他分析工作。

3,Logistic回歸的優缺點

  優點:計算代碼不多,易於理解和實現,計算代價不高,速度快,存儲資源低

  缺點:容易欠擬合,分類精度可能不高

  適用數據類型:數值型和標稱型數據

4,基於Logistic回歸和Sigmoid函數的分類

   我們想要的函數應該是:能接受所有的輸入,然后預測出類型。例如,在兩個類的情況下,上述函數輸出0或1。該函數稱為海維賽德階躍函數(Heaviside step function),或者直接稱為單位階躍函數。然而,海維賽德階躍函數的問題在於:該函數在跳躍點上從0瞬間跳躍到1,這個瞬間跳躍過程有時很難處理。幸好,另一個函數也有類似的性質(可以輸出0或者1),且數學上更易處理,這就是Sigmoid函數。Sigmoid函數具體的計算公式如下:

  自變量取值為任意實數,值域[0, 1]

  圖5-1給出了Sigmoid函數在不同坐標尺度下的兩條曲線圖。當x為0時,Sigmoid函數值為0.5。隨着x的增大,對應的Sigmoid值將逼近於1;而隨着x的減少,Sigmoid值將逼近於0.如果橫坐標刻度足夠大,Sigmoid函數看起來很像一個階躍函數。

  解釋Sigmoid函數:將任意的輸入映射到了 [0, 1]區間,我們在線性回歸中可以得到一個預測值,再將該值映射到 Sigmoid函數中這樣就完成了由值到概率的轉換,也就是分類任務

  因此,為了實現Logistic回歸分類器,我們可以在每個特征上都乘以一個回歸系數,然后把所有的結果值相加,將這個總和帶入Sigmoid函數中,進而得到一個范圍在0~1之間的數值。任何大於0.5的數據被分入1類,小於0.5即被歸入0類,所以,Logistic回歸也可以被看成是一種概率估計。

  確定了分類器的函數形式之后,現在的問題變成了:最佳回歸系數是多少?如何確定其大小

5,基於最優化方法的最佳回歸系數確定

  Sigmoid函數的輸入記為z,由下面公式得到:

  如果采用向量的寫法,上述公式可以寫成  z = wTx  ,它表示將這兩個數值向量對應元素相乘,然后全部加起來即得到z值。

  其中的向量x是分類器的輸入數據,向量w也就是我們要找到的最佳參數(系數),從而使得分類器盡可能的准確,為了尋找該最佳參數,需要用到最優化理論的一些知識。

  然后再看看我們的Logistic回歸模型的公式:

  這里假設 W>0,Y與X各維度疊加的圖形關系,如下圖所示(x為了方便取1維):

  下面首先學習梯度上升的最優化方法,我們將學習到如何使用該方法求得數據集的最佳參數,接下來,展示如何繪制梯度上升法產生的決策邊界圖,該圖將梯度上升法的分類效果可視化的呈現出來,最后我們將學習隨機梯度上升算法,以及如何對其進行修改以獲得很好地結果。

    可能我們最常聽到的是梯度下降算法,它與這里的梯度上升算法是一樣的,只是公式中的
加法需要變成減法,梯度上升算法用來求函數的最大值,而梯度下降算法是用來求函數的最小值

6,梯度上升法 

  梯度上升法基於的思想是:要找到某函數的最大值,最好的方法是沿着該函數的梯度方向探尋,如果梯度記為,則函數 f(x,y) 的梯度由下面式子表示:

  這個梯度意味着要沿着x的方向移動,沿着y方向移動,其中函數f(x,y)必須要在待計算的點上有定義並且可微,一個具體的函數例子見圖5-2:

   上圖中的梯度上升算法沿梯度方向移動了一步,可以看出,梯度算子總是指向函數值增長最快的方向。這里所說的移動方向,而未提到移動量的大小。該量值稱為步長,記為。用向量來表示的話,梯度算法的迭代公式如下:

  該公式將一直被迭代執行,直至達到某個停止條件為止,比如迭代次數達到某個指定值或算法達到某個可以允許的誤差范圍。

  基於上面的內容,我們來看一個Logistic回歸分類器的應用例子,從圖5-3可以看到我們采用的數據集。

梯度上升法的公式推導(LR 損失函數)

  在LR中,應用極大似然估計法估計模型參數,由於Sigmoid函數的特性,我們可以做如下的假設:

  上式即為在已知樣本X和參數θ的情況下。樣本X屬性正類(y=1)和負類(y=0)的條件概率,將兩個公式合並成一個,如下:

  假定樣本與樣本之間相互獨立,那么整個樣本集生成的概率即為所有樣本生成概率的乘積(也就是n個獨立樣本出現的似然函數如下):

  為了簡化問題,我們對整個表達式求對數(即為LR 損失函數)

  滿足似然函數(θ)的最大的θ值即時我們需要求解的模型。

  那么梯度上升法就像爬坡一樣,一點一點逼近極值,而上升這個動作用數學公式表達即為:

  其中,α 為步長。

  回到Logistic回歸問題,我們同樣對函數求偏導。

對這個公式進行分解,先看:

  我們可以看到,對函數求偏導,分解為三部分,然后我們對這三部分分布求導。

其中:

再由:

可得:

接下來:

最后:

綜合三部分即得到:

  如果上面鏈式分解不好理解的話,可以看下面直接求導(結果是一樣的):

   注意上面是將梯度上升求最大值,轉換為梯度下降了,本質沒變。

因此梯度迭代公式為:

 

 如果為梯度下降,我們注意符號的變化,如下:

7,訓練算法:使用梯度上升找到最佳參數

   上圖有100個樣本點,每個點包含兩個數值型特征:X1和X2,在此數據集上,我們將通過使用梯度上升法找到最佳回歸系數,也就是擬合出Logistic回歸模型的最佳參數。

  所以我們的目標:建立分類器,求解出theta參數

  設定閾值,根據閾值判斷結果

  梯度上升法的偽代碼如下:

每個回歸系數初始化為1

重復R次:
    
    計算整個數據集的梯度

    使用alpha * gradient 更新回歸系數的向量

    返回回歸系數

  testSet.txt的文件內容請去我的GitHub下載。地址:https://github.com/LeBron-Jian/MachineLearningNote

  下面具體實現梯度上升算法的代碼:

#_*_coding:utf-8_*_
from numpy import *

# 讀取數據
def loadDataSet(filename):
    '''
        對於testSet.txt,每行前兩個值分別是X1和X2,第三個值數據對應的類別標簽
        而且為了設置方便,該函數還將X0的值設置為1.0
        :return:
        '''
    dataMat = []
    labelMat = []
    fr = open(filename)
    for line in fr.readlines():
        lineArr = line.strip().split()
        dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])
        labelMat.append(int(lineArr[2]))
    return dataMat,labelMat

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

def gradAscent(dataMatIn,classLabels):
    '''
        :param dataMatIn: 是一個2維Numpy數組,每列分別代表每個不同的特征
        每行則代表每個訓練樣本。
        :param classLabels: 是類別標簽,是一個1*100的行向量,為了便於矩陣運算,需要將行向量
        轉換為列向量,就是矩陣的轉置,再將其賦值與labelMat。
        :return:
        '''
    dataMatrix = mat(dataMatIn)
    labelMat = mat(classLabels).transpose()
    # labelMat = mat(classLabels).T
    m,n = shape(dataMatrix)
    # alpha是向目標移動的步長
    alpha = 0.001
    # 迭代次數
    maxCycles = 500
    weights = ones((n,1))
    for k in range(maxCycles):
        h = sigmoid(dataMatrix*weights)
        error = (labelMat-h)
        weights = weights + alpha*dataMatrix.transpose()*error
    return weights

  注意上面代碼(為什么將x0的值設置為1呢,主要是我們要使用矩陣進行運算,可以看下圖比較明顯)

  測試結果如下:

if __name__  == '__main__':
    filename = 'testSet.txt'
    dataArr,labelMat = loadDataSet(filename)
    weights_res = gradAscent(dataArr,labelMat)
    print(weights_res)
    
'''
[[ 4.12414349]
 [ 0.48007329]
 [-0.6168482 ]]
 '''

  上面已經解出了一組回歸系數,它確定了不同類別數據之間的分割線,那么怎樣畫出該分割線,從而使得優化的過程便於理解呢?下面代碼來解決這個問題。

  畫出數據集和Logistic回歸最佳擬合直線的函數代碼:

def plotBestFit(wei):
    import matplotlib.pyplot as plt
    weights = wei.getA()
    dataMat,labelMat = loadDataSet(filename)
    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()

  輸出的結果和代碼如下圖所示:

if __name__  == '__main__':
    filename = 'testSet.txt'
    dataArr,labelMat = loadDataSet(filename)
    weights_res = gradAscent(dataArr,labelMat)
    print(weights_res)
    plotBestFit(weights_res)

梯度上升算法在500次迭代后得到的Logistic回歸最佳擬合直線

  這個分類結果相當不錯,從圖上看,只錯分了四個點。但是,盡管例子簡單且數據集很小,這個方法卻需要大量的計算(300次乘法),因此下一節將對算法稍作改進,從而使它可以用在真實數據集上。

8,訓練算法:隨機梯度上升

   梯度上升算法在每次更新回歸系數時都需遍歷整個數據集,該方法在處理100個左右的數據集尚可,但是若有數十億樣本和成千上萬的特征,那么該方法的計算復雜度就太高了。一種改進方法是一次僅用一個樣本點來更新回歸系數,該方法稱為隨機梯度上升算法。由於可以在新樣本到來時對分類器進行增量式更新,因而隨機梯度上升算法是一個在線學習算法,與“在線學習”相對應,一次處理所有數據被稱作是“批處理”。

  隨機梯度上升算法可以寫成如下的偽代碼:

所有回歸系數初始化為1

對數據集中每個樣本
    
    計算該樣本的梯度

    使用alpha*gradient 更新回歸系數值

返回回歸系數值

  以下是隨機梯度上升算法的實現代碼:

# 隨機梯度上升算法
def stocGradAscent0(dataMatrix,classLabels):
    m,n = shape(dataMatrix)
    alpha = 0.01
    weights = ones(n)
    for i in range(m):
        h = sigmoid(sum(dataMatrix[i]*weights))
        error = classLabels[i] - h
        weights = weights + alpha*error*dataMatrix[i]
    return weights

  實現的代碼如下:

if __name__  == '__main__':
    filename = 'testSet.txt'
    dataArr,labelMat = loadDataSet(filename)
    weights_res = stocGradAscent0(array(dataArr),labelMat)
    print(weights_res)
    plotBestFit(weights_res)

  

圖5-5  隨機梯度上升算法在上述數據集上的執行結果,最佳擬合直線並非最佳分類線

  可以看出,擬合出來的直線效果還不錯,但並不像,但是不像上面那個完美,這里的分類器錯分了三分之一的樣本。

  直接比較結果兩個結果是不公平的,后者的結果是在整個數據集上迭代了500次才得到的。一個判斷優化算法優劣的可靠方法是看它是否收斂,也就是說參數是否達到了穩定值,是否還會不斷地變化?對此,我們在上面的的隨機梯度算法上做了些修改,使其在整個數據集上運行200次,最終繪制的三個回歸系數的變化情況如下圖所示:

  上圖展示了隨機梯度上升算法在200次迭代過程中回歸系數的變化情況,其中的系數2,也就是圖5-5中的X2只經過了50次迭代就達到了穩定值,但系數1和0則需要更多次的迭代。另外值得注意的是,在大的波動停止后,還有一些小的周期性波動。不難理解,產生這種現象的原因是存在一些不能正確分類的樣本點(數據集並非現象可分),在每次迭代時會引發系數的劇烈改變。我們期望算法能避免來回波動,從而收斂到某個值。另外,收斂速度也需要加快。

  改進的隨機梯度上升算法代碼如下:

# 改進的隨機梯度上升算法
def stocGradAscent1(dataMatrix,classLabels,numIter=150):
    m,n = shape(dataMatrix)
    weights = ones(n)
    for j in range(numIter):
        dataIndex = list(range(m))
        for i in range(m):
            alpha = 4/(1.0+j+i)+0.01
            randIndex = int(random.uniform(0,len(dataIndex)))
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
            error = classLabels[randIndex] - h
            weights = weights + alpha *error*dataMatrix[randIndex]
            del(dataIndex[randIndex])
    return weights

  上述代碼大體上與之前的隨機梯度上升算法一致,修改了兩處,一處是alpha在每次迭代的時候都會調整,這會環節之前的數據波動或者高頻波動。另外,雖然alpha會隨着迭代次數不斷減少,但永遠不會減少到0,。必須這樣做的原因是為了保證在多次迭代之后新數據仍然具有一定的影響。如果要處理的問題是動態變化的,那么可以適當增加常數項,來確保新的值獲得更大的回歸系數。另外一點值得注意的是,在降低alpha的函數中,alpha每次減少1/(j+1),其中j是迭代次數,i是樣本點的下標,這樣當j<max(i)的時候,alpha就不是嚴格下降的,避免參數的嚴格下降也常見於模擬退火算法等其他優化算法中。程序的第二個改進地方就是通過隨機選取樣本來更新回歸系數,這種方法將減少周期性的波動,並且改進的算法還增加一個迭代次數作為第三個參數,如果該參數沒有給定的話,算法將默認迭代150次,如果給定,那么算法將按照新的參數值進行迭代。

  與stocGradAscent1()類似,下圖顯示了每次迭代時各個回歸系數的變化情況。

  比較5-7和5-6可以看到兩點不同,第一點是,圖5-7中的系數沒有像5-6里那樣出現周期性的波動,這歸功於stocGradAscent1()里的樣本隨機選擇機制,第二點是5-7的水平軸比5-6的短了很多,這是由於stocGradAscent1()可以收斂的更快,這次我們僅僅對數據集做了20次遍歷,之前是500次。

  下面看看在同一個數據集上的分類效果,將程序運行可以看到:

if __name__  == '__main__':
    filename = 'testSet.txt'
    dataArr,labelMat = loadDataSet(filename)
    weights_res = stocGradAscent1(array(dataArr),labelMat)
    print(weights_res)
    plotBestFit(weights_res)

  該分割線達到了與GradientAscent()差不多的效果,但是所使用的計算量更少。

  默認的迭代次數是150次,但是我們通過stocGradAscent()的第三個參數來對此進行修改,例如:

weights_res = stocGradAscent1(array(dataArr),labelMat,500)

  迄今為止我們分析了回歸系數的變化情況,但是還沒有達到目的,,即完成具體的分類任務,下面我們將使用隨機梯度上升算法來解決病馬的生死預測問題。,

9,示例:從疝氣病症預測病馬的死亡率

  本次將使用logistic回歸來預測患有疝氣的馬的存活問題,這里的數據包含了368個樣本和28個特征。(疝氣指的是馬胃腸痛的術語,然而這種病不一定源於馬 的腸胃問題,其他問題也可以引發疝氣)該數據集包含了醫院檢測馬疝氣病的一些指標,有的指標比較主觀,有的指標難以預測,例如馬的疼痛級別

使用Logistic回歸估計馬疝氣病的死亡率的流程

(1)收集數據:使用給定數據文件

(2)准備數據:用Python解析文本文件並填充缺失值

(3)分析數據:可視化並觀察數據

(4)訓練算法:使用優化算法,找到最佳的系數

(5)測試算法:為了量化回歸的結果,需要觀察錯誤率,根據錯誤率決定是否回退到訓練階段,通過改變迭代的次數和步長等參數來得到更好的回歸系數

(6)使用算法:實現一個簡單的命令行程序來收集馬的症狀並輸出預測結果

准備數據:處理數據中的缺失值

  數據中的缺失值是個非常棘手的問題,有很少文獻都致力於解決這個問題。那么,數據缺少究竟帶來了什么問題?假設有100個樣本和20個特征,這些數據都是機器收集回來的,若機器上的某個傳感器損壞導致一個特征無效時該怎么辦?此時是否要扔掉整個數據?這種情況下,另外19個特征怎么辦?它們是否還可用?答案是肯定的。因為有時候數據相當昂貴,扔掉和重新獲取都是不可取的,所以必須采取一些方法來解決這個問題。

  下面給出了一些可選的做法:

  • 使用可用特征的均值來填補缺失值;
  • 使用特征值來填充缺失值,如-1
  • 忽略有缺失值的樣本
  • 使用相似樣本的均值填補缺失值
  • 使用另外的機器學習算法預測缺失值

   現在,我們對下一節要用的數據集進行預處理,使其可以順利地使用分類算法。在預處理階段需要做兩件事:第一,所有的缺失值必須用一個實數值來替換,因為我們使用的Numpy數據類型不允許包含缺失值。這里選擇實數0來替換,因為我們使用的Numpy數據類型不允許包括缺失值,這里選擇實數0來替換所有缺失值,恰好能適應於Logistic回歸。這樣做的直覺在於,我們需要的是一個在更新時不會影響系數的值,回歸系數的更新公式如下:

weights = weights + alpha *error*dataMatrix[randIndex]

  如果dataMatirx的某特征對應值為0,那么該特征的系數將不做更新,即:

weights = weights

  另外,由於sigmoid(0) = 0.5 ,即對它結果的預測不具有任何傾向性,因此上述做法也不會對誤差項造成任何影響,基於上述原因,將缺失值用0代替既可以保留現有數據,也不需要對優化算法進行修改,此外,數據集中的特征值一般不取0,因此在某種意義上說它也滿足“特殊值”這個要求。

  預處理中做的第二件事,是如果在測試數據集中發現了一條數據的類別標簽以及缺失,那么我們的簡單做法是將該條數據丟棄。這是因為類別標簽與特征不同,很難確定采用某個合適的值來替換,采用Logistic回歸進行預處理之后保存兩個文件,horseColicTest.txt和horseColicTraining.txt。如果想對於原始數據和預處理之后的數據做個比較。

  我們有一個“干凈”可用的數據集和一個不錯的優化算法,下面將這些部分融合在一起訓練出一個分類器,然后利用該分類器來預測病馬的生死問題。

  horseColicTest.txt的數據 horseColicTraining.txt的數據  請去我的GitHub下載:https://github.com/LeBron-Jian/MachineLearningNote

 完整代碼:

#-*_coding:utf-8_*_
import math
from numpy import *

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

# 改進的隨機梯度上升算法
def stocGradAscent1(dataMatrix,classLabels,numIter=150):
    m,n = shape(dataMatrix)
    weights = ones(n)
    for j in range(numIter):
        dataIndex = list(range(m))
        for i in range(m):
            alpha = 4/(1.0+j+i)+0.01
            randIndex = int(random.uniform(0,len(dataIndex)))
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
            error = classLabels[randIndex] - h
            weights = weights + alpha *error*dataMatrix[randIndex]
            del(dataIndex[randIndex])
    return weights


# Logistic回歸分類函數
def  classifyVetor(inX,weights):
    prob = sigmoid(sum(inX*weights))
    if prob > 0.5:
        return 1.0
    else:
        return 0.0

def colicTest(filetrain,filetest):
    frTrain = open(filetrain)
    frTest = open(filetest)
    trainingSet = []
    trainingLabeles = []
    for line in frTrain.readlines():
        currLine = line.strip().split('\t')
        lineArr = []
        for i in range(21):
            lineArr.append(float(currLine[i]))
        trainingSet.append(lineArr)
        trainingLabeles.append(float(currLine[21]))
    trainWeights = stocGradAscent1(array(trainingSet),trainingLabeles,500)
    errorCount = 0
    numTestVec = 0
    for line in frTest.readlines():
        numTestVec += 1.0
        currLine = line.strip().split('\t')
        lineArr =[]
        for i in range(21):
            lineArr.append(float(currLine[i]))
        if int(classifyVetor(array(lineArr),trainWeights)) != int(currLine[21]):
            errorCount += 1
    errorRate = (float(errorCount)/numTestVec)
    print('the error rate of this test is  %f'%errorRate)
    return errorRate

def multTest(filetrain,filetest):
    numTests = 6
    errorSum = 1.0
    for k in range(numTests):
        errorSum += colicTest(filetrain,filetest)
    print('after %d iterations  the average error rate is %f'%(numTests,errorSum/float(numTests)))


if __name__ == '__main__':
    filetrain = 'horseColicTraining.txt'
    filetest = 'horseColicTest.txt'
    multTest(filetrain,filetest)

  運行結果:

the error rate of this test is  0.373134
the error rate of this test is  0.358209
the error rate of this test is  0.402985
the error rate of this test is  0.432836
the error rate of this test is  0.462687
the error rate of this test is  0.343284
after 6 iterations  the average error rate is 0.562189

  從結果來看,6次迭代之后的平均錯誤率為0.56,事實上,這個結果還不錯,因為數據有30%的數據缺失,當然如果調整colicTest()中的迭代次數和stochGradAscent1()中的步長,平均錯誤率就可以降到20%左右。

10,使用梯度下降求解邏輯回歸

  這里再補充一個示例,我們將建立一個邏輯回歸模型來預測一個學生是否被大學錄取。假設你是一個大學習的管理員,你想根據兩次考試的結果來決定每個申請人的錄取機會。你又以前的申請人的歷史數據,你可以用它作為邏輯回歸的訓練集。對於每一個培訓例子,你有兩個考試的申請人的分數和錄取決定。為了做到這一點,我們將建立一個分類模型,根據考試成績估計入學概率。

  數據和代碼請參考我的GitHub:https://github.com/LeBron-Jian/MachineLearningNote

  首先,拿到數據,並查看部分數據:

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

path = 'LogiReg_data.txt'
pdData = pd.read_csv(path, header=None, names=['Exam 1', 'Exam 2', 'Admitted'])
print(pdData.head(10))

   部分數據如下:

   拿到數據后,我們將其按照結果分為錄取和不被錄取的集合,並展示:

positive = pdData[pdData['Admitted'] == 1]
negative = pdData[pdData['Admitted'] == 0]
fig, ax = plt.subplots(figsize=(5, 5))
ax.scatter(positive['Exam 1'], positive['Exam 2'], s=30, c='b', marker='o', label='Admitted')
ax.scatter(negative['Exam 1'], negative['Exam 2'], s=30, c='r', marker='x', label='Not Admitted')
ax.legend()
ax.set_xlabel('Exam 1 Score')
ax.set_ylabel('Exam 2 Score')

   圖示如下:

   下面使用邏輯回歸來建立分類器,我們需要求解出三個參數 theta0,  theta1,  theta2。 然后設定閾值,根據閾值判斷錄取結果。

  要完成的模塊函數:

  • Sigmoid:映射到概率的函數
  • model:返回預測結果值
  • cost:根據參數計算損失
  • gradient:計算每個參數的梯度方法
  • descent:進行參數更新
  • accuracy:計算精度

  首先,完成Sigmoid函數:

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

# 如果想查看Sigmoid函數,可以運行下面代碼
nums = np.arange(-10, 10, step=1)
fig, ax = plt.subplots(figsize=(12, 4))
ax.plot(nums, sigmoid(nums), 'r')

   然后對數據集進行划分,分為訓練數據和標簽,其中我們添加了一列,並且初始化了theta值

# 增加一行1  在第0列
pdData.insert(0, 'Ones', 1)
print(pdData.head())
# 設置X和Y  即訓練數據和訓練變量
orig_data = pdData.values
cols = orig_data.shape[1]  # (100, 4)
X = orig_data[:, 0:cols-1]   # 等價於 np.matrix(X.values)
Y = orig_data[:, cols-1:cols]  # 等價於 np.matrox(data.iloc[:, 3:4].value

# 傳遞numpu矩陣 並初始化參數矩陣 theta
theta = np.zeros([1, 3])  # [[0. 0. 0.]]

   我們可以看一下添加了一列后,總體的數據:

   其他的數據均可以自己打印。

  下面是損失函數:

   代碼整理如下:

def cost(X, y, theta):
    left = np.multiply(-y, np.log(model(X, theta)))
    right = np.multiply(1-y, np.log(1 - model(X, theta)))
    return np.sum(left - right) / (len(X))

   下面是計算梯度,梯度函數如下(計算梯度就是我們對theta求偏導):

   代碼可寫成下面方式:

def gradient(X, y, theta):
    grad = np.zeros(theta.shape)
    error = (model(X, theta)-y).ravel()
    for j in range(len(theta.ravel())):  # for each parmeter
        term = np.multiply(error, X[:, j])
        grad[0, j] = np.sum(term) / len(X)
    return grad

   下面比較三種不同梯度下降方法,三種梯度下降分別是 批量梯度下降,隨機梯度下降,小批量梯度下降。如下圖:

 

   我們的方法為1,根據迭代次數停止  2,根據損失函數停止  3,根據梯度下降停止

  (注意:我們每次進行迭代更新的時候,都會對數據進行洗牌,一把情況數據都會有某種規律)

# 比較三種不同梯度下降方法
STOP_ITER = 0
STOP_COST = 1
STOP_GRAD = 2

def stopCriterion(type, value, threshold):
    # 設定三種不同的停止策略
    if type == STOP_ITER:
        return value > threshold
    elif type == STOP_COST:
        return abs(value[-1]-value[-2]) < threshold
    elif type == STOP_GRAD:
        return np.linalg.norm(value) < threshold

# 洗牌
def shuffleData(data):
    np.random.shuffle(data)
    cols = data.shape[1]
    X = data[:, 0:cols-1]
    y = data[:, cols-1:]
    return X, y

def descent(data, theta, batchSize, stopType, thresh, alpha):
    # 梯度下降求解
    init_time = time.time()
    i = 0  # 迭代次數
    k = 0  # batch
    X, y = shuffleData(data)
    grad = np.zeros(theta.shape)  # 計算的梯度
    costs = [cost(X, y, theta)]  # 損失值

    while True:
        grad = gradient(X[k:k+batchSize], y[k:k+batchSize], theta)
        k += batchSize  # 取batch數量個數據
        if k>=n:
            k = 0
            X, y = shuffleData(data)  # 重新洗牌
        theta = theta = alpha*grad  # 參數更新
        costs.append(cost(X, y, theta))  # 計算新的損失
        i += 1

        if stopType == STOP_ITER:
            value = i
        elif stopType == STOP_COST:
            value = costs
        elif stopType == STOP_GRAD:
            value = grad
        if stopCriterion(stopType, value, thresh):
            break
    return theta, i-1, costs, grad, time.time()-init_time

   下面我們定義一個功能性的函數,此函數的作用就是執行一次梯度下降,並完成更新,選擇停止策略。其中主要代碼就第一句。。。

 

def runExpe(data, theta, batchSize, stopType, thresh, alpha):
    theta, iter, costs, grad, dur = descent(data, theta, batchSize, stopType, thresh, alpha)
    name = 'Original' if (data[:, 1] > 2).sum() > 1 else 'Scaled'
    name += 'data - learning rate: {} -'.format(alpha)
    if batchSize == n:
        strDescType = 'Gradient'
    elif batchSize == 1:
        strDescType = 'Stochastic'
    else:
        strDescType = 'Min-batch({})'.format(batchSize)
    name += strDescType + 'descent - Stop: '
    if stopType == STOP_ITER:
        strStop = "{} iterations".format(thresh)
    elif stopType == STOP_COST:
        strStop = "costs change < {}".format(thresh)
    else:
        strStop = 'gradient norm < {}'.format(thresh)
    name += strStop
    print("***{}\nTheta: {} - Iter:{}-Last cost: {:03.2f}--Duration: {:03.2f}s".format(
        name, theta, iter, costs[-1], dur
    ))
    fig, ax = plt.subplots(figsize=(12, 4))
    ax.plot(np.arange(len(costs)), costs, 'r')
    ax.set_xlabel('Iterations')
    ax.set_ylabel('Cost')
    ax.set_title(name.upper() + '- Error vs. Iteration')
    return theta

   下面按照不同的停止策略畫圖展示:

1,設定迭代次數  5000

2,設定閾值1e-6,差不多需要 110000次迭代

3,設置閾值 0.05,差不多需要 40000次迭代

if __name__ == '__main__':
    # 設定迭代次數   選擇的梯度下降方法是基於所有樣本的
    n = 100
    # runExpe(orig_data, theta, n, STOP_ITER, thresh=5000, alpha=0.000001)
    # 根據損失值停止    設定閾值 1e-6 差不多需要110 000 次迭代
    # runExpe(orig_data, theta, n, STOP_COST, thresh=0.000001, alpha=0.001)
    # 根據梯度變換停止  設定閾值 0.05,差不多需要 40000次迭代
    # runExpe(orig_data, theta, n, STOP_GRAD, thresh=0.05, alpha=0.001)
    # 對比不同的梯度下降方法  stochastic descent  這個不收斂
    # runExpe(orig_data, theta, 1, STOP_ITER, thresh=5000, alpha=0.001)
    # 我們將學習率調小,提高thresh再試試,發現收斂了
    # runExpe(orig_data, theta, 1, STOP_ITER, thresh=15000, alpha=0.000001)
    # Mine-batch descent  這個也不收斂 ,我們將alpha調小 0.000001  發現收斂了
    runExpe(orig_data, theta, 16, STOP_ITER, thresh=15000, alpha=0.001)

   浮動仍然比較大,我們來嘗試下對數據進行標准化,將數據按其屬下(按列進行)減去其均值,然后除以其方差。最后得到的結果是,對每個屬性/每列來說所有數據都聚集在 0 附加,方差值為1。

    scaled_data = orig_data.copy()
    scaled_data[:, 1:3] = pp.scale(orig_data[:, 1:3])
    runExpe(scaled_data, theta, n, STOP_ITER, thresh=5000, alpha=0.01)

   下面兩幅圖,我們分別將學習率設置為 0.001  0.01

 ·  我們發現,這里已經降到了0.3以下,所以說做數據預處理非常重要。

    scaled_data = orig_data.copy()
    scaled_data[:, 1:3] = pp.scale(orig_data[:, 1:3])
    # runExpe(scaled_data, theta, n, STOP_ITER, thresh=5000, alpha=0.01)
    runExpe(scaled_data, theta, n, STOP_GRAD, thresh=0.02, alpha=0.01)

   我們修改一下,根據梯度變換停止,發現更多的迭代次數會使得損失下降的更多

runExpe(scaled_data, theta, 1, STOP_GRAD, thresh=0.002/5, alpha=0.01)

   隨機梯度下降更快,但是我們需要迭代的次數也需要更多,所以還是用 batch的比較合適。

   精度

# 設置閾值
def predict(x, theta):
    return [1 if x >= 0.5 else 0 for x in model(X, theta)]


if __name__ == '__main__':
    # 設定迭代次數   選擇的梯度下降方法是基於所有樣本的
    n = 100
    scaled_data = orig_data.copy()
    scaled_X = scaled_data[:, :3]
    y = scaled_data[:, 3]
    predictions = predict(scaled_X, theta)
    correct = [1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0 for a, b in zip(predictions, y)]
    accuracy = (sum(map(int, correct)) % len(correct))
    print('accuracy = {0}%'.format(accuracy))
    print(len(correct))

   不知道為什么,我求出來的精度不高。。。。百分之60。

11,總結

  Logistic回歸的目的是尋找一個非線性函數Sigmoid的最佳擬合參數,求解過程可以由最優化算法來完成。在最優化算法中,最常用的就是梯度上升算法,而梯度上升算法又可以簡化為隨機梯度上升算法。

  隨機梯度上升算法與梯度上升算法的效果相當,但是占用更少的計算資源。此外,隨機梯度上升是一個在線算法,它可以在新數據到來時就完成參數更新,不需要重新讀取整個數據集來進行批處理運算。

  機器學習的一個重要問題就是如何處理缺失數據,這個問題沒有標准答案,取決於實際應用中的需求。

12,推廣

  前面也說了Logistic回歸模型主要用於二分類,那么下面說一下多分類問題中的推廣——softmax回歸。

  softmax與Logistic回歸的主要區別就是,Logistic處理二分類問題,只有一組權重參數θ,而softmax處理多分類問題,如果有k個類別,那么softmax就有k組權值參數。每組權值對應一種分類,通過k組權值求解出樣本數據對應每個類別的概率,最后取概率最大的類別作為該數據的分類結果,它的概率函數為:

  softmax函數經常用於神經網絡的最后一層,用於對神經網絡已經處理好的特征進行分類。

  總結來說:邏輯回歸真的真的很好用!

基於Sklearn構建Logistic回歸分類器 

  下面讓我們看一下Sklearn的Logistic回歸分類器

  英文的Sklearn文檔地址:請點擊我

  sklearn.linear_model模塊提供了很多模型供我們使用,比如Logistic回歸、Lasso回歸、貝葉斯脊回歸等,可見需要學習的東西還有很多很多。本次,我們使用LogisticRegressioin。

1,LogisticRegression

  讓我們先看一下LogisticRegression這個函數,一共有14個參數

參數說明如下:

  penalty:懲罰項,str類型,可選參數為l1和l2,默認為l2。用於指定懲罰項中使用的規范。newton-cg、sag和lbfgs求解算法只支持L2規范。L1G規范假設的是模型的參數滿足拉普拉斯分布,L2假設的模型參數滿足高斯分布,所謂的范式就是加上對參數的約束,使得模型更不會過擬合(overfit),但是如果要說是不是加了約束就會好,這個沒有人能回答,只能說,加約束的情況下,理論上應該可以獲得泛化能力更強的結果。
  dual:對偶或原始方法,bool類型,默認為False。對偶方法只用在求解線性多核(liblinear)的L2懲罰項上。當樣本數量>樣本特征的時候,dual通常設置為False。
  tol:停止求解的標准,float類型,默認為1e-4。就是求解到多少的時候,停止,認為已經求出最優解。
  c:正則化系數λ的倒數,float類型,默認為1.0。必須是正浮點型數。像SVM一樣,越小的數值表示越強的正則化。
  fit_intercept:是否存在截距或偏差,bool類型,默認為True。如果使用中心化的數據,可以考慮設置為False,不考慮截距。注意這里是考慮,一般還是要考慮截距
  intercept_scaling:僅在正則化項為”liblinear”,且fit_intercept設置為True時有用。float類型,默認為1。
  class_weight:用於標示分類模型中各種類型的權重,可以是一個字典或者balanced字符串,默認為不輸入,也就是不考慮權重,即為None。如果選擇輸入的話,可以選擇balanced讓類庫自己計算類型權重,或者自己輸入各個類型的權重。舉個例子,比如對於0,1的二元模型,我們可以定義class_weight={0:0.9,1:0.1},這樣類型0的權重為90%,而類型1的權重為10%。如果class_weight選擇balanced,那么類庫會根據訓練樣本量來計算權重。某種類型樣本量越多,則權重越低,樣本量越少,則權重越高。當class_weight為balanced時,類權重計算方法如下:n_samples / (n_classes * np.bincount(y))。n_samples為樣本數,n_classes為類別數量,np.bincount(y)會輸出每個類的樣本數,例如y=[1,0,0,1,1],則np.bincount(y)=[2,3]。
    那么class_weight有什么作用呢? 在分類模型中,我們經常會遇到兩類問題:

  • 1,第一種是誤分類的代價很高。比如對合法用戶和非法用戶進行分類,將非法用戶分類為合法用戶的代價很高,我們寧願將合法用戶分類為非法用戶,這時可以人工再甄別,但是卻不願將非法用戶分類為合法用戶。這時,我們可以適當提高非法用戶的權重。
  • 2,第二種是樣本是高度失衡的,比如我們有合法用戶和非法用戶的二元樣本數據10000條,里面合法用戶有9995條,非法用戶只有5條,如果我們不考慮權重,則我們可以將所有的測試集都預測為合法用戶,這樣預測准確率理論上有99.95%,但是卻沒有任何意義。這時,我們可以選擇balanced,讓類庫自動提高非法用戶樣本的權重。提高了某種分類的權重,相比不考慮權重,會有更多的樣本分類划分到高權重的類別,從而可以解決上面兩類問題。

  random_state:隨機數種子,int類型,可選參數,默認為無,僅在正則化優化算法為sag,liblinear時有用。
  solver:優化算法選擇參數,只有五個可選參數,即newton-cg,lbfgs,liblinear,sag,saga。默認為liblinear。solver參數決定了我們對邏輯回歸損失函數的優化方法,有四種算法可以選擇,分別是:

  • liblinear:使用了開源的liblinear庫實現,內部使用了坐標軸下降法來迭代優化損失函數。
  • lbfgs:擬牛頓法的一種,利用損失函數二階導數矩陣即海森矩陣來迭代優化損失函數。
  • newton-cg:也是牛頓法家族的一種,利用損失函數二階導數矩陣即海森矩陣來迭代優化損失函數。
  • sag:即隨機平均梯度下降,是梯度下降法的變種,和普通梯度下降法的區別是每次迭代僅僅用一部分的樣本來計算梯度,適合於樣本數據多的時候。
  • saga:線性收斂的隨機優化算法的的變重。

總結:

  • liblinear適用於小數據集,而sag和saga適用於大數據集因為速度更快。
  • 對於多分類問題,只有newton-cg,sag,saga和lbfgs能夠處理多項損失,而liblinear受限於一對剩余(OvR)。啥意思,就是用liblinear的時候,如果是多分類問題,得先把一種類別作為一個類別,剩余的所有類別作為另外一個類別。一次類推,遍歷所有類別,進行分類。
  • newton-cg,sag和lbfgs這三種優化算法時都需要損失函數的一階或者二階連續導數,因此不能用於沒有連續導數的L1正則化,只能用於L2正則化。而liblinear和saga通吃L1正則化和L2正則化。
  • 同時,sag每次僅僅使用了部分樣本進行梯度迭代,所以當樣本量少的時候不要選擇它,而如果樣本量非常大,比如大於10萬,sag是第一選擇。但是sag不能用於L1正則化,所以當你有大量的樣本,又需要L1正則化的話就要自己做取舍了。要么通過對樣本采樣來降低樣本量,要么回到L2正則化。
  • 從上面的描述,大家可能覺得,既然newton-cg, lbfgs和sag這么多限制,如果不是大樣本,我們選擇liblinear不就行了嘛!錯,因為liblinear也有自己的弱點!我們知道,邏輯回歸有二元邏輯回歸和多元邏輯回歸。對於多元邏輯回歸常見的有one-vs-rest(OvR)和many-vs-many(MvM)兩種。而MvM一般比OvR分類相對准確一些。郁悶的是liblinear只支持OvR,不支持MvM,這樣如果我們需要相對精確的多元邏輯回歸時,就不能選擇liblinear了。也意味着如果我們需要相對精確的多元邏輯回歸不能使用L1正則化了。

  max_iter:算法收斂最大迭代次數,int類型,默認為10。僅在正則化優化算法為newton-cg, sag和lbfgs才有用,算法收斂的最大迭代次數。
  multi_class:分類方式選擇參數,str類型,可選參數為ovr和multinomial,默認為ovr。ovr即前面提到的one-vs-rest(OvR),而multinomial即前面提到的many-vs-many(MvM)。如果是二元邏輯回歸,ovr和multinomial並沒有任何區別,區別主要在多元邏輯回歸上。

  • OvR和MvM有什么不同?
  1. OvR的思想很簡單,無論你是多少元邏輯回歸,我們都可以看做二元邏輯回歸。具體做法是,對於第K類的分類決策,我們把所有第K類的樣本作為正例,除了第K類樣本以外的所有樣本都作為負例,然后在上面做二元邏輯回歸,得到第K類的分類模型。其他類的分類模型獲得以此類推。
  2. 而MvM則相對復雜,這里舉MvM的特例one-vs-one(OvO)作講解。如果模型有T類,我們每次在所有的T類樣本里面選擇兩類樣本出來,不妨記為T1類和T2類,把所有的輸出為T1和T2的樣本放在一起,把T1作為正例,T2作為負例,進行二元邏輯回歸,得到模型參數。我們一共需要T(T-1)/2次分類。
  3. 可以看出OvR相對簡單,但分類效果相對略差(這里指大多數樣本分布情況,某些樣本分布下OvR可能更好)。而MvM分類相對精確,但是分類速度沒有OvR快。如果選擇了ovr,則4種損失函數的優化方法liblinear,newton-cg,lbfgs和sag都可以選擇。但是如果選擇了multinomial,則只能選擇newton-cg, lbfgs和sag了。

  verbose:日志冗長度,int類型。默認為0。就是不輸出訓練過程,1的時候偶爾輸出結果,大於1,對於每個子模型都輸出。
  warm_start:熱啟動參數,bool類型。默認為False。如果為True,則下一次訓練是以追加樹的形式進行(重新使用上一次的調用作為初始化)。
  n_jobs:並行數。int類型,默認為1。1的時候,用CPU的一個內核運行程序,2的時候,用CPU的2個內核運行程序。為-1的時候,用所有CPU的內核運行程序。

除此之外,LogisticRegression也有一些方法供我們使用:

有一些方法和MultinomialNB的方法都是類似的,因此不再累述。

 2,編寫代碼

  了解了基本知識點,我們就可以編寫Sklearn分類器的代碼了,代碼如下:

#_*_ codingLutf-8_*_

from sklearn.linear_model import LogisticRegression

def colicSklearn(filetrain,filetest):
    frTrain = open(filetrain)
    frTest = open(filetest)
    trainingSet = []
    trainingLabels = []
    testSet = []
    testLabels = []
    for line in frTrain.readlines():
        currLine = line.strip().split('\t')
        lineArr = []
        for i in range(len(currLine)-1):
            lineArr.append(float(currLine[i]))
        trainingSet.append(lineArr)
        trainingLabels.append(float(currLine[-1]))
    for line in frTest.readlines():
        currLine = line.strip().split('\t')
        lineArr = []
        for i in range(len(currLine)-1):
            lineArr.append(float(currLine[i]))
        testSet.append(lineArr)
        testLabels.append(float(currLine[-1]))
    classifier = LogisticRegression(solver='liblinear',max_iter=20).fit(trainingSet,trainingLabels)
    test_accurcy = classifier.score(testSet,testLabels)*100
    print("正確率為%s%%"%test_accurcy)

if __name__ == '__main__':
    filetrain = 'horseColicTraining.txt'
    filetest = 'horseColicTest.txt'
    colicSklearn(filetrain,filetest)

  執行結果如下:

正確率為73.13432835820896%

  可以看到,正確率又搞了,更改solver參數,比如設置為sag,使用隨機平均梯度下降算法,看一看效果,你會發現:正確率高了,但是發出了警告,發出警告是因為算法還沒有收斂,更改迭代次數即可。

 當我們迭代到5000的時候,就不會報錯了,所以說,對於小數據集,sag算法需要迭代上千次才收斂,而liblinear只需要不到10次。

  所以,我們需要根據數據集情況,選擇最優化算法。

3,sklearn利用LR模型進行三分類的原理及其代碼

  首先,LR將線性模型利用Sigmoid函數進一步做了非線性映射。將分類超平面兩側的正負樣本點通過壓縮函數轉化成了以 0.5 為分類的兩類:類別0 和類別1。

  在上面我們講述了線性邊界與LR分布函數(即Sigmoid函數)的映射對應關系;同樣,對於非線性判定邊界Sigmoid函數也使用,如下:

  然后,再去考慮如何去獲得模型參數,就是判定邊界的參數怎么獲得,這里是利用MLE進行求解的,具體求解過程可以百度。

  代碼如下:

#_*_coding:utf-8_*_
import numpy as np
import matplotlib.pyplot as plt

def sigmoid(z):
    return 1.0/(1.0 + np.exp(-z))

z = np.arange(-7, 7, 0.1)
print(z)
phi_z = sigmoid(z)

def cost_1(z):
    return -np.log(sigmoid(z))

def cost_0(z):
    return -np.log(1-sigmoid(z))


z = np.arange(-10, 10, 0.1)
phi_z = sigmoid(z)

c1 = [cost_1(x) for x in z]
plt.plot(phi_z, c1, label='J(w) if y=1')

c0 = [cost_0(x) for x in z]
plt.plot(phi_z, c0, linestyle='--', label='J(w) if y=0')

plt.ylim(0.0, 5.1)
plt.xlim([0, 1])
plt.xlabel('$\phi$(z)')
plt.ylabel('J(w)')
plt.legend(loc='best')
plt.tight_layout()
# plt.savefig('./figures/log_cost.png', dpi=300)
plt.show()

  結果如下:

  最后,考慮LR進行三分類(多分類)時,是特征的線性組合和Sigmoid函數符合的函數進行概率計算和分類的。

  代碼:

from IPython.display import Image

# Added version check for recent scikit-learn 0.18 checks
from distutils.version import LooseVersion as Version
from sklearn import __version__ as sklearn_version

from sklearn import datasets
import numpy as np

iris = datasets.load_iris()
# http://scikit-learn.org/stable/auto_examples/datasets/plot_iris_dataset.html
X = iris.data[:, [2, 3]]
print(X.shape)
y = iris.target  # 取species列,類別

if Version(sklearn_version) < '0.18':
    from sklearn.cross_validation import train_test_split
else:
    from sklearn.model_selection import train_test_split

# train_test_split方法分割數據集
X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.3, random_state=0)

from sklearn.preprocessing import StandardScaler

sc = StandardScaler()  # 初始化一個對象sc去對數據集作變換
sc.fit(X_train)  # 用對象去擬合數據集X_train,並且存下來擬合參數
X_train_std = sc.transform(X_train)
X_test_std = sc.transform(X_test)

from sklearn.linear_model import LogisticRegression


def sigmoid(z):
    return 1.0 / (1.0 + np.exp(-z))


lr = LogisticRegression(C=1000.0, random_state=0)
# http://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LogisticRegression.html#sklearn.linear_model.LogisticRegression
lr.fit(X_train_std, y_train)

# 計算該預測實例點屬於各類的概率
lr.predict_proba(X_test_std[0, :].reshape(1, -1))
# Output:array([[  2.05743774e-11,   6.31620264e-02,   9.36837974e-01]])

# 驗證predict_proba的作用
c = lr.predict_proba(X_test_std[0, :].reshape(1, -1))
# c[0, 0] + c[0, 1] + c[0, 2]
# Output:0.99999999999999989

# 查看lr模型的特征系數
lr = LogisticRegression(C=1000.0, random_state=0)
# http://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LogisticRegression.html#sklearn.linear_model.LogisticRegression
lr.fit(X_train_std, y_train)
print(lr.coef_)
# Output:[[-7.34015187 -6.64685581]
#        [ 2.54373335 -2.3421979 ]
#        [ 9.46617627  6.44380858]]

# 驗證predict_proba工作原理
Zz = np.dot(lr.coef_, X_test_std[0, :].T) + lr.intercept_
np.array(sigmoid(Zz)) / sum(np.array(sigmoid(Zz)))
# Output:array([  2.05743774e-11,   6.31620264e-02,   9.36837974e-01])
# 此結果就是預測實例點各類的概率

  

 

4,使用Sklearn的Logistic回歸算法計算鳶尾花

  對鳶尾花進行分類代碼:

from sklearn import datasets
from numpy import *
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split


def colicSklearn():
    iris = datasets.load_iris()
    X = iris.data
    Y = iris.target
    trainingSet,testSet,trainingLabels,testLabels = train_test_split(X,Y,test_size=0.25,random_state=40)
    classifier = LogisticRegression(solver='sag', max_iter=5000).fit(trainingSet, trainingLabels)
    test_accurcy = classifier.score(testSet, testLabels) * 100
    print("正確率為%s%%" % test_accurcy)

if __name__  == '__main__':
    colicSklearn()

  結果:

正確率為100.0%

  

5,使用Sklearn的Logistic回歸算法預測線性回歸函數

 代碼:

#_*_coding:utf-8_*_
'''
下面這個例子,從數據產生,到數據提取,數據標准化
模型訓練和評估來說明各個API的調用過程
'''
print(__doc__)
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
# import matplotlib as mpl
import pylab as mpl

# 設置字符集,防止中文亂碼
mpl.rcParams['font.sans-serif'] = [u'simHei']
mpl.rcParams['axes.unicode_minus'] = False

# 定義目標函數通過改函數產生對應的y
# y = 1 *x[0] + 2 *x[1] + ... (n+1) *x[n]
def l_model(x):
    params = np.arange(1, x.shape[-1] +1)
    y = np.sum(params *x) + np.random.randn(1) * 0.1
    return y

# 定義數據集
x = pd.DataFrame(np.random.rand(500,6))
# print(x)
y = x.apply(lambda x_rows:pd.Series(l_model(x_rows)),axis=1)
# print(y)

# 划分數據集
x_train, x_test, y_train, y_test = train_test_split(x,y,test_size=0.3,random_state=2)

# 數據標准化
ss = StandardScaler()
x_train_s = ss.fit_transform(x_train)
x_test_s = ss.fit_transform(x_test)

# 輸出下元數據的標准差和平均數
print(ss.scale_)
print(ss.mean_)

# 訓練模型
lr = LinearRegression()
lr.fit(x_train_s , y_train)
# 訓練后的輸入端模型系數,如果label有兩個,即y值有兩列,那么是一個2D的array
print(lr.coef_)
# 截距
print(lr.intercept_)

# 用模型預測,並計算得分
y_predict = lr.predict(x_test_s)
test_accuracy = lr.score(x_test_s, y_test)

print("正確率為%s%%" % test_accuracy)

# 預測值和實際值畫圖比較
t = np.arange(len(x_test_s))
# 建一個畫布,facecolor是背景色
plt.figure(facecolor='W')
plt.plot(t, y_test, 'r-', linewidth= 2, label = '真實值')
plt.plot(t, y_predict, 'b-', linewidth= 1, label = '預測值')
# 顯示圖例,設置圖例的位置
plt.legend(loc= 'upper left')
plt.title("線性回歸預測真實值之間的關系", fontsize = 20)
# 加網格
plt.grid(b = True)
plt.show()

  結果展示:

下面這個例子,從數據產生,到數據提取,數據標准化
模型訓練和評估來說明各個API的調用過程


[0.28299035 0.3138177  0.29522932 0.30125841 0.28402977 0.28989605]
[0.53282546 0.51045318 0.48925819 0.5202345  0.4989613  0.54335117]
[[0.29202621 0.5787407  0.87280326 1.12786483 1.50043119 1.73271735]]
[10.38441274]

正確率為0.9706613427158242%

  

完整代碼及其數據,請移步小編的GitHub

  傳送門:請點擊我

  如果點擊有誤:https://github.com/LeBron-Jian/MachineLearningNote

 

 

參考:https://www.cnblogs.com/bonelee/p/7253508.html

https://blog.csdn.net/c406495762/article/details/77851973


免責聲明!

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



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