在火山图上标记基因(方法二)

测试开头

今天是生信星球陪你的第414天

大神一句话，菜鸟跑半年。我不是大神，但我可以缩短你走弯路的半年~

就像歌儿唱的那样，如果你不知道该往哪儿走，就留在这学点生信好不好~

这里有豆豆和花花的学习历程，从新手到进阶，生信路上有你有我！

关于火山图加标签的需求，有几种方法来实现，昨天已经写过方法一 (点击查看)，似乎只能处理较少的数据，应用有所限制，今天介绍方法二。

示例数据

方法一的示例数据是data.Rdata，方法二三的示例数据是test.Rdata。我将数据打包放在了“生信星球”公众号后台，回复“火山图”即可获得。你解压后双击文件夹里的volcano.Rproj，复制粘贴运行本文代码即可。

代码来源

以下代码出自R数据科学笔记第21章，原书第312页：

best_in_class <- mpg %>%
  group_by(class) %>%
  filter(row_number(desc(hwy)) == 1)

ggplot(mpg, aes(displ, hwy)) +
  geom_point(aes(color = class)) +
  geom_point(size = 3, shape = 1, data = best_in_class) +
  ggrepel::geom_label_repel(
    aes(label = model),
    data = best_in_class
  )

这个方法与数据量大小无关。

端详代码找思路

1.从原来数据中挑选了一部分，生成新数据
2.用新数据作图，向原数据做的点图上叠加两个图层，一个空心点图，一个geom_label_repel。

step1：先把火山图画出

load("test.Rdata")
p <- ggplot(data = test, 
            aes(x = logFC, 
                y = `-log10(P.value)`)) +
  geom_point(alpha=0.4, size=3.5, 
             aes(color=change)) +
  scale_color_manual(values=c("blue", "grey","red"))+
  geom_vline(xintercept=c(-1,1),lty=4,col="black",lwd=0.8) +
  geom_hline(yintercept = -log10(0.01),lty=4,col="black",lwd=0.8) +
  theme_bw()
p

step2：生成用于添加图层的新数据

⭐重点在这里

新数据框的内容是你想要标记的基因，这里根据logFC和Pvalue的大小来筛选，可以自定义阈值来调整要显示的基因的数量：

for_label <- test %>% 
  filter(abs(logFC) >4& `-log10(P.value)`> -log10(0.000001))

step3:新图层叠加到原图上去

p +
  geom_point(size = 3, shape = 1, data = for_label) +
  ggrepel::geom_label_repel(
    aes(label = symbol),
    data = for_label,
    color="black"
  )