Y轴为log10刻度的柱状图stat_summary均值显示错误求助

问题描述 #

要绘制微生物计数的柱状图（误差棒为标准差）并叠加散点图，数据包含Steel/Copper两种材料、T0/T+1h/T+2h三个时间点。因部分原始数据为0，已替换为1以便在log轴上显示，但绘制后Copper组全部、Steel组部分的柱状图均值显示异常，推测与Y轴log10转换有关。

问题根源 #

当前代码中stat_summary(fun=mean)计算的是原始数据的算术均值，但微生物计数这类右偏分布的数据，在log10刻度轴上应该使用几何均值才合理：

• 算术均值会被大量小值（替换后的1）拉低，在log轴上显示位置远低于真实数据的集中趋势
• 基于原始数据计算的标准差，在log轴上同样不适用，应该基于log转换后的数据计算标准差再转换回原始尺度

修正方案 #

步骤1：提前计算每组的几何统计量 #

按Material分组，计算几何均值、log转换后的标准差，再推导误差棒的上下限：

• 几何均值：exp(mean(log(CFU)))
• 误差棒上限：exp(mean(log(CFU)) + sd(log(CFU)))
• 误差棒下限：exp(mean(log(CFU)) - sd(log(CFU)))，同时确保下限不小于1（匹配替换后的最小数值）

步骤2：重构ggplot代码 #

用提前计算好的统计量绘制柱状图和误差棒，散点图保留原始数据的抖动显示。

完整修正代码 #

library(ggplot2)
library(dplyr)
library(scales)

# 原始数据
Material<-c("Copper T0", "Copper T0","Copper T0","Copper T0","Copper T0","Copper T0","Copper T0","Copper T0","Copper T0",
             "Copper T+1h","Copper T+1h","Copper T+1h","Copper T+1h","Copper T+1h","Copper T+1h","Copper T+1h","Copper T+1h","Copper T+1h",
             "Copper T+2h","Copper T+2h","Copper T+2h","Copper T+2h","Copper T+2h","Copper T+2h","Copper T+2h","Copper T+2h","Copper T+2h",
             "Steel T0", "Steel T0","Steel T0","Steel T0","Steel T0","Steel T0","Steel T0","Steel T0","Steel T0",
             "Steel T+1h", "Steel T+1h","Steel T+1h","Steel T+1h","Steel T+1h","Steel T+1h","Steel T+1h","Steel T+1h","Steel T+1h",
             "Steel T+2h","Steel T+2h","Steel T+2h","Steel T+2h","Steel T+2h","Steel T+2h","Steel T+2h","Steel T+2h","Steel T+2h")
CFU<-c(1,1,1,1,1,1,1,1,12500,
       1,1,1,1,1,1,1,1,12500,
       1,1,1,1,1,1,1,1,12500,
       262500,137500,437500,112500,37500,62500,250000,225000,50000,
       50000,112500,37500,62500,75000,1,12500,87500,1,
       1,1,25000,12500,1,1,137500,150000,112500)
data3<-data.frame(Material,CFU)

# 计算每组的几何统计量
summary_stats <- data3 %>%
  group_by(Material) %>%
  summarise(
    geom_mean = exp(mean(log(CFU))),
    log_sd = sd(log(CFU)),
    ymin = pmax(exp(mean(log(CFU)) - log_sd), 1),  # 下限不小于1
    ymax = exp(mean(log(CFU)) + log_sd)
  ) %>%
  mutate(Material = factor(Material, levels = c("Steel T0", "Steel T+1h", "Steel T+2h", "Copper T0", "Copper T+1h", "Copper T+2h"))) %>%
  arrange(Material)

# 绘图
ggplot() +
  # 柱状图：用提前计算的几何均值
  geom_col(data = summary_stats, aes(x = Material, y = geom_mean, fill = Material), 
           color = "black", show.legend = FALSE) +
  # 误差棒：用几何标准差的上下限
  geom_errorbar(data = summary_stats, aes(x = Material, ymin = ymin, ymax = ymax), 
                width = 0.2) +
  # 散点图：原始数据抖动显示
  geom_point(data = data3, aes(x = factor(Material, levels = levels(summary_stats$Material)), y = CFU), 
             position = position_jitter(width = 0.2), color = "black", show.legend = FALSE) +
  # Y轴log10刻度设置
  scale_y_continuous(trans = 'log10', expand = c(0, 0.0001),
                     breaks = trans_breaks("log10", function(x) 10^x, n = 8),
                     labels = trans_format("log10", math_format(10^.x))) +
  # 主题设置
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        panel.background = element_blank(), axis.line = element_line(colour = "darkgrey"),
        axis.title.x = element_blank(), axis.ticks.x = element_blank()) +
  # 填充色设置
  scale_fill_manual(values = c("sienna3", "sienna3", "sienna3", 
                               "lightslategrey", "lightslategrey", "lightslategrey")) +
  # Y轴标签
  ylab(bquote(CFU.cm^-2))

关键说明 #

• 用dplyr分组计算统计量，确保每组的均值和误差棒符合log刻度下的分布特征
• 误差棒下限设置为不小于1，避免出现小于替换值的不合理数值
• Y轴标签修正为CFU（原代码误写为PFU），如果确实是PFU可改回

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