Skip to content

Instantly share code, notes, and snippets.

# cmtucker/caRd2015.RSecret Last active Dec 17, 2015

The EEB & Flow's annual holiday caRd, 2015
 ##2015 Holiday caRd from the EEB & Flow, by C. Tucker ##If you don't have these libraries, please install before running... #I didn't test this on different OS/versions of R this year, so if you run into problems refer to the blog post instead :-) ##Req'd libraries require(ape) require(rworldmap) require(FD) # function modified from on phytools (& dependencies) and rworldmap (& dependencies) (written by Liam J. Revell 2013) phylo.to.map2<-function(tree,coords){ # open & size a new plot par(mai=c(0.1,0.1,0.1,0.1)) map<-getMap(resolution="low") plot(map,ylim=c(-30,100)) # rescale tree so it fits in the upper half of the plot # with enough space for labels sh<-max(strwidth(tree\$tip.label))/((par()\$usr-par()\$usr)*(par()\$usr-par()\$usr)) tree\$edge.length<-tree\$edge.length/max(nodeHeights(tree))*(90-sh) n<-length(tree\$tip.label) # reorder cladewise to assign tip positions cw<-reorder(tree,"cladewise") x<-vector(length=n+cw\$Nnode) x[cw\$edge[cw\$edge[,2]<=n,2]]<-0:(n-1)/(n-1)*360-180 # reorder pruningwise for post-order traversal pw<-reorder(tree,"pruningwise") nn<-unique(pw\$edge[,1]) # compute horizontal position of each edge for(i in 1:length(nn)){ xx<-x[pw\$edge[which(pw\$edge[,1]==nn[i]),2]] x[nn[i]]<-mean(range(xx)) } # compute start & end points of each edge Y<-180-nodeHeights(cw) # plot vertical edges for(i in 1:nrow(Y)) lines((rep(x[cw\$edge[i,2]],2)+10),(Y[i,]+10),lwd=2,lend=2, col="forestgreen") # plot horizontal relationships for(i in 1:tree\$Nnode+n) lines((range(x[cw\$edge[which(cw\$edge[,1]==i),2]])+10),(Y[which(cw\$edge[,1]==i),1]+10),lwd=2,lend=2, col="forestgreen") # plot tip labels lab_short <- c("Coca \ncola Santa","Mall \nSanta","Salvation \nArmy.Santa","Kris.Kringle", "Santa \nClaus","Pere.Noel","Sinterklaas","Saint \nNicholas", "Father \nChristmas", "Ded \nMoroz") for(i in 1:n) text(x[i],Y[which(cw\$edge[,2]==i),2],lab_short[i],srt=0, col="red", cex=0.5) coords<-coords[tree\$tip.label,2:1] points(coords,pch=16,cex=1,col="red") for(i in 1:n) lines(c(x[i],coords[i,1]),c(Y[which(cw\$edge[,2]==i),2],coords[i,2]),col="red",lty="dashed") } #also from phytools nodeHeights <- function (tree) { if (attr(tree, "order") != "cladewise" || is.null(attr(tree, "order"))) t <- reorder(tree) else t <- tree root <- length(t\$tip.label) + 1 X <- matrix(NA, nrow(t\$edge), 2) for (i in 1:nrow(t\$edge)) { if (t\$edge[i, 1] == root) { X[i, 1] <- 0 X[i, 2] <- t\$edge.length[i] } else { X[i, 1] <- X[match(t\$edge[i, 1], t\$edge[, 2]), 2] X[i, 2] <- X[i, 1] + t\$edge.length[i] } } if (attr(tree, "order") != "cladewise" || is.null(attr(tree, "order"))) o <- apply(matrix(tree\$edge[, 2]), 1, function(x, y) which(x == y), y = t\$edge[, 2]) else o <- 1:nrow(t\$edge) return(X[o, ]) } ######################## ##Make tree edge <- matrix(c(11,12,12,13,13,14,14,9,14,15,15,6,15,7,13,8, 12, 10,11,16,16,17,17,18,18,1,18,2,17,3,16,19,19,4,19,5), nrow=18, ncol=2, byrow=TRUE) edge.length <- c(0.1, 0.4, 0.1, 0.5, 0.3, 0.2, 0.2, 0.6, 1.0, 0.2, 0.3, 0.3, 0.3, 0.3, 0.6, 0.1, 0.8, 0.8) tip.label <- c("Coca.cola.Santa", "Mall.Santa", "Salvation.Army.Santa", "Kriss.Kringle", "Santa.Claus", "Pere.Noel", "Father.Christmas", "Sinterklaas", "Saint.Nicholas", "Ded.Moroz") Nnode <- 9 xmas.tree <- list(edge, edge.length, Nnode, tip.label) class(xmas.tree) <- "phylo" names(xmas.tree) <- c("edge", "edge.length", "Nnode", "tip.label") ####Start #taxa exploration #functional ecology heft <- c("fat", "fat", "fat", "fat", "fat", "thin", "fat","thin", "thin", "thin") transport <- c("reindeer", "reindeer", "reindeer", "reindeer", "reindeer", "donkey", "foot", "horse", "foot", "foot") first.appearance <- c(1900, 1900, 1890, 1800, 1700, 1400, 1400, 400, 300, 1937) santatraits <- cbind.data.frame(heft, transport, first.appearance) rownames(santatraits) <- xmas.tree\$tip.label ###Card plotting #Make device with set size, works on multiple OS dev.new <- function(width = 5, height = 5){ platform <- sessionInfo()\$platform if (grepl("linux",platform)) { x11(width=width, height=height) }else{ if (grepl("pc",platform)) { windows(width=width, height=height) }else{ if (grepl("w32",platform)) { windows(width=width, height=height) }else{ if (grepl("apple", platform)){ quartz(width=width, height=height)} }}}} dev.new(10, 6) par(mfrow=c(1,2), bg="darkseagreen1", mai=c(0.1,0.1,0.1,0.1), oma=c(0.1,0.1,0.1,0.1)) plot(xmas.tree, type = "c", FALSE, edge.color="darkgreen", edge.lty=1, edge.width=18, label.offset = 1, direction="downward", font=3, tip.color="darkred") #reconstruct ancestral state cc <- ace(transport, xmas.tree, type="discrete") co2 <- c("yellow", "gold", "darkorange", "red") nodelabels(pie = cc\$lik.anc, piecol = co2, cex = c(1.5, rep(1, 8))) #Plot traits (fatness) against Santa co1 <- c("blue", "purple") tiplabels(pch = 12, col = co1[as.factor(heft)], cex = 3.5, adj=c(0.5, 0), lwd=2) #Let's see the transportation mode trait too: tiplabels(pch = 8, col = co2[as.factor(transport)], cex = 2, adj=c(0.5, 0), lwd=2) #where do santa taxa live? locales <- cbind(c(40.7, 40.7, 40.7, 70.93, 70.93, 47, 51.4, 52.2, 39.8, 55.7), c( -73.93, -73.93, -73.93, -95.4, -95.4, 2.7, -0.85, 5.7, 33.3, 37.6)) rownames(locales) <- xmas.tree\$tip.label phylo.to.map2(xmas.tree, locales) text(x=-20, y=-120, "Happy Holidays from \nthe EEB & Flow", col="darkred", cex=2, font=4)
to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.