Friday, July 25, 2014

Update to locate.yeti permitting the missing tip to attach below the root node

I just updated the phytools function locate.yeti (a method to attach a missing tip to an ultrametric base tree using continuous character data) so that the missing leaf is now permitted to attach below the root node.

When I first started working on this, I initially mistakenly thought that I would have to add one additional parameter to when attaching a new leaf to the root node - that being the length of the edge leading to the tip. Actually, this is not the case. Whereas for a leaf attached to an edge we have to optimize the position along the edge from whence the lineage splits; in the case of a leaf attached to the root we only have to optimize its total length. Since our tree is invariably ultrametric - the lengths of the root stem and the terminal edge can be found by simply midpoint rooting our phylogeny.

The new version of locate.yeti is here; however since locate.yeti uses phytools internal functions it is probably best & easiest to simply install the latest non-CRAN phytools version.

Here's a quick demo using a phylogeny in which the true position of the missing lineage is sister to our ultrametric base tree:

> packageVersion("phytools")
[1] ‘0.4.25’
> tree<-pbtree(n=26,tip.label=LETTERS,scale=1)
> tt<-bind.tip(tree,where=length(tree$tip.label)+1,
edge.length=1.2*max(nodeHeights(tree)),
tip="Missing-lineage")
> tt<-midpoint.root(tt)
> ## this is the full true tree
> plotTree(tt)
> ## simulate on the full tree
> X<-fastBM(tt,nsim=10)
> ## place our missing taxon
> mltree<-locate.yeti(tree,X,method="exhaustive")
Optimizing the phylogenetic position of Missing-lineage. Please wait....
Done.
> mltree$logL
[1] -221.0057
> plotTree(mltree)

Although it's pretty obvious that we have the correct tree in this case, we can also quantify it:

> require(phangorn)
Loading required package: phangorn
> RF.dist(tt,mltree)
[1] 0
> require(Rphylip)
Loading required package: Rphylip
> Rtreedist(tt,trees2=mltree,quiet=TRUE)
          2,1
1,1 0.0862447

That's pretty much it for now.

No comments:

Post a Comment