Tuesday, August 9, 2016

Co-phylogenetic plotting for trees with polytomies

One issue with the phytools function for co-phylogenetic plotting, cophylo, is how it performs when a tree has one or multiple multifurcations. Because the function uses tipRotate internally, which itself depends on rotate and a greedy algorithm to attempt to maximize tip matching, polytomies are not handled particularly well.

Just to demonstrate, the following is a co-phylogenetic plot between two trees that are identical, except in tip rotation, but that have polytomies:

library(phytools)
t1<-read.tree(text="(((A1,A2),(B1,B2,B3),C,D),E,F);")
t1$edge.length<-runif(n=nrow(t1$edge))
## perform some random rotations at the nodes with polytomies
t2<-sample(rotate.multi(t1,13),1)[[1]]
t2<-sample(rotate.multi(t2,11),1)[[1]]
t2<-sample(rotate.multi(t2,10),1)[[1]]
par(mfrow=c(1,2))
plotTree(t1)
plotTree(t2)

plot of chunk unnamed-chunk-1

And to confirm this observation:

all.equal.phylo(t1,t2)
## [1] TRUE

Unfortunately, cophylo doesn't do as well as we might hope here:

obj<-cophylo(t1,t2,print=TRUE,methods=c("pre","post"))
## Rotating nodes to optimize matching...
## objective: 118
## objective: 88
## objective: 86
## objective: 86
## objective: 86
## objective: 86
## objective: 86
## objective: 86
## objective: 86
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## objective: 62
## Done.
plot(obj,link.type="curved")

plot of chunk unnamed-chunk-3

I just pushed an update to cophylo which adds the option rotate.multi - a logical argument which tells the function whether or not to use the new phytools function rotate.multi to attempt to optimize the node rotation around polytomies.

Note that the optimization routine is still greedy, so it is probably not guaranteed to find the best solution. However, in my experience to date, it tends to - particularly using methods=c("pre","post"), indicating to conduct a pre- and post-order tree traversal.

For instance:

obj<-cophylo(t1,t2,print=TRUE,methods=c("pre","post"),
    rotate.multi=TRUE)
## Rotating nodes to optimize matching...
## objective: 188
## objective: 158
## objective: 158
## objective: 158
## objective: 156
## objective: 156
## objective: 126
## objective: 2
## objective: 2
## objective: 2
## objective: 2
## objective: 2
## objective: 2
## objective: 2
## objective: 2
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## objective: 0
## Done.
plot(obj,link.type="curved")

plot of chunk unnamed-chunk-4

That's it.

This update can be obtained by installing from GitHub; and I would be more than happy to hear how it works!

No comments:

Post a Comment