多元线性回归

文基于MT5753 统计建模复习笔记

谢绝任何形式的转载

第二章 多元线性回归

1. 建模

多元线性回归是建立在一元简单线性回归的基础上的，只不过多引入了一些变量。yi是第i个响应变量，β为协方差，x1i是第i个第一个解释变量，以此类推,xki就是第i个第k个解释变量，e为误差项，同样需要服从均值为0的正态分布。同样地，这个公式可以通过线性代数写成更统一的矩阵相乘的形式：Y=βX+e.

这里我们依然使用最小二乘来对β进行估计。最小二乘有一个特性，叫做“最佳线性无偏估计” （Best Linear Unbiased Predictors, BLUE):

（1）估计值是无偏的。例如，估计值的分布是集中于真值的。

（2）估计值是一致的。根据大数定理，如果我们增加样本数量，估计值应该越接近于真值。

（3）估计值是目前最佳的——因为它是无偏的。

最佳线性无偏预测

BLUP best linear unbiased prediction

等同于固定效应（fixed effects）的最佳线性无偏估计

对于得到的估计参数β，我们需要用下式置信区间来评估（任何不用置信区间的估计值都是耍流氓）。注意：双尾所以α除以2，自由度为n-p-1,p为参数数量，依旧不包括β0；在样本数量超过30,000时，t分布于正态分布基本相等。

对于置信区间包含0的参数，这个参数基本就是废了（但我们依然可以得出结论：这个模型解释性不强）。例如估计出值为1的β1，其置信区间为[-1.5, 2]，它可以是下图中的任何一种可能。像图3这样，根本就是不存在线性关系了。

置信区间评估完参数之后，一个重要环节就是查看参数的显著性。这里有一个误区就是觉得p值非常小就说明解释变量和响应变量无关。其实有可能是解释变量对响应变量的影响比较小。每个变量的影响都是有对应的权重的，显著性检验只是在看这个权重的大小。

在模型中，我们常常会引入年份、月份等因子型变量，线性回归会自动将因子型变量的第一个选为基底baseline。例如，我们有2009-2017年9个因子，建模时，2009年会自动变成baseline（baseline的选择在代码里是可以调整的），其他8个因子的显著性和参数估计都是根据2009这个baseline进行计算的。如果我们改变了baseline，将2010作为baseline，那么其他因子的参数估计和显著性也会随之改变。在模型中，baseline是不会出现的。例如我们理解的模型应该是y = β0+β1Year(2009)+β2Year(2010)+...，其实它是y =β0+β1Year(2010)+β2Year(2011)+...

不同于一元模型，多元模型我们可以引入交叉项：interaction terms。比如在现有模型基础上，我们想知道x1和x2对于响应变量有没有交叉影响，就可以引入x1:x2这样的交互项。用公式表达就是：

y=β0+β1x1+β2x2+β3x1:x2

2. 模型选择

建立好模型之后，我们需要考虑以下几个问题：

（1）多重共线性 collinearity

（2）变量自动选择 automated variable selection

（1）多重共线性

在线性回归模型中，解释变量之间存在着高度相关或精确相关的关系，这使得模型估计失真或者难以准确估计。想象一个y = β0+β1x1+β2x2+e的模型，如下图所示。如果x1和x2之间没有共线性，那么这个三维坐标轴无论如何旋转，y的值都不会受到影响——因为x1和x2轴是独立垂直的。如果两者存在共线性，那么模型会随着旋转而改变。有兴趣的同学可以拿自己的手臂模拟试试。

那么如何探查共线性呢，一个简单的判断是大致可以从变量意义上判断，比如x1是年份，x2是月份，他们都是时间相关，就是极有可能出现共线性的。还有就是当我们在模型中加入一个变量时，模型估计参数变化明显，就可能说明这个变量和模型中已存在的变量有共线性。

正式的判断我们可以用方差膨胀因子VIF来判断。VIF越大，共线性越严重。这里没有统一的标准VIF应该多大表示共线性严重。一般认为：当 0 =100，存在严重共线性。

方差膨胀因子 Variance inflation factors (VIFs)

方差膨胀因子是指解释变量之间存在多重共线性时的方差与不存在时的方差之比容忍度的倒数。

——摘自百度百科《方差膨胀因子》【1】

R代码：

library(car)

vif(modelName)

最后一点，注意共线性和自相关不是一个东西。

（2）变量自动选择 automated variable selection

如何处理共线性呢？一个最简单的方法是删掉这个变量，但是这个方法显然是不科学的。还有就是使用automated variable seleciton. 我们分为三个方法：

a. forward selection 前向选择

从只引入截距开始，每次加一个变量，然后通过一些标准判断这个变量的去留。

b. backward elimination 后向消除

一次引入所有变量，然后开始逐个删除，通过一些标准判断这个变量的去留。

c. stepwise selection 逐步选择

结合前向选择和后向消除。当我们使用前向选择选择了两个变量时，下一步不是继续添加变量了，而是使用后向消除的检验方法来检查当前模型的变量是否显著，不显著的删掉。（因为向前选择第二步可能使得之前的变量不再显著了，所以要在检验一下整个模型）。它的停止规则也很简单：在某一阶段，模型中的每个变量都是显著的，没在模型中的每个变量都是不显著的。

上面我们说到这三个方法都是通过一些标准来判断模型变量的引入，这个标准有很多种，我们使用AIC或BIC。

AIC：最小信息准则，计算公式如下。log-likelihood value是通过模型拟然估计求得的值，p依旧是参数的个数。

AIC 越小越好，但是对于1234和1238这样的AIC值，并没有什么实质性的区别，并不能说明1234的模型比1238好，我们应该结合其他衡量标准来判断。

AIC衡量的是相对模型，而不是绝对模型。意思就是，AIC只会在现有的模型里选择它觉得最好的，而不会给我们提供全局最优解。如下图所示，AIC只是在有偏bias和高方差 variance之间帮我们找到一个参数para的平衡点。

R代码：

AIC(modelName)

step(modelName, direction = "both")

step 的default 判断标准是AIC

在出现在第一行时，代表模型选择结束

除了AIC，我们还有BIC：贝叶斯信息准则。L和p同AIC。BIC依旧是越小越好。

在R代码方面，step使用BIC需要加上“k=log(nrow(y))”（AIC的k=2，为默认值）。

对于AIC和BIC，用一个就足够了，两者没有太大区别。但是由于AIC可能会存在过拟合问题，所以我个人更倾向于BIC。

3. 参数检验ANOVA

在线性回归模型中，F检验（ANOVA）十分重要。其原假设H0为：所有估计参数都为0（除截距）；H1为至少有一个参数不为0（除截距）。ESS为模型无法解释的内容（例如残差，误差等），RSS为模型可以解释的内容（如回归项等）。p是全模型full model的参数个数，q是削减模型 reduced model的参数个数。

——摘自5753统计建模课件，侵删

R代码：

anova(Model1,Model2)

或

library(car)

Anova(model1, model2)

3. 模型评估

（1）adj-R square

不同于之前的R方，这里我们使用调整后的R方（adjusted R square）。好处就是比R方更稳定。当我们每次引入变量，R方都会不同程度的增加，这样过拟合等问题都没办法体现，调整后的R方可以在一定程度上回避这个问题。

（2）有影响的值 influential points

大致分为三种，如下图所示。蓝点是观测值，红点为观测值中不正常的y，绿点为不正常的x,黑点是二者结合体。

探查这些特殊值，我们可以用The Cook's Distance. S是标准差，p是除去截距的参数数量，yk是当除去第k个观察值而建成的模型对于第k个观察值的预测。有点拗口，附上原话吧："yk is the fitted values for the t-th observation when the kth observation is omitted from the model fitting process"——摘自MT5753 统计建模。

下图就是一个实例。可以看出来有三个点是远远脱离大部分散点的区域的，他们可能就是influential points.

——摘自MT5753 统计建模，侵删

（3）模型假设

同一元简单线性模型一样。

a. 线性性

对于检查线性性（无论有无interactions），最好使用局部残差（partial residuals）。r为原本的残差，这里我觉得应该是使用标准化残差。后缀项是对应的协方差和解释变量。

相比原来的残差图，它有这么几个好处：

i. 它的斜率直接反应了我们的估计参数

ii. 散点图能反应变量x的重要性

iii. 极大残差值可以比较容易被发现

不同于原本的拟合值-残差图，这里我们引入的线性量，所以如果符合线性性，拟合值-局部残差图应该是线性的。

b. 同方差

同一元，Breusch-Pagan test. 对于同方差的模型，不应该出现均值-方差关系（mean-variance relationship）。

正态分布和独立性同一元，就不再赘述。

R代码：

model1

model2

.表示保留原有的东西，第一个.是保留y,第二个点是保留x1

update也可以用于删除变量

y~x1*x2 = y~ x1+x2+x1:x2

参考文献

【1】https://baike.baidu.com/item/%E6%96%B9%E5%B7%AE%E8%86%A8%E8%83%80%E5%9B%A0%E5%AD%90/4993652?fr=aladdin