Tuesday, August 6, 2013

New version of mergeMappedStates; some tutorials in progress

I just posted a new version of the function mergeMappedStates that can "merge" a set of states with only one entry into a new state (basically, change the name of state mapped on the tree). This might be useful for (for example) converting a tree with mapped states a, b, and c into one with two states: a & not-a (coded as 0 & 1, perhaps) for use with densityMap.

This function version is available in a new phytools build (phytools 0.3-21).

Let's try and use the function to do exactly that:

> require(phytools)
Loading required package: phytools
> packageVersion("phytools")
[1] ‘0.3.21’
> ## simulate tree & character data
> tree<-pbtree(n=100,scale=1)
> Q<-matrix(c(-2,1,1,1,-2,1,1,1,-2),3,3)
> colnames(Q)<-rownames(Q)<-letters[1:3]
> x<-sim.history(tree,Q)$states
> ## simulate stochastic histories
> mtrees<-make.simmap(tree,x,nsim=80,model="ER")
make.simmap is sampling character histories conditioned on the transition matrix
Q =
          a        b        c
a -2.433741  1.216870  1.216870
b  1.216870 -2.433741  1.216870
c  1.216870  1.216870 -2.433741
(estimated using likelihood);
and (mean) root node prior probabilities
pi =
        a        b        c
0.3333333 0.3333333 0.3333333
> ## now let's merge states to make a "a"/"not-a" set of trees
> statea<-mergeMappedStates(mtrees,c("b","c"),"0")
> statea<-mergeMappedStates(statea,"a","1")
> densityMap(statea,type="fan")
sorry - this might take a while; please be patient

I'm also posting a series of tutorials that will form the basis of the computer exercises I will lead this week at the NESCent Academy Evolutionary Quantitative Genetics Workshop in Durham NC this week. I have put draft version of three of these up so far. Check them out here. These tutorials were all created using knitr and have not yet been proofed thoroughly.

No comments:

Post a Comment