I have aligned sequences in class of phyDat
. Here is a minimal example:
library(seqinr)
seqs <- DNAMultipleAlignment( c(a1 = "------ATGTTCATTAACCGCTGACTATTCTCAACCA",
a2 = "------ATGTTCACCGACCGCTGACTATTCTCTACAA",
a3 = "---GTGACCTTCATCAACCGATGATTATTCTCAA---",
a4 = "---TCAGTCGTCACCAGGCGTTG-CAGGACCCGAC--",
a5 = "ATGGGGGTCTTCCTCA-TCGCCGTCGCCGCGT-----"))
library(phangorn)
phyDat_seqs <- as.phyDat(seqs)
I followed the manual of phangorn package to construct a phylogenic tree. (Below, I wrote a simplified code. You don't need to check manual, you can run the code directly.)
dm <- dist.ml(phyDat_seqs)
treeNJ <- NJ(dm)
fit <- pml(treeNJ, data=phyDat_seqs)
fitGTR <- optim.pml(fit, model="GTR")
bs <- bootstrap.pml(fitGTR, bs=5, optNni=TRUE)
And I can simply plot this tree with plotBS
:
plotBS(fitGTR$tree, bs, type="p")
There is no problem so far. My question starts from this point.
I want to customize this tree by using ggtree
package but I couldn't figure out how to do that. When I use ggtree(fitGTR$tree)
, it gives me a different tree.
plotBS
somehow uses bs
values to change tree branching. If you compare plot(fitGTR$tree, type="p")
and plotBS(fitGTR$tree, bs, type="p")
difference can easily be seen. ggtree()
and plot()
construct the same trees, but not the tree that I want to customize.
I've tried different approaches such as saving tree with write.tree
or write.nexus
functions and calling the file again with ggtree
package functions but I can't save the plotBS
version of the tree. How can I get an output of plotBS
tree?
Thank you in advance.
Problem was occuring because I misread the documentation. I did a mistake in this line:
plotBS(fitGTR$tree, bs, type="p")
Tree must be constructed through midpoint()
function. Therefore, when I changed the line as this:
plotBS(midpoint(fitGTR$tree), bs, type="p")
no problem occurs. I can use ggtree
and I can save/write the tree data without a problem.
Even if I can't figure out what midpoint()
function does to tree actually (you can comment about this), problem is solved.