I'm running some tests on phytools right now to try and put a new version on CRAN, but one of them required I use the function lambdaTree in the geiger package. geiger is not a dependency of phytools, so I wanted to be able to do this (ideally) without loading geiger. Here is a simplified lambdaTree (now deprecated into transform.phylo):
lambdaTree<-function(tree,lambda){
ii<-which(tree$edge[,2]>length(tree$tip.label))
H1<-nodeHeights(tree)
tree$edge.length[ii]<-lambda*tree$edge.length[ii]
H2<-nodeHeights(tree)
tree$edge.length[-ii]<-tree$edge.length[-ii]+ H1[-ii,2]-H2[-ii,2]
tree
}
ii<-which(tree$edge[,2]>length(tree$tip.label))
H1<-nodeHeights(tree)
tree$edge.length[ii]<-lambda*tree$edge.length[ii]
H2<-nodeHeights(tree)
tree$edge.length[-ii]<-tree$edge.length[-ii]+ H1[-ii,2]-H2[-ii,2]
tree
}
Here's a quick check:
> library(geiger)
> library(phytools)
> tree<-pbtree(n=26,tip.label=LETTERS[26:1])
> plotTree(transform(tree,model="lambda",lambda=0.5), mar=c(0,1,4,1))
> title(main="a) transform.phylo",adj=0)
> plotTree(lambdaTree(tree,0.5),mar=c(0,1,4,1))
> title(main="b) simple lambdaTree",adj=0)
> library(phytools)
> tree<-pbtree(n=26,tip.label=LETTERS[26:1])
> plotTree(transform(tree,model="lambda",lambda=0.5), mar=c(0,1,4,1))
> title(main="a) transform.phylo",adj=0)
> plotTree(lambdaTree(tree,0.5),mar=c(0,1,4,1))
> title(main="b) simple lambdaTree",adj=0)
No comments:
Post a Comment
Note: due to the very large amount of spam, all comments are now automatically submitted for moderation.