我仍在学习如何使用 r,并且在制作基于 Tukey 事后测试的效果图时遇到了一些麻烦。我有一个包含 3 个菌落的数据集(和A
,见下文)。每个菌落分为 2 个处理:和。我必须运行 GLM 来测试每个菌落的治疗组之间的预期寿命 last.scan 是否存在差异。我将“治疗”作为固定因素,将“菌落”作为随机因素。B
C
control
DWV
我使用了以下代码:
fit_life = lmer(last.scan~treatment + (1|colony), data = data)
Anova(fit_life, type = 3)
我还需要 95% 的置信限,但我不确定我是否使用了正确的命令:
confint(lsmeans(fit_life, list(~treatment)), type = "response")
以下命令为我提供了治疗的 Tukey 事后测试的结果,但我想要一个 Tukey 测试,它可以让我了解每个群体的治疗之间的差异。我怎样才能做到这一点?
mcp = glht(fit_life, linfct = mcp(treatment = "Tukey"))
summary(mcp)
# DWV treatment significantly changes life expectancy (z = -9.734, p = < 2e-16)
有了事后测试的结果,我想制作一个像这里显示的效果图,但我不知道使用什么命令。
以下是数据集的示例:
> head(data, 30)
RFID colony treatment last.scan
1 A0 01 03 C0 00 C0 20 01 A Control 24
2 A0 01 03 C0 00 C0 20 0C A Control 21
3 A0 01 03 C0 00 C0 20 1D A Control 19
4 A0 01 03 C0 00 C0 20 1E A Control 18
5 A0 01 03 C0 00 C0 20 1F A Control 31
6 A0 01 03 C0 00 C0 20 21 A Control 19
7 A0 01 03 C0 00 C0 20 2F A Control 18
8 A0 01 03 C0 00 C0 20 37 A Control 16
9 A0 01 03 C0 00 C0 20 5E A Control 23
10 A0 01 03 C0 00 C0 20 79 A Control 19
11 A0 01 03 C0 00 C0 20 7F A Control 17
12 A0 01 03 C0 00 C0 20 8C A Control 24
13 A0 01 03 C0 00 C0 20 92 A Control 26
14 A0 01 03 C0 00 C0 20 95 A Control 19
15 A0 01 03 C0 00 C0 20 98 A Control 21
16 A0 01 03 C0 00 C0 20 B8 A Control 21
17 A0 01 03 C0 00 C0 20 B9 A Control 20
18 A0 01 03 C0 00 C0 20 D5 A Control 17
19 A0 01 03 C0 00 C0 20 D9 A Control 27
20 A0 01 03 C0 00 C0 20 E4 A Control 26
21 A0 01 03 C0 00 C0 20 FE A Control 31
22 A0 01 03 C0 00 C0 3A 02 A DWV 11
23 A0 01 03 C0 00 C0 3A 0C A DWV 22
24 A0 01 03 C0 00 C0 3A 0D A DWV 12
25 A0 01 03 C0 00 C0 3A 11 A DWV 19
26 A0 01 03 C0 00 C0 3A 14 A DWV 21
27 A0 01 03 C0 00 C0 3A 1D A DWV 24
28 A0 01 03 C0 00 C0 3A 24 A DWV 9
29 A0 01 03 C0 00 C0 3A 2A A DWV 16
30 A0 01 03 C0 00 C0 3A 2C A DWV 23