Saturday, May 18, 2013

Function to merge mapped states

I just posted a new function (mergeMappedStates) to merge mapped states on a phylogeny with a discrete character map. This is pretty easy. We just do two passes through all the edges of the tree. In the first pass we rename any mapped state in the old states with the new merged state:

rr<-function(map,oo,nn){
  for(i in 1:length(map)) if(names(map)[i]%in%oo)
    names(map)[i]<-nn
  map
}
maps<-lapply(maps,rr,oo=old.states,nn=new.state)
In the second pass, we join any adjacent map elements in the same (new) merged state:
mm<-function(map){
  if(length(map)>1){
    new.map<-vector()
    j<-1
    new.map[j]<-map[1]
    names(new.map)[j]<-names(map)[1]
    for(i in 2:length(map)){
      if(names(map)[i]==names(map)[i-1]){
        new.map[j]<-map[i]+new.map[j]
        names(new.map)[j]<-names(map)[i]
      } else {
        j<-j+1
        new.map[j]<-map[i]
        names(new.map)[j]<-names(map)[i]
      }
    }
    map<-new.map
  }
  map
}
maps<-lapply(maps,mm)

Here's a quick demo of how it works:

> tree<-pbtree(n=100,scale=1)
> Q
   a  b  c
a -2  1  1
b  1 -2  1
c  1  1 -2
> tree<-sim.history(tree,Q)
> plotSimmap(tree,lwd=3,ftype="off",pts=F)
no colors provided. using the following legend:
      a        b        c
 "black"    "red" "green3"
> merged<-mergeMappedStates(tree,c("a","c"),"ac")
> plotSimmap(merged,lwd=3,ftype="off",pts=F)
no colors provided. using the following legend:
    ac      b
"black"  "red"

Cool - seems to work.

No comments:

Post a Comment