计划中的对比
计划中的对比
提问于 2020-01-19 18:15:28
我想要计算一个特定子集的计划对比使用to，但有困难的编码这些。
在我的样本数据集中，我有两个条件，"drugA“和"drugB”。有6只动物，每只动物的体重在每种药物的影响下被测量了3次。
id <- rep(c("A","B","C","D","E","F"),6)
drug <- c(rep(c("drugA"), 18), rep(c("drugB"), 18))
time <- rep(rep(1:3, each = 6),2)
value <- c(rnorm(6, 1, 0.4), rnorm(6, 3, 0.5), rnorm(6, 6, 0.8), rnorm(6, 1.1, 0.4), rnorm(6, 0.8, 0.2), rnorm(6, 1, 0.6))
df <- data.frame(id,drug, time, value)
df$id <- as.factor(df$id)
df$drug <- as.factor(df$drug)
df$time <- as.factor(df$time)
stats <- lmer(value ~ drug*time + drug + time + (1|id), data = df)
summary(stats)
emm <- emmeans(stats, list(pairwise ~ drug + time), adjust = "tukey")
emm但是，我只想来计算以下对比：
DrugA，time1与DrugB，time1
DrugA，time2与DrugB，time2
DrugA，time3与DrugB，time3
DrugA，time1和time2
DrugA，time2和time3
DrugB，time1和time2
DrugB，time2和time3
我要如何编码这些对比呢？非常感谢你的建议。
回答 2
Stack Overflow用户
回答已采纳
发布于 2020-01-22 08:23:48
这里的线索是看看结果网格的to。前两栏是形成对比的基础，每一行代表一个因素的组合。
emm <- emmeans(stats, list(~ drug + time)) # not used afterwards, but to check result Grid
con <- list(
DrugA1_DrugB1 = c(1,-1,0,0,0,0),
DrugA2_DrugB2 = c(0,0,1,-1,0,0),
DrugA3_DrugB3 = c(0,0,0,0,-1,1),
DrugA1_DrugA2 = c(1,0,-1,0,0,0),
DrugA2_DrugA3 = c(0,0,1,0,-1,0),
DrugB1_DrugB2 = c(0,1,0,-1,0,0),
DrugB2_DrugB3 = c(0,0,0,1,0,-1)
)下面给出的正是这些对比。
emm <- emmeans(stats, list(~ drug + time), contr = con, adjust = "mvt")Stack Overflow用户
发布于 2020-06-24 10:42:43
您可以使用contrast()：
library(lme4)
library(emmeans)
id <- rep(c("A","B","C","D","E","F"),6)
drug <- c(rep(c("drugA"), 18), rep(c("drugB"), 18))
time <- rep(rep(1:3, each = 6),2)
value <- c(rnorm(6, 1, 0.4), rnorm(6, 3, 0.5), rnorm(6, 6, 0.8), rnorm(6, 1.1, 0.4), rnorm(6, 0.8, 0.2), rnorm(6, 1, 0.6))
df <- data.frame(id,drug, time, value)
df$id <- as.factor(df$id)
df$drug <- as.factor(df$drug)
df$time <- as.factor(df$time)
stats <- lmer(value ~ drug*time + drug + time + (1|id), data = df)
emm <- emmeans(stats, ~ drug + time)
emm
#> drug time emmean SE df lower.CL upper.CL
#> drugA 1 1.16 0.187 28.8 0.778 1.55
#> drugB 1 1.05 0.187 28.8 0.666 1.43
#> drugA 2 3.30 0.187 28.8 2.917 3.68
#> drugB 2 0.84 0.187 28.8 0.457 1.22
#> drugA 3 5.99 0.187 28.8 5.602 6.37
#> drugB 3 1.30 0.187 28.8 0.920 1.69
#>
#> Degrees-of-freedom method: kenward-roger
#> Confidence level used: 0.95
contrast(emm, method = "pairwise", by = "time")
#> time = 1:
#> contrast estimate SE df t.ratio p.value
#> drugA - drugB 0.113 0.253 25 0.446 0.6593
#>
#> time = 2:
#> contrast estimate SE df t.ratio p.value
#> drugA - drugB 2.461 0.253 25 9.745 <.0001
#>
#> time = 3:
#> contrast estimate SE df t.ratio p.value
#> drugA - drugB 4.683 0.253 25 18.545 <.0001
#>
#> Degrees-of-freedom method: kenward-roger
contrast(emm, method = "consec", by = "drug")
#> drug = drugA:
#> contrast estimate SE df t.ratio p.value
#> 2 - 1 2.139 0.253 25 8.471 <.0001
#> 3 - 2 2.685 0.253 25 10.634 <.0001
#>
#> drug = drugB:
#> contrast estimate SE df t.ratio p.value
#> 2 - 1 -0.209 0.253 25 -0.828 0.6244
#> 3 - 2 0.463 0.253 25 1.834 0.1383
#>
#> Degrees-of-freedom method: kenward-roger
#> P value adjustment: mvt method for 2 tests页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/59813002
复制相关文章
相似问题
腾讯云开发者
Copyright © 2013 - 2026 Tencent Cloud. All Rights Reserved. 腾讯云 版权所有
深圳市腾讯计算机系统有限公司 ICP备案/许可证号：粤B2-20090059 
粤公网安备44030502008569号
腾讯云计算（北京）有限责任公司 京ICP证150476号 | 京ICP备11018762号