[ ML In Action ] Logistic Regression

Preface :
以下為 Machine Learning in Action 第四章內容, 相關代碼可以在 這裡下載 (Ch05).
這一章要提的 Machine Learning 是跟 Optimization algorithms 相關的 Logistic Regression. 而底下是有關這個演算法的描述 :

Finding the best fit is similar to regression, and in this method it's how we train our classifier. We'll use optimization algorithms to find these best-fit parameters. This best-fit stuff is where the name regression comes from.

簡單來說就是我們希望找出一條方程式構成的 hyperplane , 並透過它將標示於座標上不同 class 的 data point group 可以區隔開來. 底下為此 Classification 的簡介 :

General approach to logistic regression
1. Collect: Any method
2. Prepare: Numeric values are needed for a distance calculation. A structured data format is best.
3. Analyze: Any method
4. Train: We will spend most of the time training, where we try to find optional coefficients to classify our data.
5. Test: Classification is quick and easy once the training step is done.
6. Use:
This application needs to get some input data and output structured numeric values. Next, the application applies the simple regress calculation on this input data and determines which class the input data should belong to.

Classification with logistic regression :
在開始了解 logistic regression 先要來知道一些常識. 底下為使用 logistic regression 的說明 :

Logistic regression
Pros: Computationally inexpensive, easy to implement, knowledge representation easy to interpret.
Cons: Prone to underfitting, may have low accuracy
Works with: Numeric values, nominal values.

在 logistic regression 中, 我們透過一個函式稱為 sigmoid 來幫我們界定某個 data 的類別, 而該函式的長相如下 :



而上述方程式的 t 實際上為 logistic regression 中 features 值的合. 透過 sigmoid 計算出來的值如果大於 0.5 則歸類為類別 C1 ; 反之則歸類為類別 C2. 而 features 通常會用下面的方程式代表 :


而 x0 到 xn 的值為 input data ; 而我們的目的即是找出最佳的 coefficients w0~wn, 讓計算出來的 t 值套入sigmoid 得到的類別與預期的類別能有最大的相似度.

- Gradient ascent
第一個要接觸的 optimization algorithm 是 Gradient ascent. 詳細的演算法過程不在這裡介紹, 而這邊只說明概念. 考慮下面的式子 :


這個函式的概念就是透過不斷的計算預期值與運算值的差(預期值 - 運算值) :

1. 如果是負值 (預期值 < 運算值) 則原本的 weight 將會扣掉這個差與設定的 α 乘積, 如此 weight 將變小, 則下次計算 sigmoid() 的運算值將會變小 (更接近預期值) 而逐漸收斂.
2. 如果是正值 (預期值 > 運算值) 則原本的 weight 將會加上這個差與設定的 α 乘積, 如此 weight 將變大, 則下次計算 sigmoid() 的運算值將會變大 (更接近預期值) 而逐漸收斂.

- Train: Using gradient ascent to find the best parameters
首先我們先來撰寫載入 dataset , 計算 sigmoid 與 plot dataset 的函數 :
- Listing 5.1 (logRegres.py)

view plaincopy to clipboardprint?

1. # -*- coding: utf-8 -*-  
2. from numpy import *  
3.   
4. def plotDataSet(dataArr, labelArr):  
5.     import matplotlib.pyplot as plt  
6.     xcord1 = []; ycord1 = [];  
7.     xcord2 = []; ycord2 = [];  
8.     dataSize = len(dataArr)  
9.     for i in range(dataSize):  
10.         if labelArr[i] == 1:  
11.             xcord1.append(dataArr[i][1]); ycord1.append(dataArr[i][2])  
12.         else:  
13.             xcord2.append(dataArr[i][1]); ycord2.append(dataArr[i][2])  
14.     fig = plt.figure()  
15.     ax = fig.add_subplot(111)  
16.     ax.scatter(xcord1, ycord1, s=30, c='red', marker='s')  
17.     ax.scatter(xcord2, ycord2, s=30, c='green')  
18.     plt.show()  
19.   
20. def loadDataSet():  
21.     dataMat = []; labelMat = []  
22.     fr = open('testSet.txt')  
23.     for line in fr.readlines():  
24.         lineArr = line.strip().split()  
25.         dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])  
26.         labelMat.append(int(lineArr[2]))  
27.     return dataMat,labelMat  
28.   
29. def sigmoid(inX):  
30.     return 1.0/(1+exp(-inX))  

接著可以如下畫出 data set 的分布 :

>>> import logRegres
>>> dataArr, labelMat = logRegres.loadDataSet()
>>> logRegres.plotDataSet(dataArr, labelMat)

執行後會出現下圖 :


接著是 Gradient ascent algorithm 的實作 :

view plaincopy to clipboardprint?

1. def gradAscent(dataMatIn, classLabels, alpha=0.001, maxCycles=500):  
2.     dataMatrix = mat(dataMatIn)             # 轉 list 為 matrix  
3.     labelMat = mat(classLabels).transpose() # 轉 class label list 為 matrix, 並做 transpose  
4.     m,n = shape(dataMatrix)                 # m 為 data set size ; n 為 features size  
5.     weights = ones((n,1))                   # 初始 weight matrix 值為 1.  
6.     for k in range(maxCycles):              # for loop 進行 Gradient ascent algorithm 的計算.  
7.         h = sigmoid(dataMatrix*weights)     # 計算 t=w0x0 + w1x1 + ... + wnxn. 並將 t 帶入 sigmoid.  
8.         error = (labelMat - h)              # 計算 error > exp=1,0 - trt=1,0 = 0 ;  
9.                                             #              exp=1 - trt=0 = 1 (加大 weight) ; exp=0 - trt=1 = -1 縮小 weight)  
10.         weights = weights + alpha * dataMatrix.transpose()* error # 重新計算 weight matrix  
11.     return weights  

接著可以如下測試 :

>>> import logRegres
>>> dataArr, labelsArr = logRegres.loadDataSet()
>>> logRegres.gradAscent(dataArr, labelsArr)
matrix([[ 4.12414349],
[ 0.48007329],
[-0.6168482 ]])

接著有了 weight matrix, 我們便可以計算 sigmoid(), 因為最佳界定類別 C1 與 類別 C2 在 t=0 的位置上故 :

(x1->x ; x2->y)

我們找到分開類別 C1 與 類別 C2 的 hyperplane 的方程式, 接著便可以透過下面函數畫出 logistic regression best-fin line :

view plaincopy to clipboardprint?

1. def plotBestFit(weights):  
2.     import matplotlib.pyplot as plt  
3.     dataMat,labelMat=loadDataSet()  
4.     dataArr = array(dataMat)  
5.     n = shape(dataArr)[0]  
6.     xcord1 = []; ycord1 = []  
7.     xcord2 = []; ycord2 = []  
8.     for i in range(n):  
9.         if int(labelMat[i])== 1:  
10.             xcord1.append(dataArr[i,1]); ycord1.append(dataArr[i,2])  
11.         else:  
12.             xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])  
13.     fig = plt.figure()  
14.     ax = fig.add_subplot(111)  
15.     ax.scatter(xcord1, ycord1, s=30, c='red', marker='s')  
16.     ax.scatter(xcord2, ycord2, s=30, c='green')  
17.     x = arange(-3.0, 3.0, 0.1)  
18.     y = (-weights[0]-weights[1]*x)/weights[2]  
19.     ax.plot(x, y)  
20.     plt.xlabel('X1'); plt.ylabel('X2');  
21.     plt.show()  

接著可以如下測試

>>> dataArr, labelsArr = logRegres.loadDataSet()
>>> weights = logRegres.gradAscent(dataArr, labelsArr)
>>> weights.getA()
array([[ 4.12414349],
[ 0.48007329],
[-0.6168482 ]])
>>> logRegres.plotBestFit(weights.getA())

產生座標圖如下 :


- Train: stochastic gradient ascent
在使用 Gradient ascent algorithm 可以發現 loop 中的 data set 的每一筆 record 都會被拿來計算 sigmoid(). 這裡的 data set 只有 100 筆所以感覺不出來影響, 在處理大量的 data set (billions), 就會有 performance concern. 因此有另一個替代的 Stochastic gradient ascent algorithm, 一次只拿 data set 的一筆 record 來重新計算 weight matrix, 而透過它當有新的 data set 進來, 我們可以利用已經算好的 weight matrix 繼續往下運算而不用重新計算舊的 data set. 底下為其 pseudo code :


接著底下為其代碼實作 :

view plaincopy to clipboardprint?

1. def stocGradAscent0(dataMatrix, classLabels):  
2.     m,n = shape(dataMatrix)  
3.     alpha = 0.01  
4.     weights = ones(n)   #initialize to all ones  
5.     for i in range(m):  
6.         h = sigmoid(sum(dataMatrix[i]*weights))  
7.         error = classLabels[i] - h  
8.         weights = weights + alpha * error * dataMatrix[i]  
9.     return weights  

但 α(alpha) 值應該隨著 loop 過程中而變小, 因為 error 會逐漸收斂, 故底下為優化過後的代碼 :

view plaincopy to clipboardprint?

1. def stocGradAscent1(dataMatrix, classLabels, numIter=150):  
2.     m,n = shape(dataMatrix)  
3.     weights = ones(n)   #initialize to all ones  
4.     for j in range(numIter):  
5.         dataIndex = range(m)  
6.         for i in range(m):  
7.             alpha = 4/(1.0+j+i)+0.0001    #apha decreases with iteration, does not  
8.             randIndex = int(random.uniform(0,len(dataIndex)))#go to 0 because of the constant  
9.             h = sigmoid(sum(dataMatrix[randIndex]*weights))  
10.             error = classLabels[randIndex] - h  
11.             weights = weights + alpha * error * dataMatrix[randIndex]  
12.             del(dataIndex[randIndex])  
13.     return weights