如何在 R 中使用 lmer 设置自定义对比

机器算法验证 r lme4-nlme 对比
2022-03-15 17:01:37

我在 R 中使用 lmer 来检查条件 ( cond) 对某些结果的影响。这是一些组成的数据，其中 s 是主题标识符a，b和c是条件。

library("tidyr")
library("dplyr")
set.seed(123)
temp <- data.frame(s = paste0("S", 1:30), 
                   a = rnorm(30, -2, 1), 
                   b = rnorm(30, -3, 1), 
                   c = rnorm(30, -4, 1))

我想比较

  1. 水平a到水平的平均值b和c
  2. 水平b到水平c。

我的问题是，如何设置对比以使截距反映三个条件的平均值，而两个计算的估计值直接反映 1. 和 2. 中定义的差异？

我试过了

c1 <- cbind(c(-0.5, 0.25, 0.25), c(0, -0.5, 0.5))
gather(temp, cond, result, a, b, c) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = c1))

哪里cond2似乎可以，但cond1不是。

以下如何解释这些自定义对比？，我尝试使用广义逆，但这些估计也没有意义。

c2 <- t(ginv(c1))
gather(temp, cond, result, a, b, c) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = c2))

我也尝试过 Helmert 对比，但方法仍然不匹配。

gather(temp, cond, result, a, b, c) %>%
  mutate(cond = factor(cond, levels = c("c", "b", "a"))) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = contr.helmert))

这样做的正确方法是什么？

1个回答

对于以下步骤，我们需要长格式的数据框。数据框dat包含因变量result、分类预测变量cond（水平：a、b和c）和随机因子s。

library(tidyr)
dat <- gather(temp, cond, result, a, b, c)

在下文中，我将说明两种方法来创建与您要比较的条件相对应的对比矩阵：

  1. a−b+c2
  2. b−c

自定义对比

该矩阵mat对应于水平差异。

mat <- rbind(c(1, -0.5, -0.5),     # a vs. (b + c) / 2
             c(0, 1, -1))          # b vs. c

为了创建实际的对比矩阵，我们用ginv(from MASS) 计算广义逆矩阵。

library(MASS)
cMat <- ginv(mat)
#            [,1]          [,2]
# [1,]  0.6666667 -7.130169e-17
# [2,] -0.3333333  5.000000e-01
# [3,] -0.3333333 -5.000000e-01

此对比矩阵cMat可用于lmer.

library(lme4)
res <- lmer(result ~ cond + (1|s), data = dat, 
            contrasts = list(cond = cMat))
coef(summary(res))    
#              Estimate Std. Error    t value
# (Intercept) -2.948115  0.0946025 -31.163182
# cond1        1.351517  0.2006822   6.734612
# cond2        1.153918  0.2317279   4.979625

如您所见，固定效应估计对应于上面指定的差异。此外，截距代表整体平均值。

赫尔默特对比contr.helmert

您还可以使用内置contr.helmert函数来创建对比度矩阵。

cHelmert <- contr.helmert(3)
#   [,1] [,2]
# 1   -1   -1
# 2    1   -1
# 3    0    2

但是，该顺序与您在问题中指定的顺序不对应。因此，我们必须颠倒列和行的顺序。第一列对应于bvs. a，第二列对应于 vs.和c的均值。ba

cHelmert2 <- cHelmert[c(3:1), 2:1]
#   [,1] [,2]
# 3    2    0
# 2   -1    1
# 1   -1   -1

将对比矩阵cHelmert2与进行比较cMat。您会注意到这些列是另一个矩阵的缩放版本。

结果lmer是： 

library(lme4)
res2 <- lmer(result ~ cond + (1|s), data = dat, 
             contrasts = list(cond = cHelmert2))
coef(summary(res2))    
#               Estimate Std. Error    t value
# (Intercept) -2.9481150 0.09460250 -31.163182
# cond1        0.4505056 0.06689407   6.734612
# cond2        0.5769590 0.11586393   4.979625

此对比矩阵允许进行与自定义对比矩阵相同的比较。但是，由于矩阵中的值不同，固定效应系数也不同。不出所料，t-值相同。

其它你可能感兴趣的问题
上一篇为什么离散度的度量不如中心度那么直观？ 下一篇R：增强型迪基富勒 (ADF) 测试