我试图从简单的配对中排除对比。
虚构数据：
library(emmeans)
library(lme4)
set.seed(1234)
dat <- data.frame(
dv = c(
# Add in t1 ctrl
rlnorm(meanlog=0.8, sdlog=0.1, n=1000), #Long Healthy
rlnorm(meanlog=0.7, sdlog=0.1, n=1000), #Lat Healthy
rlnorm(meanlog=0.6, sdlog=0.1, n=500), #Long Damaged
rlnorm(meanlog=0.5, sdlog=0.1, n=500), #Lat Damaged
# Add in t2 ctrl
rlnorm(meanlog=0.7, sdlog=0.1, n=1000), #Long Healthy
rlnorm(meanlog=0.6, sdlog=0.1, n=1000), #Lat Healthy
rlnorm(meanlog=0.5, sdlog=0.1, n=500), #Long Damaged
rlnorm(meanlog=0.4, sdlog=0.1, n=500), #Lat Damaged
# Add in t1 trt
rlnorm(meanlog=0.8, sdlog=0.1, n=1000), #Long Healthy
rlnorm(meanlog=0.7, sdlog=0.1, n=1000), #Lat Healthy
rlnorm(meanlog=0.6, sdlog=0.15, n=500), #Long Damaged
rlnorm(meanlog=0.5, sdlog=0.15, n=500), #Lat Damaged
# Add in t2 trt
rlnorm(meanlog=0.7, sdlog=0.1, n=1000), #Long Healthy
rlnorm(meanlog=0.6, sdlog=0.1, n=1000), #Lat Healthy
rlnorm(meanlog=0.65, sdlog=0.15, n=500),#Long Damaged
rlnorm(meanlog=0.55, sdlog=0.15, n=500) #Lat Damaged
),
id=c(rep(c("subj_1", "subj_2"), times=c(6000, 6000))),
intervention=c(rep(c("ctrl", "trt"), times=c(6000, 6000))),
timepoint=c(rep(rep(c("t1", "t2"), times=c(3000, 3000)),2)),
direction=c(rep(rep(c("long", "lat", "long", "lat"), times=c(1000, 1000, 500, 500)),4)),
region=c(rep(rep(c("healthy", "damaged"), times=c(2000, 1000)),4))
)
型号：
lmm_1 <- lmer(log(dv) ~ intervention*timepoint*region + direction + (1|id), data=dat)
emmeans参考网格：
lmm_1_emm <- emmeans(lmm_1,
pairwise ~ intervention*region*timepoint,
type = "response",
bias.adj=TRUE,
sigma=sigma(lmm_1))
简单配对：
pairs(regrid(lmm_1_emm,
bias.adjust = TRUE,
sigma = sigma(lmm_1)),
simple = "each",
combine = TRUE)
I bet you wanted to call this with just object[[1]] - use '[[]]' or which' if I'm wrong.
See '? emm_list' for more information
region timepoint intervention contrast estimate SE df z.ratio p.value
damaged t1 . ctrl - trt -0.00263 0.06162 Inf -0.043 1.0000
healthy t1 . ctrl - trt -0.00554 0.07468 Inf -0.074 1.0000
damaged t2 . ctrl - trt -0.26043 0.06057 Inf -4.300 0.0002
healthy t2 . ctrl - trt -0.00388 0.06743 Inf -0.057 1.0000
. t1 ctrl damaged - healthy -0.37968 0.01228 Inf -30.928 <.0001
. t1 trt damaged - healthy -0.38259 0.01234 Inf -31.001 <.0001
. t2 ctrl damaged - healthy -0.33741 0.01099 Inf -30.691 <.0001
. t2 trt damaged - healthy -0.08085 0.00814 Inf -9.932 <.0001
damaged . ctrl t1 - t2 0.16378 0.00907 Inf 18.057 <.0001
damaged . trt t1 - t2 -0.09403 0.00906 Inf -10.382 <.0001
healthy . ctrl t1 - t2 0.20605 0.00863 Inf 23.868 <.0001
healthy . trt t1 - t2 0.20771 0.00867 Inf 23.957 <.0001
Results are averaged over some or all of the levels of: direction
Degrees-of-freedom method: inherited from asymptotic when re-gridding
P value adjustment: bonferroni method for 12 tests
你将如何排除4"；受损-健康"；对比？我可以单独计算简单的效果，并省略"；区域"；，但我喜欢"simple="；每个"&"combine=TRUE"选项在所有测试的p值校正中滚动。
我参考了文档，试图与以下代码进行对比，但失败了：
contrast(lmm_1_emm, method = list("ctrl damaged t2" - "trt damaged t2") )
Error in "ctrl damaged t2" - "trt damaged t2" :
non-numeric argument to binary operator
p.s：什么是"我敢打赌，你只想用对象[[1]-如果我错了，请使用"[[]]"或"which"。看到了吗？emm_list’以获取更多信息"；意思是
首先，对于自定义对比度，最好在emmeans()中的公式中而不是有左侧。只需获得你想要的手段，然后单独进行对比，例如
emm <- emmeans(lmm_1,
~ intervention * region * timepoint,
type = "response",
bias.adj = TRUE,
sigma = sigma(lmm_1))
其次，我看不出文档中的任何内容是如何建议contrast(lmm_1_emm, method = list("ctrl damaged t2" - "trt damaged t2") )的。method参数必须是命名列表或对比度族的名称。
我们可以通过列出对象(例如(来查看手段，以及它们的呈现顺序
emm
在本例中，将列出8种方式。然后，您可以通过将1和-1放在列表顺序的位置，将0放在其他位置来选择要比较的，例如
contrast(emm, list(
`t d 1 - c h 2` = c(0, 1, -1, 0, 0, 0, 0, 0),
`t h 1 - t d 1` = c(0, -1, 0, 1, 0, 0, 0, 0)) )
这可能相当乏味，因此在可能的情况下，最好尽可能使用by变量和命名对比族。例如，要分别比较每个区域的连续interventionxtimepoint组合，请执行
contrast(emm, "consec", by = "region")