使用ggtree绘制系统发育树

绘制系统发育树的方法

一、基本概念：nwk格式文件

nwk格式是一种存储系统发育树的格式，更具体的信息可以参考https://evolution.genetics.washington.edu/phylip/newicktree.html。在这一类型的文件中，系统发育树的节点和分支由圆括号的层次关系所表示。

一颗系统发育树的表达式以分号结尾，用圆括号和逗号标记节点和分支。节点名称可以为空，也可以是一串可打印字符，但是不能包含空格、逗号、分号、圆括号和方括号。

一种拓展的表示方法是在节点名称后加分支长度，二者以冒号隔开。例如

(B:6.0,(A:5.0,C:3.0,E:4.0):5.0,D:11.0);

另外，由于分支和节点除了拓扑关系以外，不存在其他的排序问题，因此对于同一棵数，可以有多种表示方法。例如，下面五条表达式，所表达的是同一棵树。

(A,(B,C),D);    
(A,(C,B),D);    
(D,(C,B),A);    
(D,A,(C,B));    
((C,B),A,D);

二、使用ggtree绘制系统发育树

ggtree是一个基于ggplot的R包，可以在R中对系统发育树进行可视化。

具体使用说明可以参考bioconductor上的主页和treedata-book文档进行了解。

安装方法：

BiocManager::install("ggtree")

一个简单的例子：

library("ggtree")
mytree <- read.tree("ggtree_test.tree")

## file content of "ggtree_test.tree":
# ((seq0:0.01222,((((seq1:0.13069,((seq3:0.09087,seq8:0.14246):0.01517,seq13:0.09038):0.01931):0.00424,seq14:0.14715):0.01204,(seq6:0.10536,seq10:0.18353):0.00463):0.00946,((seq2:0.14280,seq11:0.09053):0.01233,(seq7:0.09201,seq15:0.14132):0.02656):0.00825):0.00566):0.00167,seq4:0.10935,((seq5:0.10304,seq12:0.16363):0.01016,seq9:0.05095):0.00454);

svg("ggtree_test.svg")
myfig <- ggtree(mytree,layout = "rectangular") + geom_tiplab() + geom_point(color='firebrick')
myfig # show plot
dev.off()

上面的代码读取了ggtree_test.tree这个文件，并使用ggtree()这个函数实现了可视化。为了后期修图方便，此处我们将其保存为svg矢量图格式。得到的结果如下图所示。

ggtree包依赖于ggplot，因此一些针对ggplot的图片参数调整方法同样适用于ggtree。例如，我们可以通过下面的代码，得到一棵更加美观的树。

myfig2 <- ggtree(mytree,layout = "circular")+
	geom_tiplab()+
	geom_highlight(node = 10,fill="blue ",alpha=0.5)+
	geom_highlight(node=9,fill=" orange", alpha=0.5)+
	geom_highlight(node=8,fill=" red ")+
	geom_cladelabel(node=3,label="Virtual creature 1",offset=0.048,barsize = 4,color="blue")+
	geom_cladelabel(node=4,label="Virtual creature 2",offset=0.048,barsize = 4,color="green",)+
	geom_cladelabel(node=5,label="Virtual creature 3",offset=0.048,barsize = 4,color="red")
myfig2