初级算法_2线性回归

理论部分

1 线性回归概念

1.1 线性回归的一般形式

有数据集 { ( x 1 , y 1 ) , ( x 2 , y 2 ) , . . . , ( x n , y n ) } \{(x_1,y_1),(x_2,y_2),...,(x_n,y_n)\} {(x1​,y1​),(x2​,y2​),...,(xn​,yn​)},其中,$x_i=(x_{i1};x_{i2};x_{i3};...;x_{id}),y_i\in R$
其中n表示变量的数量，d表示每个变量的维度。
可以用以下函数来描述y和x之间的关系：
$$\begin{array}{l}f(x)={\theta }_0+{\theta }_1{x}_1+{\theta }_2{x}_2+...+{\theta }_d{x}_d\\ =\sum _{i=0}^d{\theta }_i{x}_i\end{array}$$
如何来确定$\theta$的值，使得$f(x)$尽可能接近y的值呢？均方误差是回归中常用的性能度量，即：
$$J(\theta )=\frac{1}{2}\sum _{j=1}^n(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$$

• 损失函数(Loss Function)：度量单样本预测的错误程度，损失函数值越小，模型就越好。
• 代价函数(Cost Function)：度量全部样本集的平均误差。
• 目标函数(Object Function)：代价函数和正则化函数，最终要优化的函数。

1.2 线性回归的优化方法

1.2.1 梯度下降法

设定初始参数$\theta$,不断迭代，使得$J(\theta )$最小化：
$${\theta }_j:={\theta }_j-\alpha \frac{\partial J(\theta )}{\partial \theta }$$
即：
$${\theta }_j={\theta }_j+\alpha \sum _{i=1}^n(y^{(i)}-f_{\theta }(x{)}^{(i)})x_j^{(i)}$$
注：下标j表示第j个参数，上标i表示第i个数据点。
向量表示：
$$\theta =\theta +\alpha \sum _{i=1}^n(y^{(i)}-f_{\theta }(x{)}^{(i)})x^{(i)}$$
梯度下降法的缺陷：如果函数为非凸函数，有可能找到的并非全局最优值，而是局部最优值。

1.2.2 最小二乘法矩阵求解

结果为：
$$\theta =(X^{T}X{)}^{(-1)}{X}^{T}Y$$

1.2.3 牛顿法

$$f(\theta {)}^{\prime }=\frac{f(\theta )}{\Delta },\Delta ={\theta }_0-{\theta }_1$$$$可求得，{\theta }_1={\theta }_0-\frac{f({\theta }_0)}{f({\theta }_0{)}^{\prime }}$$ 重复迭代，可以让逼近取到$f(\theta )$的最小值
当我们对损失函数$l(\theta )$进行优化的时候，实际上是想要取到$l^{\prime }(\theta )$的最小值，因此迭代公式为：
$$\theta :=\theta -\frac{l^{\prime }(\theta )}{l^{\prime \prime }(\theta )}$$ $$当\theta 是向量值的时候，\theta :=\theta -H^{-1}{\Delta }_{\theta }l(\theta )$$ 其中，${\Delta }_{\theta }l(\theta )$是$l(\theta )$对${\theta }_i$的偏导数，$H$是$J(\theta )$的海森矩阵，
$$H_{ij}=\frac{{\partial }^{2}l(\theta )}{\partial {\theta }_i\partial {\theta }_j}$$

1.2.4 拟牛顿法

拟牛顿法的思路是用一个矩阵替代计算复杂的海森矩阵H，因此要找到符合H性质的矩阵
要求得海森矩阵符合的条件，同样对泰勒公式求导$f^{\prime }(x)=f^{\prime }(x_0)+f^{\prime \prime }(x_0)x-f^{\prime \prime }(x_0)x_0$
令$x=x_1$，即迭代后的值，代入可得：
$$f^{\prime }(x_1)=f^{\prime }(x_0)+f^{\prime \prime }(x_0)x_1-f^{\prime \prime }(x_0)x_0$$更一般的，
$$f^{\prime }(x_{k+1})=f^{\prime }(x_k)+f^{\prime \prime }(x_k)x_{k+1}-f^{\prime \prime }(x_k)x_k$$ $$f^{\prime }(x_{k+1})-f^{\prime }(x_k)=f^{\prime \prime }(x_k)(x_{k+1}-x_k)=H(x_{k+1}-x_k)$$ $x_k$为第k个迭代值

1.3 线性回归的评价指标

均方误差(MSE):$\frac{1}{m}{\sum }_{i=1}^m(y^{(i)}-{\hat{y}}^{(i)}{)}^2$
均方根误差(RMSE)：$\sqrt{MSE}=\sqrt{\frac{1}{m}{\sum }_{i=1}^m(y^{(i)}-{\hat{y}}^{(i)}{)}^2}$
平均绝对误差(MAE)：$\frac{1}{m}\sum^{m}_{i=1} | (y^{(i)} - \hat y^{(i)} | $
但以上评价指标都无法消除量纲不一致而导致的误差值差别大的问题，最常用的指标是$R^2$,可以避免量纲不一致问题
$$R^2:=1-\frac{{\sum }_{i=1}^m(y^{(i)}-{\hat{y}}^{(i)}{)}^2}{{\sum }_{i=1}^m(\bar{y}-{\hat{y}}^{(i)}{)}^2}=1-\frac{\frac{1}{m}{\sum }_{i=1}^m(y^{(i)}-{\hat{y}}^{(i)}{)}^2}{\frac{1}{m}{\sum }_{i=1}^m(\bar{y}-{\hat{y}}^{(i)}{)}^2}=1-\frac{MSE}{VAR}$$ 我们可以把$R^2$理解为，回归模型可以成功解释的数据方差部分在数据固有方差中所占的比例，$R^2$越接近1，表示可解释力度越大，模型拟合的效果越好。

2 代码实现

• 生成随机数据

#生成数据
import numpy as np
#生成相同随机数
np.random.seed(1234)
x = np.random.rand(500,3)
#构建映射关系，模拟真实的数据待预测值,映射关系为y = 4.2 + 5.7*x1 + 10.8*x2，可自行设置值进行尝试
y = x.dot(np.array([4.2,5.7,10.8]))

2.1 调用sklearn的线性回归模型训练数据

import numpy as np
from sklearn.linear_model import LinearRegression
import matplotlib.pyplot as plt
%matplotlib inline

# 调用模型
lr = LinearRegression(fit_intercept=True)
# 训练模型
lr.fit(x,y)
print("估计的参数值为：%s" %(lr.coef_))
# 计算R平方
print('R2:%s' %(lr.score(x,y)))
# 任意设定变量，预测目标值
x_test = np.array([2,4,5]).reshape(1,-1)
y_hat = lr.predict(x_test)
print("预测值为: %s" %(y_hat))

结果：

估计的参数值为：[ 4.2  5.7 10.8]
R2:1.0
预测值为: [85.2]

2.2 最小二乘法的矩阵求解

class LR_LS():
    def __init__(self):
        self.w = None      
    def fit(self, X, y):
        # 最小二乘法矩阵求解
        #============================= show me your code =======================
        self.w = np.dot(np.dot(np.linalg.inv(np.dot(X.T,X)),X.T),y)
        #============================= show me your code =======================
    def predict(self, X):
        # 用已经拟合的参数值预测新自变量
        #============================= show me your code =======================
        y_pred = X.dot(self_w)
        #============================= show me your code =======================
        return y_pred

if __name__ == "__main__":
    lr_ls = LR_LS()
    lr_ls.fit(x,y)
    print("估计的参数值：%s" %(lr_ls.w))
    x_test = np.array([2,4,5]).reshape(1,-1)
    print("预测值为: %s" %(lr_ls.predict(x_test)))

结果：

估计的参数值：[ 4.2  5.7 10.8]
预测值为: [85.2]

2.3 梯度下降法

class LR_GD():
    def __init__(self):
        self.w = None     
    def fit(self,X,y,alpha=0.02,loss = 1e-10): # 设定步长为0.002,判断是否收敛的条件为1e-10
        [m,d] = np.shape(X) #自变量的维度
        self.w = np.zeros((d)) #将参数的初始值定为0
        tol = 1e5
        #============================= show me your code =======================
        while tol > loss:
            for i in range(d):
                temp = y - X.dot(self.w)
                self.w[i] = self.w[i] + alpha *np.sum(temp * X[:,i])/m
            tol = np.sum(np.abs(y -  X.dot(self.w)))
        #============================= show me your code =======================
    def predict(self, X):
        # 用已经拟合的参数值预测新自变量
        y_pred = X.dot(self.w)
        return y_pred  

if __name__ == "__main__":
    lr_gd = LR_GD()
    lr_gd.fit(x,y)
    print("估计的参数值为：%s" %(lr_gd.w))
    x_test = np.array([2,4,5]).reshape(1,-1)
    print("预测值为：%s" %(lr_gd.predict(x_test)))

结果：

估计的参数值为：[ 4.2  5.7 10.8]
预测值为：[85.2]