线性回归学习笔记

基本概念

线性回归试图学习一个线性模型来尽可能预测给定输入数据的输出值。假定我们的数据集的形式为\(D=\{\left(x_1,y_1\right),\left(x_2,y_2\right),\left(x_n,y_n\right)\}\)包含了n个样本,线性回归的目的就是对于给定的数据\(x_i=\left(x_{i1};x_{i2};\cdots;x_{im}\right)\)(包含m个特征)尽可能做到能够输出与\(y_i\)相较接近的值。

我们可以将线性回归的预测函数表示为\(h\left(x\right)=\theta_0+\theta_1x_1+\theta_2x_2+\cdots+\theta_mx_m\),其中\(\theta_i\)称之为参数,用于参数化\(X\to Y\)的线性函数映射空间。在补充上\(x_0=1\)后,这里可以将式子进一步表示为\(h(x)=\sum_{i=0}^m\theta_ix_i=\theta^Tx\)

我们的目标其实就是求解参数\(\theta\),使得计算出的\(h(x)\)与数据集的\(y\)尽量接近。为了描述所谓接近,这里定义一个关于\(\theta\)的函数,使其能够描述\(h(x_i)\)与相应\(y_i\)间的差距,通常称作成本函数(cost function)

\[J(\theta)=\frac{1}{2}\sum_{i=1}^n\left(h\left(x_i\right)-y_i\right)^2 \]

这里除以2的原因是为了后续求导的时候能够抵消从而方便运算。

求解线性回归

这里介绍两种常见的求解方式:梯度下降法和最小二乘法。

梯度下降法(Gradient Descent)

算法为初始化某个\(\theta\)值,随后不断修正\(\theta\)值使得\(J(\theta)\)变小,最终收敛到使得\(J(\theta)\)取到最小值。梯度下降法通过每次更新\(\theta_j:=\theta_j-\alpha\frac{\partial}{\partial\theta_j}J(\theta)\quad (j=0,\cdots,m)\)使得\(J(\theta)\)收敛于最小值。其中\(\alpha\)称为学习率(learning rate)。

\[\begin{align} \frac{\partial}{\partial\theta_j}J(\theta)&=\frac{\partial}{\partial\theta_j}\frac{1}{2}\sum_{i=1}^n\left(h\left(x_i\right)-y_i\right)^2 \\&=\sum_{i=1}^n \left(h(x_i)-y_i\right)\cdot\frac{\partial}{\partial\theta_j}\left(h(x_i)-y_i\right) \\&=\sum_{i=1}^n \left(h(x_i)-y_i\right)\cdot\frac{\partial}{\partial\theta_j}\left(\sum_{k=0}^m\theta_kx_{ik}-y_i\right) \\&=\sum_{i=1}^n \left(h(x_i)-y_i\right)x_{ij} \end{align}\]

批量梯度下降法(Batch Gradient Descent)

这里引出批量梯度下降法,按下式更新参数\(\theta\)直到收敛:

\[\theta_j:=\theta_j-\alpha\sum_{i=1}^n \left(h(x_i)-y_i\right)x_{ij}\quad (j=0,\cdots,m) \]

可以看出批量梯度下降法每次更新参数都需要遍历数据集中所有的样本。

随机梯度下降法(Stochastic Gradient Descent)

随机梯度下降法也称为增量梯度下降(Incremental Gradient Descent),顾名思义,每有一个新的数据,则调整一次参数\(\theta\),而不需要遍历完整的数据集。
也即对于数据集中的数据按下式进行更新:

\[\quad \theta_j:=\theta_j-\alpha\left(h(x_i)-y_i\right)x_{ij}\quad (j=0,\cdots,m) \]

随机梯度下降有可能无法收敛于最小值点,\(\theta\)可能会在\(J(\theta)\)最小值附近变化,但通常来说都足够接近。并且在实践中会慢慢减小\(\alpha\)值来使得参数收敛于全局最小值。在数据量比较大的情况下,随机梯度下降法会是一个不错的选择。

最小二乘法

这里介绍用最小二乘法的矩阵解法来做线性回归的参数估计。
设计一个\(n\times (m+1)\)的矩阵:

\[X= \left[ \begin{matrix} 1 & x_{11} & x_{12} & \cdots & x_{1m}\\ 1 & x_{21} & x_{22} & \cdots & x_{2m}\\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & x_{n1} & x_{n2} & \cdots & x_{nm} \end{matrix} \right] \]

定义目标矩阵Y:

\[Y= \left[ \begin{matrix} y_1 \\ \vdots \\ y_n \end{matrix} \right] \]

这时有

\[X\theta-Y= \left[ \begin{matrix} h(x_1)-y_1 \\ \vdots \\ h(x_n)-y_n \end{matrix} \right] \]

那么我们可以将\(J(\theta)\)表示为

\[J(\theta)=\frac{1}{2}\sum_{i=1}^n\left(h\left(x_i\right)-y_i\right)^2=\frac{1}{2}(X\theta-Y)^T(X\theta-Y) \]

对其进行求导可以得到

\[\begin{align} \frac{\partial}{\partial\theta} J(\theta) &=\frac{\partial}{\partial\theta} \frac{1}{2}(X\theta-Y)^T(X\theta-Y) \\&=\frac{1}{2} \frac{\partial}{\partial\theta} (\theta^TX^TX\theta-\theta^TX^TY-Y^TX\theta+Y^TY) \\&=\frac{1}{2} (X^TX\theta+X^TX\theta-2X^TY) \\&=X^TX\theta-X^TY \end{align}\]

令式(8)取零,可得

\[\theta=(X^TX)^{-1}X^TY \]

实践

接下来通过实践来运用一下线性回归。这里使用很经典的波士顿房价数据,在sklearn的dataset中已经包含了这份数据。波士顿房价数据描述中定义了数据的格式,共有14列属性,其中前13列用于描述房屋的各种信息,而第14列为房价的中位数,也是我们需要预测的值。

属性名 解释 类型
CRIM 该镇的人均犯罪率 连续值
ZN 占地面积超过25,000平方呎的住宅用地比例 连续值
INDUS 非零售商业用地比例 连续值
CHAS 是否邻近 Charles River 离散值,1=邻近;0=不邻近
NOX 一氧化氮浓度 连续值
RM 每栋房屋的平均客房数 连续值
AGE 1940年之前建成的自用单位比例 连续值
DIS 到波士顿5个就业中心的加权距离 连续值
RAD 到径向公路的可达性指数 连续值
TAX 全值财产税率 连续值
PTRATIO 学生与教师的比例 连续值
B 1000(BK - 0.63)^2,其中BK为黑人占比 连续值
LSTAT 低收入人群占比 连续值
MEDV 同类房屋价格的中位数 连续值

由于主要是以实践线性回归为目的,所以这里不对数据作一些特征工程处理。

下面是我写的一个分别用sklearn中的LinearRegression和SGDRegressor对这份波士顿房价数据进行线性回归预测的代码。

from sklearn import datasets
from sklearn.metrics import mean_squared_error
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import SGDRegressor
import pandas as pd
import matplotlib.pyplot as plt

boston = datasets.load_boston()
data = boston.data
target = boston.target
# 切分数据,80%数据用作训练,剩下的用于测试
X_train, X_test, y_train, y_test = train_test_split(data, target, test_size=0.2, random_state=1)

# 正规化数据
ss = StandardScaler()
X_train = ss.fit_transform(X_train)
X_test = ss.transform(X_test)

# 用LinearRegression训练切分完的数据
lr = LinearRegression()
lr.fit(X_train, y_train)
lr_y_predict = lr.predict(X_test)

# 用SGDRegressor训练切分完的数据
sgd = SGDRegressor()
sgd.fit(X_train, y_train)
sgd_y_predict = sgd.predict(X_test)

# 比较LinearRegression和SGDRegressor的R2与MSE
compares = pd.DataFrame(
    {'Regressor': ['LinearRegression', 'SGDRegressor'],
     'Score': [lr.score(X_test, y_test), sgd.score(X_test, y_test)],
     'MSE': [mean_squared_error(y_test, lr_y_predict), mean_squared_error(y_test, sgd_y_predict)],
     },
    columns=['Regressor', 'Score', 'MSE'])
print(compares)

# 作图
fig, axes = plt.subplots(1, 2)
combine = [(lr_y_predict, 'LinearRegression Result'), (sgd_y_predict, 'SGDRegressor Result')]
for ax, predict in zip(axes, combine):
    result, title = predict
    ax.set_title(title)
    ax.scatter(y_test, result, c='b')
    ax.plot([result.min(), result.max()], [result.min(), result.max()], 'r--')

plt.show()

可以看到比较结果为:
Regressor Score MSE
0 LinearRegression 0.763481 23.374563
1 SGDRegressor 0.750579 24.649675

最后作出的图中展示了模型最终的预测效果。

参考

posted @ 2017-10-17 19:36  活在夢裡  阅读(286)  评论(0)    收藏  举报