Thursday, June 27, 2013

Fix for bind.tip

At this year's Evolution meeting it was reported to me that bind.tip has some weird behavior if you try and attach a tip to the end of a terminal edge. One might want to do this if, say, identical haplotypes were removed for inference - and then one wanted to add them back in for plotting.

Here's a demo of the problem:

> tree<-pbtree(n=20)
> plotTree(tree,setEnv=TRUE,offset=0.7)
setEnv=TRUE is experimental. please be patient with bugs
> tiplabels()
> aa<-bind.tip(tree,"t21",where=1,position=0)
> plotTree(aa,setEnv=TRUE,offset=0.7)
setEnv=TRUE is experimental. please be patient with bugs
> tiplabels()

Holy cow! Weird, right? Instead of binding our new tip zero-length below the tip node, bind.tip is replacing the tip label of our target terminus! Note that this will not happen even if we add the the new tip just below the target tip.

In this case, bind.tip actually just wraps around the ape function bind.tree which attaches two trees together.

I have built a work-around into the latest version of bind.tip (source code now here). This is also in a new build of phytools (phytools 0.2-89). Let's verify that this works:

> require(phytools)
Loading required package: phytools
> packageVersion("phytools")
[1] ‘0.2.89’
> bb<-bind.tip(tree,"t21",where=1,position=0)
> plotTree(bb,setEnv=TRUE,offset=0.7)
setEnv=TRUE is experimental. please be patient with bugs
> tiplabels()
> bb$edge.length
[1] 0.49852021 0.90133001 0.63759396 0.36295139 0.00000000 0.00000000 ....

This fix is a bit of a hack - but at least it works.

No comments:

Post a Comment