機器學習——Logistic迴歸

參考《機器學習實戰》

利用Logistic迴歸進行分類的主要思想：

根據現有數據對分類邊界線創建迴歸公式，以此進行分類。

分類藉助的Sigmoid函數：

Sigmoid函數圖：

Sigmoid函數的做用：

將全部特徵都乘上一個迴歸係數，而後將全部結果值相加，將這個總和代入Sigmoid函數中，進而獲得一個0～1之間的數值。任何大於0.5的數據被分1類，小於0.5分入0類。

綜上，Sigmoid的輸入能夠記爲z：

因此向量w便是咱們要經過最優化方法找的係數。

w向量的求解：

1）、梯度上升法（思想：沿函數梯度方向找函數最大值）

梯度上升法僞代碼：

更新w係數的細節實現代碼：

要注意的是，做者在這裏dataMatrix的特徵矢量維度比實際特徵矢量多了一維。做者使用的數據是二維[x1,x2]，而程序中增長了一維[x0=1,x1,x2].奇怪的是x0加在了最前面的位置，而不是最後的位置。此外，上圖中畫紅線處的公式做者未給出由來，網上搜索了下，找到一篇博文，寫得還不錯。這裏帖上點簡要概述：

具體過程以下：（參考:http://blog.csdn.net/yangliuy/article/details/18504921?reload）

參數機率方程：

其中x爲訓練特徵，y爲x對應的類，θ爲待估計參數

利用上式中y只取0或1的特色，恰好能夠表示爲:

似然函數：（這裏就是Logistic Regression的目標函數，原書中並未指明，因此若是不網上找logistic的資料區先前學過機器學習，很沒法理解書中的作法的）

對數似然函數：

因此極大似然估計：

從而獲得梯度上升法的遞歸公式：

這裏也就是上面的圖中，我畫紅線處公式的由來了。

這裏再上傳下本身寫的代碼（未優化的logistic算法），代碼中的數據來源還是《機器學習實戰》一書提供的數據：

#-*- coding:cp936 -*-
import numpy as np
import matplotlib.pyplot as plt

class Log_REG():
    def __init__(self):
        self._closed=False
        
    def loadData(self, dataFile='testSet.txt'):
        f_file = open(dataFile)
        lines = f_file.readlines()
        line_data = lines[0].strip().split()
        self.num_feature = len(line_data) - 1
        self.xData = np.zeros((len(lines), self.num_feature + 1))
        self.label = np.zeros((len(lines), 1))
        self.num_label = len(lines)
        line_cnt = 0
        for iLine in lines:
            line_data = iLine.strip().split()  
            for i in range(self.num_feature):
                self.xData[line_cnt][i] = float(line_data[i])
            self.xData[line_cnt][self.num_feature] = 1
            self.label[line_cnt] = float(line_data[-1])
            line_cnt+=1
    
    def _sigmoid(self, z):
        return 1.0 / (1 + np.exp(-z))
    

    def gradAscendClass(self):
        maxIter = 500
        self.omiga = np.ones((1, self.num_feature+1))
        xDataMat = np.matrix(self.xData)
        alpha = 0.01
        self.omiga_record=[]
        for i in range(maxIter):
            h = self._sigmoid(self.omiga * xDataMat.transpose())  # 矩陣乘
            error = self.label - h.transpose()
            self.omiga = self.omiga + alpha * (xDataMat.transpose()*error).transpose() 
            self.omiga_record.append(self.omiga)
            if np.sum(np.abs(error)) < self.num_label * 0.05:
                print  "error very low",i
                break
    
    def stochasticGradAscend(self):
        pass
#         maxIter = 150
#         self.omiga = np.ones((1,self.num_feature+1))
#         for 


    def plotResult(self):
        self._close()
        if self.num_feature != 2:
            print "Only plot data with 2 features!"
            return
        label0x = []
        label0y = []
        label1x = []
        label1y = []
        for i in range(self.num_label):
            if int(self.label[i]) == 1:
                label1x.append(self.xData[i][0])
                label1y.append(self.xData[i][1])
            else:
                label0x.append(self.xData[i][0])
                label0y.append(self.xData[i][1])
        fig = plt.figure()
        ax = fig.add_subplot(111)
        ax.scatter(label0x, label0y, c='b',marker='o')
        ax.scatter(label1x, label1y, c='r',marker='s')
        
        minx = min(min(label0x),min(label1x))
        maxx = max(max(label0x),max(label1x))
        wx = np.arange(minx,maxx,0.1)
        wy = (-self.omiga[0,2]-self.omiga[0,0]*wx)/self.omiga[0,1] 
        ax.plot(wx,wy)
        
    
    def plotIteration(self):
        self._close()
        iterTimes = len(self.omiga_record)
        w0=[i[0][0,0] for i in self.omiga_record]
        w1=[i[0][0,1] for i in self.omiga_record]
        w2=[i[0][0,2] for i in self.omiga_record]
        fig = plt.figure()
        ax1 = fig.add_subplot(3,1,1)
        ax1.plot(range(iterTimes),w0,c='b')#,marker='*')
        plt.xlabel('w0')
        ax2 = fig.add_subplot(3,1,2)
        ax2.plot(range(iterTimes),w1,c='r')#,marker='s')
        plt.xlabel('w1')
        ax3 = fig.add_subplot(3,1,3)
        ax3.plot(range(iterTimes),w2,c='g')#,marker='o')
        plt.xlabel('w2')

    def show(self):
        plt.show()

    def _close(self):
        pass
            
        
        
if __name__ =='__main__':
    testclass = Log_REG()
    testclass.loadData()
    testclass.gradAscendClass()
    testclass.plotResult()
    
    testclass.plotIteration()
    testclass.show()

顯示結果：

分類結果

分類參數收斂結果

梯度上升（或降低）算法的改進：

當數據量很大時，上述梯度上升算法每次迭代都要對全部數據進行處理，會形成計算量異常龐大。解決的方法是引入隨機梯度的思想。

隨機梯度降低的基本原理是：不直接計算梯度的精確值，而是用梯度的無偏估計g(w)來代替梯度：

實際操做時，隨機地選取單個數據而非整個數據集參與迭代，詳細的原理推導可參見：http://www.52ml.net/2024.html

改進的隨機梯度上升法：

def stochasticGradAscend2(self):
        maxIter = 150
        self.omiga = np.ones((1,self.num_feature+1))
        self.omiga_record=[]
        
        for j in range(maxIter):
            randRange = range(self.xData.shape[0])
            for i in range(self.xData.shape[0]):
                alpha = 4/(1.0+i+j)+0.01
                randIndex  = int(random.uniform(0,len(randRange)-1))
                index = randRange[randIndex]
                h = self._sigmoid(np.matrix(self.omiga)[0]*np.matrix(self.xData[index,:]).transpose())
                error = self.label[index]-h
                
                self.omiga  = self.omiga+alpha*error*self.xData[index,:]
                self.omiga_record.append(np.matrix(self.omiga))
                del(randRange[randIndex])

從上圖能夠看出，改進後的隨機梯度算法收斂很快，上圖只對全部數據作150次迭代。

參考文獻：

http://blog.csdn.net/yangliuy/article/details/18504921?reload

隨機梯度：http://www.52ml.net/2024.html