原文链接：http://tecdat.cn/?p=3015
介绍
首先，请留神，围绕多层次模型的术语有很大的不一致性。例如，多层次模型自身可能被称为分层线性模型、随机效应模型、多层次模型、随机截距模型、随机斜率模型或汇合模型。依据不同的学科、应用的软件和学术文献，许多这些术语可能指的是雷同的个别建模策略。
在本文中，我将试图通过演示如何在 R 中拟合多层次模型，并试图将模型拟合过程与无关这些模型的罕用术语分割起来，为用户提供一个多层次模型的指南。
读入数据
多层次模型适宜于一种非凡的数据结构，即单位嵌套在组内（个别是 5 个以上的组），咱们想对数据的组构造进行建模。
## id extro open agree social class school
## 1 1 63.69 43.43 38.03 75.06 d IV
## 2 2 69.48 46.87 31.49 98.13 a VI
## 3 3 79.74 32.27 40.21 116.34 d VI
## 4 4 62.97 44.41 30.51 90.47 c IV
## 5 5 64.25 36.86 37.44 98.52 d IV
## 6 6 50.97 46.26 38.83 75.22 d I在这里，咱们有对于嵌套在班级和学校内的科目的内向性的数据。
在咱们开始之前，让咱们先理解一下数据的构造。
## 'data.frame': 1200 obs. of 7 variables:
## $ id : int 1 2 3 4 5 6 7 8 9 10 ...
## $ extro : num 63.7 69.5 79.7 63 64.2 ...
## $ open : num 43.4 46.9 32.3 44.4 36.9 ...
## $ agree : num 38 31.5 40.2 30.5 37.4 ...
## $ social: num 75.1 98.1 116.3 90.5 98.5 ...
## $ class : Factor w/ 4 levels "a","b","c","d": 4 1 4 3 4 4 4 4 1 2 ...
## $ school: Factor w/ 6 levels "I","II","III",..: 4 6 6 4 4 1 3 4 3 1 ...这里咱们看到咱们有两个可能的分组变量 – 班级和学校。让咱们进一步探讨一下它们。
##
## a b c d
## 300 300 300 300
##
## I II III IV V VI
## 200 200 200 200 200 200
##
## I II III IV V VI
## a 50 50 50 50 50 50
## b 50 50 50 50 50 50
## c 50 50 50 50 50 50
## d 50 50 50 50 50 50这是一个齐全均衡的数据集。让咱们来绘制一下数据。能够摸索学校和班级之间的变量关系。
plot(extro ~ open + social + agree | school + class
在这里咱们看到，在班级外部有显著的分层，咱们也看到社会变量与凋谢和同意变量有强烈的区别。咱们还看到，A 班和 D 班在最低和最高段的散布显著更多。接下来让咱们按学校绘制数据。
按学校划分，咱们看到学生的分组很严密，但学校一和学校六比其余学校显示出显著的分散性。在学校之间，咱们的预测指标的模式与班级之间的模式雷同。
这里咱们能够看到，学校和班级仿佛亲密辨别了咱们的预测因素和内向性之间的关系。
摸索 merMod 对象
咱们对咱们的嵌套数据拟合了一系列的随机截距模型。咱们将更深刻地钻研咱们拟合这个模型时产生的 lmerMod 对象，以理解如何在 R 中应用混合效应模型。咱们首先拟合上面这个按班级分组的根本例子。
## \[1\] "lmerMod"
## attr(,"package")
## \[1\] "lme4"首先，咱们看到 MLexamp1 当初是一个 lmerMod 类的 R 对象。这个 lmerMod 对象是一个 S4 类，为了摸索其构造，咱们应用 slotNames。
## \[1\] "resp" "Gp" "call" "frame" "flist" "cnms" "lower"
## \[8\] "theta" "beta" "u" "devcomp" "pp" "optinfo"在 lmerMod 对象中，咱们看到了一些咱们可能心愿摸索的对象。比如说。
## lmer(formula = extro ~ open + agree + social + (1 | school),
## data = lmm.data)
## \[1\] 59.116514 0.009751 0.027788 -0.002151
## \[1\] "data.frame"
## extro open agree social school
## 1 63.69 43.43 38.03 75.06 IV
## 2 69.48 46.87 31.49 98.13 VI
## 3 79.74 32.27 40.21 116.34 VI
## 4 62.97 44.41 30.51 90.47 IV
## 5 64.25 36.86 37.44 98.52 IV
## 6 50.97 46.26 38.83 75.22 ImerMod 对象有许多可用的办法 – 在此无奈一一列举。然而，咱们将在上面的列表中介绍一些比拟常见的办法。
## \[1\] anova.merMod* as.function.merMod* coef.merMod*
## \[4\] confint.merMod deviance.merMod* df.residual.merMod*
## \[7\] drop1.merMod* extractAIC.merMod* extractDIC.merMod*
## \[10\] family.merMod* fitted.merMod* fixef.merMod*
## \[13\] formula.merMod* isGLMM.merMod* isLMM.merMod*
## \[16\] isNLMM.merMod* isREML.merMod* logLik.merMod*
## \[19\] model.frame.merMod* model.matrix.merMod* ngrps.merMod*
## \[22\] nobs.merMod* plot.merMod* predict.merMod*
## \[25\] print.merMod* profile.merMod* qqmath.merMod*
## \[28\] ranef.merMod* refit.merMod* refitML.merMod*
## \[31\] residuals.merMod* sigma.merMod* simulate.merMod*
## \[34\] summary.merMod* terms.merMod* update.merMod*
## \[37\] VarCorr.merMod* vcov.merMod weights.merMod*
##
## Non-visible functions are asterisked一个常见的需要是从 merMod 对象中提取固定效应。fixef 提取一个固定效应的命名数字向量，这很不便。
## (Intercept) open agree social
## 59.116514 0.009751 0.027788 -0.002151理解这些参数的 p 值或统计显着性：
如果你想理解这些参数的 P 值或统计学意义，首先通过运行 “mcmcsamp “ 查阅 lme4 帮忙，理解各种办法。包中内置的一个不便的办法是。
## 0.5 % 99.5 %
## .sig01 4.91840 23.88758
## .sigma 2.53287 2.81456
## (Intercept) 46.27751 71.95610
## open -0.02465 0.04415
## agree -0.01164 0.06721
## social -0.01493 0.01063从这里咱们能够首先看到咱们的固定成果参数重叠 0 示意没有成果。咱们还能够看到.sig01，这是咱们对随机效应变动的预计。这表明咱们的群体之间的影响很小，要么失去更准确的预计的群体太少。
另一个常见的需要是提取残余标准误差，这是计算成果大小所必须的。要取得残差规范误的命名向量：
## \[1\] 2.671例如，教育钻研中的常见做法是通过将固定效应参数除以残差标准误差来将固定成果标准化为“成果大小”：
## (Intercept) open agree social
## 22.1336705 0.0036508 0.0104042 -0.0008055从这一点，咱们能够看到，咱们对开放性，同意性和社交性的预测因素在预测内向性方面简直毫无用处。让咱们把注意力转向下一个随机效应。
摸索组变动和随机成果
您很可能适宜混合成果模型，因为您间接对模型中的组级变动感兴趣。目前还不分明如何从后果中摸索这种群体程度的变动summary.merMod。咱们从这个输入失去的是组效应的方差和标准偏差。
## $school
## (Intercept)
## I -14.091
## II -6.183
## III -1.971
## IV 1.966
## V 6.331
## VI 13.948运行 ranef 函数能够失去每个学校的截距，但没有太多的额定信息 – 例如这些估计值的精确度。咱们须要一些额定的命令。
## \[1\] "ranef.mer"
## , , 1
##
## \[,1\]
## \[1,\] 0.03565
##
## , , 2
##
## \[,1\]
## \[1,\] 0.03565
##
## , , 3
##
## \[,1\]
## \[1,\] 0.03565
##
## , , 4
##
## \[,1\]
## \[1,\] 0.03565
##
## , , 5
##
## \[,1\]
## \[1,\] 0.03565
##
## , , 6
##
## \[,1\]
## \[1,\] 0.03565ranef.mer 对象是一个列表，其中蕴含每个组程度的数据框。数据框蕴含每个组的随机效应（这里咱们只有每个学校的截距）。lme4 提供随机效应的条件方差时，它贮存在这些数据框的一个属性中，作为方差 - 协方差矩阵的一个列表。
这种构造确实很简单，但它很弱小，因为它容许嵌套的、分组的和跨级别的随机效应。
``````
## $school
## (Intercept)
## I -14.091
## II -6.183
## III -1.971
## IV 1.966
## V 6.331
## VI 13.948
##
## with conditional variances for "school"
## $school
这个图形显示了随机效应项的点阵图，也被称为毛毛虫图。在这里，你能够分明地看到每个学校对外向性的影响，以及它们的标准误差，以帮忙确定随机效应之间的区别。
应用模仿和图来摸索随机效应
一种常见的计量经济学办法是创立所谓的群体程度项的教训贝叶斯预计。然而，对于什么是随机效应项的适当标准误差，甚至对于如何统一地定义教训贝叶斯估计值，目前还没有太多的共识。在 R 中，有一些额定的办法能够取得随机效应的估计值。
## X1 X2 mean median sd
## 1 I (Intercept) -14.053 -14.093 3.990
## 2 II (Intercept) -6.141 -6.122 3.988
## 3 III (Intercept) -1.934 -1.987 3.981
## 4 IV (Intercept) 2.016 2.051 3.993
## 5 V (Intercept) 6.378 6.326 3.981
## 6 VI (Intercept) 13.992 14.011 3.987REsim 函数为每个学校返回程度名称 X1、预计名称 X2、估计值的平均值、中位数和估计值的标准差。
另一个不便的函数能够帮忙咱们绘制这些后果，看看它们与 dotplot 的后果相比如何。
ymin = "ymin")) +
geom\_pointrange() + theme\_dpi() +
theme(panel.grid.major = element\_blank(), panel.grid.minor = element\_blank(),
axis.text.x = element\_blank(), axis.ticks.x = element\_blank()) +
geom_hline(yintercept = 0, color = I("red"), size = I(1.1))
}
依据你的数据收集形式和你的钻研问题，预计这些效应大小的其余办法是可能的。
咱们能够测试是否蕴含随机效应来进步模型的拟合度，咱们能够应用基于模仿的似然比测试来评估额定的随机效应项的 p 值。
##
## simulated finite sample distribution of LRT. (p-value based on
## 10000 simulated values)
##
## data:
## LRT = 2958, p-value < 2.2e-16这里程序收回了一个正告，因为咱们最后是用 REML 而不是齐全最大似然来拟合模型。但 lme4 中的 refitML 函数容许咱们应用齐全最大似然法轻松地从新拟合咱们的模型。
##
## simulated finite sample distribution of LRT. (p-value based on
## 10000 simulated values)
##
## data:
## LRT = 2958, p-value < 2.2e-16在这里，咱们能够看到，即便每个独自的群体的影响可能是实质性的小和 / 或不准确的测量，但纳入咱们的分组变量是显著的。这对于了解模型是很重要的。
随机效应有什么意义？
如何解释咱们的随机效应的实质性影响？在应用多层次构造的察看工作中，这往往是至关重要的，理解分组对个体察看的影响。为了做到这一点，咱们抉择了 12 个随机案例，而后模仿他们的预测值，如果他们被安顿在 6 所学校中的每一所，他们的 extro 预测值。请留神，这是一个非常简单的模仿，只是应用了固定效应的平均值和随机效应的条件模式，并没有进行复制或抽样。
qplot(school, yhat, data = simDat) + facet\_wrap(~case) + theme\_dpi()
这张图通知咱们，在每个图中，代表一个案例，不同的学校有很大的变动。因而，把每个学生移到不同的学校，对外向性得分有很大的影响。然而，每个案例在每个学校都有差别吗？
+ facet\_wrap(~school) + theme\_dpi() +
theme(axis.text.x = element_blank())
在这里，咱们能够分明地看到，在每所学校内，案例绝对雷同，表明群体效应大于个体效应。
这些图在以实质性的形式展现群体和个体效应的绝对重要性方面很有用。甚至还能够做更多的事件来使图表更有信息量，比方提到后果的总变异性，也能够看一下挪动群体使每个察看值离其实在值的间隔。
论断
R 为解决混合效应模型提供了一个十分弱小的面向对象的工具集。了解模型拟合和置信区间须要一些钻研。
参考文献
1. 基于 R 语言的 lmer 混合线性回归模型
2.R 语言用 Rshiny 摸索 lme4 狭义线性混合模型（GLMM）和线性混合模型（LMM）
3.R 语言线性混合效应模型实战案例
4.R 语言线性混合效应模型实战案例 2
5.R 语言线性混合效应模型实战案例
6. 线性混合效应模型 Linear Mixed-Effects Models 的局部折叠 Gibbs 采样
7.R 语言 LME4 混合效应模型钻研老师的受欢迎水平
8.R 语言中基于混合数据抽样 (MIDAS) 回归的 HAR-RV 模型预测 GDP 增长
9. 应用 SAS，Stata，HLM，R，SPSS 和 Mplus 的分层线性模型 HLM