1、梯度下降-矩陣形式

上篇文章介紹了一元線性回歸，包括Python實作和sklearn實作的實體、對比，以及一些問題點，詳情可以看這里：
鏈接: 手寫演算法-Python代碼實作一元線性回歸

里面封裝的one_variable_linear()類只適用于一元線性回歸，
本篇文章修改代碼，推廣至多元線性回歸，并介紹2種更簡潔的方法，
先給大家復習一下矩陣的基本知識：





轉置矩陣：

損失函式可表示為：

可以求得：矩陣形式下，偏導的運算式是：


下面附上我的推導證明程序（剛寫的）：


有了上述運算式，我們修改上次的代碼如下：

class linear():
    def __init__(self):
        pass
    
    #梯度下降法迭代訓練模型引數,x為特征資料，y為標簽資料，a為學習率，epochs為迭代次數
    def fit(self,x,y,a,epochs):  
        #計算總資料量
        m=x.shape[0]
        #給x添加偏置項
        X = np.concatenate((np.ones((m,1)),x),axis=1)
        #計算總特征數
        n = X.shape[1]
        #初始化W的值,要變成矩陣形式
        W=np.mat(np.ones((n,1)))
        #X轉為矩陣形式
        xMat = np.mat(X)
        #y轉為矩陣形式，這步非常重要,且要是m x 1的維度格式
        yMat =np.mat(y.reshape(-1,1))
        #回圈epochs次
        for i in range(epochs):
            W=W-a*xMat.T*(xMat*W-yMat)
        return W
    def predict(self,x,w):  #這里的x也要加偏置，訓練時x是什么維度的資料，預測也應該保持一樣
        return np.dot(x,w)

依然用上次的測驗資料集，2個代碼比較如下：

import numpy as np
import pandas as pd
from sklearn import datasets   #sklearn生成資料集都在這里
from matplotlib import pyplot as plt


#生成一個特征的回歸資料集
x,y=datasets.make_regression(n_features=1,noise=15,random_state=2020)  
plt.scatter(x,y)
plt.show()




class one_variable_linear():
    #初始化引數，k為斜率，b為截距，a為學習率，n為迭代次數
    def __init__(self,k,b,a,n):
        self.k =k 
        self.b=b
        self.a=a
        self.n = n
     
    #梯度下降法迭代訓練模型引數
    def fit(self,x,y):
        #計算總資料量
        m=len(x)
        #回圈n次
        for i in range(self.n):
            b_grad=0
            k_grad=0
            #計算梯度的總和再求平均
            for j in range(m):
                b_grad += (1/m)*((self.k*x[j]+self.b)-y[j])
                k_grad += (1/m)*((self.k*x[j]+self.b)-y[j])*x[j]

            #更新k,b
            self.b=self.b-(self.a*b_grad)
            self.k=self.k-(self.a*k_grad)

            #每迭代10次，就輸出一次影像
            if i%10==0:
                print('迭代{0}'.format(i)+'次')
                plt.plot(x,y,'b.')
                plt.plot(x,self.k*x+self.b,'r')
                plt.show()
        self.params= {'k':self.k,'b':self.b}
        #輸出系數
        return self.params
    
    #預測函式
    def predict(self,x):
        y_pred =self.k * x + self.b
        return y_pred

lr=one_variable_linear(k=1,b=1,a=0.1,n=60)
lr.fit(x,y)

舊代碼得到的引數如上，
新代碼得到的引數如下：

model = linear()
w = model.fit(x,y,a=0.1,epochs=50)
print(w)

我的天，這是什么鬼，怎么和上面的差的這么多（其實是我故意的），明顯這個是不正確的，模型完全沒有收斂，問題在哪里？
我們細想一下，正常的梯度下降法，前面是帶m的，而矩陣形式，我們直接約掉了m,相當于學習率就被放大了m倍，所以這里學習率a應該設定為a/m=0.001，這樣a就相等了，迭代次數也相等，冥冥中感覺這次的系數應該是一樣的才對，再跑一下代碼：

w = model.fit(x,y,a=0.001,epochs=50)
print(w)

哈哈，完全一樣，破案了，這里也解釋了，之前說的，為什么損失函式前面1/2m這個值，其實對模型的引數沒有影響，
但是你的學習率要選擇得對，不能可能無法收斂，

現在這個類就是Python線性回歸代碼的一般式了，設定合理的學習率和迭代次數，就會得到不錯的結果，

2、標準方程法

下面來介紹第二種方法，標準方程法，
有了前面的鋪墊，這里就很容易理解了，損失函式：

因為這是一個凸函式，因此一定有極小值，根據最小二乘法的原理，我們要對這個損失函式對θ向量求導取0，結果如下式：
這個推導程序中也可以看到，1/2m對最終的系數沒有影響，可以直接被約掉，
撰寫標準方程法代碼如下：

class normal():
    def __init__(self):
        pass

    def fit(self,x,y):
        m=x.shape[0]
        X = np.concatenate((np.ones((m,1)),x),axis=1)
        xMat=np.mat(X)
        yMat =np.mat(y.reshape(-1,1))

        xTx=xMat.T*xMat
        #xTx.I為xTx的逆矩陣
        ws=xTx.I*xMat.T*yMat
        return ws
       
model =normal()
model.fit(x,y)    

求出來的引數為：

這里要注意：XTX的逆矩陣不是什么時候都可以求得出來的，以下情況求不到XTX的逆矩陣：
1、特征資料高度線性相關；
2、n >>m，即特征數量大于樣本數量，此時為非滿秩矩陣；

sklearn實作對比標準方程法

from sklearn.linear_model import LinearRegression
LR=LinearRegression()
LR.fit(x,y)
LR.intercept_,LR.coef_

和撰寫的標準方程法得到的引數一模一樣，這里回答了之前說過為什么梯度下降法得到的引數和sklearn里面得到的引數不一樣的問題，也說明了sklearn中封裝的是標準方程法，畢竟真的簡單！

下篇介紹非線性回歸，當資料表現為非線性時，該怎么處理，