Re:從零開始的機器學習 - Machine Learning(一) 線性迴歸

樑王發表於2017-12-21

從我對整個職業生涯的規劃出發，我不僅想做一些高質量的應用(軟體工程的角度),還想做一些激動人心的應用，所以我希望能在機器學習的方向走，儘管我在大學粗淺的學了些皮毛，但如果要把機器學習作為職業發展的話這些還差得遠，所以我開始寫了這個系列的文章。

我希望通過這個系列能對機器學習領域的知識點做一個總結，所以對於Machine Learning這個部分，我的目標是寫出能讓高中生看得懂。

前言

這篇文章的主角是線性迴歸，也就是LR(Linear Regression)。不過高中生肯定不知道迴歸是什麼吧？我現在前言裡面簡單介紹一下。

迴歸(Regression)

迴歸(Regression)問題是機器學習裡面很大的一塊。

統計學中，迴歸分析（regression analysis)指的是確定兩種或兩種以上變數間相互依賴的定量關係的一種統計分析方法。迴歸分析按照涉及的變數的多少，分為一元迴歸和多元迴歸分析；按照因變數的多少，可分為簡單迴歸分析和多重回歸分析；按照自變數和因變數之間的關係型別，可分為線性迴歸分析和非線性迴歸分析。

不玩定義，直接上例子，設想一下這樣的場景：如果我現在有了一組(身高, 體重)的資料集，可以令身高為x值，體重為y值將其畫在紙上，就和下圖一樣。

那麼現在問題來了，假設我身高為1.8m，能不能預測出我的體重？這就是迴歸要解決的問題，我希望得到一個函式，我輸入我的身高之後他能夠幫我預測出我的體重。

在這個例子裡，我假設身高和體重的關係是線性關係，所以假設模型是y=kx+b這種一次函式，為了確定變數k和b，我需要利用之前的資料去學習出這兩個引數。

儘管上面的例子誤差很大。因為體重的函式肯定不單單由身高這一個特徵(feature)決定，而且即使只有身高一個特徵，最佳函式也很可能不是一次函式。

迴歸問題相關概念

通常，這類預測問題可以用迴歸(regression)模型進行解決，迴歸模型定義了輸入與輸出的關係，輸入即現有知識，而輸出則為預測。

一個預測問題在迴歸模型下的解決步驟為：

積累知識：我們將儲備的知識稱之為訓練集(Training Set)，畢竟預測是需要過去的資料作為參考的，這個很好理解。
學習：當我們有了資料的時候，我們就需要去學習這些資料，為什麼說機器學習智慧？因為當我告訴計算機我的模型是線性(一次函式)或其他型別的函式然後把資料直接倒進去之後它就能返回給我最後的函式(引數都訓練好了)。
預測：學習完成後，當接受了新的資料（輸入）後，我們就能通過學習階段獲得的對應關係來預測輸出。

kaggle上有個給萌新練習的比賽(competition)，關於泰坦尼克號的。

主要內容就是給出幾千個人的個人資訊(性別、年齡，船艙登記，登船口等資訊)以及他們是否存活的資料，然後給出一些測試資料，即上面說的個人資訊，讓你去預測他們是否會存活下來。

有興趣可以瞭解一下: kaggle-titanic

Coursera - 史丹佛 Machine Learning

正文

這篇部落格主要講的是線性迴歸(Linear regression)，經過前言後大家也知道了，迴歸裡面用到的函式有多種多樣的，這個需要開發者自己去選擇，這次先介紹最簡單的線性迴歸(Linear regression)。

線性迴歸LR(Linear Regression)

從數學上來說，給定由d個屬性描述的示例x=(x1;x2;...;xd),其中xi是x在第i個屬性上的取值，線性模型試圖學得一個通過屬性的線性組合來進行預測的函式，即

一般用向量形式

什麼，你說你是高中生不知道這個T是什麼？去看看矩陣或向量的轉置吧。當然你還得去看看矩陣的乘法，畢竟這裡有個1行n列的矩陣和一個n行一列的矩陣相乘

好了，現在我們知道了線性迴歸模型的基本形式了，那麼我們現在的任務就是要去學習出這個w向量和b這些引數的值，有了這個值我們才能去做預測。

一般來說，我們會給w和b一個初始值，然後我們去修正這些值來讓這些值符合預期，那麼我們怎麼去修正這些值呢？我們需要一個損失函式(loss function)，這個函式用來指明我的預測值和訓練資料實際值的差別。

那麼這個損失函式(loss function)是什麼，以及我們怎麼用它來修正我們的引數w和b，看下面。

梯度下降法

這裡我只講一個梯度下降法，如果以後有需要我會再回來補充。

不知道大家知不知道啟發式搜尋(Heuristically Search)？

啟發式搜尋(Heuristically Search)又稱為有資訊搜尋(Informed Search)，它是利用問題擁有的啟發資訊來引導搜尋，達到減少搜尋範圍、降低問題複雜度的目的，這種利用啟發資訊的搜尋過程稱為啟發式搜尋。——百度百科

舉個例子，大一大二在大學裡面的時候我是有單車的，晚上下課的時候在停車場裡面經常忘記自己的車停在哪，要找半天才能找得到，這種就叫盲目搜尋，廣度優先搜尋(BFS)和寬度優先搜尋(DFS)都是盲目搜尋。

那麼如果我手裡有個黑科技，能夠告訴我我現在距離我的愛車多遠，那麼我是不是可以以此為依據，往距離逐漸減少的方向去搜尋呢？這就是啟發式搜尋(Heuristically Search)，Astar(A*)尋路演算法等就是啟發式搜尋。

啟發式搜尋與機器學習有一些概念是共通的，那麼對於機器學習來說，我也同樣需要一個和“距離我愛車多遠”這樣一個指標來判斷我現在引數距離我最優引數有“多遠”。我們可以統一的把這種東西叫做損失函式(loss function)

損失函式(Loss Function)

我們給之前例子裡的函式一個名字，叫假設函式(hypothesis function)，意為預估值函式。損失函式則是用來衡量假設函式(hypothesis function)的準確性，具體衡量指標有很多，這裡我們和吳恩達教程裡面一樣採用平方差的方式計算。

其中函式J就是損失函式(為什麼是J我好好奇)，函式h是假設函式，可以看到後面是平方差的平均值然後除以了個2。

公式裡面的2其實就是為了後面求偏導的時候可以把分數消掉，沒什麼太大的實際意義。

梯度下降演算法(Gradient descent algorithm)

既然我們現在知道了當前怎麼評價當前引數的好壞，那麼我怎麼去修正引數讓引數更好(損失函式的最小值)呢？

高中生都知道，在一元函式裡面，導數的幾何意義就是函式變化增加最快的方向。梯度其實類似，也是類似的，說白了就是一個向量，表示上升最快的方向。

梯度、偏導數部分的補充大家可以自己去看高數書或者網上的一些資料。

那麼我們就可以得到一個修正的公式，我們迭代這個公式許多次來修正引數。

這裡是減號是因為梯度表示上升最快，所以用的是負梯度。

然後其中其中α表示學習速率(learning rate)，這個值越大每次修正的就越多，不過這個不是越高越好，如果太高了可能會一直在最低點“擺動”而無法收斂。也有的使用可變的學習速率，一開始設定較高，接近最低點的過程中逐漸降低。

接下來我們看一下求導之後的結果

這裡有一點值得注意的是，在這個梯度下降中我們每次迭代都使用了所有的m個訓練資料，這種也叫批量梯度下降法Batch Gradient Descent，BGD

這樣每次迭代都將對m個樣本進行計算，計算量大。所以有些優化方案，有興趣的可以去看一下

現在我們知道了如何去修正引數了，但我們實際上修正之後得到的是損失函式(loss function)的極小值而不一定是最小值

由於起始點(起始引數)的不同，最後可能得到的並不是全域性的最優解(損失函式最小)。我說一下從西瓜書上看到的幾個優化方法。

以多組不同的引數值初始化，反正就是找多個起點選最好的結果。
使用“模擬退火”(Simulated Annealing)技術，模擬退火在每一步都以一定概率接受比當前解更差的結果，從而有助於“跳出”區域性極小。在每步迭代過程中，接受“次優解”的概率會隨著時間的推移而逐步降低，從而保證演算法的穩定。
使用隨機梯度下降，它在計算梯度的時候加入了隨機因素，所以即使陷入了區域性極小值，它計算出來的梯度仍然可能不為0，這樣就有機會跳出區域性極小繼續搜尋。

實踐

說了那麼多理論，是時候寫一些程式碼了。我打算使用Python來做一下史丹佛Machine Learning課程裡面關於線性迴歸的練習。pdf以及資料都可以在我的GitHub庫上下載到

環境

如果你不想被配環境煩死的話，我真的推薦裝Anaconda，除此之外要說的就是我用的都是Python3.x。

背景

在這個練習中，我們要用簡單線性迴歸實現預測食物卡車的利潤。我們現在已經有了(城市人數，城市利潤)這樣的許多對資料，現在我要做的就是用線性迴歸模型並訓練出引數來預測我如果給另一個城市(城市人數)，那麼卡車的利潤是多少。

程式碼及註釋

# 參考http://www.johnwittenauer.net/machine-learning-exercises-in-python-part-1/
import os
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# 計算損失，用了向量化程式設計而不是for迴圈
def computeLoss(X, y, theta):  
    inner = np.power(((X * theta.T) - y), 2)
    return np.sum(inner) / (2 * len(X))

# 梯度下降部分
def gradientDescent(X, y, theta, alpha, iters):  
    temp = np.matrix(np.zeros(theta.shape))
    parameters = int(theta.ravel().shape[1])
    cost = np.zeros(iters)

    for i in range(iters):
        
        error = (X * theta.T) - y
        for j in range(parameters):
            term = np.multiply(error, X[:,j])
            temp[0,j] = theta[0,j] - ((alpha / len(X)) * np.sum(term))

        theta = temp
        cost[i] = computeLoss(X, y, theta)

    return theta, cost

# 讀入訓練資料
# windows使用者路徑可能需要修改下，後期有時間可能會做統一
def loadData(path):
    trainingData = pd.read_csv(path, header=None, names=['Population', 'Profit'])

    trainingData.head()

    trainingData.describe()

    trainingData.plot(kind='scatter', x='Population', y='Profit', figsize=(12,8))
    plt.show()
    return trainingData

trainingData = loadData(os.getcwd() + '/../data/ex1data1.txt')

# 在資料集前插入一列Ones作為常數係數，也就是y=k*x+b*1這種形式
trainingData.insert(0, 'Ones', 1)

# 將輸入X以及輸出y從資料集中分割
cols = trainingData.shape[1]  
X = trainingData.iloc[:,0:cols-1]  
y = trainingData.iloc[:,cols-1:cols]  

# 把pandas的DataFrames轉換成numpy的矩陣
X = np.matrix(X.values)  
y = np.matrix(y.values)  
# 初始化引數為全0的，當然也可以初始化成其他的
theta = np.matrix(np.array([0,0]))  

# 各向量的維度
X.shape, theta.shape, y.shape  

# 初始損失函式值
computeLoss(X, y, theta)   # 32.07，後面可以看看訓練完後的損失函式值

# 設定學習速率以及迭代次數
alpha = 0.01  
iters = 2000

# 使用梯度下降得到模型引數
theta_fin, loss = gradientDescent(X, y, theta, alpha, iters)  
theta_fin

# 計算訓練後的引數的損失值
computeLoss(X, y, theta_fin)  # 4.47

# 為了畫線用的，畫出訓練好後的直線
x = np.linspace(trainingData.Population.min(), trainingData.Population.max(), 100)  
f = theta_fin[0, 0] + (theta_fin[0, 1] * x)

fig, ax = plt.subplots(figsize=(12,8))  
ax.plot(x, f, 'r', label='Prediction')  
ax.scatter(trainingData.Population, trainingData.Profit, label='Traning Data')  
ax.legend(loc=2)  
ax.set_xlabel('Population')  
ax.set_ylabel('Profit')  
ax.set_title('Predicted Profit vs. Population Size')  
plt.show()

# 損失隨著迭代次數的變化
fig, ax = plt.subplots(figsize=(12,8))  
ax.plot(np.arange(iters), loss, 'r')  
ax.set_xlabel('Iterations')  
ax.set_ylabel('Loss')  
ax.set_title('Error vs. Training Epoch')  
plt.show()
複製程式碼