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)
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)
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")
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)
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)
}
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.
Hi Liam,
ReplyDeleteHow to plot a "contMap" objecting showing aligned all tip labels?
Rafa B.