首页 > 解决方案 > 如何创建带有基因 ID 的火山图?

问题描述

我从来没有创建过火山图,我不太确定我是否正确地绘制了我需要的东西。INJ 是我的数据文件,文件的列是 GeneID、Treatment 1、Treatment 2、log2FoldChange、pvalue、padj、Gene...

我想创建一个火山图,将治疗 1 和治疗 2 与标记到每个点的基因进行比较。有人可以就如何做到这一点或如何开始编写代码提供一些建议。谢谢你。

在此处输入图像描述

标签: rggplot2

解决方案


这是一个非常基本的设置,数据的形状与我期望的大致相同,其中填充了一些虚拟数据。

这些是您至少需要的数据类型,您需要记录 2 倍的变化、一个(FDR 校正的)p 值和标签的 ID。

library(ggplot2)
library(scales)
library(ggrepel)

df <- data.frame(
  log2FoldChange = log2(rcauchy(1000, 50) / rcauchy(1000, 50)),
  padj = runif(1000),
  gene_id = paste0("ENSG00000", 1:1000)
)
#> Warning in data.frame(log2FoldChange = log2(rcauchy(1000, 50)/rcauchy(1000, :
#> NaNs produced

接下来,您通常会决定如何展示您的基因。我倾向于指出 htey 是否重要,如果它们重要,它们会朝哪个方向发展。

df$category <- with(df, ifelse(padj > 0.05 | is.na(padj), "n.s.",
                               ifelse(log2FoldChange > 0, "Up", "Down")))

此外,我建议您不要绘制所有基因的标签,而可能只是绘制最强的效应大小。我在这里将截止值设置为绝对值 1,但您应该使这个值与您的数据相得益彰。

df$labels <- with(df, ifelse(abs(log2FoldChange) > 1 & padj < 0.05, gene_id, ""))

接下来,您只需绘制 log2FoldChange 与 -log10 FDR 校正的 p 值。对于标签,我建议你使用 ggrepel 包。

ggplot(df, aes(log2FoldChange, log10(padj))) +
  geom_point(aes(colour = category)) +
  geom_text_repel(aes(label = labels)) +
  scale_y_continuous(trans = "reverse",
                     labels = math_format()) +
  scale_x_continuous(name = expression("Log"[2]*" Fold Change"),
                     limits = function(x){c(-1, 1) * max(abs(x))})
#> Warning: Removed 13 rows containing missing values (geom_point).
#> Warning: Removed 13 rows containing missing values (geom_text_repel).

reprex 包(v0.3.0)于 2020 年 7 月 24 日创建


推荐阅读