Tuesday, July 23, 2013

Plotting the variance among stochastic character maps on the tree

The phytools function densityMap can be used to visualize the posterior probability of a binary character on the nodes & edges on the tree obtained via stochastic character mapping. By default, this function uses a continuous color gradient from blue to red (with low posterior probability for state 1 in blue); thus purple colors indicate high uncertainty about the state of the character mapped on the tree.

I thought it might be interesting to visualize the variance (i.e., uncertainty) of our ancestral reconstruction directly. We can compute the variance in our posterior probability (which is just a proportion from the posterior sample) as p(1-p). Here is what this looks like:

First, simulate a tree & some data, and then perform stochastic mapping × 100:

> tree<-pbtree(n=30,scale=2)
> Q<-matrix(c(-1,1,1,-1),2,2)
> rownames(Q)<-colnames(Q)<-0:1
> x<-sim.history(tree,Q)$states
> mtrees<-make.simmap(tree,x,nsim=100)

Now we're ready for my somewhat hacky solution for plotting the variance among stochastic maps on the edges & nodes of the tree. Note that this requires the latest phytools build, phytools 0.3-13.

# collapse each prob with the same variance into each other
for(i in 0:499){
 if(sum(chk)==2) tt<-mergeMappedStates(tt,oo,nn)
 else if(chk[2])
  for(j in 1:length(tt$maps))
# compute the variances
# this is a hack
XX$cols<-setNames(as.vector(rbind(grey(v/0.25), grey(v/0.25))),as.vector(rbind(0:500,0:500)))
# plot!

OK, now let's compare the normal contMap plot with our "uncertainty" plot. First, here is the contMap plot:

Now here is the variance plot:

That's it.

No comments:

Post a Comment