查看原文
其他

bootstrap分段标记

2017-05-05 Y叔 biobabble

想要分段标记bootstrap,比如0-70, 70-90, 90-100之类的。这个对于ggtree来说,太简单。
首先我们要读树,用treeio的read.newick,它和read.tree的不同之处是,加入了参数node.label,当node.label存的不是label,而是bootstrap等数字型的时候,你可以传入node.label=’support’,这样它会把node label解析为support value,另存为树注释数据,而不是和tip label混在一起。label变量只能存为字符型,因为和tip label混了,而字符型会让数字操作稍微复杂点。所以我们要用read.newick。

tree <- read.newick("RMI.phy_phyml_tree_rooted_labeled.txt", node.label='support')
root <- getRoot(tree@phylo)  
p <- ggtree(tree, color="black", size=1.5, linetype=1,  ladderize=TRUE) + ggtitle(label="Figure A") + geom_tiplab(size=4.5, hjust = -0.060, fontface="bold") +  xlim(0, 0.09)

首先用ggtree画出树的结构,用ggtitle加title,用geom_tiplab加tip label,出图如下:

然后我们就可以用geom_point2来打点了,我们要指定只打内部节点,但不打root node,因为root和外部节点都是没有值的。然后这个support值是数字型,我们想按区间分段,这太容易,用cut来切就可以了。这里的值是做了1000次的bootstrap。于是得到下面的图:

p <- p + geom_point2(aes(subset=!isTip & node != root, fill=cut(support, c(0, 700, 900, 1000))), shape=21, size=4) + theme_tree(legend.position=c(0.2, 0.3))

已然是近乎完美,但我们仍需要修图,想用简单的颜色来填充,低值用白色,中间值用灰色,而高值用黑色。 label想除以1000,换成百分比,legend的名字想换成bootstrap百分值,legend的label想到数学公式表达区间,这些都可以在scale_fill_manual中实现:

p + scale_fill_manual(values=c("white", "grey", "black"), guide='legend', name='Bootstrap Percentage(BP)', breaks=c('(900,1e+03]', '(700,900]', '(0,700]'), labels=expression(BP>=90,70 <= BP * " < 90", BP < 70))

华丽丽的图就产生了:


您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存