作 者：崔家華 編 輯：李文臣

Python版本： Python3.x

執行平臺 ： Windows

IDE ： Sublime text3

一、前言

本系列文章對梯度上升演算法和改進的隨機梯度上升演算法進行了對比，總結了各自的優缺點，並對sklearn.linear_model.LogisticRegression進行了詳細介紹。

二、改進的隨機梯度上升演算法

梯度上升演算法在每次更新迴歸係數(最優引數)時，都需要遍歷整個資料集。可以看一下我們之前寫的梯度上升演算法：

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.01 
      #移動步長,也就是學習速率,控制更新的幅度。
    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(),weights_array 
           #將矩陣轉換為陣列，返回權重陣列
 

假設，我們使用的資料集一共有100個樣本。那麼，dataMatrix就是一個100*3的矩陣。每次計算h的時候，都要計算dataMatrix*weights這個矩陣乘法運算，要進行100*3次乘法運算和100*2次加法運算。同理，更新迴歸係數(最優引數)weights時，也需要用到整個資料集，要進行矩陣乘法運算。總而言之，該方法處理100個左右的資料集時尚可，但如果有數十億樣本和成千上萬的特徵，那麼該方法的計算複雜度就太高了。因此，需要對演算法進行改進，我們每次更新迴歸係數(最優引數)的時候，能不能不用所有樣本呢？一次只用一個樣本點去更新迴歸係數(最優引數)？這樣就可以有效減少計算量了，這種方法就叫做隨機梯度上升演算法。

1、隨機梯度上升演算法

讓我們直接看程式碼：

def stocGradAscent1(dataMatrix, classLabels, numIter=150):
    m,n =np.shape(dataMatrix)
            #返回dataMatrix的大小,m為行數,n為列數。
    weights = np.ones(n)  #引數初始化
    for j in range(numIter): 
        dataIndex = list(range(m)
        for i in range(m): 
               #降低alpha的大小，每次減小1/(j+i)。
            randIndex = int(random.uniform(0,len(dataIndex))) 
                #隨機選取樣本
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
                #選擇隨機選取的一個樣本，計算h
            error = classLabels[randIndex] - h 
                #計算誤差
            weights = weights + alpha 
                   * error * dataMatrix[randIndex] 
               #更新迴歸係數
            del(dataIndex[randIndex]) 
                     #刪除已經使用的樣本
    return weights #返回
 

該演算法第一個改進之處在於，alpha在每次迭代的時候都會調整，並且，雖然alpha會隨著迭代次數不斷減小，但永遠不會減小到0，因為這裡還存在一個常數項。必須這樣做的原因是為了保證在多次迭代之後新資料仍然具有一定的影響。如果需要處理的問題是動態變化的，那麼可以適當加大上述常數項，來確保新的值獲得更大的迴歸係數。另一點值得注意的是，在降低alpha的函式中，alpha每次減少1/(j+i)，其中j是迭代次數，i是樣本點的下標。第二個改進的地方在於跟新迴歸係數(最優引數)時，只使用一個樣本點，並且選擇的樣本點是隨機的，每次迭代不使用已經用過的樣本點。這樣的方法，就有效地減少了計算量，並保證了迴歸效果。

編寫程式碼如下，看下改進的隨機梯度上升演算法分類效果如何：

# -*- coding:UTF-8 -*-from matplotlib.font_manager import FontPropertiesimport matplotlib.pyplot as pltimport numpy as npimport random"""函式說明:載入資料Parameters:    無Returns:    dataMat - 資料列表    labelMat - 標籤列表Author:    Jack CuiBlog:    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 CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-28"""def sigmoid(inX):    return 1.0 / (1 + np.exp(-inX))"""函式說明:繪製資料集Parameters:    weights - 權重引數陣列Returns:    無Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    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()"""函式說明:改進的隨機梯度上升演算法Parameters:    dataMatrix - 資料陣列    classLabels - 資料標籤    numIter - 迭代次數Returns:    weights - 求得的迴歸係數陣列(最優引數)Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-31"""def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n= np.shape(dataMatrix)
           #返回dataMatrix的大小m為行數,n為列數。
    weights = np.ones(n)       #引數初始化
    for j in range(numIter):
        dataIndex = list(range(m))
        for i in range(m):            alpha = 4/(1.0+j+i)+0.01                #降低alpha的大小，每次減小1/(j+i)。            randIndex = int(random.uniform(0,len(dataIndex))) 
            #隨機選取樣本
            h = sigmoid(sum(dataMatrix[randIndex]*weights)) 
            #選擇隨機選取的一個樣本，計算h
            error = classLabels[randIndex] - h #計算誤差
            weights=weights+alpha*error*dataMatrix[randIndex] 
             #更新迴歸係數
            del(dataIndex[randIndex])#刪除已經使用的樣本    return weights    #返回if __name__ == '__main__':    dataMat, labelMat = loadDataSet()
    weights = stocGradAscent1(np.array(dataMat), labelMat)
    plotBestFit(weights)

程式碼執行結果：

2、迴歸係數與迭代次數的關係

可以看到分類效果也是不錯的。不過，從這個分類結果中，我們不好看出迭代次數和迴歸係數的關係，也就不能直觀的看到每個迴歸方法的收斂情況。因此，我們編寫程式，繪製出迴歸係數和迭代次數的關係曲線：

# -*- coding:UTF-8 -*-from matplotlib.font_manager import FontPropertiesimport matplotlib.pyplot as pltimport numpy as npimport random"""函式說明:載入資料Parameters:    無Returns:    dataMat - 資料列表    labelMat - 標籤列表Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    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 CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    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() - 求得的權重陣列(最優引數)    weights_array - 每次更新的迴歸係數Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    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.01 
        #移動步長,也就是學習速率,控制更新的幅度。
    maxCycles = 500   #最大迭代次數
    weights = np.ones((n,1))
    weights_array = np.array([])
    for k in range(maxCycles):
        h = sigmoid(dataMatrix * weights)
                    #梯度上升向量化公式
        error = labelMat - h
        weights = weights + alpha 
                   * dataMatrix.transpose() * error
        weights_array = np.append(weights_array,weights)
    weights_array = weights_array.reshape(maxCycles,n)    return weights.getA(),weights_array 
                 #將矩陣轉換為陣列,並返回"""函式說明:改進的隨機梯度上升演算法Parameters:    dataMatrix - 資料陣列    classLabels - 資料標籤    numIter - 迭代次數Returns:    weights - 求得的迴歸係數陣列(最優引數)    weights_array - 每次更新的迴歸係數Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-31"""def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n= np.shape(dataMatrix)
       #返回dataMatrix的大小,m為行數,n為列數。
    weights = np.ones(n)    #引數初始化
    weights_array = np.array([]) 
               #儲存每次更新的迴歸係數
    for j in range(numIter):                                      
        dataIndex = list(range(m))
        for i in range(m):           
            alpha = 4/(1.0+j+i)+0.01
             #降低alpha的大小,每次減小1/(j+i)。
            randIndex =int(random.uniform(0,len(dataIndex)))
                  #隨機選取樣本
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
               #選擇隨機選取的一個樣本，計算h
            error = classLabels[randIndex] - h
                 #計算誤差
            weights = weights + alpha * error 
           * dataMatrix[randIndex] #更新迴歸係數
            weights_array = np.append(weights_array,
             weights,axis=0) #添加回歸係數到陣列中
            del(dataIndex[randIndex]) 
                 #刪除已經使用的樣本
    weights_array = weights_array.reshape(numIter*m,n) 
                 #改變維度    return weights,weights_array      #返回"""函式說明:繪製迴歸係數與迭代次數的關係Parameters:    weights_array1 - 迴歸係數陣列1    weights_array2 - 迴歸係數陣列2Returns:    無Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-30"""def plotWeights(weights_array1,weights_array2):    #設定漢字格式
    font =
    FontProperties(fname=r"c:windowsfontssimsun.ttc",size=14)
    #將fig畫布分隔成1行1列,不共享x軸和y軸
    #fig畫布的大小為（13,8）當nrow=3,nclos=2時
    #代表fig畫布被分為個區域
    #axs[0][0]表示第一行第一列

    fig, axs = 
      plt.subplots(nrows=3, ncols=2,sharex=False,
                  sharey=False, figsize=(20,10))
    x1 = np.arange(0, len(weights_array1), 1)
    #繪製w0與迭代次數的關係
    axs[0][0].plot(x1,weights_array1[:,0])
    axs0_title_text = 
       axs[0][0].set_title(u'梯度上升演算法：
         迴歸係數與迭代次數關係',FontProperties=font)
    axs0_ylabel_text = axs[0][0].
             set_ylabel(u'W0', FontProperties=font)
    plt.setp(axs0_title_text, 
                size=20,weight='bold',color='black') 
    plt.setp(axs0_ylabel_text,
                  size=20,weight='bold',color='black')
    #繪製w1與迭代次數的關係
    axs[1][0].plot(x1,weights_array1[:,1])
    axs1_ylabel_text = axs[1][0].
            set_ylabel(u'W1',FontProperties=font)
    plt.setp(axs1_ylabel_text,
                size=20,weight='bold',color='black')
                        #繪製w2與迭代次數的關係
    axs[2][0].plot(x1,weights_array1[:,2])
    axs2_xlabel_text = axs[2][0].
        set_xlabel(u'迭代次數',FontProperties=font)
    axs2_ylabel_text = axs[2][0].
         set_ylabel(u'W1', FontProperties=font)
    plt.setp(axs2_xlabel_text,
              size=20,weight='bold',color='black') 
    plt.setp(axs2_ylabel_text,
              size=20,weight='bold',color='black')


    x2 = np.arange(0, len(weights_array2), 1)
    #繪製w0與迭代次數的關係
    axs[0][1].plot(x2,weights_array2[:,0])
    axs0_title_text = 
       axs[0][1].set_title(u'改進的隨機梯度上升演算法：
         迴歸係數與迭代次數關係',FontProperties=font)
    axs0_ylabel_text = axs[0][1].
             set_ylabel(u'W0',FontProperties=font)
    plt.setp(axs0_title_text,
            size=20, weight='bold',color='black') 
    plt.setp(axs0_ylabel_text,
            size=20,weight='bold',color='black')
    #繪製w1與迭代次數的關係
    axs[1][1].plot(x2,weights_array2[:,1])
    axs1_ylabel_text = axs[1][1].
             set_ylabel(u'W1',FontProperties=font)
    plt.setp(axs1_ylabel_text,
            size=20,weight='bold',color='black')
    #繪製w2與迭代次數的關係
    axs[2][1].plot(x2,weights_array2[:,2])
    axs2_xlabel_text = axs[2][1].
       set_xlabel(u'迭代次數', FontProperties=font)
    axs2_ylabel_text = axs[2][1].
          set_ylabel(u'W1', FontProperties=font)
    plt.setp(axs2_xlabel_text,
            size=20,weight='bold',color='black') 
    plt.setp(axs2_ylabel_text,
            size=20,weight='bold',color='black')

    plt.show()

if __name__ == '__main__':
    dataMat, labelMat = loadDataSet()           
    weights1,weights_array1 =
       stocGradAscent1(np.array(dataMat), labelMat)
    weights2,weights_array2 = 
            gradAscent(dataMat, labelMat)
    plotWeights(weights_array1, weights_array2)

執行結果如下：

由於改進的隨機梯度上升演算法，隨機選取樣本點，所以每次的執行結果是不同的。但是大體趨勢是一樣的。我們改進的隨機梯度上升演算法收斂效果更好。為什麼這麼說呢？讓我們分析一下。我們一共有100個樣本點，改進的隨機梯度上升演算法迭代次數為150。而上圖顯示15000次迭代次數的原因是，使用一次樣本就更新一下回歸係數。因此，迭代150次，相當於更新迴歸係數150*100=15000次。簡而言之，迭代150次，更新1.5萬次迴歸引數。從上圖左側的改進隨機梯度上升演算法迴歸效果中可以看出，其實在更新2000次迴歸係數的時候，已經收斂了。相當於遍歷整個資料集20次的時候，迴歸係數已收斂。訓練已完成。

再讓我們看看上圖右側的梯度上升演算法迴歸效果，梯度上升演算法每次更新迴歸係數都要遍歷整個資料集。從圖中可以看出，當迭代次數為300多次的時候，迴歸係數才收斂。湊個整，就當它在遍歷整個資料集300次的時候已經收斂好了。

沒有對比就沒有傷害，改進的隨機梯度上升演算法，在遍歷資料集的第20次開始收斂。而梯度上升演算法，在遍歷資料集的第300次才開始收斂。想像一下，大量資料的情況下，誰更牛逼？

本系列篇章：

Logistic迴歸實戰篇之預測病馬死亡率(一)

Logistic迴歸實戰篇之預測病馬死亡率(二)

Logistic迴歸實戰篇之預測病馬死亡率(三)