应用机器学习（八）：线性模型

回归分析

在统计建模中，回归分析( Regression analysis ) 是用来刻画变量之间的统计关系的一种统计技术。当一个变量( 因变量 dependent variable )受另外一些变量( 自变量 independent variables or predictors )的强烈影响时，适宜用回归的方法。例如，”回归”的思想和方法最早由英国著名的统计学家F. Galton 和他的学生 K. Pearson 提出。他们在研究父母身高与其子女身高的遗传问题时，以每对夫妇的平均身高作为自变量 x, 而取他们的一个成年孩子的身高作为 y, 共观察了1,078对夫妇，建立了回归方程 y=33.73+0.516x, 从此产生了回归分析的方法。

回归模型

回归模型的一般形式

设因变量 $y$, 自变量向量 $\mathbf{x}=(x_1, x_2, \dots, x_p)$, 则刻画 $y$ 与 $\mathbf{x}$关系的回归模型的一般形式为

$$\begin{align*} y=f(\mathbf{x})+\varepsilon \tag{1} \end{align*}$$

其中， $\varepsilon$ 为随机误差，它表示除了 $\mathbf{x}$ 外的其它随机干扰因素。正是因为 $\varepsilon$ 的存在，给定 $\mathbf{x}$ 的值，不能唯一确定 $y$, 称这种变量关系为统计关系。随机变量 $\varepsilon$ 满足基本假设 $E(\varepsilon)=0$, 即零均值，它表示没有系统误差。

线性回归

线性回归( Linear regression ), 指的是 $y$ 与 $\mathbf{x}$ 之间是线性关系，即

$$\begin{align*} y=\beta_0+\beta_1 x_1 + \beta_2 x_2 + \dots + \beta_p x_p +\varepsilon \tag{2} \end{align*}$$
其中的参数 $\beta_0, \beta_1, \dots, \beta_p$ 称为回归系数。特别地， $p=1$, 即一个自变量时，称模型为一元线性回归。

假设有 $n$ 次观测数据 $( y_i; x_{i1}, x_{i2}, \dots, x_{ip} ),\, i=1,2,\dots, n$.
记$\mathbf{y}=(y_1, y_2, \dots, y_n)'$, $\mathbf{\beta}=(\beta_0, \beta_1, \dots, \beta_p)'$,
$\mathbf{\varepsilon}=(\epsilon_1, \epsilon_2, \dots, \epsilon_p)'$,

$$\begin{equation*} \mathbf{X}=\left( \begin{array}{ccccc} 1 & x_{11} & x_{12} & \dots & x_{1p}\\ 1 & x_{21} & x_{22} & \dots & x_{2p}\\ \vdots & \vdots & \vdots & & \vdots\\ 1 & x_{n1} & x_{n2} & \dots & x_{np}\\ \end{array} \right) \end{equation*}$$

则观测数据的线性回归模型可以表示为

$$\begin{align*} y_i=\beta_0+\beta_1 x_{i1} + \beta_2 x_{i2} + \dots + \beta_p x_{ip} +\varepsilon_i, \, i=1,2,\dots,n \tag{3} \end{align*}$$
等价的矩阵形式为
$$\begin{align*} \mathbf{y}=\mathbf{X}\beta +\varepsilon \tag{4} \end{align*}$$

线性回归模型的基本假设

1. 自变量 $\mathbf{x}$ 是确定性变量，不是随机变量，因变量 $y$ 是随机变量。要求
   设计矩阵 $\mathbf{X}$ 是列满秩的，即 $rank(\mathbf{X})=p+1<n$, 这表明 $\mathbf{X}$
   的列向量之间线性无关;

2. 随机误差零均值同方差，即 $E(\varepsilon)=0,\, cov(\varepsilon, \varepsilon)=\sigma^2 E_n$.
   这个假设常称为高斯-马尔科夫条件( Gauss-Markov assumption );

3. 通常，假定随机误差向量 $\varepsilon$ 服从正态分布，即 $\varepsilon\sim N(0, \sigma^2 E_n)$,
   因此，有 $\mathbf{y}\sim N(\mathbf{X}\beta, \sigma^2 E_n)$.

回归系数的估计

• 普通最小二乘估计

普通最小二乘估计( ordinary least square estimation, OLSE )的思想是，考虑 $n$ 次观测的 $y_i$ 与均值 $E(y_i)=\beta_0+\beta_1 x_{i1} + \beta_2 x_{i2} + \dots + \beta_p x_{ip}$
的离差平方和

$$\begin{align*} Q(\beta)=\sum\limits_{i=1}^n (y_i - E(y_i))^2=\sum\limits_{i=1}^n (y_i - \beta_0-\beta_1 x_{i1} - \beta_2 x_{i2} - \dots - \beta_p x_{ip})^2 \tag{5} \end{align*}$$

所谓最小二乘估计，就是找回归系数$\hat{\beta}=(\hat{\beta_0}, \hat{\beta_1}, \dots, \hat{\beta_p})'$,
使得

$$\begin{align*} Q(\hat{\beta})=\mathop{\min}_{\beta}Q(\beta) \tag{6} \end{align*}$$

注意到， 函数 $Q$ 关于 $\beta_0, \beta_1, \dots, \beta_p$ 是非负可导的二次函数，故存在最小值，且在导数为0的点处。这样，解得 $\beta$ 最小二乘估计为

$$\begin{align*} \hat{\beta}=(\mathbf{X}'\mathbf{X})^{-1} \mathbf{X}' \mathbf{y} \tag{7} \end{align*}$$
$\mathbf{y}$ 的回归值向量

$$\begin{align*} \hat{\mathbf{y}}=\mathbf{X}\hat{\beta} \tag{8} \end{align*}$$

• 最大似然估计

最大似然估计( maximum likelihood estimation ), 是在正态假设下，$\mathbf{y}\sim N(\mathbf{X}\beta, \sigma^2 E_n)$. 因此， $\mathbf{y}$ 的似然函数为

$$\begin{align*} L(\beta, \sigma^2)=(2\pi)^{-\frac{n}{2}}(\sigma^2)^{-\frac{n}{2}}exp(-\dfrac{1}{2\sigma^2} (\mathbf{y}-\mathbf{X}\beta)'(\mathbf{y}-\mathbf{X}\beta)) \tag{9} \end{align*}$$
等价地转换为对数似然
$$\begin{align*} \ln L(\beta, \sigma^2)=-\dfrac{n}{2}\ln 2\pi -\dfrac{n}{2}\ln \sigma^2 -\dfrac{1}{2\sigma^2} (\mathbf{y}-\mathbf{X}\beta)'(\mathbf{y}-\mathbf{X}\beta) \tag{10} \end{align*}$$
注意，这里的(对数)似然函数为 $\beta, \,\sigma^2$ 的函数，固定 $\sigma^2$, 求 $\ln L$ 的最大值，等价于求 $(\mathbf{y}-\mathbf{X}\beta)'(\mathbf{y}-\mathbf{X}\beta)$ 的最小值，而这时离差平方和 $Q$, 这样，回归系数向量 $\beta$ 的最大似然估计与最小二乘估计是等价的，下面，固定 $\hat{\beta}$, 得到误差项方差 $\sigma^2$ 为

$$\begin{align*} \hat{\sigma}^2=\dfrac{1}{n}(\mathbf{y}-\mathbf{X}\hat{\beta})'(\mathbf{y}-\mathbf{X}\hat{\beta}) \tag{11} \end{align*}$$
定义残差向量
$$\begin{align*} \mathbf{e}=\mathbf{y}-\hat{\mathbf{y}}=\mathbf{y}-\mathbf{X}\hat{\beta} \end{align*}$$

残差 $\mathbf{e}$ 可以看作是随机误差 $\varepsilon$ 的估计。那么， $\sigma^2$ 的最大似然估计还可以记为
$\hat{\sigma}^2=\dfrac{1}{n}\mathbf{e}'\mathbf{e}$.

回归方程的显著性检验

确定回归系数的估计后，在使用回归方程作预测或估计前，还需要对回归方程进行显著性检验，即，检验自变量和因变量之间是否有线性关系。下面介绍两种统计检验方法，一是回归方程的 F 检验，另一个是回归系数的 t 检验。同时介绍度量回归拟合数据程度的拟合优度。

• 回归方程的 F-检验

考虑所有自变量 $x_1, x_2, \dots, x_p$ 作为一个整体是否与 $y$ 有线性关系，为此提出假设

$$\begin{align*} H_0: \beta=0 \qquad H_1: \beta \ne 0 \end{align*}$$

定义平方和

$$\begin{align*} & SST=\sum\limits_{i=1}^n (y_i - \bar{y})^2\\ & SSR=\sum\limits_{i=1}^n (\hat{y}_i -\bar{y})^2\\ & SSE=\sum\limits_{i=1}^n (y_i - \hat{y}_i)^2\\ \end{align*}$$

$SST$ 称离差平方和，它表示因变量的波动，即不确定性；
$SSR$ 称回归平方和，它表示回归方程的波动；
$SSE$ 称残差平方和，它表示随机波动。
可以证明平方和分解式

$$\begin{align*} SST=SSR+SSE \tag{12} \end{align*}$$

令检验统计量

$$\begin{align*} F=\dfrac{SSR/p}{SSE/n-p-1} \tag{13} \end{align*}$$
$H_0$ 为真时， $F\sim F(p, n-p-1)$, 因此，给定水平 $\alpha$, 拒绝域的形式 $W=\{F> F_{\alpha}(p, n-p-1)\}$. 当 $F> F_{\alpha}(p, n-p-1)$ 时，拒绝 $H_0$, $y$ 与 $x_1, x_2, \dots, x_p$ 有线性关系，即，回归方程是显著的；当 $F\le F_{\alpha}(p, n-p-1)$时，则认为没有线性关系，回归方程不显著。

• 回归系数的t-检验

回归方程显著，并不意味着每个自变量都与因变量有线性关系。因此，还需要对每个自变量是否与因变量有线性关系作显著性检验。为此提出假设

$$\begin{align*} H_{0j}: \beta_j = 0 \qquad H_{1j}: \beta_j \ne 0,\qquad j=1,2,\dots,p. \end{align*}$$

在正态假设下，可以证明 $\hat{\beta}\sim N(\beta, \sigma^2 (\mathbf{X}'\mathbf{X})^{-1})$, 记 $(\mathbf{X}'\mathbf{X})^{-1}=(c_{ij}),\, i,j=0,1,2,\dots,p$. 构造检验统计量 $t_j = \dfrac{\hat{\beta}_j}{\sqrt{c_{jj}}\hat{\sigma}}$,
$\hat{\sigma}=\sqrt{\dfrac{1}{n-p-1}SSE}$. 在 $H_{0j}$ 下， $t_j\sim t(n-p-1)$, 因此，在水平 $\alpha$ 下，拒绝域
$W=\{|t_j|\ge t_{\frac{\alpha}{2}}\}$. 当 $|t_j|\ge t_{\frac{\alpha}{2}}$ 时，拒绝 $H_{0j}$, 认为自变量 $x_j$ 与 $y$ 有线性关系。否则，可以认为没有线性关系。

• 拟合优度

拟合优度( Goodness of fit ), 用于考察回归方程对观测数据的拟合程度，在回归中称决定系数，定义为

$$\begin{align*} R^2=\dfrac{SSR}{SST}=1-\dfrac{SSE}{SST} \tag{14} \end{align*}$$

$R^2$ 在区间 $[0, 1]$ 内取值，$R^2$ 越接近1，表明回归拟合的效果越好；相反，越接近0，表明拟合的越差。与 F-检验比， $R^2$ 可以更清楚地反映回归拟合的效果，但不能作为严格的显著性检验。可以证明，当观测数 n=2 时，不论变量间有没有相关性， $R^2=1$.

逐步回归法

在建立一个实际问题的回归模型时，首先遇到的问题是如何确定自变量。如果遗漏了某些重要的自变量，回归效果不好；如果包括过多的自变量，这些变量可能产生共线性，回归模型的稳定性较差。首先，介绍自变量选择的重要准则，赤池信息量准则。

自变量的选择

赤池信息量准则( Akaike information criterion, AIC ), 是根据最大似然估计的原理提出的一种模型选择准则。对一般情况，设模型的似然函数为 $L(\mathbf{\theta}, \mathbf{x})$, $\theta$ 的维数为 $p$, $\mathbf{x}$ 为随机样本，则 AIC 定义为

$$\begin{align*} AIC=-2\ln L(\hat{\theta}, \mathbf{x}) + 2 p \tag{15} \end{align*}$$

其中， $\hat{\theta}$ 是 $\theta$ 的最大似然估计。似然函数越大的估计量越好， $AIC$ 是似然函数的对数乘 -2 再加上惩罚因子 $2p$ , 因而选择使 $AIC$ 达到最小的模型是”最优”模型。
下面讨论把 $AIC$ 用于回归模型的选择。假定回归模型的随机误差 $\varepsilon\sim N(0, \sigma^2)$,
那么，最大对数似然

$$\begin{align*} \ln L_{max} =-\dfrac{n}{2}\ln 2\pi -\dfrac{n}{2}\ln \hat{\sigma}^2 - \dfrac{1}{2} \hat{\sigma}^2 SSE \tag{16} \end{align*}$$

将 $\hat{\sigma}^2 = \dfrac{1}{n} SSE$ 代入，得

$$\begin{align*} \ln L_{max} =-\dfrac{n}{2}\ln 2\pi -\dfrac{n}{2}\ln \dfrac{SSE}{n} - \dfrac{n}{2} \tag{17} \end{align*}$$

代入式(15)，得

$$\begin{align*} AIC = n\ln SSE + 2p \tag{18} \end{align*}$$

在回归建模过程中，对每一个回归的自变量子集，计算 $AIC$, 其中的 $AIC$ 最小的模型为”最优”的回归模型。

逐步回归法

逐步回归的基本思想是”有进有出”。在自变量的选择过程中，将变量一个一个引入，每引入一个自变量后，对已选入的变量进行逐个显著性检验。当原先引入的变量由于后面变量的引入而变得不再显著时，要将其从回归模型中剔除。引入一个变量或从回归模型中剔除一个变量，为逐步回归的一步，每一步都要进行 F 检验，以确保每次引入新的变量之前回归模型中只包含显著的变量。这个过程反复进行，直到既无显著的自变量选入，也无不显著的自变量剔除为止。这样，保证了最后所得的回归模型是最优的。在逐步回归过程，要求引入自变量的检验水平 $\alpha_{entry}$, 小于剔除自变量的检验水平 $\alpha_{removal}$, 否则，可能产生”死循环”。

一个实例：带交互项的线性回归

我们以 R 包 datasets 的数据集 mtcars 为例，以 mpg (每加仑汽油行驶的英里数)为因变量，
hp (马力)、 wt (汽车重量)为自变量，使用函数 lm() 建立线性回归模型。

fit <- lm(mpg ~ hp + wt, data = mtcars)
summary(fit)

从 Pr(>|t|) 栏可以看到，马力和车重都是显著的，回归拟合优度 R-squared 和 adjusted R-squared 都超过了0.8，说明线性程度较高。进一步考虑，马力和车重是否存在交互效应，为此，建立带交互项的回归模型。

fit <- lm(mpg ~ hp + wt + hp:wt, data = mtcars)
summary(fit)

从 Pr(>|t|) 栏可以看到，马力和车重的交互项是显著的，这说明响应变量 mpg 与其中一个变量的关系依赖于另一个变量的值。即，每加仑汽油行驶的英里数与汽车马力的关系依车重不同而不同。而且，通过与不带交互项的模型对比发现，带交互项的回归模型，拟合优度 R-squared 有所提高，说明带交互项的回归模型拟合数据更佳。预测模型为 mpg = 49.81 - 0.12hp - 8.22wt + 0.03hp·wt

使用 effects 包中的 effect() 函数，你可以用图形展示交互项的结果。

library(effects)
plot(effect("hp:wt", fit,
     xlevels=list(wt=c(2.2, 3.2, 4.2))),
     multiline = TRUE)

从结果图可以很清晰地看出，随着车重的增加（2.2、3.2、4.2），马力与每加仑汽油行驶英里数的关系减弱了。特别地，当 wt = 4.2 时，直线几乎是水平的，表明随着 hp 的增加，mpg 不会发生改变。

阅读更多精彩内容，请关注微信公众号“统计学习与大数据”!