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(fit$tree,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(fitJC$tree), 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))