A phytools user, Bruno Vilela,
recently reported the bug that the phytools function contMap
can fail for very short internal edges. This tends to be the case when
an edge has a length that is less than 1/res
× the tree
height, in which res
is an argument of contMap
specifying the graininess of the color gradient to be plotted along edges.
[Actually, I believe I already fixed this bug except for edges descended from the root node. This is because for those edges the function will attempt to create a list of mappings along edges, but create a vector instead. Silly!]
I have previously suggested collapsing very short internal edges or increasing
res
; however Bruno did one better by providing an incredibly
simple solution, which was to create a list of mappings along edges before
iterating across the edges of the tree to compute the mappings. Perhaps
obvious in retrospect, please keep in mind that I wrote the first version of
this function in 2012 when I was quite an R programming novice. (I'm not sure
what I'd call myself now. Slightly more than an a novice? I guess in an era
when
Donald
Trump is running for present I'm a world-class R programming expert!)
Here's a quick demo of the problem, along with a demonstration that Bruno's fix does indeed restore function:
library(phytools)
tree
##
## Phylogenetic tree with 26 tips and 25 internal nodes.
##
## Tip labels:
## t5, t26, t18, t17, t19, t22, ...
##
## Rooted; includes branch lengths.
x
## t5 t26 t18 t17 t19 t22
## 0.20063042 0.90298749 2.28555569 2.11331867 1.38012423 0.01110355
## t15 t12 t6 t3 t2 t24
## 0.31029532 -0.31270146 0.36006396 0.86328160 -0.16281977 -1.88665386
## t4 t9 t25 t11 t8 t21
## -2.40507772 -2.69645998 -1.57485299 1.30619095 1.96062919 -0.91760839
## t1 t14 t7 t16 t20 t23
## -1.45220994 -2.58959187 -1.44919866 0.93419977 0.19530914 0.06853757
## t13 t10
## 1.79467390 0.29857041
plotTree(tree) ## short node from the root
obj<-contMap(tree,x,plot=FALSE)
## Error in tree$maps[[i]] <- XX[, 2] - XX[, 1]: more elements supplied than there are to replace
## increase res
obj<-contMap(tree,x,plot=FALSE,res=1000)
plot(obj)
Now let's load the new version from source off GitHub:
source("https://raw.githubusercontent.com/liamrevell/phytools/master/R/contMap.R")
source("https://raw.githubusercontent.com/liamrevell/phytools/master/R/densityMap.R")
source("https://raw.githubusercontent.com/liamrevell/phytools/master/R/read.simmap.R")
obj<-contMap(tree,x,plot=FALSE)
plot(obj)
This fix also allows us to set the value of res
to any
arbitrarily low value. Although the resultant plot is grainier, this might
nonetheless be useful for users working with very large trees where
computation is an issue.
contMap(tree,x,res=100)
contMap(tree,x,res=10)
Thanks Bruno!
No comments:
Post a Comment
Note: due to the very large amount of spam, all comments are now automatically submitted for moderation.