在 R 中使用混合效应模型分析纵向数据答案

【问题标题】：Analyze longitudinal data with a mixed effects model in R在 R 中使用混合效应模型分析纵向数据
【发布时间】：2016-10-13 15:16:30
【问题描述】：

我尝试使用混合效应模型（lme4 包）分析 R 中的一些模拟纵向数据。

模拟数据：25 名受试者必须在 5 个连续时间点执行 2 项任务。

#Simulate longitudinal data
N    <- 25
t    <- 5
x    <- rep(1:t,N)

#task1
beta1 <- 4
e1    <- rnorm(N*t, mean = 0, sd = 1.5)
y1    <- 1 + x * beta1 + e1

#task2
beta2 <- 1.5
e2    <- rnorm(N*t, mean = 0, sd = 1)
y2    <- 1 + x * beta2 + e2

data1 <- data.frame(id=factor(rep(1:N, each=t)), day = x, y = y1, task=rep(c("task1"),length(y1)))
data2 <- data.frame(id=factor(rep(1:N, each=t)), day = x, y = y2, task=rep(c("task2"),length(y2)))
data <- rbind(data1, data2)

问题1：如何分析受试者如何学习每项任务？

library(lme4)
m1 <- lmer(y ~ day + (1 | id), data=data1)
summary(m1)

...
Fixed effects:
            Estimate Std. Error       df t value Pr(>|t|)    
(Intercept)   1.2757     0.3561 123.0000   3.582 0.000489 ***
day           3.9299     0.1074 123.0000  36.603  < 2e-16 ***

使用ranef(m1)，我得到了每个科目的随机截距，我认为这反映了每个科目在第 1 天时的基线值。但我不明白如何判断个人如何学习任务，或者是否科目在学习任务的方式上有所不同。

问题2：如何分析task1和task2的受试者学习方式是否不同。

【问题讨论】：

标签： r lme4 mixed-models

【解决方案1】：

我扩展了您的示例以简要回答您的问题，但我建议您阅读 Snijders & Bosker (2012) 的第 15 章或 Singer & Willet (2003) 的书以获得更好的解释。天在您的模型中被视为连续变量，因为您有面板数据（即每个人都在同一天进行测量）并且除了指示不同的测量场合之外，天没有任何意义，最好将天视为一个因素（即使用虚拟变量）。但是，现在我将继续您的示例

您的第一个模型（我认为您希望 data1 的数据插入）给出了固定的线性斜率（即平均斜率，任务没有差异，个体之间没有差异）。固定截距是 day 为 0 时的性能，没有任何意义，因此您可能需要考虑将 day 的效果居中以获得更好的解释（或者实际上使用假人）。随机效应使个体偏离该截距，在您的示例中估计方差为 0.00，因此个体在起始位置上几乎没有区别。

m1 <- lmer(y ~ day + (1 | id), data=data)
summary(m1)
Random effects:
 Groups   Name        Variance Std.Dev.
 id       (Intercept)  0.00    0.000   
 Residual             18.54    4.306   
Number of obs: 250, groups:  id, 25

我们可以通过添加与任务的交互来扩展此模型。这意味着 task1 和 task2 的固定斜率不同，我相信可以回答问题 2（您也可以使用 update() 更新您的模型）

m2 <- lmer(y ~ day*task + (1|id), data = data)
summary(m2)

该模型中日的影响是您的参考类别（task1）的固定斜率，交互作用是 task1 和 task2 的斜率之间的差异。任务的固定效果是截距的差异。模型拟合可以通过偏差检验来评估，请阅读 Snijders & Boskers (2012) 了解 ML 和 REML 估计的解释。

anova(m1,m2)

为了给个体的成长添加随机效应，我们可以再次更新模型，这回答了问题 1

m3 <- lmer(y ~ day*task + (day|id), data = data)
summary(m3)
ranef(m3)

随机效应表示斜率和截距的个体偏差。模型摘要中包含随机效应分布的摘要（与 m1 相同）。最后，我认为您可以在日常任务交互中添加随机效应，以评估个人在任务 1 和任务 2 上的绩效增长是否存在差异。但这在很大程度上取决于您的数据和以前模型的性能。

m4 <- lmer(y ~ day*task + (day*task|id), data = data)
summary(m4)
ranef(m4)

希望这会有所帮助。我推荐的书当然应该。两者都提供了出色的示例和理论解释（遗憾的是没有 R 示例）。如果您决定使用固定场合模型（由假人表示的一天的影响），nlme 包提供了很好的选项来控制随机效应的协方差结构。 Pinheiro & Bates (2000) 提供了该软件包的良好文档。

【讨论】：