分享

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

 生物_医药_科研 2019-07-20

 今天是生信星球陪你的第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'
  )

加号连接两句代码就实现了图层的叠加。
如果对ggplot2不了解,请看R数据科学第一章和第21章,是专门讲图的,我写过对应的笔记:来自小洁的《R数据科学》目录。但21章是整本书的错误重灾区,我的笔记有改正后的代码。

    本站是提供个人知识管理的网络存储空间,所有内容均由用户发布,不代表本站观点。请注意甄别内容中的联系方式、诱导购买等信息,谨防诈骗。如发现有害或侵权内容,请点击一键举报。
    转藏 分享 献花(0

    0条评论

    发表

    请遵守用户 评论公约

    类似文章 更多