Tuesday, April 4, 2017

Showing tip labels for some, but not all, of the tips in a plotted tree

A phytools user recently contacted me to ask me how to plot a "contMap" objecting showing some but not all tip labels. This is totally possible (of course!), but a bit trickier than it might seem. The following is a basic demo - then I may illustrate some more “advanced” options.

library(phytools)
tree
## 
## Phylogenetic tree with 26 tips and 25 internal nodes.
## 
## Tip labels:
##  P._kgfqhxm, B._dhwfteyp, J._lnrpb, X._kgbdui, A._mgcj, F._rxaikrligq, ...
## 
## Rooted; includes branch lengths.
x
##    P._kgfqhxm   B._dhwfteyp      J._lnrpb     X._kgbdui       A._mgcj 
##   1.221447299   4.999731654   3.324608498   4.245941152   1.867780375 
## F._rxaikrligq Y._fmejxidzgj  Q._hurzjgrgt   C._hurhmpyj N._inqhwrfkkg 
##   0.788551106   3.944202674   1.948093627   0.263497722  -0.418167671 
##       V._ohvc    O._dzohmrm   E._ifdmucjo  T._mplzwkmew       Z._bwuy 
##   2.253341739   3.153452128   2.249732267  -0.208776645   0.711507433 
##  K._iejohjylu     D._qpbdfz       H._hqes       U._olon  S._slfeoiwaq 
##   0.407903549   0.581238603   0.824611185  -0.219487212   0.490048186 
##     I._mzgvaq      W._nxtrp  R._sombrvlvo     G._zvjklv  L._wnjhcuwki 
##  -1.475669703  -0.221141585   0.447807409  -0.006821951   0.261657895 
## M._kkgayweykn 
##   0.789239770
obj<-contMap(tree,x)

plot of chunk unnamed-chunk-1

lastPP<-get("last_plot.phylo",envir=.PlotPhyloEnv)
xlim<-lastPP$x.lim
ylim<-lastPP$y.lim
plot(obj,xlim=xlim,ylim=ylim,ftype="off")
lastPP<-get("last_plot.phylo",envir=.PlotPhyloEnv)
tips ## tips to label
##  [1] "S._slfeoiwaq"  "N._inqhwrfkkg" "U._olon"       "R._sombrvlvo" 
##  [5] "G._zvjklv"     "E._ifdmucjo"   "Q._hurzjgrgt"  "L._wnjhcuwki" 
##  [9] "B._dhwfteyp"   "T._mplzwkmew"
xpos<-lastPP$xx[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]
ypos<-lastPP$yy[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]
text(xpos,ypos,gsub("_"," ",tips),pos=4,font=3)

plot of chunk unnamed-chunk-1

Sometimes, however, we might want to label only a subset of the tips because all the labels don't fit (at the desired font size), and we'd like to space them out. For this we could try to use linking lines as follows:

plot(obj,xlim=xlim,ylim=ylim,ftype="off")
lastPP<-get("last_plot.phylo",envir=.PlotPhyloEnv)
xpos<-lastPP$xx[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]+strwidth("i")
ypos<-lastPP$yy[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]
xmax<-rep(max(lastPP$xx),length(tips))+strwidth(" ")
ylab<-seq(1,Ntip(tree),by=(Ntip(tree)-1)/(length(tips)-1))
ylab<-ylab[rank(ypos)]
text(xmax,ylab,gsub("_"," ",tips),pos=4,font=3,cex=1.2,
    offset=0)
segments(xpos,ypos,xmax,ylab,lty="dotted",col="blue")

plot of chunk unnamed-chunk-2

Of course, here we have the problem that the linking lines cross the edges of the tree. We could try to substitute sigmoidal linking lines as follows:

plot(obj,xlim=1.2*xlim,ylim=ylim,ftype="off")
lastPP<-get("last_plot.phylo",envir=.PlotPhyloEnv)
xpos<-lastPP$xx[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]+strwidth("i")
ypos<-lastPP$yy[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]
xmax<-rep(max(lastPP$xx),length(tips))+0.2/1.2*xlim[2]
ylab<-seq(1,Ntip(tree),by=(Ntip(tree)-1)/(length(tips)-1))
ylab<-ylab[rank(ypos)]
text(xmax,ylab,gsub("_"," ",tips),pos=4,font=3,cex=1.2,
    offset=0)
for(i in 1:length(tips))
    phytools:::drawCurve(c(xpos[i],xmax[i]),c(ypos[i],ylab[i]),
        scale=0.2,lty="dotted",col="blue",lwd=2)

plot of chunk unnamed-chunk-3

This looks really good in this case; however in general it might still result in line crossing. A final alternative that I offer is as follows:

plot(obj,xlim=1.2*xlim,ylim=ylim,ftype="off")
lastPP<-get("last_plot.phylo",envir=.PlotPhyloEnv)
xpos<-lastPP$xx[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]+strwidth("i")
ypos<-lastPP$yy[sapply(tips,function(x,y) which(y==x),
    y=obj$tree$tip.label)]
xmax<-rep(max(lastPP$xx),length(tips))+0.2/1.2*xlim[2]
ylab<-seq(1,Ntip(tree),by=(Ntip(tree)-1)/(length(tips)-1))
ylab<-ylab[rank(ypos)]
text(xmax,ylab,gsub("_"," ",tips),pos=4,font=3,cex=1.2,
    offset=0)
tipmax<-max(nodeHeights(obj$tree))
for(i in 1:length(tips)){
    ff<-strwidth("W")
    segments(xpos[i],ypos[i],tipmax,ypos[i],lty="dotted",
        col="blue",lwd=2)
    segments(tipmax,ypos[i],tipmax+ff,ylab[i],lty="dotted",
        col="blue",lwd=2)
    segments(tipmax+ff,ylab[i],xmax[i],ylab[i],lty="dotted",
        col="blue",lwd=2)
}

plot of chunk unnamed-chunk-4

Neat.

The data for this exercise were simulated as follows:

foo<-function() paste(sample(letters,sample(4:10,1),replace=TRUE),
    collapse="")
tip.labels<-paste(LETTERS,"._",replicate(26,foo()),sep="")
tree<-rtree(n=26,tip.label=tip.labels)
x<-fastBM(tree)
tips<-sample(tree$tip.label,10)

That's it.

No comments:

Post a Comment