Commit 045253f9 by Jannes Höpke

### Upload New File

parent f2b04835
 ### Script by: ### GLOR, Rich (2008). "R Tip: Labeling Trees w/ Posterior Probability and Bootstrap Support". ### http://treethinkers.blogspot.de/2008/10/labeling-trees-posterior-probability.html ### (access date: 4th August 2015) ### Modified by: ### HÖPKE, Jannes #Set working directory setwd("~/Documents/Uni Oldenburg/6. Semester/Bachelorarbeitsmodul/Phylogenetische Bäume in R/trnLF_2") #Fist we need to open some necessary libraries library(ape) library(geiger) #Some plot options #define cex-size for tiplabel tcex<-0.7 #define cex-size for support scex<-0.6 #define cex-size for tree line lcex<-1.5 #Read in trees read.nexus("trnLF_best.tre")->bestTree #FigTree output (rooted) read.nexus("trnLF_mb_consensus.tre")->bayesTree #SumTrees.py output read.nexus("trnLF_garli_consensus.tre")->bootTree #SumTrees.py output #Setting plot options opar<-par() par(mar=c(0,0,0,0),xpd=T) #Reroot bestTree plot(bestTree,cex=tcex) library(phytools) #Adjust bestTree via "rotate" nodelabels() rrbestTree<-rotate(bestTree, 59, polytom = c(1,2)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 60, polytom = c(1,2)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 61, polytom = c(1,2)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 67, polytom = c(1,7)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 67, polytom = c(2,10)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 67, polytom = c(3,11)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 67, polytom = c(4,6)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 67, polytom = c(5,7)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 70, polytom = c(1,2)) plot(rrbestTree,show.node.label=F,cex=tcex) nodelabels() rrbestTree<-rotate(rrbestTree, 72, polytom = c(3,4)) plot(rrbestTree,show.node.label=F,cex=tcex) #Create new labels for bestTree bestLabels<-read.csv("bestTree_tips.csv", header=T, sep=";") str(bestLabels) attach(bestLabels) bestLabels<-data.frame(lapply(bestLabels, as.character), stringsAsFactors=FALSE) #converts all colums from a dataframe into characters detach(bestLabels) attach(bestLabels) str(bestLabels) bestTips <- mixedFontLabel(new1," ",new2,", ",new3," ",new4, italic = 1:3, parenthesis = 7, sep="") #Reroot bootTree plot(bootTree,cex=tcex) nodelabels() tiplabels() library(phytools) i<-63 rbootTree <- reroot(bootTree, i, position = 0.5 * bootTree\$edge.length[which(bootTree\$edge[,2] == i)]) plot(rbootTree,show.node.label=T,cex=tcex) #Adjust bootTree via "rotate" nodelabels() rrbootTree<-rotate(rbootTree, 64, polytom = c(1,2)) plot(rrbootTree,show.node.label=T,cex=tcex) nodelabels() rrbootTree<-rotate(rrbootTree, 72, polytom = c(5,9)) plot(rrbootTree,show.node.label=T,cex=tcex) nodelabels() rrbootTree<-rotate(rrbootTree, 66, polytom = c(1,8)) plot(rrbootTree,show.node.label=T,cex=tcex) #Delete "Root" labeling rrbootTree\$node.label[which(rrbootTree\$node.label=="Root")]<-"" plot(rrbootTree,show.node.label=T,cex=tcex) #Reroot bayesTree plot(bayesTree, show.node.label=T,cex=tcex) nodelabels() i<-63 rbayesTree <- reroot(bayesTree, i, position = 0.5 * bayesTree\$edge.length[which(bayesTree\$edge[,2] == i)]) plot(rbayesTree,show.node.label = T,cex=tcex) #Adjust bootTree via "rotate" nodelabels() rrbayesTree<-rotate(rbayesTree, 64, polytom = c(1,2)) plot(rrbayesTree,show.node.label=T,cex=tcex) nodelabels() rrbayesTree<-rotate(rrbayesTree, 73, polytom = c(1,3)) plot(rrbayesTree,show.node.label=T,cex=tcex) nodelabels() rrbayesTree<-rotate(rrbayesTree, 68, polytom = c(1,2)) plot(rrbayesTree,show.node.label=T,cex=tcex) #Delete "Root" labeling rrbayesTree\$node.label[which(rrbayesTree\$node.label=="Root")]<-"" plot(rrbayesTree,show.node.label=T,cex=tcex) #The getAllSubTrees function below is a necessary subfunction that atomizes a tree into #each individual subclade and was provided compliments of Luke Harmon. getAllSubtrees<-function(phy, minSize=2) { res<-list() count=1 ntip<-length(phy\$tip.label) for(i in 1:phy\$Nnode) { l<-tips(phy, ntip+i) bt<-match(phy\$tip.label, l) if(sum(is.na(bt))==0) { st<-phy } else st<-drop.tip(phy, phy\$tip.label[is.na(bt)]) if(length(st\$tip.label)>=minSize) { res[[count]]<-st count<-count+1 } } res } #The plotBootBayes function below plots both posterior probability and bootstrap values #on each node of the consensus tree obtained from your Maximum-Likelihood analysis. #Bootstrap values will appear above and to the left of the node they support, whereas #Bayesian posterior probabilies will appear below and to the left of the node. tiff("trnLF_bestTree.tiff", width = 20, height = 25, units = 'cm', res = 300) par(mar=c(2,0,0,8),xpd=T) getAllSubtrees(rrbayesTree)->bayesSub getAllSubtrees(rrbootTree)->bootSub getAllSubtrees(rrbestTree)->bestSub bayesList<-matrix("",Nnode(rrbestTree),1) bootList<-matrix("",Nnode(rrbestTree),1) #The commands below compare all the subclades in the Bayes tree to all the subclades #in the bootstrap tree, and vice versa, and identifies all those clades that are #identical. for(i in 1:Nnode(rrbestTree)) { for(j in 1:Nnode(rrbayesTree)) { match(bestSub[[i]]\$tip.label[order(bestSub[[i]]\$tip.label)], bayesSub[[j]]\$tip.label[order(bayesSub[[j]]\$tip.label)])->shared match(bayesSub[[j]]\$tip.label[order(bayesSub[[j]]\$tip.label)], bestSub[[i]]\$tip.label[order(bestSub[[i]]\$tip.label)])->shared2 if(sum(is.na(c(shared,shared2)))==0) { rrbayesTree\$node.label[j]->bayesList[i] }}} for(i in 1:Nnode(rrbestTree)) { for(j in 1:Nnode(rrbootTree)) { match(bestSub[[i]]\$tip.label[order(bestSub[[i]]\$tip.label)], bootSub[[j]]\$tip.label[order(bootSub[[j]]\$tip.label)])->shared match(bootSub[[j]]\$tip.label[order(bootSub[[j]]\$tip.label)], bestSub[[i]]\$tip.label[order(bestSub[[i]]\$tip.label)])->shared2 if(sum(is.na(c(shared,shared2)))==0) { rrbootTree\$node.label[j]->bootList[i] }}} plot(rrbestTree, cex=tcex, edge.width =lcex, show.tip.label=F) #Plots your Maximum-Likelihood consensus tree tiplabels(bestTips, adj = -0.02, frame = "none", col = "black", bg="white",cex=tcex) nodelabels(bootList, adj=c(1.2, -0.3), frame="n", cex=scex, font=1) nodelabels(bayesList, adj=c(1.12, 1.3), frame="n", cex=scex, font=1) add.scale.bar(x = 0, y = -1, length = NULL, cex=tcex) dev.off() # --> graphical output! par(opar) \ No newline at end of file
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment