treeio
treeio copied to clipboard
Feature Request: Import MCMCTree Result Tree (PAML)
Thank you for the excellent packages!
I was wondering if it would be possible to implement functionality that would enable the importing of MCMCTree 'nexus' result trees, and their subsequent use in ggtree with geom_range and other functions to visualize ranges on divergence time estimates. I have attached one such tree that comes from the examples provided with MCMCTree:
#NEXUS
BEGIN TREES;
UTREE 1 = ((((human: 0.068414, (chimpanzee: 0.027796, bonobo: 0.027796) [&95%={0.023087, 0.0332708}]: 0.040618) [&95%={0.0600119, 0.0786158}]: 0.024982, gorilla: 0.093396) [&95%={0.0810593, 0.105526}]: 0.052890, (orangutan: 0.054156, sumatran: 0.054156) [&95%={0.044843, 0.0635921}]: 0.092129) [&95%={0.12543, 0.160626}]: 0.010679, gibbon: 0.156965) [&95%={0.132252, 0.180134}];
END;
Hello, you can use read.beast( ) function, Just make two changes.
- remove the letter "U" from the word "UTREE" in your MCMCtree file.
- replace "95%" with "reltime_95%_CI"
and when you use the following code, your 95%_CI bar center would be the mean value of the range.
> mcmctree <- read.beast("FigTree.tre") > ggtree(mcmctree)+ geom_range(range='reltime_0.95_CI', color='red', size=3, alpha=0.5)
You can manually add a parameter reltime (the estimated value, i.e. the node time) in the square brackets, this would help you anchor the right 0.95_CI. Then use the following code.
> ggtree(mcmctree)+ geom_range(range='reltime_0.95_CI', color='red', size=3, alpha=0.5, center='reltime')
You can compare the Beast format (instance files that come with the treeio packages) with the MCMCtree format.
C:\Users\...\Documents\R\win-library\4.1\treeio\extdata\BEAST\beast_mcc.tree
and you will understand what I'm saying.
Good luck.
The following is my code.
library(treeio)
mcmctree <- read.beast("FigTree.tre")
library(ggplot2)
library(ggtree)
p <- ggtree(mcmctree, size=1, mrsd='00-01-01') + # 翻转时间轴
geom_tiplab(as_ylab=TRUE) # 将物种名定义为y轴的标签
insert_minor <- function(major_labs, n_minor) { # 生成label序列:"-200" "" "" "" "" "-175" ""
labs <- c(sapply(major_labs, function(x) c(x, rep("", n_minor))))
labs[1:(length(labs)-n_minor)]}
p + theme_tree2(plot.margin=margin(16, 16, 30, 16)) + # 设置绘图区域(上,右,下,左)
geom_range(range='reltime_0.95_CI', color='red', # 添加置信区间色条
size=3, alpha=0.5, center='reltime') + # reltime 定义了锚点
geom_nodelab(aes(label=reltime*100), # 添加每个节点的估计时间,并扩大100倍
hjust=1.2, vjust=-0.5, # 调整位置
color="red", # 设置节点标签的颜色
size=5) + # 设置节点标签的大小
scale_x_continuous(limits=c(-2,0), # 设置时间轴范围,单位:百万年
breaks=seq(-2,0,0.10), # 每0.10 为间隔(实际上是10 个百万年)
labels=insert_minor(seq(-200,0,20),1)) + # 标签直接扩大100倍
theme(axis.line=element_line(size=1, linetype="solid"), # 时间轴加粗
axis.text.x=element_text(size=15, vjust=-0.8), # 时间轴的字体大小
axis.text.y=element_text(size=20, color='black', face='italic'), # 设置物种名的字体大小及斜体
axis.ticks.x=element_line(size=1.2, lineend=10)) +
geom_segment(x=-0.158273, y=0, xend=-0.158273, yend=6, # 添加一条竖直的虚线
lty=2, size=1.5, color='lightblue')
@xiangpin we need to implement a read.mcmctree
function.
ok
@xiangpin we need to implement a
read.mcmctree
function.
Thanks for the new read.mcmctree feature. McMcTree does not provide the reltime parameter in square brackets, thus the confidence interval cannot be anchored to the node using the estimated node time directly (with the center=reltime in geom_range( ) ). Node time is obtained by adding the branch times together. Manual editing the tree file is time consuming and error prone for a big tree. Could you please help to deal with this?