在R中绘制带Tukey HSD置信区间的均值图及版本适配问题

针对你在R 3.5.2、ggplot2 3.1.0、dplyr 0.7.8环境下，因版本参数差异导致Tukey置信区间绘制失败的问题，以下是兼容低版本的解决方案，实现Statgraphics风格的均值对比图：

1. 数据准备与ANOVA分析 #

先完成单因素方差分析，确认因子对响应变量的显著影响：

# 加载依赖包（确保已安装对应版本）
library(ggplot2)
library(dplyr)
library(multcomp)

# 模拟含poison因子的示例数据集（可替换为你的真实数据）
set.seed(123)
data <- data.frame(
  poison = factor(rep(c("A", "B", "C", "D"), each = 10)),
  response = c(rnorm(10, 5, 1), rnorm(10, 7, 1), rnorm(10, 6, 1), rnorm(10, 8, 1))
)

# 执行单因素ANOVA
anova_model <- aov(response ~ poison, data = data)
summary(anova_model)

2. 计算Tukey HSD置信区间并整理数据 #

老版本ggplot2的stat_summary对自定义置信区间函数支持有限，因此直接提前计算各组均值与Tukey置信区间，避免参数兼容问题：

# 用multcomp获取各组Tukey置信区间
mc <- glht(anova_model, linfct = mcp(poison = "Tukey"))
ci <- confint(mc)
tukey_ci_df <- data.frame(
  poison = rownames(ci$confint),
  lower = ci$confint[, "lwr"],
  upper = ci$confint[, "upr"]
)

# 提取组均值并合并数据
group_means <- model.tables(anova_model, type = "means")$tables$poison
plot_data <- merge(
  data.frame(poison = names(group_means), mean = as.numeric(group_means)),
  tukey_ci_df,
  by = "poison"
)

# 标记最优水平（以响应值最高组为例）
plot_data <- plot_data %>%
  mutate(is_optimal = mean == max(mean))

3. 添加显著性字母标记 #

为直观呈现组间差异，用multcompView生成显著性字母：

library(multcompView)
tukey_result <- TukeyHSD(anova_model)$poison
tukey_letters <- multcompLetters4(anova_model, tukey_result)
plot_data <- merge(
  plot_data,
  data.frame(poison = names(tukey_letters$poison$Letters), 
             letters = tukey_letters$poison$Letters),
  by = "poison"
)

4. 绘制Statgraphics风格均值图 #

用基础ggplot图层实现目标样式，彻底规避版本参数冲突：

ggplot(plot_data, aes(x = poison, y = mean)) +
  # 绘制Tukey置信区间
  geom_errorbar(aes(ymin = lower, ymax = upper), width = 0.2, color = "black") +
  # 绘制均值点，最优水平高亮放大
  geom_point(aes(color = is_optimal, size = is_optimal), shape = 19) +
  # 添加显著性字母
  geom_text(aes(y = upper + 0.2, label = letters), size = 4) +
  # 样式配置，贴近Statgraphics风格
  scale_color_manual(values = c("FALSE" = "black", "TRUE" = "red")) +
  scale_size_manual(values = c("FALSE" = 3, "TRUE" = 5)) +
  labs(title = "Poison因子均值与Tukey HSD置信区间",
       x = "Poison类型",
       y = "响应值") +
  theme_bw() +
  theme(legend.position = "none",
        plot.title = element_text(hjust = 0.5),
        axis.title = element_text(size = 12),
        axis.text = element_text(size = 10))

关键版本兼容说明 #

• 放弃使用stat_summary(fun.data = ...)的方式，直接用预计算的置信区间数据绘制，彻底解决老版本参数被忽略、默认调用mean_se()的问题。
• 所有数据处理函数（dplyr::mutate、merge）均兼容0.7.8版本，model.tables在R 3.5.2中返回结构稳定。
• 若multcompView安装失败，可手动根据TukeyHSD的两两比较结果标记字母：无显著差异的组分配相同字母。

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