R语言绘制带t检验显著性标记的多组频次统计柱状图方法

实现逻辑 #

• 先统计Phenotype列的样本计数、占比，按分子分型大类拆分簇分组，补全无样本的表型空值
• 逐组执行t检验，按p值自动生成显著性标记：p<0.001标记***，p<0.01标记**，p<0.05标记*，其余统一标记ns
• 用ggplot2绘图，直接映射显著性标注的位置、文本，比基础barplot调整效率更高

完整代码 #

# 加载依赖包，未安装先执行 install.packages(c("ggplot2","dplyr","ggsignif"))
library(ggplot2)
library(dplyr)
library(ggsignif)

# 1. 数据统计整理
pheno_stat <- as.data.frame(table(ttcluster_dataset$Phenotype)) |> 
  rename(Phenotype = Var1, Count = Freq) |> 
  mutate(Percentage = Count/sum(Count)*100)

# 补全当前数据中计数为0的表型（Neural (Cluster 2)）
full_pheno <- data.frame(
  Phenotype = c("Proneural (Cluster 1)", "Proneural (Cluster 2)", 
                "Neural (Cluster 1)", "Neural (Cluster 2)", 
                "Classical (Cluster 1)", "Classical (Cluster 2)", 
                "Mesenchymal (Cluster 1)", "Mesenchymal (Cluster 2)")
)
pheno_stat <- left_join(full_pheno, pheno_stat, by = "Phenotype") |> 
  mutate(
    Count = ifelse(is.na(Count), 0, Count),
    Percentage = ifelse(is.na(Percentage), 0, Percentage),
    Pheno_Type = gsub(" \\(Cluster.*", "", Phenotype),
    Cluster = gsub(".*\\((Cluster \\d)\\)", "\\1", Phenotype)
  )

# 2. 批量t检验生成显著性结果
sig_res <- data.frame(
  Pheno_Type = unique(pheno_stat$Pheno_Type),
  start = seq(0.8, 6.8, 2),
  end = seq(1.2, 7.2, 2)
)
sig_res$p <- sapply(sig_res$Pheno_Type, \(x){
  t.test(Percentage ~ Cluster, data = subset(pheno_stat, Pheno_Type == x))$p.value
})
sig_res <- sig_res |> 
  mutate(
    label = case_when(
      p < 0.001 ~ "***",
      p < 0.01 ~ "**",
      p < 0.05 ~ "*",
      TRUE ~ "ns"
    ),
    y = sapply(Pheno_Type, \(x){
      max(subset(pheno_stat, Pheno_Type == x)$Percentage) + 2
    })
  )

# 3. 绘图
# 自定义配色，可替换为原有clustercolor向量
clustercolor <- rep(c("#E64B35","#4DBBD5","#00A087","#3C5488"), each = 2)

ggplot(pheno_stat, aes(x = Phenotype, y = Percentage, fill = Phenotype)) +
  geom_col(colour = "black", width = 0.7) +
  scale_fill_manual(values = clustercolor) +
  geom_signif(data = sig_res,
              aes(xmin = start, xmax = end, annotations = label, y_position = y),
              textsize = 5, tip_length = 0.01, manual = T) +
  labs(
    x = "Phenotypes",
    y = "Percentage of Samples expressed",
    title = "Sample wise Phenotype distribution"
  ) +
  theme_bw() +
  theme(
    axis.text.x = element_text(size = 8, angle = 45, hjust = 1),
    legend.position = "none",
    plot.title = element_text(hjust = 0.5)
  )

参数调整说明 #

• 如果比较组不是同分型下的两个cluster，直接修改sig_res中的分组规则、x轴起止坐标即可
• 显著性标注的高度通过sig_res中y列的偏移值调整，避免和柱顶重叠
• 如果不需要展示计数为0的表型，删除补全full_pheno的关联步骤即可
• t检验默认使用方差不齐的Welch校正版本，如果需要普通 Student's t检验，在t.test函数内加参数var.equal = TRUE即可

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