ggplot2中geom_errorbar()使用求助：二次模型绘图添加误差棒遇阻

嘿，我懂你现在的困扰——已经画出了分性别（Sexo）的二次拟合曲线和散点，但不知道怎么把误差棒加上去对吧？别慌，咱们一步步来搞定这个问题。

首先得明确：geom_errorbar()必须要有明确的ymin和ymax值，所以咱们得先把包含这两个值的数据准备好。根据你的需求，我分两种最常见的场景来给你说明：

场景1：给原始数据的分组统计（均值±误差）加误差棒 #

如果你想展示每个性别下，不同density水平对应的size均值±标准误/标准差，可以先对数据做分组汇总：

# 加载需要的包
library(tidyverse)

# 按Sexo和density分组，计算size的均值、标准误
summary_df <- df %>%
  group_by(Sexo, density) %>%
  summarise(
    mean_size = mean(size, na.rm = TRUE),
    se_size = sd(size, na.rm = TRUE)/sqrt(n()), # 计算标准误，换成sd(size)就是标准差
    .groups = "drop"
  )

然后把汇总数据和原数据结合到ggplot里，用geom_errorbar()绘制误差棒：

ggplot() +
  # 先画原始散点
  geom_point(data = df, aes(x = density, y = size, col = Sexo)) +
  # 画二次拟合曲线
  geom_smooth(data = df, aes(x = density, y = size, col = Sexo), 
              method = lm, formula = y ~ x + I(x^2), size = 1) +
  # 画分组统计的误差棒
  geom_errorbar(data = summary_df, 
                aes(x = density, y = mean_size, col = Sexo, 
                    ymin = mean_size - se_size, ymax = mean_size + se_size),
                width = 0.1) # width控制误差棒横向宽度，可根据数据调整

场景2：给模型拟合的置信区间加误差棒 #

如果你想展示二次模型拟合值的置信区间（和geom_smooth的置信区间对应），可以用predict()函数从模型中提取拟合值和置信区间：

首先按Sexo分组拟合模型，生成预测数据：

# 按Sexo分组生成预测数据
pred_df <- df %>%
  group_by(Sexo) %>%
  group_modify(function(data, key) {
    # 拟合二次模型
    model <- lm(size ~ poly(density, 2, raw = TRUE), data = data)
    # 生成包含density序列的预测数据框
    new_data <- tibble(density = seq(min(data$density), max(data$density), length.out = 100))
    # 预测拟合值和95%置信区间
    preds <- predict(model, newdata = new_data, interval = "confidence")
    # 合并数据
    cbind(new_data, preds, Sexo = key$Sexo)
  }) %>%
  ungroup()

然后用预测数据绘制误差棒和拟合曲线：

ggplot() +
  geom_point(data = df, aes(x = density, y = size, col = Sexo)) +
  # 用预测数据画拟合曲线
  geom_line(data = pred_df, aes(x = density, y = fit, col = Sexo), size = 1) +
  # 画模型置信区间的误差棒
  geom_errorbar(data = pred_df, 
                aes(x = density, y = fit, col = Sexo, 
                    ymin = lwr, ymax = upr),
                width = 0.1)

实用小提示 #

• 如果你的density是连续变量，场景1里直接按density分组可能会因为每个density值的个体太少，导致误差棒过于密集。这时候可以用cut_width()或cut_number()把density分成若干区间后再汇总：

  summary_df <- df %>%
    group_by(Sexo, density_bin = cut_width(density, width = 0.5)) %>% # 按宽度0.5分箱
    mutate(density_mid = mean(density, na.rm = TRUE)) %>% # 取每个箱的中点作为x轴位置
    summarise(
      mean_size = mean(size, na.rm = TRUE),
      se_size = sd(size, na.rm = TRUE)/sqrt(n()),
      density = first(density_mid), # 用中点作为x值
      .groups = "drop"
    )
  

• 误差棒的width参数可以根据你的x轴刻度灵活调整，让图表更美观。

内容的提问来源于stack exchange，提问作者Marcos Patrício