Comparing phylogeny in R

So I want to compare the phylogeny created using two methods for example Maximum likelihood and maximum parsimony.Is there any way to compare the two phylogeny ?

I did read about phangorn but not sure if its the right R library for comparative analysis.

Any suggestion or helped would be highly appreciated

My data file

library(phangorn)
library(phytools)
library(dendextend)

data <- read.dna("abhi_seq/clean_dup_align_fast.fas", format = "fasta")
data
dat <- as.phyDat(data)

dm <- dist.ml(dat)
treeUPGMA <- upgma(dm)
treeNJ <- NJ(dm)

layout(matrix(c(1,2), 2, 1), height=c(1,2))
par(mar = c(0,0,2,0)+ 0.1)
plot(treeUPGMA, main="UPGMA")
plot(treeNJ, "phylogram", main="NJ")

dev.off()

parsimony(treeUPGMA, dat)
parsimony(treeNJ,dat)

tr.mp <- optim.parsimony(treeNJ, dat)

#tr.ml = optim.pml(treeNJ, dat)

fit <- pml(treeNJ, dat)
fit <- optim.pml(fit, rearrangement="NNI")

fit.ini <- pml(treeNJ, dat)
fit.ini

fit <- optim.pml(fit.ini, optNni=TRUE, optBf=TRUE, optQ=TRUE, optGamma=TRUE)
fit

tr.ml <- root(fittree,1) tr.mp.ultra<-force.ultrametric(tr.mp) tr.ml.ultra<-force.ultrametric(tr.ml) is.ultrametric(tr.mp.ultra) is.binary.tree(tr.mp.ultra) is.rooted(tr.mp.ultra) dd.ml.ultra<-as.dendrogram(tr.mp.ultra)  Error in ape::as.hclust.phylo(object) : the tree is not rooted I ran into this error this error which says trees is not rooted Comparing phylogeny code updated ** Working thanks to ***thomas duge de bernonville* thomas for putting codes together and fixing the errors library(dendextend) library(seqinr) library(phytools) library(phangorn) a<-read.alignment("abhi_seq/clean_dup_align_fast.fas", format="fasta") a.phydat<-as.phyDat(a) dist.a.phydat<-dist.dna(as.DNAbin(a.phydat)) upgma.a<-upgma(dist.a.phydat) parsimony(upgma.a,a.phydat) pars.a <- optim.parsimony(upgma.a, a.phydat) pars.a<-acctran(pars.a, a.phydat) pars.a.rooted<-root(pars.a, outgroup="AAA64460", resolve.root=T) pars.a.rooted.dd<-as.dendrogram(force.ultrametric(pars.a.rooted)) mt <- modelTest(a.phydat, tree=upgma.a,multicore = TRUE,mc.cores=10) #ml.a = pml(upgma.a,a.phydat) #fitJC <- optim.pml(ml.a, TRUE) #ml.a.rooted<-root(midpoint(fitJCtree), outgroup="AAA64460", resolve.root=T)
#ml.a.rooted.dd<-as.dendrogram(force.ultrametric(ml.a.rooted))

#########################################################################

######################################################################33

#mt <- modelTest(dat, tree=tree, multicore=TRUE)
mt[order(mt$$AICc),]# choose best model from the table according to AICc bestmodel <- mt$$Model[which.min(mt$AICc)] env = attr(mt, "env") fitStart = eval(get(bestmodel, env), env) fit = optim.pml(fitStart, rearrangement = "stochastic",optGamma=TRUE, optInv=TRUE, model="GTR")#tree bs=bootstrap.pml(fit, bs=25, optNni=TRUE, multicore=TRUE) ml.a.rooted<-root(midpoint(fit$tree), outgroup="AAA64460", resolve.root=T)
ml.a.rooted.dd<-as.dendrogram(force.ultrametric(ml.a.rooted))


4 Answers

phangorn is a really powerful package for phylogenies. But to compare the trees, I think you may convert them into dendrograms and calculate a correlation measure, such as the Fowlkes-Mallows Index or a distance measure such as Baker’s Gamma Index. These can be easily computed using the dendextend R package (https://cran.r-project.org/web/packages/dendextend/vignettes/dendextend.html#correlation-measures).

• Yes, you can build the two trees and then convert them into dendrograms. Be careful, trees should be ultrametric for that. You can convert them using the force.ultrametric function from phytools package followed by as.dendrogram from dendextend package. May 25, 2020 at 19:56
• I think you can still use phangorn to construct the phylogenies.You just have to read your alignment with the read.alignment function of seqinr package, then convert it into a phyDat object, then calculate a distance with dist.ml. From that point, you may construct ML or parsimony trees (cran.r-project.org/web/packages/phangorn/vignettes/Trees.pdf). Once you've obtained the trees, convert them into dendrograms and calculate a correlation as described above. May 26, 2020 at 4:51
• from your phyDat object, you may estimate a distance eg d<-dist.ml(dat), then construct a rooted starting tree with UPGMA, tr.upgma<-upgma(d) then from this starting tree construct either ML or MP following the phangorn vignette. Once you get the bootstrapped trees, set them as ultrametric tr.ml.ultra<-force.ultrametric(tr.ml) and convert them as dendrogram dd.ml.ultra<-as.dendrogram(tr.ml.ultra). Please load the appropriate libraries to do that, and to have a look at the parameters that can be tuned at each step. May 26, 2020 at 7:13
• I think you should use the UPGMA tree as a starting tree, both for the MP and ML, because I think it produces rooted trees, in contrast to NJ. May 26, 2020 at 12:08
• May 26, 2020 at 12:55

Another 2 tree compare tools: ggtree and blatic 3

• The tree building is usually used IQ-Tree, MrBayes, RAxML, and beast to get a robust maximum likelihood or maximum-credibility tree. This process is time-consuming and computational resource consuming. If I want to use the dendextend to compare the phylogenetic trees which are built by the methods mentioned above rather than in R, can I do it with dendextend? May 29, 2020 at 1:29

Bootstrapping I agree correlation is one way in, but the classic way to do this is via bootstrapping, i.e. resampling the alignment with replacement and making a consensus phylogeny of of 100 to 1000 replicates. Values exceeding 75 to 80% define node robustness and you can calculate incongruence from there.

Phangorn certainly does bootstrapping, but the efficiency of the algorhithm for maximum likelihood (ML) I don't know (might take a very long time), just think how long one tree takes and times it by 100 (or even 1000).

If you just performed a correlation measure there is a risk the differences are not robust, but if you've got loads of taxa you could bootstrap collapse any value less than the robustness threshold into a polytomy and then perform a correlation method. Its a bit complicated, but it would work.

Bootstrapping parsimony isn't hard nor computationally expensive, bootstrapping maximum likelihood sometimes requires specific efficient algorhithms notably RAxML.

Non-R The other way in is just dump the trees into Dendroscope3, okay its not R but will work. Dendroscrope is hard for a non-tree person to read, however it will flag the differences for a specialist. Dendroscope will make a network tree for the topological differences, in this case between methods, i.e. non-bifurcating. You really then usually have to present the differences as a "mirror tree", i.e. parsimony and ML side by side, because usually non-tree people wonder what a Dendroscope output means. Technically there is nothing wrong with describing incongruence via networks.

Basically, the heart-beat is bootstrapping.

Generally speaking parsimony will give the same basic output as ML, except when you have some taxa undergoing rapid evolution against their sister group. This does happen in pathogen evolution and parsimony falls over into a phenomina known as 'long-branch attraction'.

• The advantage of the correlation measure will be the metrics which provide a more objective evaluation. Michael's recommandation with bootstrapping is really good, especially using collapsed branches. May 25, 2020 at 11:49

If you only have to compare two phylogenies (as opposed to more), I would suggest that an alternative to using R would be using iqtree . One of its parameters (-rf) can calculate the Robinson-Foulds metric between two trees. (you can also use iqtree to calculate the phylogenies in the first place)

It is command-line based, but can be easily installed with conda. Having anaconda installed in general makes bioinformatics work easier, I feel.

Alternatively, differences in trees can be visualized (without numbers) in Dendroscope, or Splitstree. They have different strategies for it. With Dendroscope it is easier to compare side-by-side and, if having precisely the same taxon names, to draw lines between taxa in different trees. Splitstree can accept an input of multiple trees (File>Tools>Load Multiple trees) and build a network (which can be a way to visualize the consensus of the trees).

• Hi @Laura, SplitsTree is very cool, but by default will not work here because classically it will identify incongruence within the sequence itself (not just between trees). There may be options to disable this (if so please post them), but by default thats what it will attempt.
– M__
May 25, 2020 at 11:27
• @Michael Thank you for pointing that out. I meant Splitstree starting from trees, not from sequence (File>Tools>Load Multiple trees). I've used it before for visualizing bootstraping results, but technically it should work on fewer trees too? May 25, 2020 at 13:22