如何用R语言agricolae包结合给定数据绘制含Tukey/LSD/Duncan分组标记的图表？

1. 准备数据（示例/替换为你的实际数据） #

这里用模拟的农业试验数据演示，你可以直接替换为自己的数据集：

set.seed(123)
# 构造示例：5种处理，每种处理5个重复的产量数据
data <- data.frame(
  Treatment = rep(c("T1", "T2", "T3", "T4", "T5"), each = 5),
  Yield = c(rnorm(5, 22, 1.5), rnorm(5, 28, 1.5), rnorm(5, 25, 1.5), rnorm(5, 20, 1.5), rnorm(5, 30, 1.5))
)

2. 单因素方差分析 #

先验证处理间差异是否显著，这是多重比较的前提：

anova_model <- aov(Yield ~ Treatment, data = data)
summary(anova_model)

3. 执行三种多重比较检验 #

使用agricolae包完成LSD、Duncan、Tukey检验：

library(agricolae)

# LSD检验
lsd_result <- LSD.test(anova_model, "Treatment", p.adj = "none")
# Duncan检验
duncan_result <- duncan.test(anova_model, "Treatment", p.adj = "none")
# Tukey HSD检验
tukey_result <- HSD.test(anova_model, "Treatment", p.adj = "none")

4. 合并检验结果与均值统计 #

提取三种检验的分组字母，和均值、标准差合并成绘图用数据集：

library(dplyr)

# 提取各检验的分组字母
lsd_groups <- lsd_result$groups %>% 
  rownames_to_column("Treatment") %>% 
  select(Treatment, LSD = groups)
duncan_groups <- duncan_result$groups %>% 
  rownames_to_column("Treatment") %>% 
  select(Treatment, Duncan = groups)
tukey_groups <- tukey_result$groups %>% 
  rownames_to_column("Treatment") %>% 
  select(Treatment, Tukey = groups)

# 计算均值、标准差并合并分组信息
plot_data <- aggregate(Yield ~ Treatment, data, mean) %>% 
  rename(Mean_Yield = Yield) %>% 
  left_join(aggregate(Yield ~ Treatment, data, sd) %>% rename(Std_Yield = Yield), by = "Treatment") %>% 
  left_join(lsd_groups, by = "Treatment") %>% 
  left_join(duncan_groups, by = "Treatment") %>% 
  left_join(tukey_groups, by = "Treatment")

5. 绘制带分组字母的图表 #

结合ggplot2绘制柱状图，添加三种检验的分组字母：

library(ggplot2)

ggplot(plot_data, aes(x = Treatment, y = Mean_Yield)) +
  geom_col(fill = "#2563eb", alpha = 0.8) +
  # 添加误差棒
  geom_errorbar(aes(ymin = Mean_Yield - Std_Yield, ymax = Mean_Yield + Std_Yield),
                width = 0.2, color = "black") +
  # 依次添加三种检验的分组字母（可根据数据调整y轴偏移量）
  geom_text(aes(label = Tukey), y = Mean_Yield + Std_Yield + 0.4, size = 4, fontface = "bold") +
  geom_text(aes(label = LSD), y = Mean_Yield + Std_Yield + 0.9, size = 4, fontface = "bold", color = "#dc2626") +
  geom_text(aes(label = Duncan), y = Mean_Yield + Std_Yield + 1.4, size = 4, fontface = "bold", color = "#059669") +
  labs(title = "不同处理的作物产量对比",
       x = "处理类型",
       y = "平均产量",
       caption = "分组字母：黑色=Tukey 红色=LSD 绿色=Duncan") +
  theme_bw() +
  theme(plot.title = element_text(hjust = 0.5, size = 14, face = "bold"),
        axis.title = element_text(size = 12),
        axis.text = element_text(size = 10))

注意事项 #

• 若你的实际数据格式不同，只需调整data数据框的结构即可。
• 若方差分析结果不显著，无需进行后续多重比较。
• 分组字母的y轴偏移量可根据数据的数值范围灵活调整，避免重叠。

内容的提问来源于stack exchange，提问作者Loop Vinyl