用python做逻辑回归梯度上升_机器学习实例---4.1、Logistic回归基础篇之梯度上升算法...
一 前言
本文從Logistic回歸的原理開始講起，補充了書上省略的數學推導。本文可能會略顯枯燥，理論居多，Sklearn實戰內容會放在下一篇文章。自己慢慢推導完公式，還是蠻開心的一件事。
二 Logistic回歸與梯度上升算法
Logistic回歸是眾多分類算法中的一員。通常，Logistic回歸用于二分類問題，例如預測明天是否會下雨。當然它也可以用于多分類問題，不過為了簡單起見，本文暫先討論二分類問題。首先，讓我們來了解一下，什么是Logistic回歸。
1 Logistic回歸
假設現在有一些數據點，我們利用一條直線對這些點進行擬合(該線稱為最佳擬合直線)，這個擬合過程就稱作為回歸，如下圖所示：
Logistic回歸一種二分類算法，它利用的是Sigmoid函數閾值在[0,1]這個特性。Logistic回歸進行分類的主要思想是：根據現有數據對分類邊界線建立回歸公式，以此進行分類。其實，Logistic本質上是一個基于條件概率的判別模型(Discriminative Model)。
所以要想了解Logistic回歸，我們必須先看一看Sigmoid函數 ，我們也可以稱它為Logistic函數。它的公式如下：
整合成一個公式，就變成了如下公式：
下面這張圖片，為我們展示了Sigmoid函數的樣子。
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函數的特性，我們可以做出如下的假設：
上式即為在已知樣本x和參數θ的情況下，樣本x屬性正樣本(y=1)和負樣本(y=0)的條件概率。理想狀態下，根據上述公式，求出各個點的概率均為1，也就是完全分類都正確。但是考慮到實際情況，樣本點的概率越接近于1，其分類效果越好。比如一個樣本屬于正樣本的概率為0.51，那么我們就可以說明這個樣本屬于正樣本。另一個樣本屬于正樣本的概率為0.99，那么我們也可以說明這個樣本屬于正樣本。但是顯然，第二個樣本概率更高，更具說服力。我們可以把上述兩個概率公式合二為一：
合并出來的Cost，我們稱之為代價函數(Cost Function)。當y等于1時，(1-y)項(第二項)為0；當y等于0時，y項(第一項)為0。為了簡化問題，我們對整個表達式求對數，(將指數問題對數化是處理數學問題常見的方法)：
這個代價函數，是對于一個樣本而言的。給定一個樣本，我們就可以通過這個代價函數求出，樣本所屬類別的概率，而這個概率越大越好，所以也就是求解這個代價函數的最大值。既然概率出來了，那么最大似然估計也該出場了。假定樣本與樣本之間相互獨立，那么整個樣本集生成的概率即為所有樣本生成概率的乘積，再將公式對數化，便可得到如下公式：
其中，m為樣本的總數，y(i)表示第i個樣本的類別，x(i)表示第i個樣本，需要注意的是θ是多維向量，x(i)也是多維向量。
綜上所述，滿足J(θ)的最大的θ值即是我們需要求解的模型。
怎么求解使J(θ)最大的θ值呢？因為是求最大值，所以我們需要使用梯度上升算法。如果面對的問題是求解使J(θ)最小的θ值，那么我們就需要使用梯度下降算法。面對我們這個問題，如果使J(θ) := -J(θ)，那么問題就從求極大值轉換成求極小值了，使用的算法就從梯度上升算法變成了梯度下降算法，它們的思想都是相同的，學會其一，就也會了另一個。本文使用梯度上升算法進行求解。
2 梯度上升算法
說了半天，梯度上升算法又是啥？J(θ)太復雜，我們先看個簡單的求極大值的例子。一個看了就會想到高中生活的函數：
來吧，做高中題。這個函數的極值怎么求？顯然這個函數開口向下，存在極大值，它的函數圖像為：
求極值，先求函數的導數：
令導數為0，可求出x=2即取得函數f(x)的極大值。極大值等于f(2)=4
但是真實環境中的函數不會像上面這么簡單，就算求出了函數的導數，也很難精確計算出函數的極值。此時我們就可以用迭代的方法來做。就像爬坡一樣，一點一點逼近極值。這種尋找最佳擬合參數的方法，就是最優化算法。爬坡這個動作用數學公式表達即為：
其中，α為步長，也就是學習速率，控制更新的幅度。效果如下圖所示：
比如從(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()
代碼運行結果如下：
結果很顯然，已經非常接近我們的真實極值2了。這一過程，就是梯度上升算法。那么同理，J(θ)這個函數的極值，也可以這么求解。公式可以這么寫：
由上小節可知J(θ)為：
sigmoid函數為：
那么，現在我只要求出J(θ)的偏導，就可以利用梯度上升算法，求解J(θ)的極大值了。
那么現在開始求解J(θ)對θ的偏導，求解如下(數學推導)：
其中：
再由：
可得：
接下來，就剩下第三部分：
綜上所述：
因此，梯度上升迭代公式為：
知道了，梯度上升迭代公式，我們就可以自己編寫代碼，計算最佳擬合參數了。
三 Python3實戰
1 數據準備
這就是一個簡單的數據集，沒什么實際意義。讓我們先從這個簡單的數據集開始學習。先看下數據集有哪些數據：
這個數據有兩維特征，因此可以將數據在一個二維平面上展示出來。我們可以將第一列數據(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()
運行結果如下：
從上圖可以看出數據的分布情況。假設Sigmoid函數的輸入記為z，那么z=w0x0 + w1x1 + w2x2，即可將數據分割開。其中，x0為全是1的向量，x1為數據集的第一列數據，x2為數據集的第二列數據。另z=0，則0=w0 + w1x1 + w2x2。橫坐標為x1，縱坐標為x2。這個方程未知的參數為w0，w1，w2，也就是我們需要求的回歸系數(最優參數)。
2 訓練算法
在編寫代碼之前，讓我們回顧下梯度上升迭代公式：
將上述公式矢量化：
根據矢量化的公式，編寫代碼如下：
# -*- coding:UTF-8 -*-
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 #返回
"""
函數說明:sigmoid函數
Parameters:
inX - 數據
Returns:
sigmoid函數
Author:
Jack Cui
Blog:
http://blog.csdn.net/c406495762
Zhihu:
https://www.zhihu.com/people/Jack--Cui/
Modify:
2017-08-28
"""
def sigmoid(inX):
return 1.0 / (1 + np.exp(-inX))
"""
函數說明:梯度上升算法
Parameters:
dataMatIn - 數據集
classLabels - 數據標簽
Returns:
weights.getA() - 求得的權重數組(最優參數)
Author:
Jack Cui
Blog:
http://blog.csdn.net/c406495762
Zhihu:
https://www.zhihu.com/people/Jack--Cui/
Modify:
2017-08-28
"""
def gradAscent(dataMatIn, classLabels):
dataMatrix = np.mat(dataMatIn) #轉換成numpy的mat
labelMat = np.mat(classLabels).transpose() #轉換成numpy的mat,并進行轉置
m, n = np.shape(dataMatrix) #返回dataMatrix的大小。m為行數,n為列數。
alpha = 0.001 #移動步長,也就是學習速率,控制更新的幅度。
maxCycles = 500 #最大迭代次數
weights = np.ones((n,1))
for k in range(maxCycles):
h = sigmoid(dataMatrix * weights) #梯度上升矢量化公式
error = labelMat - h
weights = weights + alpha * dataMatrix.transpose() * error
return weights.getA() #將矩陣轉換為數組，返回權重數組
if __name__ == '__main__':
dataMat, labelMat = loadDataSet()
print(gradAscent(dataMat, labelMat))
運行結果如圖所示：
可以看出，我們已經求解出回歸系數[w0,w1,w2]。
通過求解出的參數，我們就可以確定不同類別數據之間的分隔線，畫出決策邊界。
3 繪制決策邊界
我們已經解出了一組回歸系數，它確定了不同類別數據之間的分隔線。現在開始繪制這個分隔線，編寫代碼如下：
# -*- 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 #返回
"""
函數說明:sigmoid函數
Parameters:
inX - 數據
Returns:
sigmoid函數
Author:
Jack Cui
Blog:
http://blog.csdn.net/c406495762
Zhihu:
https://www.zhihu.com/people/Jack--Cui/
Modify:
2017-08-28
"""
def sigmoid(inX):
return 1.0 / (1 + np.exp(-inX))
"""
函數說明:梯度上升算法
Parameters:
dataMatIn - 數據集
classLabels - 數據標簽
Returns:
weights.getA() - 求得的權重數組(最優參數)
Author:
Jack Cui
Blog:
http://blog.csdn.net/c406495762
Zhihu:
https://www.zhihu.com/people/Jack--Cui/
Modify:
2017-08-28
"""
def gradAscent(dataMatIn, classLabels):
dataMatrix = np.mat(dataMatIn) #轉換成numpy的mat
labelMat = np.mat(classLabels).transpose() #轉換成numpy的mat,并進行轉置
m, n = np.shape(dataMatrix) #返回dataMatrix的大小。m為行數,n為列數。
alpha = 0.001 #移動步長,也就是學習速率,控制更新的幅度。
maxCycles = 500 #最大迭代次數
weights = np.ones((n,1))
for k in range(maxCycles):
h = sigmoid(dataMatrix * weights) #梯度上升矢量化公式
error = labelMat - h
weights = weights + alpha * dataMatrix.transpose() * error
return weights.getA() #將矩陣轉換為數組，返回權重數組
"""
函數說明:繪制數據集
Parameters:
weights - 權重參數數組
Returns:
無
Author:
Jack Cui
Blog:
http://blog.csdn.net/c406495762
Zhihu:
https://www.zhihu.com/people/Jack--Cui/
Modify:
2017-08-30
"""
def plotBestFit(weights):
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) #繪制負樣本
x = np.arange(-3.0, 3.0, 0.1)
y = (-weights[0] - weights[1] * x) / weights[2]
ax.plot(x, y)
plt.title('BestFit') #繪制title
plt.xlabel('X1'); plt.ylabel('X2') #繪制label
plt.show()
if __name__ == '__main__':
dataMat, labelMat = loadDataSet()
weights = gradAscent(dataMat, labelMat)
plotBestFit(weights)
運行結果如下：
這個分類結果相當不錯，從上圖可以看出，只分錯了幾個點而已。但是，盡管例子簡單切數據集很小，但是這個方法卻需要大量的計算(300次乘法)。因此下篇文章將對改算法稍作改進，從而減少計算量，使其可以應用于大數據集上。
四 總結
Logistic回歸的一般過程：
收集數據：采用任意方法收集數據。
準備數據：由于需要進行距離計算，因此要求數據類型為數值型。另外，結構化數據格式則最佳。
分析數據：采用任意方法對數據進行分析。
訓練算法：大部分時間將用于訓練，訓練的目的是為了找到最佳的分類回歸系數。
測試算法：一旦訓練步驟完成，分類將會很快。
使用算法：首先，我們需要輸入一些數據，并將其轉換成對應的結構化數值；接著，基于訓練好的回歸系數，就可以對這些數值進行簡單的回歸計算，判定它們屬于哪個類別；在這之后，我們就可以在輸出的類別上做一些其他分析工作。
其他：
Logistic回歸的目的是尋找一個非線性函數Sigmoid的最佳擬合參數，求解過程可以由最優化算法完成。
本文講述了Logistic回歸原理以及數學推導過程。
下篇文章將講解Logistic回歸的改進以及Sklearn實戰內容。
如有問題，請留言。如有錯誤，還望指正，謝謝！
總結
以上是生活随笔為你收集整理的用python做逻辑回归梯度上升_机器学习实例---4.1、Logistic回归基础篇之梯度上升算法...的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： python程序的基本结构知识点总结_P
- 下一篇： python内置字符串处理函数_Pyth