Hi Nick, seems ladderize is getting popular these days, second reply on the topic today (https://stat.ethz.ch/pipermail/r-sig-phylo/2015-October/004352.html). ladderize silently assume that trees are in "cladewise" order.
just run tr_wFossils = reorder(tr_wFossils) before ladderize, should fix your problem. The read.tree(write.tree) trick also returns a "cladewise" tree. Regards, Klaus On Sun, Oct 4, 2015 at 4:55 AM, Nick Matzke <[email protected]> wrote: > Hi all, > > I have been hitting intermittent problems using trees generated by > ape::rphylo. Here is a reproducible example. > > > ############################ > library(ape) > sessionInfo() > > # Simulate a tree with e.g. 5 species > nspecies = 5 > set.seed(123465) > tr_wFossils = rphylo(n=nspecies, birth=0.3, death=0.25, T0=50, > fossils=TRUE, eps=1e-06) > > # Ladderize tree with fossils > ltr = ladderize(tr_wFossils, right=FALSE) > ltr$edge > round(ltr$edge.length, 4) > write.tree(phy=ltr, file="") > > tr_wFossils = read.tree(file="", text=write.tree(phy=ltr, file="") ) > tr_wFossils$tip.label = paste0("sp", 1:length(tr_wFossils$tip.label)) > plot(tr_wFossils); axisPhylo(); title("Sim Tree w Fossils") > ############################ > > Output: > ############################ > > library(ape) > > sessionInfo() > R version 3.2.2 (2015-08-14) > Platform: x86_64-apple-darwin13.4.0 (64-bit) > Running under: OS X 10.10.5 (Yosemite) > > locale: > [1] en_AU.UTF-8/en_AU.UTF-8/en_AU.UTF-8/C/en_AU.UTF-8/en_AU.UTF-8 > > attached base packages: > [1] stats graphics grDevices utils datasets methods base > > other attached packages: > [1] ape_3.3 > > loaded via a namespace (and not attached): > [1] nlme_3.1-121 grid_3.2.2 lattice_0.20-33 > > > > # Simulate a tree with e.g. 5 species > > nspecies = 5 > > set.seed(123465) > > tr_wFossils = rphylo(n=nspecies, birth=0.3, death=0.25, T0=50, > fossils=TRUE, eps=1e-06) > > > > # Ladderize tree with fossils > > ltr = ladderize(tr_wFossils, right=FALSE) > > ltr$edge > [,1] [,2] > [1,] 7 11 > [2,] 11 2 > [3,] 8 5 > [4,] 8 9 > [5,] 9 6 > [6,] 9 10 > [7,] 10 4 > [8,] 10 3 > > round(ltr$edge.length, 4) > [1] 2.6862 0.0670 1.2078 0.1028 0.4985 0.0915 1.0135 1.0135 > > write.tree(phy=ltr, file="") > [1] "((t2:0.06701971429):2.686176768);" > > > > tr_wFossils = read.tree(file="", text=write.tree(phy=ltr, file="") ) > Warning message: > In FUN(X[[i]], ...) : NAs introduced by coercion > > tr_wFossils$tip.label = paste0("sp", 1:length(tr_wFossils$tip.label)) > > plot(tr_wFossils); axisPhylo(); title("Sim Tree w Fossils") > NULL > Warning message: > In plot.phylo(tr_wFossils) : found less than 2 tips in the tree > Error in get("last_plot.phylo", envir = .PlotPhyloEnv) : > object 'last_plot.phylo' not found > ############################ > > > > It looks like a problem with the edge matrix generated by rphylo - perhaps > older simulation code. The workaround is Yea Olde Writeth To And Readeth > >From Newick trick: > > ########################### > # Simulate a tree with e.g. 5 species > nspecies = 5 > set.seed(123465) > tr_wFossils = rphylo(n=nspecies, birth=0.3, death=0.25, T0=50, > fossils=TRUE, eps=1e-06) > tr_wFossils = read.tree(file="", text=write.tree(phy=tr_wFossils, file="") > ) > # Ladderize tree with fossils > ltr = ladderize(tr_wFossils, right=FALSE) > ltr$edge > round(ltr$edge.length, 4) > write.tree(phy=ltr, file="") > > tr_wFossils = read.tree(file="", text=write.tree(phy=ltr, file="") ) > tr_wFossils$tip.label = paste0("sp", 1:length(tr_wFossils$tip.label)) > plot(tr_wFossils); axisPhylo(); title("Sim Tree w Fossils") > ########################### > > > > Now there are no issues: > ########################### > > # Simulate a tree with e.g. 5 species > > nspecies = 5 > > set.seed(123465) > > tr_wFossils = rphylo(n=nspecies, birth=0.3, death=0.25, T0=50, > fossils=TRUE, eps=1e-06) > > tr_wFossils = read.tree(file="", text=write.tree(phy=tr_wFossils, > file="") ) > > # Ladderize tree with fossils > > ltr = ladderize(tr_wFossils, right=FALSE) > > ltr$edge > [,1] [,2] > [1,] 7 8 > [2,] 8 1 > [3,] 8 2 > [4,] 7 9 > [5,] 9 6 > [6,] 9 10 > [7,] 10 3 > [8,] 10 11 > [9,] 11 4 > [10,] 11 5 > > round(ltr$edge.length, 4) > [1] 2.6862 0.0670 0.0670 1.5454 1.2078 0.1028 0.4985 0.0915 1.0135 1.0135 > > write.tree(phy=ltr, file="") > [1] > > "((t2:0.06701971429,t1:0.06701971429):2.686176768,(t5:1.207768321,(t6:0.4984853178,(t4:1.013489985,t3:1.013489985):0.09146602338):0.1028123128):1.545428161);" > > > > tr_wFossils = read.tree(file="", text=write.tree(phy=ltr, file="") ) > > tr_wFossils$tip.label = paste0("sp", 1:length(tr_wFossils$tip.label)) > > plot(tr_wFossils); axisPhylo(); title("Sim Tree w Fossils") > ########################### > > Cheers, > Nick > > [[alternative HTML version deleted]] > > _______________________________________________ > R-sig-phylo mailing list - [email protected] > https://stat.ethz.ch/mailman/listinfo/r-sig-phylo > Searchable archive at > http://www.mail-archive.com/[email protected]/ > -- Klaus Schliep Postdoctoral Fellow Revell Lab, University of Massachusetts Boston [[alternative HTML version deleted]] _______________________________________________ R-sig-phylo mailing list - [email protected] https://stat.ethz.ch/mailman/listinfo/r-sig-phylo Searchable archive at http://www.mail-archive.com/[email protected]/
