线性回归模型（OLS）3

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

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

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

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

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

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

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

1. 模型评估

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

y=Xβ+ϵ\begin{equation} \mathbf{y} = \mathbf{X} \boldsymbol{\beta} + \boldsymbol{\epsilon} \tag{3.1} \end{equation}y=Xβ+ϵ(3.1)

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

β^=arg⁡min⁡β∥y−Xβ∥2\begin{align*} \hat{\boldsymbol{\beta}} &= \arg \min_{\boldsymbol{\beta}} \| \mathbf{y} - \mathbf{X}\boldsymbol{\beta} \|^2 \tag{3.2} \end{align*}β^=argβmin∥y−Xβ∥2(3.2)

其解析解为：

β^=(XTX)−1XTy\begin{equation} \hat{\boldsymbol{\beta}} = (\mathbf{X}^\mathsf{T} \mathbf{X})^{-1} \mathbf{X}^\mathsf{T} \mathbf{y} \tag{3.3} \end{equation}β^=(XTX)−1XTy(3.3)

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

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

1.1. 拟合优度的计算

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

R2=ESSTSS=1−RSSTSS\begin{equation} R^2 = \frac{\mathrm{ESS}}{\mathrm{TSS}} = 1 - \frac{\mathrm{RSS}}{\mathrm{TSS}} \tag{3.4} \end{equation}R2=TSSESS=1−TSSRSS(3.4)

其中

ESS=∑i=1n(y^i−yˉ)2=∥y^−yˉ∥2\begin{equation} \mathrm{ESS} = \sum_{i=1}^n (\hat{y}_i - \bar{y})^2 = \| \hat{\mathbf{y}} - \bar{y} \|^2 \tag{3.5} \end{equation}ESS=i=1∑n(y^i−yˉ)2=∥y^−yˉ∥2(3.5)

RSS=∑i=1n(yi−y^i)2=∥y−y^∥2\begin{equation} \mathrm{RSS} = \sum_{i=1}^n (y_i - \hat{y}_i)^2 = \| \mathbf{y} - \hat{\mathbf{y}} \|^2 \tag{3.6} \end{equation}RSS=i=1∑n(yi−y^i)2=∥y−y^∥2(3.6)

TSS=∑i=1n(yi−yˉ)2=∥y−yˉ∥2\begin{equation} \mathrm{TSS} = \sum_{i=1}^n (y_i - \bar{y})^2 = \| \mathbf{y} - \bar{y} \|^2 \tag{3.7} \end{equation}TSS=i=1∑n(yi−yˉ)2=∥y−yˉ∥2(3.7)

上面的 yˉ=1n∑i=1nyi\displaystyle \bar{y} = \frac{1}{n} \sum_{i=1}^{n}y_iyˉ=n1i=1∑nyi。

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

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

我们会在下文中证明

TSS=ESS+RSS\begin{equation} \mathrm{TSS} = \mathrm{ESS} + \mathrm{RSS} \tag{3.8} \end{equation}TSS=ESS+RSS(3.8)

此外，值得注意的是：

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

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

adj.R2=1−RSSn−p−1TSSn−1=1−(1−R2)n−1n−p−1\begin{align*} \mathrm{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} \tag{3.9} \end{align*}adj.R2=1−n−1TSSn−p−1RSS=1−(1−R2)n−p−1n−1(3.9)

其中 nnn 是观测值的数量；而 ppp 是自变量的数量。调整 R2R^2R2 并不总是随着自变量数量的增多而增大。

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

1.2. 模型的 FFF 检验

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

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

E[ϵ∣X]=0\begin{equation} \mathbb{E}[\boldsymbol{\epsilon}|\mathbf{X}] = \boldsymbol{0} \tag{3.10} \end{equation}E[ϵ∣X]=0(3.10)

var[ϵ∣X]=σ2In\begin{equation} \text{var}[\boldsymbol{\epsilon}|\mathbf{X}] = \sigma^2 \mathbf{I}_n \tag{3.11} \end{equation}var[ϵ∣X]=σ2In(3.11)

ϵ∼N(0,σ2In)\begin{equation} \boldsymbol{\epsilon} \sim \mathcal{N}(\boldsymbol{0}, \sigma^2 \mathbf{I}_n) \tag{3.12} \end{equation}ϵ∼N(0,σ2In)(3.12)

rank(X)=p+1\begin{equation} \text{rank}(\mathbf{X}) = p + 1 \tag{3.13} \end{equation}rank(X)=p+1(3.13)

FFF 检验

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

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

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

Test Model:y=Xβ+ϵ\begin{equation} \text{Test Model} \, \text{:} \quad \mathbf{y} = \mathbf{X} \boldsymbol{\beta} + \boldsymbol{\epsilon} \tag{3.14} \end{equation}Test Model:y=Xβ+ϵ(3.14)

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

Test Model:yi=β0+β1xi1+β2xi2+⋯+βpxip+ϵi\begin{equation} \text{Test Model} \, \text{:} \quad y_i = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip} + \epsilon_i \tag{3.15} \end{equation}Test Model:yi=β0+β1xi1+β2xi2+⋯+βpxip+ϵi(3.15)

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

H0:β1=β2=⋯=βp=0\begin{equation} \mathrm{H}_0 \, \text{:} \quad \beta_1 = \beta_2 = \cdots = \beta_p = 0 \tag{3.16} \end{equation}H0:β1=β2=⋯=βp=0(3.16)

Null Model:yi=β0+ϵi\begin{equation} \text{Null Model} \, \text{:} \quad y_i = \beta_0 + \epsilon_i \tag{3.17} \end{equation}Null Model:yi=β0+ϵi(3.17)

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

Null Model:y=β0+ϵ\begin{equation} \text{Null Model} \, \text{:} \quad \mathbf{y} = \beta_0 + \boldsymbol{\epsilon} \tag{3.18} \end{equation}Null Model:y=β0+ϵ(3.18)

而备择假设 H1\mathrm{H}_1H1 为自变量整体上对模型的拟合效果有“显著”的贡献，不应该将所有自变量从模型中移除；其对应的模型记为扩展模型（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+β1xi1+β2xi2+⋯+βpxip+ϵ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。此外，在进行 FFF 检验之前，我们添加了式(3.10)到式(3.13)所描述的几个条件；那么式(3.18)和式(3.21)都有 ϵ∼N(0,σ2In)\boldsymbol{\epsilon} \sim \mathcal{N}(\boldsymbol{0}, \sigma^2 \mathbf{I}_n)ϵ∼N(0,σ2In)。

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

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

RSS(null)=TSS(ext)\begin{equation} \mathrm{RSS}^{\text{(null)}} = \mathrm{TSS}^{\text{(ext)}} \tag{3.22} \end{equation}RSS(null)=TSS(ext)(3.22)

于是，

ESS(ext)=TSS(ext)−RSS(ext)=RSS(null)−RSS(ext)\begin{align*} \mathrm{ESS}^{\text{(ext)}} &= \mathrm{TSS}^{\text{(ext)}} - \mathrm{RSS}^{\text{(ext)}} \\ &= \mathrm{RSS}^{\text{(null)}} - \mathrm{RSS}^{\text{(ext)}} \tag{3.23} \end{align*}ESS(ext)=TSS(ext)−RSS(ext)=RSS(null)−RSS(ext)(3.23)

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

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

ESS(ext)/pRSS(ext)/(n−p−1)∼Fp,n−p−1\begin{equation} \frac{ \mathrm{ESS}^{\text{(ext)}} / p} { \mathrm{RSS}^{\text{(ext)}} / (n-p-1) } \sim F_{p,n-p-1} \tag{3.24} \end{equation}RSS(ext)/(n−p−1)ESS(ext)/p∼Fp,n−p−1(3.24)

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

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

1.3. 回归系数的 ttt 检验

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

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

H0:βj=0\begin{equation} \mathrm{H}_0 \, \text{:} \quad \beta_j = 0 \tag{3.25} \end{equation}H0:βj=0(3.25)

注意这里 ttt 检验在检验回归系数是否为 0，也可以用类似的方法对偏移项 β0\beta_0β0 是否为 0 进行检验。

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

H1:βj≠0\begin{equation} \mathrm{H}_1 \, \text{:} \quad \beta_j \ne 0 \tag{3.26} \end{equation}H1:βj=0(3.26)

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

Tj=β^js2(XTX)jj−1∼tn−p−1\begin{equation} T_j = \frac{\hat{\beta}_j}{\sqrt{s^2(\mathbf{X}^\mathsf{T}\mathbf{X})^{-1}_{jj}}} \sim t_{n-p-1} \tag{3.27} \end{equation}Tj=s2(XTX)jj−1β^j∼tn−p−1(3.27)

其中

s2=(y−Xβ^)T(y−Xβ^)n−p−1\begin{equation} s^2 = \frac{(\mathbf{y} - \mathbf{X}\hat{\boldsymbol{\beta}})^{\mathsf{T}} (\mathbf{y} - \mathbf{X}\hat{\boldsymbol{\beta}}) }{n-p-1} \tag{3.28} \end{equation}s2=n−p−1(y−Xβ^)T(y−Xβ^)(3.28)

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

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

置信区间

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

(β^j−Tn−p−1−1(1+γ2)s.e.^(β^j),β^j+Tn−p−1−1(1+γ2)s.e.^(β^j))\begin{equation} \left( \hat{\beta}_j - \mathcal{T}^{-1}_{n-p-1} \left( \frac{1+\gamma}{2} \right) \, \widehat{\mathrm{s.e.}}(\hat{\beta}_j) \, , \ \hat{\beta}_j + \mathcal{T}^{-1}_{n-p-1} \left( \frac{1+\gamma}{2} \right) \, \widehat{\mathrm{s.e.}}(\hat{\beta}_j) \right) \tag{3.29} \end{equation}(β^j−Tn−p−1−1(21+γ)s.e.(β^j), β^j+Tn−p−1−1(21+γ)s.e.(β^j))(3.29)

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

s.e.^(β^j)=s2(XTX)jj−1\begin{equation} \widehat{\mathrm{s.e.}}(\hat{\beta}_j) = \sqrt{s^2(\mathbf{X}^\mathsf{T}\mathbf{X})^{-1}_{jj}} \tag{3.30} \end{equation}s.e.(β^j)=s2(XTX)jj−1(3.30)

Tm(⋅)\mathcal{T}_{m}(\cdot)Tm(⋅) 是自由度为 mmm 的 ttt 分布的累积分布函数；而 Tm−1()\mathcal{T}^{-1}_{m}()Tm−1() 是 Tm(⋅)\mathcal{T}_{m}(\cdot)Tm(⋅) 对应的分位数（quantile）函数。

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

2. 示例：mtcars数据集