Problem with bind.tree()

Today I ran into a slight problem with the {ape} function bind.tree(). I posted this issue to the R-sig-phylo email listserve. Emmanuel Paradis, the author of {ape}, is normally excellent at responding to queries like this - but in the meantime I thought I would also describe the issue in this forum.

I said in my listserve post: "For some reason, when I bind a tree with a root.edge to a tip, the length of the edge leading to that tip is attached not only to the root.edge of the bound subtree, but also to another tip of the tree." Actually, this is imprecise. In fact, the edge length from the edge leading to the target tip is for some reason moved to another tip in the tree.

The easiest way to illustrate this problem is with an example.

First, take the following tree:

text="(1:3.785,(((((8:0.13,9:0.13):0.205,6:0.335):0.154,(10:0.278,7:0.278):0.211):0.1,(4:0.386,5:0.386):0.203):2.655,(2:1.08,3:1.08):2.164):0.54);"
tree<-read.tree(text=text)
plot(tree)
edgelabels(round(tree$edge.length,2),adj=c(0.5,-0.2),frame="none")

Which produces the plot below in which I have labeled edges with their (rounded) lengths:

Now to reproduce our issue, let's cut out the clade with tips "2" and "3", and then regraft it in the same place:

node<-19 # this is the node we will extract
position<-1.0 # this is the length of the root.edge
# extract clade and attach root edge
tr1<-extract.clade(tree,node=node); tr1$root.edge<-1.0
x11(); plot(tr1,root.edge=T) # plot the extracted clade
edgelabels(round(tr1$edge.length,2),adj=c(0.5,-0.2),frame="none")

The length of the root edge is not shown, but remember we have set it to be 1.0.

# now remove tips in tr1 from tree
tr2<-drop.tip(tree,tr1$tip.label,trim.internal=FALSE)
# subtract root.edge of tr1 from tip ending in "NA"
temp<-match(which(tr2$tip.label=="NA"),tr2$edge[,2])
tr2$edge.length[temp]<-tr2$edge.length[temp]-position
x11(); plot(tr2,root.edge=T) # plot the pruned tree
edgelabels(round(tr2$edge.length,2),adj=c(0.5,-0.2),frame="none")

Here is the tree with clade (2,3) extracted.
# now bind tr1 to tr2, where it was removed
tr.bound<-bind.tree(tr2,tr1,where=which(tr2$tip.label=="NA"))
# now plot, to visualize the result
x11(); plot(tr.bound)
edgelabels(round(tr.bound$edge.length,2),adj=c(0.5,-0.2),frame="none")

So, for some reason the edge leading to "NA" before we use bind.tree(), is actually moving to taxon "5" when we reattach the subtree with species "2" and "3". Weird.

Please let me know if you can repeat this error; or, even better, if you figure out that I have made some kind of mistake in my example, above. Thanks!

In a future post (hopefully soon, if I can get this figured out) I will describe what I have been using extract.clade() and bind.tree() for today as it may be of interest to some readers.