遍历数据框绘制直方图时ggplot报Frequency未找到的问题排查

问题背景 #

我在做RNA-seq分析，需要展示不同基因的CPM（每百万转录本数）频率分布，重点关注零CPM的占比。现有数据集dat包含Gene、CPM和Frequency列，尝试循环遍历目标基因列表goi绘制直方图时，使用ggplot的geom_histogram并设置y = after_stat(Frequency / sum(Frequency))，报错提示“未找到对象Frequency”，但数据框中明确存在该列。

错误原因 #

after_stat() 仅用于访问ggplot统计变换后生成的变量（比如geom_histogram默认生成的count、density等），而你数据中的Frequency是原始数据列，不属于统计变换的输出结果，因此无法在after_stat()中引用。

另外，你的数据已经是汇总格式（每行对应一个CPM值及其出现频率），geom_histogram是针对原始观测值进行分箱计数的工具，不适合直接用于已汇总的数据。

修正方案 #

方案一：使用geom_col（推荐，适配已汇总数据） #

直接利用原始数据的Frequency列计算占比，用geom_col绘制每个CPM值的频率百分比：

# 创建空列表存储绘图对象
histogram_plots <- list()

# 遍历目标基因
for (element in goi) {
  # 筛选单个基因的数据
  element_data <- dat[dat$Gene == element, ]
  
  # 绘制柱状图
  plot <- ggplot(element_data, aes(x = CPM, y = Frequency / sum(Frequency))) +
    geom_col(fill = "cornflowerblue", color = "white") +
    labs(title = paste(element, "的CPM频率分布"),
         y = "百分比",
         x = "CPM") +
    scale_y_continuous(labels = scales::percent_format()) +
    theme_minimal()
  
  # 保存到列表
  histogram_plots[[element]] <- plot
  
  # 保存为图片
  ggsave(paste0("histogram_", gsub("-", "_", element), ".png"), plot, width = 8, height = 6)
}

方案二：使用geom_histogram（需先展开数据） #

如果需要保留直方图的分箱逻辑（合并相近CPM值为区间），先将汇总数据展开为原始观测格式（按Frequency重复CPM行），再用直方图绘制：

# 创建空列表存储绘图对象
histogram_plots <- list()

# 遍历目标基因
for (element in goi) {
  # 筛选单个基因的数据
  element_data <- dat[dat$Gene == element, ]
  # 展开数据：按Frequency重复每行
  element_data_expanded <- element_data[rep(seq(nrow(element_data)), element_data$Frequency), ]
  
  # 绘制直方图
  plot <- ggplot(element_data_expanded, aes(x = CPM)) +
    geom_histogram(aes(y = after_stat(count / sum(count))), 
                   fill = "cornflowerblue", color = "white", binwidth = 5) +
    labs(title = paste(element, "的CPM频率分布"),
         y = "百分比",
         x = "CPM") +
    scale_y_continuous(labels = scales::percent_format()) +
    theme_minimal()
  
  # 保存到列表
  histogram_plots[[element]] <- plot
  
  # 保存为图片
  ggsave(paste0("histogram_", gsub("-", "_", element), ".png"), plot, width = 8, height = 6)
}

验证说明 #

两种方案都能准确展示零CPM的频率占比：比如对于hsa-miR-15a-5p，零CPM的Frequency为30，总Frequency为该基因所有行Frequency的和，计算出的占比会直观显示在图表中。

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