Monday, August 15, 2016

Getting all pairs of sister species in a tree

Here is some code to make a matrix containing all pairs of sister species - that is, the descendants of all nodes leading to two & only leafs. It uses phangorn's Descendants function & apply family functions:

library(phytools)
library(phangorn)
tree<-rtree(n=100)
## starting here
dd<-lapply(1:tree$Nnode+Ntip(tree),function(n,t)
    Descendants(t,n)[[1]],t=tree)
nodes<-c(1:tree$Nnode+Ntip(tree))[which(sapply(dd,
    length)==2)]
sisters<-t(sapply(nodes,function(n,t) 
    t$tip.label[Descendants(t,n)[[1]]],t=tree))
rownames(sisters)<-nodes
sisters
##     [,1]  [,2]  
## 106 "t74" "t99" 
## 110 "t70" "t53" 
## 111 "t49" "t10" 
## 114 "t19" "t52" 
## 117 "t68" "t50" 
## 120 "t72" "t87" 
## 127 "t9"  "t8"  
## 131 "t75" "t46" 
## 133 "t66" "t84" 
## 137 "t86" "t81" 
## 138 "t73" "t24" 
## 141 "t1"  "t94" 
## 146 "t42" "t33" 
## 150 "t83" "t18" 
## 151 "t96" "t67" 
## 155 "t88" "t85" 
## 156 "t57" "t22" 
## 160 "t35" "t89" 
## 161 "t71" "t14" 
## 167 "t29" "t91" 
## 169 "t65" "t16" 
## 170 "t60" "t2"  
## 175 "t28" "t95" 
## 177 "t13" "t63" 
## 179 "t54" "t21" 
## 181 "t90" "t36" 
## 185 "t26" "t59" 
## 186 "t6"  "t30" 
## 187 "t3"  "t79" 
## 191 "t47" "t77" 
## 195 "t82" "t31" 
## 197 "t41" "t100"
## 198 "t11" "t76" 
## 199 "t15" "t5"

Let's cross check them:

plotTree(tree,fsize=0.6,lwd=1)
nodelabels(node=nodes,pch=21,cex=1.2,bg="red")

plot of chunk unnamed-chunk-2

That's it.

No comments:

Post a Comment