I was wondering if this function (countSimmap) could work on SIMMAP files that contain multiple mappings (like the attached). I got it to work on single mappings like the example on the blog but get errors on my other files.
This is not too hard. Here's how we can do it using sapply:
library(phytools)
## first simulate a tree & some mapped histories
tree<-pbtree(n=100,1)
Q<-matrix(c(-2,1,1,1,-2,1,1,1,-2),3,3)
rownames(Q)<-colnames(Q)<-letters[1:3]
trees<-sim.history(tree,Q,nsim=10)
## now let's run the analysis using sapply
foo<-function(tree){
XX<-countSimmap(tree)
setNames(c(XX$N,as.vector(t(XX$Tr))),c("N",
sapply(rownames(XX$Tr),paste,colnames(XX$Tr),sep=",")))
}
XX<-t(sapply(trees,foo))
## first simulate a tree & some mapped histories
tree<-pbtree(n=100,1)
Q<-matrix(c(-2,1,1,1,-2,1,1,1,-2),3,3)
rownames(Q)<-colnames(Q)<-letters[1:3]
trees<-sim.history(tree,Q,nsim=10)
## now let's run the analysis using sapply
foo<-function(tree){
XX<-countSimmap(tree)
setNames(c(XX$N,as.vector(t(XX$Tr))),c("N",
sapply(rownames(XX$Tr),paste,colnames(XX$Tr),sep=",")))
}
XX<-t(sapply(trees,foo))
And here are what the results look like for this toy example:
> XX
N c,c c,a c,b a,c a,a a,b b,c b,a b,b
[1,] 212 0 38 28 36 0 41 33 36 0
[2,] 179 0 25 36 28 0 23 31 36 0
[3,] 188 0 29 31 33 0 31 33 31 0
[4,] 184 0 36 33 30 0 30 24 31 0
[5,] 189 0 36 30 27 0 32 28 36 0
[6,] 196 0 30 31 36 0 36 28 35 0
[7,] 177 0 24 25 37 0 33 26 32 0
[8,] 175 0 30 19 42 0 38 19 27 0
[9,] 184 0 24 27 31 0 31 36 35 0
[10,] 214 0 41 47 36 0 30 32 28 0
N c,c c,a c,b a,c a,a a,b b,c b,a b,b
[1,] 212 0 38 28 36 0 41 33 36 0
[2,] 179 0 25 36 28 0 23 31 36 0
[3,] 188 0 29 31 33 0 31 33 31 0
[4,] 184 0 36 33 30 0 30 24 31 0
[5,] 189 0 36 30 27 0 32 28 36 0
[6,] 196 0 30 31 36 0 36 28 35 0
[7,] 177 0 24 25 37 0 33 26 32 0
[8,] 175 0 30 19 42 0 38 19 27 0
[9,] 184 0 24 27 31 0 31 36 35 0
[10,] 214 0 41 47 36 0 30 32 28 0
That's it.
I just added this automatically within countSimmap if the input object is of class "multiPhylo" (described here).
ReplyDelete