microbiomeViz:绘制lefse结果中Cladogram

EndNote相关资讯 | 2019-02-10 13:09

本文为R包microbiomeViz作者李陈浩原创,博文链接见文末。作者授权发布于本平台,刘永鑫对本文进行测试,有修改。

平日经常会分析shotgun宏基因组的数据,我们的pipeline使用MetaPhlAn,Kraken等profiler。这种数据经常会产生一个表格,如下

官网上相关的教程很详细,但是问题是,这个完全封闭的python程序,想要hack,还真的是挺难得。Krona可能是另一个选择,但是同样还是会有同样的问题。最近发布的R包Metacoder,画出的图个人真心不是很喜欢:

跟Y叔讨论了一下用ggtree实现像GraPhlAn那样图的可能性,得到了肯定的答复,于是开始自己造轮子。

MicrobiomeViz–千里之行,始于足下其实可以写一个简单的函数,但是还是想做一个拓展性更强的东西,所以就有了这个包(不断完善中):

让我们产生lefse调用graphlan绘制的物种树标记差异物种的Cladogram

输入数据为metaphlan2结果合并的矩阵。如何生成详见:

ID      BM_SRS    BM_SRS015374    BM_SRS015646    BM_SRS017687    BM_SRS019221    BM_SRS019329    BM_SRS020336    BM_SRS022145    BM_SRS022532     k__Bacteria     100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   100.0   k__Bacteria|p__Actinobacteria   1.33609 2.90435 0.45117 6.7964  14.08966        2.30709 7.30108 0.53534 3.57207 8.47622 7.07037 17.30722        30.62601 k__Bacteria|p__Actinobacteria|c__Actinobacteria 1.33609 2.90435 0.45117 6.7964  14.08966        2.30709 7.30108 0.53534 3.57207 8.47622 7.07037 17.30722包安装和加载# microbiomeViz需要 R 3.5 以上,依赖包安装 source("") biocLite("ggtree") devtools::install_github("lch14forever/microbiomeViz") library(microbiomeViz)物种相对丰对矩阵绘制物种树# 加载测试数据 df <- read.table("http://bailab.genetics.ac.cn/markdown/R/microbiomeViz/merged_abundance_table.txt", head=TRUE, stringsAsFactors = FALSE) ## 计算均值用于呈现结点大小 dat <- data.frame(V1=df[,1], V2=rowMeans(df[,-1]), stringsAsFactors = FALSE) # 用物种和丰度生成树骨架 tr <- parseMetaphlanTSV(dat, node.size.offset=2, node.size.scale=0.8) p <- tree.backbone(tr, size=0.5) p

# 读取需要颜色标注的差异物种列表,本质上是两列和颜色对应表 lefse_lists = data.frame(node=c('s__Haemophilus_parainfluenzae','p__Proteobacteria',                                'f__Veillonellaceae','o__Selenomonadales',                                'c__Negativicutes', 's__Streptococcus_parasanguinis',                                'p__Firmicutes','f__Streptococcaceae',                                'g__Streptococcus','o__Lactobacillales',                                'c__Bacilli','s__Streptococcus_mitis'),                         color=c(rep('darkgreen',6), rep('red','6')),                         stringsAsFactors = FALSE ) # 注释树 p <- clade.anno(p, lefse_lists, alpha=0.3) p

简单几行代码,美图大功告成。

-04-20-r-microbiomeviz_example/

系列教程:

专业技能:

一文读懂:

必备技能:

文献阅读

扩增子分析:

在线工具:

科研经验:

编程模板:

生物科普:

学习16S扩增子、宏基因组科研思路和分析实战,关注“宏基因组”