线性回归

线性回归

文章分三块，第一块线性回归模型，第二块模型策略，第三块学习算法。

线性回归模型

线性回归的模型通常定义如下：

$$h_w(x)=w\cdot x+b$$
这里， $x\in R^N，h_w(x)\in R$, $w,b$分别是权值和偏置。

模型策略

对于给定的训练集数据$T=\{(x_1,y_1),(x_2,y_2),\cdots,(x_n,y_n)\}$,如何去选择权值$w$?一个合理的想法是使得$h_w(x_i)$接近于$y_i$，至少在训练集中是如此。
因此，我们需要定义一个损失函数并最小化。
则定义如下：

$$L(w,b)=\frac{1}{2}\sum_{i=1}^N(h_w(x_i)-y_i)^2$$
对 $L(w,b)$求极小值，得到 $w$的估计值。
线性回归模型的策略即是求使得$L(w,b)$极小的权值 $w,b$。

学习算法

这里分别采用批次梯度下降算法和随机梯度下降算法，并进行比较。

批次梯度下降算法

由之前一篇文章最小二乘法的两种表示可知，对w的求导向量表示法如下：

$$\begin{align} \nabla_w(L(w))&=\nabla_w\frac{1}{2}(Xw-y)^T(Xw-y)\\ &=\frac{1}{2}\nabla_w(w^TX^TXw-w^TX^Ty-y^TXw+y^Ty)\\ &=\frac{1}{2}\nabla_wtr(w^TX^TXw-w^TX^Ty-y^TXw+y^Ty)\\ &=\frac{1}{2}(2X^TXw-2X^Ty)\\ &=X^T(Xw-y) \end{align}$$
$$w=w-\alpha \nabla_w(L(w))$$ 其中 $\alpha$为学习率，过大或过小都不会训练处合适的模型。
这个算法在每一次迭代中需要用到所用的样本，因此也称为批次梯度下降算法。尽管梯度下降很容易收敛到局部最优解，但在线性回归这种情况下它只有一个全局最优解，没有局部最优解。

import numpy as np

class LinearRegression():
    def __init__(self,eta=0.001,n_iter=10):
        self.eta=eta
        self.n_iter=n_iter
    def fit(self,X,y):
        m,n=np.shape(X)
        self.w=np.zeros(n+1).reshape(n+1,-1)
        for i in np.arange(self.n_iter):
            error=self.net_input(X) -y
            self.w[1:] -=self.eta*np.dot(X.T,error)
            self.w[0] -=self.eta*np.sum(error)
        return self

    def net_input(self,X):
        return np.dot(X,self.w[1:])+self.w[0]    
    def predict(self,X):
        return self.net_input(X)

if __name__=="__main__":
    X=np.array([[3,6],[6,12],[12,24]])
    y=np.array([3,6,12]).reshape(3,-1)
    XX=np.array([[18,36],[30,60]])
    lr=LinearRegression(eta=0.0001,n_iter=200)
    lr.fit(X,y)
    print lr.predict(XX)

随机梯度下降算法

它不是针对全部训练集进行一次$w$的更新，而是对每一个样本就更新一次。所以表达式修改如下：

$$w=w-\alpha x_i^T(x_iw-y_i)\qquad (i=1,2\cdots,m)$$
在这个算法里，我们每次更新仅仅使用了一个训练样本，因此称为随机梯度算法。在fit()中稍作修改：

import numpy as np

class LinearRegression():
    def __init__(self,eta=0.001,n_iter=10):
        self.eta=eta
        self.n_iter=n_iter
    def fit(self,X,y):
        m,n=np.shape(X)
        self.w=np.zeros(n+1).reshape(n+1,-1)
        for i in np.arange(self.n_iter):
            for j in np.arange(m):                
                error=self.net_input(X[j]) -y[j]
                self.w[1:] -=self.eta*error*(X[j].reshape(n,-1))
                self.w[0] -=self.eta*error
        return self

    def net_input(self,Xi):
        return np.dot(Xi,self.w[1:])+self.w[0]

    def predict(self,X):
        return self.net_input(X)

if __name__=="__main__":
    X=np.array([[3,6],[6,12],[12,24]])
    y=np.array([3,6,12]).reshape(3,-1)
    XX=np.array([[18,36],[30,60]])
    lr=LinearRegression(eta=0.0001,n_iter=200)
    lr.fit(X,y)
    print lr.predict(XX)

两个算法对比

• 随机梯度算法往往速度更快，但可能无法收敛到最优解，而是在其附近；
• 对于训练样本数量较大时，一个大概的最优解大多能满足我们的需求，这时随机梯度算法是较好的选择。