线性回归模型（OLS）3

本系列文章基于 R 语言中 lm 函数的输出，介绍线性回归模型的例子和原理。

本文是系列文章的第三篇。前两篇文章中：

• 《线性回归模型（OLS）1》介绍了线性回归模型的定义以及基于普通最小二乘（OLS）求解回归系数的方法。
• 《线性回归模型（OLS）2》介绍了线性回归模型常用的假设以及基于这些假设得到的一些结论。

本文将介绍线性回归模型的模型评估。包括以下 4 个小节：

1. 模型评估
2. 示例：mtcars 数据集
3. 模型推导
4. 附录代码

为求简洁，下文用《1》指代《线性回归模型（OLS）1》，而用《2》指代《线性回归模型（OLS）2》。

以下内容为免费试读部分，完整文章可到公号“生信了”付费阅读

1. 模型评估

前文《1》介绍了线性回归模型用以下关系式来描述因变量 ${\mathbf{y}}_{n\times 1}$ 与自变量 ${\mathbf{X}}_{n\times (p+1)}$ 之间的关系：

$$\begin{array}{c}\mathbf{y}=\mathbf{X}\mathbf{\beta }+\mathbf{ϵ}\end{array}\text{\hspace{1em}(3.1)}$$

其中 $\mathbf{\beta }$ 是回归系数，而 $\mathbf{ϵ}$ 是误差项。注意 $\mathbf{X}$ 一共包含 $p$ 个自变量；其第一列都是 1，是对应着偏移项 ${\beta }_0$ 的值。OLS 下回归系数的最优估计值 $\hat{\mathbf{\beta }}$ 满足下式，即

$$\begin{array}{rl}\hat{\mathbf{\beta }}& =\arg \underset{\mathbf{\beta }}{\min }\Vert \mathbf{y}-\mathbf{X}\mathbf{\beta }{\Vert }^2\end{array}\text{\hspace{1em}(3.2)}$$

其解析解为：

$$\begin{array}{c}\hat{\mathbf{\beta }}=({\mathbf{X}}^{\mathsf{T}}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathsf{T}}\mathbf{y}\end{array}\text{\hspace{1em}(3.3)}$$

在求得回归系数的估计值后，我们得到了一个“具体”的模型来描述因变量和自变量之间的关系。那接下来就是评估这个模型。在 R 语言 lm 函数的输出中，模型评估主要分为三个部分：

• 拟合优度的计算，即检查模型对数据的拟合效果有多好。
• 模型的 $F$ 检验，即检查模型是否“显著”，或模型的自变量整体是否和因变量显著相关。
• 回归系数的 $t$ 检验，即检查每个自变量是否和因变量显著相关。

1.1. 拟合优度的计算

拟合优度（goodness of fit）是用来评估模型对观测值的拟合程度的概念；在线性回归模型中，可以用决定系数 $R^2$（coefficient of determination）表示，其计算式如下：

$$\begin{array}{c}{R}^2=\frac{\mathrm{ESS}}{\mathrm{TSS}}=1-\frac{\mathrm{RSS}}{\mathrm{TSS}}\end{array}\text{\hspace{1em}(3.4)}$$

其中

$$\begin{array}{c}\mathrm{ESS}=\sum _{i=1}^n({\hat{y}}_i-\bar{y}{)}^2=\Vert \hat{\mathbf{y}}-\bar{y}{\Vert }^2\end{array}\text{\hspace{1em}(3.5)}$$

$$\begin{array}{c}\mathrm{RSS}=\sum _{i=1}^n(y_i-{\hat{y}}_i{)}^2=\Vert \mathbf{y}-\hat{\mathbf{y}}{\Vert }^2\end{array}\text{\hspace{1em}(3.6)}$$

$$\begin{array}{c}\mathrm{TSS}=\sum _{i=1}^n(y_i-\bar{y}{)}^2=\Vert \mathbf{y}-\bar{y}{\Vert }^2\end{array}\text{\hspace{1em}(3.7)}$$

上面的 $\bar{y}=\frac{1}{n}\sum _{i=1}^n{y}_i$。

$\mathrm{ESS}$（Explained Sum of Squares）是回归平方和，反映了观测值中被模型“解释”了的那部分方差；$\mathrm{RSS}$（Residual Sum of Squares）是残差平方和（也叫作剩余平方和），反映了观测值中没有被模型“解释”的那部分方差；而 $\mathrm{TSS}$（Total Sum of Squares）是观测值中所有的方差。

当 $R^2$ 越大时，代表被回归模型“解释”了的方差的比例越高，也就说明模型的拟合效果越好。

我们会在下文中证明

$$\begin{array}{c}\mathrm{TSS}=\mathrm{ESS}+\mathrm{RSS}\end{array}\text{\hspace{1em}(3.8)}$$

此外，值得注意的是：

推论(T3.1): 当模型使用的自变量增多（在原来模型的基础上增加自变量；第 3 小节有详细说明），$R^2$ 总是会不变或升高（至少不下降）。

因此，为了减少模型包含“无关”自变量使得 $R^2$ 升高的 bias，一般我们还会计算调整 $R^2$（Adjusted R-Square; $adj.R^2$）来衡量模型的拟合效果：

$$\begin{array}{rl}adj.R^2& =1-\frac{\frac{\mathrm{RSS}}{n-p-1}}{\frac{\mathrm{TSS}}{n-1}}\\ & =1-(1-R^2)\frac{n-1}{n-p-1}\end{array}\text{\hspace{1em}(3.9)}$$

其中 $n$ 是观测值的数量；而 $p$ 是自变量的数量。调整 $R^2$ 并不总是随着自变量数量的增多而增大。

另外，无论是 $R^2$ 还是调整 $R^2$ 都只是衡量模型对观测值拟合效果的指标，而不能作为模型是否“正确”的衡量指标 [ref 9]。$R^2$ 或调整 $R^2$ 都没有对模型的前提假设进行诊断。“错误”的模型可能得到较高的 $R^2$，而“正确”的模型也可能得到较低的 $R^2$。

1.2. 模型的 $F$ 检验

我们知道，统计学中仅仅求得一个点估计值往往是不够的，经常还需要计算出这个点估计值的“可靠程度”，即要进行假设检验。在求得回归系数的估计值之后，我们自然地想要对每个回归系数进行检验；一般会检验每个回归系数是否为 0，也就是去检查对应的自变量是否和因变量显著相关。我们可以利用 $F$ 检验对自变量整体进行检验或者利用 $t$ 检验对每个自变量进行单独的检验。

值得注意的是，无论是 $F$ 检验 还是 $t$ 检验，都是建立在以下几个前提条件上的（我们在前文《2》中已经介绍过这几个条件）：

$$\begin{array}{c}\mathbb{E}[\mathbf{ϵ}|\mathbf{X}]=0\end{array}\text{\hspace{1em}(3.10)}$$

$$\begin{array}{c}\text{var}[\mathbf{ϵ}|\mathbf{X}]={\sigma }^2{\mathbf{I}}_n\end{array}\text{\hspace{1em}(3.11)}$$

$$\begin{array}{c}\mathbf{ϵ}\sim \mathcal{N}(0,{\sigma }^2{\mathbf{I}}_n)\end{array}\text{\hspace{1em}(3.12)}$$

$$\begin{array}{c}\text{rank}(\mathbf{X})=p+1\end{array}\text{\hspace{1em}(3.13)}$$

$F$ 检验

在求得回归系数估计值后，我们可以首先检验模型是否“显著”，或者说检验自变量是否整体上和因变量显著相关（$F$ 检验）；更具体地，就是去检验自变量是否整体上对模型的拟合效果有“显著”的贡献。

从上一小节 1.1 中我们知道当移除模型中某些自变量后，模型的拟合效果不变或降低。如果某些自变量和因变量不显著相关，或者说这些自变量对模型的拟合没有贡献，那删除这些自变量可以使模型变得简洁，同时不”显著“降低模型的拟合效果。一个特殊情形是如果所有自变量在整体上对模型拟合没有贡献，那么可以删去所有自变量而只保留偏移项。

线性回归模型（OLS）的 $F$ 检验意在检查自变量是否整体上对模型的拟合效果有“显著”的贡献。根据式(3.1)，我们将要检验的模型记为 test model：

$$\begin{array}{c}\text{Test Model}\text{:}\mathbf{y}=\mathbf{X}\mathbf{\beta }+\mathbf{ϵ}\end{array}\text{\hspace{1em}(3.14)}$$

其中 $\mathbf{ϵ}\sim \mathcal{N}(0,{\sigma }^2{\mathbf{I}}_n)$。或者依照单个观测值表示为：

$$\begin{array}{c}\text{Test Model}\text{:}{y}_i={\beta }_0+{\beta }_1{x}_{i1}+{\beta }_2{x}_{i2}+\cdots +{\beta }_p{x}_{ip}+{ϵ}_i\end{array}\text{\hspace{1em}(3.15)}$$

$F$ 检验的零假设（${\mathrm{H}}_0$）为自变量整体上对模型的拟合效果没有“显著”的贡献，所有的自变量可以从模型中移除；其对应的模型记为零模型（null model）：

$$\begin{array}{c}{\mathrm{H}}_0\text{:}{\beta }_1={\beta }_2=\cdots ={\beta }_p=0\end{array}\text{\hspace{1em}(3.16)}$$

$$\begin{array}{c}\text{Null Model}\text{:}{y}_i={\beta }_0+{ϵ}_i\end{array}\text{\hspace{1em}(3.17)}$$

零模型用矩阵形式表示为：

$$\begin{array}{c}\text{Null Model}\text{:}\mathbf{y}={\beta }_0+\mathbf{ϵ}\end{array}\text{\hspace{1em}(3.18)}$$

而备择假设 ${\mathrm{H}}_1$ 为自变量整体上对模型的拟合效果有“显著”的贡献，不应该将所有自变量从模型中移除；其对应的模型记为扩展模型（extended model）：

H1 :∃ j∈{1,2,⋯ ,p}, s.t. βj≠0\begin{equation} \mathrm{H}_1 \, \text{:} \quad \exists \, j \in \{1,2,\cdots,p \}, \ \text{s.t.} \ \beta_j \ne 0 \tag{3.19} \end{equation}H1​:∃j∈{1,2,⋯,p}, s.t. βj​=0​(3.19)​

Extended Model :yi=β0+β1xi1+β2xi2+⋯+βpxip+ϵi;∃ j∈{1,2,⋯ ,p}, s.t. βj≠0\begin{align*} \text{Extended Model} \, \text{:} \quad &y_i = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip} + \epsilon_i; \\ &\exists \, j \in \{1,2,\cdots,p \}, \ \text{s.t.} \ \beta_j \ne 0 \tag{3.20} \end{align*}Extended Model:​yi​=β0​+β1​xi1​+β2​xi2​+⋯+βp​xip​+ϵi​;∃j∈{1,2,⋯,p}, s.t. βj​=0​(3.20)​

扩展模型用矩阵形式表示为：

Extended Model :y=Xβ+ϵ;∃ j∈{1,2,⋯ ,p}, s.t. βj≠0\begin{align*} \text{Extended Model} \, \text{:} \quad &\mathbf{y} = \mathbf{X} \boldsymbol{\beta} + \boldsymbol{\epsilon}; \\ &\exists \, j \in \{1,2,\cdots,p \}, \ \text{s.t.} \ \beta_j \ne 0 \tag{3.21} \end{align*}Extended Model:​y=Xβ+ϵ;∃j∈{1,2,⋯,p}, s.t. βj​=0​(3.21)​

简单来说，extended model 使用了 p 个自变量，其回归系数至少有一个不为 0。此外，在进行 $F$ 检验之前，我们添加了式(3.10)到式(3.13)所描述的几个条件；那么式(3.18)和式(3.21)都有 $\mathbf{ϵ}\sim \mathcal{N}(0,{\sigma }^2{\mathbf{I}}_n)$。

为了避免混淆，下文中我们分别对 null model 和 extended model 的统计量加注上标 (null) 和 (ext)。

在 $F$ 检验中，我们构建一个 $F$ 统计量 $\frac{{\mathrm{ESS}}^{\text{(ext)}}/p}{{\mathrm{RSS}}^{\text{(ext)}}/(n-p-1)}$；其中 ${\mathrm{ESS}}^{\text{(ext)}}$ 和 ${\mathrm{RSS}}^{\text{(ext)}}$ 的计算就是式(3.5)和式(3.6)。我们发现（证明过程见第 3 节）：

$$\begin{array}{c}{\mathrm{RSS}}^{\text{(null)}}={\mathrm{TSS}}^{\text{(ext)}}\end{array}\text{\hspace{1em}(3.22)}$$

于是，

$$\begin{array}{rl}{\mathrm{ESS}}^{\text{(ext)}}& ={\mathrm{TSS}}^{\text{(ext)}}-{\mathrm{RSS}}^{\text{(ext)}}\\ & ={\mathrm{RSS}}^{\text{(null)}}-{\mathrm{RSS}}^{\text{(ext)}}\end{array}\text{\hspace{1em}(3.23)}$$

所以我们可以将上述 $F$ 统计量的分子看成是两个模型的 $\mathrm{RSS}$ 的差值，从直观上这也符合 $F$ 检验的初衷。也就是说我们可以将线性回归模型的 $F$ 检验看作是比较 null model 和 extended model 之间的 $\mathrm{RSS}$ 的差异是否显著。

我们可以证明在零假设条件下，上述统计量符合自由度是 $p$ 和 的 $n-p-1$ 的 $F$ 分布（证明过程见第 3 节）：

$$\begin{array}{c}\frac{{\mathrm{ESS}}^{\text{(ext)}}/p}{{\mathrm{RSS}}^{\text{(ext)}}/(n-p-1)}\sim F_{p,n-p-1}\end{array}\text{\hspace{1em}(3.24)}$$

于是，我们可以根据式(3.24)计算出 $p$-value，从而对自变量整体是否对模型拟合有显著贡献作出统计检验。

类似上一小节讲到拟合优度时提到的，$F$ 检验只是检查模型自变量整体在拟合效果上是否显著；并不作为模型是否“正确”的指标；它并没有对模型的前提假设进行诊断。

1.3. 回归系数的 $t$ 检验

如果 $F$ 检验的结果是自变量整体没有和因变量显著相关，那么需要根据数据的特性对现有自变量作一些变换，或者重新选择一些新的自变量，然后再构建新的模型。否则，如果 $F$ 检验的结果显著，我们需要进一步对每个自变量作单独的检验，即 $t$ 检验。当然，我们需要注意 $t$ 检验中可能存在的多重检验带来的（假阳性）问题，尤其是自变量数量多的时候。

假设我们想要检验第 $j$ 个自变量是否和因变量显著相关。此时，$t$ 检验的零假设（null hypothesis）为第 $j$ 个自变量和因变量没有显著相关，即：

$$\begin{array}{c}{\mathrm{H}}_0\text{:}{\beta }_j=0\end{array}\text{\hspace{1em}(3.25)}$$

注意这里 $t$ 检验在检验回归系数是否为 0，也可以用类似的方法对偏移项 ${\beta }_0$ 是否为 0 进行检验。

$t$ 检验的备择假设为第 $j$ 个自变量和因变量显著相关，即：

$$\begin{array}{c}{\mathrm{H}}_1\text{:}{\beta }_j\ne 0\end{array}\text{\hspace{1em}(3.26)}$$

假定式(3.10)到式(3.13)所描述的几个假设条件都成立，我们已经在前文《2》中证明了，在零假设条件下，统计量 $T_j$ 符合自由度是 $n-p-1$ 的 $t$ 分布：

$$\begin{array}{c}{T}_j=\frac{{\hat{\beta }}_j}{\sqrt{s^2({\mathbf{X}}^{\mathsf{T}}\mathbf{X}{)}_{jj}^{-1}}}\sim t_{n-p-1}\end{array}\text{\hspace{1em}(3.27)}$$

其中

$$\begin{array}{c}{s}^2=\frac{(\mathbf{y}-\mathbf{X}\hat{\mathbf{\beta }}{)}^{\mathsf{T}}(\mathbf{y}-\mathbf{X}\hat{\mathbf{\beta }})}{n-p-1}\end{array}\text{\hspace{1em}(3.28)}$$

详细信息可以参见前文《2》中的式(2.11)。

于是我们可以根据式(3.27)和式(3.28)计算出各个回归系数（包括偏移项）的 $T_j$ 值以及 $p$-value，从而对各个自变量是否和因变量显著相关作出统计检验。

置信区间

顺便提一下，关于线性回归模型的 $t$ 检验，常见的操作还包括回归系数置信区间的计算。这个可以根据式(3.27)和式(3.28)推导得到。比如，${\beta }_j$ 的 $\gamma$（比如常见的95%）置信区间可以表示为：

$$\begin{array}{c}\left({\hat{\beta }}_j-{\mathcal{T}}_{n-p-1}^{-1}\left(\frac{1+\gamma }{2}\right)\widehat{s.e.}({\hat{\beta }}_j),{\hat{\beta }}_j+{\mathcal{T}}_{n-p-1}^{-1}\left(\frac{1+\gamma }{2}\right)\widehat{s.e.}({\hat{\beta }}_j)\right)\end{array}\text{\hspace{1em}(3.29)}$$

其中 ${\hat{\beta }}_j$ 的“标准差” $\widehat{s.e.}({\hat{\beta }}_j)$ 为：

$$\begin{array}{c}\widehat{s.e.}({\hat{\beta }}_j)=\sqrt{s^2({\mathbf{X}}^{\mathsf{T}}\mathbf{X}{)}_{jj}^{-1}}\end{array}\text{\hspace{1em}(3.30)}$$

${\mathcal{T}}_m(\cdot )$ 是自由度为 $m$ 的 $t$ 分布的累积分布函数；而 ${\mathcal{T}}_m^{-1}()$ 是 ${\mathcal{T}}_m(\cdot )$ 对应的分位数（quantile）函数。

下文给出一个 R 语言中线性回归模型的例子，然后介绍其拟合优度的计算、$F$ 检验和 $t$ 检验。

2. 示例：mtcars数据集