多层回归模型简介

深藏阁楼爱情的钟 2022-08-05 14:16 462阅读 0赞

多层回归模型简介

多层回归模型（Multi-level model）中有很多容易混淆的概念，因为很多概念是来源于不同的专业背景。首先让我们先罗列这些名词进行区分，再来R语言来举例。  
  
**多层回归模型**通常涉及到对同一个体进行反复测量，这样得到的数据就不再相互独立而是存在某种相关性，所以普通线性回归不再适用。当这种反复测量是在不同时点上进行时，这就称为 **面板数据分析**（panel data analysis）或者 **纵向数据分析**（longitudinal data analysis）。  
  
  
**Fixed Effect：固定效应**，也就是普通线性回归中处理的预测变量，它在模型中对响应变量的期望造成影响，其因子水平包含明确的信息，通常是在实验中加以控制的因素。  
  
**Random Effect：随机效应**，也就是一个随机变量，估计它是没有意义的。它在模型中对响应变量的方差造成影响，其因子水平包含信息不清楚，通常在实验中无法控制的因素，例如在面板数据分析中，不同的个体差异就是随机效应。  
  
**Mix effect model：混合效应模型**，模型中包含了fixed effect和random effect，根据random effect的影响，又区别为对模型截距的影响（random intercept）和对模型截距与斜率的影响（random intercept and slope）。  
  
假设一个面板数据中个体之间存在差异，那么它是一个random effect。它的影响分两种情况，一种是只影响截距，即是它与模型中其它预测变量相加，如果对不同个体分别作线性回归，那么得到回归线截距会不同，但回归线平行，此时又称为 **固定效应回归模型（Fixed Effects Regression Model）**。另一种同时影响模型截距与斜率，是指它还与其它变量相乘，那么分别得到的回归线截距和斜率都不同，此时又称为 **随机效应模型（Random Effects Regression Model）**。  
  
如果模型中不考虑random effect，也就是认为个体和时间都没有显著性差异，此时模型退化为 **混合估计模型（Pooled Regression Model）**，此时可以直接把面板数据混合在一起用普通线性回归估计参数。  
  
下面我们以faraway包中的psid数据来举例，该数据集是对美国人的收入情况进行调查所得到的，其中包括了年龄、教育、性别、时间和个体ID这几个变量，我们希望了解这些因素对收入的影响。  
  age educ sex income year person  
1  31   12   M   6000   68      1  
2  31   12   M   5300   69      1  
3  31   12   M   5200   70      1  
4  31   12   M   6900   71      1  
5  31   12   M   7500   72      1  
6  31   12   M   8000   73      1  
  
如果假设认为这些调查对象是同质的，也就是个体间没有差异性，那么可以将数据 **完全汇集（complete pooling）**到一起，直接利用lm函数进行回归。但这个混合效应模型的同质假设往往不成立，数据汇集导致过度简化。  
  
另一种思路是假设研究的异质性，将不同的个体分别进行回归，从而得到针对特定个体的估计值，这称为 **不汇集（no pooling）**。但这种方法导致每个回归所用到的样本减少，从而难以估计统计量的标准差。  
  
多层回归模型的思路是前两者的折中，所以又称为 **部分汇集（partial pooling）**。在R语言中我们使用 **lme4包**中的 **lmer函数**来完成这项工作。首先载入faraway包以便读取psid数据集，然后加载mgcv包，再将年份数据中心化以方便解释模型，最后用lmer函数进行建模。  
\-----------------  
library(faraway)  
library(lme4)  
psid$cyear <- psid$year-78   
model1=lmer(log(income) ~ cyear\*sex +age+educ+(cyear|person),psid)   
\-----------------  
lmer函数使用和lm是类似的，一般变量表示固定效应，括号内竖线右侧的person表示它是一个随机效应，它与模型中其它变量相加，而且与年份cyear变量相乘，影响其斜率。这就是一个随机效应模型。如果认为随机效应只影响模型截距，那么固定效应回归模型可以用下面的公式  
\-----------------  
model2=lmer(log(income) ~ cyear\*sex +age+educ+(1|person),psid)  
\-----------------  
为了判断哪一个模型更为合适，可以使用anova函数，从结果中观察到P值很小，判断应当使用model1  
\-----------------  
anova(model1,model2)  
  
Data: psid  
Models:  
model2: log(income) ~ cyear \* sex + age + educ + (1 | person)  
model1: log(income) ~ cyear \* sex + age + educ + (cyear | person)  
       Df    AIC    BIC  logLik  Chisq Chi Df  
model2  8 3943.9 3987.2 -1963.9                
model1 10 3805.6 3859.7 -1892.8 142.27      2  
       Pr(>Chisq)      
model2                 
model1  < 2.2e-16 \*\*\*  
\-----------------  
  
得到的模型结果还可以用各种泛型函数如summary、predict、resid进行进一步处理。当响应变量不符合正态分布假设时，还可以用 **广义多层回归进行(glmer)**建模  
  
参考资料：  
环境与生态统计--R语言的应用  
The R Book  
A\_Handbook\_of\_Statistical\_Analyses\_Using\_R  
Extending\_the\_Linear\_Model\_with\_R\_\_Generalized\_Linear\_\_Mixed\_Effects\_and\_Nonparametric\_Regression\_Model