Monday, March 6, 2023

Mapping the value of a continuous trait at only the nodes of a plotted phylogenetic tree using contMap

Inspired by some visualizations during our recent R phylogenetics workshop at UMass-Boston…

… as well as by similar plots in the recent scientific literature…

… I decided to add a simple nodes_only=TRUE option to the S3 plot.contMap method for the phytools "contMap" object class. Updated R code can be seen here.

Here’s how it works.

First, load phytools (updated from GitHub).

library(phytools)
packageVersion("phytools")
## [1] '1.5.11'

Now, load a dataset. Here, I’ll use the mammal tree & body mass data from Garland et al. (1992).

data(mammal.tree)
data(mammal.data)
lnBodyMass<-setNames(log(mammal.data$bodyMass),
	rownames(mammal.data))

We can next create a "contMap" object as follows.

Here is a typical visualization.

mammal.contMap<-contMap(mammal.tree,lnBodyMass,plot=FALSE)
mammal.contMap<-setMap(mammal.contMap,palette()[c(2,4)])
plot(mammal.contMap,fsize=c(0.7,0.9),ftype="i",
	leg.txt="log(body mass)")

plot of chunk unnamed-chunk-3

Now let’s test the nodes_only=TRUE option.

plot(mammal.contMap,nodes_only=TRUE,fsize=c(0.7,0.9),
	ftype="i",leg.txt="log(body mass)")

plot of chunk unnamed-chunk-4

So far so good.

The implementation works with trees of different types and orientation. Here is a “fan” style tree plot.

data(anoletree)
data(anole.data)
svl<-setNames(anole.data$SVL,rownames(anole.data))
anole.contMap<-contMap(anoletree,svl,plot=FALSE)
anole.contMap<-setMap(anole.contMap,viridisLite::viridis(10))
plot(anole.contMap,nodes_only=TRUE,fsize=c(0.7,0.8),
	ftype="i",leg.txt="log(SVL)",type="fan",legend=5,
	offset=2)

plot of chunk unnamed-chunk-5

For this plot, I used a phylogeny & dataset of Anolis lizards from Mahler et al. (2010).

Lastly, here’s a visualization in which I also create and add a custom legend.

data(primate.tree)
data(primate.data)
lnSkull<-setNames(log(primate.data$Skull_length),
	rownames(primate.data))
primate.contMap<-contMap(primate.tree,lnSkull,plot=FALSE)
primate.contMap<-setMap(primate.contMap,
	viridisLite::viridis(n=20,direction=-1))
plot(primate.contMap,nodes_only=TRUE,fsize=0.6,ftype="i",
	legend=FALSE,type="cladogram",nodes="weighted",
	xlim=c(-15,95))
LWD<-diff(par()$usr[1:2])/dev.size("px")[1]
lines(x=rep(-0.15*max(nodeHeights(primate.tree))+
	LWD*10/2,2),y=c(5,Ntip(primate.tree)-5))
nticks<-10
Y<-cbind(seq(5,Ntip(primate.tree)-5,
	length.out=nticks),seq(5,Ntip(primate.tree)-5,
	length.out=nticks))
X<-cbind(rep(-0.15*max(nodeHeights(primate.tree))+
	LWD*10/2,nticks),rep(-0.15*max(nodeHeights(primate.tree))+
	LWD*10/2+0.02*max(nodeHeights(primate.tree)),nticks))
for(i in 1:nrow(Y)) lines(X[i,],Y[i,])
add.color.bar(Ntip(primate.tree)-10,primate.contMap$cols,
	title="log(skull length)",lims=NULL,digits=3,
	direction="upwards",subtitle="",lwd=10,
	x=-0.15*max(nodeHeights(primate.tree)),y=5,
	prompt=FALSE)
ticks<-seq(primate.contMap$lims[1],primate.contMap$lims[2],
	length.out=nticks)
text(x=X[,2],y=Y[,2],prettyNum(exp(ticks),digits=3),pos=4,
	cex=0.8)

plot of chunk unnamed-chunk-6

In this case, the tree & data are from Kirk & Kay (2004).

That’s it!

No comments:

Post a Comment

Note: due to the very large amount of spam, all comments are now automatically submitted for moderation.