機器學習算法·回歸預測

一、回歸預測簡介

現在我們知道的回歸一詞最早是由達爾文的表兄弟Francis Galton發明的。Galton在根據上一年的豌豆種子的尺寸預測下一代豌豆種子的尺寸時首次使用了回歸預測。他在大量的對象上應用了回歸分析，包括人的身高。他注意到，如果雙親的高度比平均高度高的話，則他們的子女也傾向於比平均身高高，但尚不及雙親。孩子的身高向着平均高度回退（回歸）。Galton在多項研究上都注意到這個現象，所以盡管這個單詞跟數值預測沒有任何關系，但是仍然把這種方法稱為回歸。

二、線性回歸

要了解回歸方法，我們可以從二維平面上的回歸開始先從直觀上理解回歸的含義。例如平面中有一系列的點(x1,y1),(x2,y2)...(xn,yn)，如圖1所示。我們可以很直觀的發現這些數據有明顯的線性關系，回歸的目的就是找到這樣一條直線，使得所有數據點盡可能的與直線接近。這條直線的方程稱為回歸方程。這樣我們就能夠根據直線的回歸方程進行預測目標值的預測了。

 

圖1

當推廣到多維空間的時候同樣可以采用相同的方法，只不過找出來的方程就是一個超平面了，這與我們在支持向量機中所說的超平面是相同的。下面我們將介紹如何求出這樣的平面。

2.1線性回歸的數學理論

對於含有m維特征向量的數據X(x1,x2,x3...xn)，我們可以把回歸模型寫成如下的形式：

 

其中x1,x2...xn為n個自變量，ε為誤差項，且服從正態分布N（0,0^2）。

b1,b2...bm為這m個自變量的系數。為方便起見，我們對m組實際觀測到的樣本數據引入矩陣記號：

  

其中X為模型設計矩陣，Y和ε是隨機向量，並且有：

（En為n階單位矩陣）

所以包含n個樣本的回歸模型式可以寫成：

 

     

ε是補不可觀測到的隨機誤差向量，B是回歸系數構成的向量，系數需要待定。

2.2利用最小二乘法估計回歸參數

設β_i’為β_i(i=1,2...m)的估計值，則當β_i為β_i^’時的誤差平方和為：

 

根據極限定理，當Q對B中各參數的偏導數為0時Q取得最小即：

       

求解上式得：

 

經過整理之后得到正規方程組為：

          

由線性代數的知識可知當X為滿秩的時候，上式表達式的解為：

2.3 python代碼的實現

首先需要定義一個函數standRegress()函數，該函數輸入訓練集的特征矩陣和標簽，返回權重值。如下所示

def standRegres(xArr,yArr):#標准的回歸預測
    xMat = mat(xArr); yMat = mat(yArr).T#將數組轉化為矩陣的形式
    xTx = xMat.T*xMat
    if linalg.det(xTx) == 0.0:#判斷行列式是否為0
        print("This matrix is singular, cannot do inverse")
        return
    ws = xTx.I * (xMat.T*yMat)#計算權重值
    return ws

  為了直觀上的理解，我們仍然采用單變量數據來對回歸預測進行說明。數據集為ex0.txt。改文件中每個數據都只包含了一個特征變量和一個目標值。利用上面的函數我們可以求得回歸方程的系數，進而利用matlplotlib庫中的繪圖功能，將回歸直線和訓練集的散點圖在同一個圖像中繪制出來，結果如圖2所示。

 圖2 回歸直線

三、局部線性加權

    從上面的實例來看我們發現所有數據點都較為均勻的分布在直線兩側，但是從直觀上來看線性回歸並不能很好的擬合數據點。這樣模型預測出來的數據就有很大的偏差，不能夠得到很好的預測結果。線性回歸之所以出現這樣的問題是因為在擬合過程中求得是最小均方誤差的無偏估計，但是這樣容易導致過擬合現象。而有些方法允許在估計中出現一些偏差，從而降低預測的均方誤差。局部加權線性回歸就是這樣一個方法。

局部線性加權算法的思想就是，給待預測點附近的每一個點賦予一定的權重。然后在這個子集上基於最小均方差進行普通的回歸。因此將上述線性回歸的誤差可以改成：

 

 

我們對上式求關於B的偏導數，並令其為0。最后計算出的結果為：

 

其中W是一個權重矩陣，用於給每個數據點賦予權重。在該算法中，使用“核”來對附近點賦予權重。核的類型可以自由選擇，但最常用的就是高斯核，高斯核對應的權重為：

 

這樣就夠簡稱了一個只含有對角元素的矩陣B，並且點x與點x(i)越近，w(i,i)就會越大。高斯核函數只有一個待確定的參數k，k值的大小決定了對附近點賦予的權重大小。在x-w的曲線中，k值越大，該圖像就越平；反之則圖像月陡峭。我們預測點0.5處的附近點在不同k值下的x-w曲線如圖3所示。

 

圖3 不同k值下權重與距離的關系圖

現在我們再利用和原來一樣的數據集，來觀察局部線性加權算法的效果。分別令k=0.003,0.01,1.0，結果如圖4所示。

 

	k=1                      k=0.01                    k=0.03

圖4 不同k條件下的局部加權算法結果

    加權線性回歸由於在預測每一個點的時候都需要重新確定權重，因此計算量非常龐大，這也是加權線性回歸預測的一個缺點。因此在實際應用中，我們需要控制好參數k的大小，以免出現欠擬合和過擬合的現象。K值的合適的大小一般跟數據集本身有關，需要根據經驗來進行選擇。

代碼實現：其中定義兩個函數lwlr()和lwlrTest()前者用於計算局部線性回歸的權重值，后者用於預測給定點預測值。

def lwlr(testPoint,xArr,yArr,k=1.0):#局部加權線性回歸
    xMat = mat(xArr); yMat = mat(yArr).T
    m = shape(xMat)[0]#表示x中樣本的個數
    weights = mat(eye((m)))#初始化為單位矩陣
    for j in range(m):                      #下面兩行用於創建權重矩陣
        diffMat = testPoint - xMat[j,:]#樣本點與待預測點的距離
        weights[j,j] = exp(diffMat*diffMat.T/(-2.0*k**2))#使用的是高斯核
    xTx = xMat.T * (weights * xMat)
    if linalg.det(xTx) == 0.0:#判斷行列式是否為滿秩
        print("This matrix is singular, cannot do inverse")
        return
    ws = xTx.I * (xMat.T * (weights * yMat))
    return testPoint * ws

def lwlrTest(testArr,xArr,yArr,k=1.0):  #利用局部加權回歸預測出待測點的值
    m = shape(testArr)[0]
    yHat = zeros(m)
    for i in range(m):
        yHat[i] = lwlr(testArr[i],xArr,yArr,k)
    return yHat

四、嶺回歸

    當數據的特征比樣本點的個數還多的時候，不能夠再使用之前的線性回歸和局部加權線性回歸的方法。特征個數比樣本點還多(n>m)，也就是說輸入數據的矩陣X不是列滿秩矩陣，而非滿秩的矩陣在計算(X^TX)^-1的時候會出現錯誤。

    例如：對於矩陣X=[1,2,3;4,5,6]。它不是一個列滿秩矩陣，我們計算一下X^TX的逆，計算結果為(X^TX)^-1=0。也就是說它的行列式正好為零。在一般情況下都有相同的結果。可以證明r(AA’)=r(A),因此當A不是滿秩的時候AA’也不是滿秩，因此行列式為0.

       為了解決上述出現的問題，統計學家提出了嶺回歸的概念。簡單來說嶺回歸就是在X^TX上加上一個λE（其中E為m*m的單位矩陣）。在這種情況下，回歸系數的計算公式會變為：

 

      這里需要指明一下嶺回歸中的“嶺”的含義。嶺回歸中使用了單位矩陣乘於常量λ，我們觀察其中的單位矩陣E可以看到值1貫穿整個對角線，其余元素都是0，在0構成的平面上有一條1組成的“嶺”就是嶺回歸中的“嶺”的來由。

     我們將仍然利用上面的數據集，對嶺回歸的統計特性進行檢驗。下面是嶺回歸的python的嶺回歸實現代碼：

def ridgeRegres(xMat,yMat,lam=0.2):#用於計算回歸系數
    xTx = xMat.T*xMat
    denom = xTx + eye(shape(xMat)[1])*lam
    if linalg.det(denom) == 0.0:
        print("This matrix is singular, cannot do inverse")
        return
    ws = denom.I * (xMat.T*yMat)#I表示逆
    return ws

def ridgeTest(xArr,yArr):
    xMat = mat(xArr); yMat=mat(yArr).T
    yMean = mean(yMat,0)   #計算y的均值
    yMat=yMat - yMean     #y變成與均值之間的距離
    #數據的標准化
    xMeans = mean(xMat,0)   #計算x的均值
    xVar = var(xMat,0)      #計算x的方差
    xMat = (xMat - xMeans)/xVar#標准化后的x
    numTestPts = 30#計算30個不同lambda下的權重矩陣
    wMat = zeros((numTestPts,shape(xMat)[1]))
    for i in range(numTestPts):
        ws = ridgeRegres(xMat,yMat,exp(i-10))
        wMat[i,:]=ws.T
    return wMat

 

    其中第一函數ridgeRes()實現了對給定的lambda下的嶺回歸求解，默認為lambda為0.2。另外為了使用嶺回歸和縮減技術，首先需要對特征進行標准化處理。第二個函數ridgeTest()中首先對數據進行了標准化處理，使得每一維的特征具有相同的重要性。然后可以輸出30個不同lambda下的權重矩陣。然后繪制出嶺回歸圖如圖5所示。對於嶺回歸參數的確定我們可以采用嶺跡法，即在嶺跡圖中個系數都比較穩定的地方所取的lambda值。也可以使用GCV廣義交叉方法。

 

圖5 嶺回歸系數變化圖

 

    嶺回歸最先用來處理特征數多於樣本數的情況，現在也用於在估計中加入偏差，從而得到更好的估計。這里引入了λ來限制所有w之和，通過引入該懲罰項，能夠減少不重要的參數，這種技術就做縮減。嶺回歸是對最小二乘回歸的一個補充，它損失了無偏性，但換回了高的數值穩定性。

 

參考文獻：

【1】機器學習實戰

【1】百度百科.嶺回歸

【2】百度文庫.嶺回歸和lasso

https://wenku.baidu.com/view/a1d929bd50e2524de4187e2a.html