一 線性模型

• 給定由n個屬性描述的列向量\(f(\mathbf{x})={(x^{(1)};x^{(2)};...;x^{(n)})}\)，其中 \(x^{(j)}\)是\(\textbf{x}\)在第\(j\)個屬性的取值。線性模型即為通過對屬性進行線性組合的函數，即
  \[f(\mathbf{x})=w_0+w_1x^{(1)}+...+w_nx^{(n)}\]
  寫成向量形式如下：
  \[f(\textbf{x})=\mathbf{w}^\mathrm{T}\mathbf{x}\]
  其中列向量\(\mathbf{w}=(w_0;w_1;...;w_n)\)，列向量\(\mathbf{x}=(1;x^{(1)};...;x^{(n)})\)
  。列向量\(\mathbf{w}\)確定後，模型隨之確定。
  線性模型形式簡單，易於建模；直觀表達了各屬性在預測中的重要性，因此具有很好的可解釋性；

二 線性模型求解

• 對於給定的數據集\(\mathbf{D}=\left \{ (\mathbf{x_1},y_1),(\mathbf{x_2},y_2),...,(\mathbf{x_m},y_m)\right \}\),其中\(\mathbf{x_i}=(x_i^{(1)};...;x_i^{(n)})\)，\(y_i\)為第\(i\)個實例的實際值。“線性回歸”試圖學得一個線性模型以盡可能準確地預測實例的輸出值，使之接近實際值。
• 關鍵問題是如何衡量兩者之間的誤差。這裏采用均方誤差作為性能度量，即利用最小二乘法來進行參數估計。
  \[\underset{\mathbf{w}}{argmin}\mathbf{\mathit{J}}(\mathbf{w})=\frac{1}{2m}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)^2\]
  
  實際上在高斯分布假設前提下，用極大似然函數來進行參數估計，可以得出上述目標,推導過程如下。
  根據中心極限定理，認為誤差項 $\mathbf\xi $服從均值為零的高斯分布
  \[P(\xi_i)=\frac{1}{\sqrt{2\pi }\sigma }exp(-\frac{\xi_i^2}{2\sigma ^2})\]
  \[P(y_i|\mathbf{x_i};\mathbf{w})=\frac{1}{\sqrt{2\pi }\sigma }exp(-\frac{(y_i-\mathbf{x}_i)^2}{2\sigma ^2})\]
  由上可得似然函數，
  \[\mathbf{L}(\mathbf{w})=\prod_{i=1}^{n}P(y_i|\mathbf{x_i};\mathbf{w})\]
  
  取對數，得對數似然函數
  \[\mathbf{L}(\mathbf{w})=nlog\frac{1}{\sqrt{2\pi }\sigma }-\frac{1}{2\sigma ^2}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)^2\]
  對上式取極大值等價於下式取極小值
  \[\mathbf{\mathit{J}}(\mathbf{w})=\frac{1}{2m}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)^2\]
  推導完畢。
• 模型求解方法：矩陣直接求解和梯度下降法。
  1 .正規方程法
  對於數據集\(D\)中的每個實例組成一個矩陣，矩陣形式如下：\[\mathbf{X}=\begin{pmatrix} 1 & x_1^{(1)} & ... & x_1^{(n)}\\ 1& x_2^{(1)} & ... & x_2^{(n)}\\ .& . & ... & .\\ 1& x_m^{(1)} & . & x_m^{(n)} \end{pmatrix}=\begin{pmatrix} 1 &\mathbf{x}_1^T \\ 1 &\mathbf{x}_2^T\\ .& .\\ 1& \mathbf{x}_m^T \end{pmatrix}\]
  對應的實際值寫成列向量形式\(\mathbf{y}=(y_1;y_2;...;y_m)\),則有
  \[\mathbf{\hat{w}}^*=\underset{\hat{\mathbf{w}}}{argmin}(\mathbf{y}-\mathbf{X\hat{w}})^T(\mathbf{y}-\mathbf{X\hat{w}})\]上式argmin後面部分對\(\hat{\mathbf{w}}\)求導，另之等於零，\[\mathbf{\hat{w}}^*=(\mathbf{X^TX})^{-1}\mathbf{X^T}\mathbf{y}\]令\(\mathbf{\hat{x}}_i=(1;\mathbf{x}_i)\)從而得到線性模型\(f(\mathbf{\hat{x}}_i)=\mathbf{\hat{x}}_i^T\mathbf{\hat{w}}^*\),或者\(f(\mathbf{\hat{x}}_i)=\mathbf{\hat{w}}^{*T}\mathbf{\hat{x}}_i\)。
  但是，現實情況中\(\mathbf{X^TX}往往不可逆\)，通常原因有兩種，一是高度共線性；二是數據特征過多而訓練數據較少，此時可以通過正則化來解決。
  2 .梯度下降法
  梯度下降是一種常用的一階優化方法，是求解無約束優化問題的經典方法之一。對於連續可微函數上某一點，有各個方向導數，沿梯度方向的方向導數達到最大值，也就是說,梯度的方向是函數在這點增長最快的方向。
  因此，我們可以得到如下結論：函數在某點的梯度是這樣一個向量，它的方向與取得最大方向導數的方向一致，而它的模為方向導數的最大值。
  所以我們可以沿反梯度方向不斷一步一步叠代，得到局部極小點。當目標函數為凸函數時，局部極小點就是全局最小點，此時梯度下降法可確保收斂到全局最優解。
  將損失函數對列向量\(\mathbf{w}\)求導，得到\(w_j\)的偏導：
  \[\frac{\partial \mathbf{J(w)}}{\partial w_j}=\frac{\partial }{\partial w_j}\frac{1}{2m}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)^2=\frac{1}{m}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)\mathbf{x}_i^{(j)},j=0,1,2,...,n\]然後對各個分量都以下面形式更新\(w_j\)：\[w_j=w_j-\alpha\frac{1}{m} \sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)\mathbf{x}_i^{(j)}\]有公式可以看出對於每一個分量進行一次叠代時計算了所有訓練樣本數據，這種稱為批量梯度下降。因此在數據量很大的時候，每次叠代都要遍歷訓練集一遍，開銷會很大。
  為改善上述情況，可以在每次叠代僅選擇一個訓練樣本去計算代價函數的梯度，然後更新參數。即使是大規模數據集，隨機梯度下降法也會很快收斂。這種方法稱為隨機梯度下降。此時有，\[w_j=w_j-\alpha (f(\mathbf{x}_i)-y_i)\mathbf{x}_i^{(j)}\]
  關於隨機梯度下降、批量梯度下降與小批量隨機梯度下降的特征不再贅述。

三 代碼實現

class MyLinearRegression():

    def __init__(self, n_iterations=10000, learning_rate=0.0001, regularization=None, gradient=True):
        ‘‘‘初始化。是否正則化及L1L2的選擇；選用梯度下降法還是正規方程法。梯度下降學習率以及叠代次數‘‘‘
        self.n_iterations = n_iterations
        self.learning_rate = learning_rate
        self.gradient = gradient
        if regularization == None:
            self.regularization = lambda x: 0
            self.regularization.grad = lambda x: 0
        else:
            self.regularization = regularization

    def initialize_weights(self, n_features):
        ‘‘‘初始化權重.初始化模型參數,加入w0‘‘‘
        limit = np.sqrt(1 / n_features)
        w = np.random.uniform(-limit, limit, (n_features, 1))              #二維數組，n行一列。
        b = 0
        self.w = np.insert(w, 0, b, axis=0)                                

    def fit(self,X,y,):
        ‘‘‘進行擬合‘‘‘
        m_samples, n_features = X.shape                                      # !!!
        self.initialize_weights(n_features)
        X = np.insert(X, 0, 1, axis=1)                                      #二維數組，每行前面加上元素1
        y = np.reshape(y, (m_samples, 1))                                    #二維數組，m 行一列
        self.training_errors = []
        if self.gradient == True:                                            #批量梯度下降
            for i in range(self.n_iterations):
                y_pred = X.dot(self.w)
                loss = np.mean(0.5 * (y_pred - y) ** 2)/m_samples + self.regularization(self.w)  # 矩陣運算
                ‘‘‘mean()函數功能：求取均值
                經常操作的參數為axis，以m * n矩陣舉例：
                axis 不設置值，對 m*n 個數求均值，返回一個實數
                axis = 0：壓縮行，對各列求均值，返回 1* n 矩陣
                axis =1 ：壓縮列，對各行求均值，返回 m *1 矩陣
                np.mean(X,axis=0或者1,keepdims=True)
                ‘‘‘
                self.training_errors.append(loss)
                w_grad = X.T.dot(y_pred - y) + self.regularization.grad(self.w)  # (y_pred - y).T.dot(X)，計算梯度
                self.w = self.w - self.learning_rate * w_grad  # 更新權值w
        else:
            # 正規方程
            X = np.matrix(X)
            y = np.matrix(y)
            X_T_X = X.T.dot(X)
            X_T_X_I_X_T = X_T_X.I.dot(X.T)
            X_T_X_I_X_T_X_T_y = X_T_X_I_X_T.dot(y)
            self.w = X_T_X_I_X_T_X_T_y

    def predict(self, X):
        X = np.insert(X, 0, 1, axis=1)
        y_pred = X.dot(self.w)
        return y_pred
‘‘‘以二元為例，進行擬合‘‘‘
lr = MyLinearRegression()
X = np.array([[1,2],[2,4],[50,3],[23,59],[10,45],[10,61]])
y = np.array([3,6,53,82,55,71])
lr.fit(X,y)
y_test = lr.predict(np.array([[1,40],[2,6]]))
print(y_test)

L1正則化下的損失函數\[\mathbf{\mathit{J}}(\mathbf{w})=\frac{1}{2m}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)^2+\lambda \left \| \mathbf{w} \right \|_{1}\]其中\(\left \| \mathbf{w} \right \|_{1}=\sum w_j\)
L2正則化下的損失函數\[\mathbf{\mathit{J}}(\mathbf{w})=\frac{1}{2m}\sum_{i=1}^{m}(f(\mathbf{x}_i)-y_i)^2+\lambda \left \| \mathbf{w} \right \|_{2}^2\]其中\(\left \| W \right \|_{2}=\sqrt{\sum w_j^2}\)
加入L1正則化後稱為Lasso回歸，加入L2正則化稱為Ridge回歸，其中 \(\lambda\)為模型的超參數。
L1 regularizer ： 它的優良性質是能產生稀疏性，導致 W 中許多項變成零。 稀疏的解除了計算量上的好處之外，更重要的是更具有“可解釋性”。
L2 regularizer ：使得模型的解偏向於 norm 較小的 W，通過限制 W 的 norm 的大小實現了對模型空間的限制，從而在一定程度上避免了 overfitting 。不過 ridge regression 並不具有產生稀疏解的能力，得到的系數 仍然需要數據中的所有特征才能計算預測結果，從計算
量上來說並沒有得到改觀。

關於正則化參考[https://www.jianshu.com/p/a47c46153326]

機器學習：線性回歸——理論與代碼實現（基於正規方程與梯度下降）