Scripts to analyze data using TreeMix. This pipeline runs TreeMix with bootstrapping, helps choose number of migration events and creates a consensus tree. It plots the maximum likelihood tree with bootstrap values, drift and residuals and calculates statistics for every migration event, such as migration support, standard error and p-values.
Then I continue with the next step using one of the unique ML trees. But the plot does not show bootstrap values. Am I doing something wrong? What I do is below:
_pdf("TreeMix_output.pdf")
treemix.bootstrap("m3_3m_finalrun_1", out.file = "tmp", #stem of TreeMix files (as above) + number of run with highest likelihood and unique topology
phylip.file = "Consensus.newick", #consensus tree in newick format (from the bootstrap procedure generated with PHYLIP)
nboot = 300, fill = TRUE, #nboot is the number of bootstraps used
pop.color.file = "col.txt", #specify colours with a tab delimited pop.color.file - first column is pop name, second the colour
boot.legend.location = "topright")
treemix.drift(in.file = "m3_3m_finalrun_1", #pairwise matrix for drift estimates with specified order of pops
pop.order.color.file = "poporder.txt") +
title("Drift")
"In dcast(data = DiffX, formula = B1 ~ B2, fill = 0, value.var = "DELTA") :
The dcast generic in data.table has been passed a data.frame and will attempt to redirect to the reshape2::dcast; please note that reshape2 is deprecated, and this redirection is now deprecated as well. Please do this redirection yourself like reshape2::dcast(DiffX). In the next version, this warning will become an error."
2nd from plot_resid.R:
"Warning message:
In xy.coords(x, y, xlabel, ylabel, log) : NAs introduced by coercion"
Hello,
I come accross a problem in each step. This time about showing the bootstrap values on the tree.
After running the step 3, I use maxLL and cfTrees functions to find unique trees and the likelihood values. Then I get these results:
_Input trees: 100 Best-known ML trees: 100
Number of duplicates in the Best-known ML trees: 98 Duplicate trees #: 2 3 4 5 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
Number of unique ML tree(s): 2 Unique ML tree(s) #: m3_3m_finalrun_1.treeout.gz m3_3m_finalrun_6.treeout.gz
Pairwise topological distance ( method = PH85 ) between unique trees: tree1 tree6 tree1 0
tree6 0 0
Job done._
Then I continue with the next step using one of the unique ML trees. But the plot does not show bootstrap values. Am I doing something wrong? What I do is below:
_pdf("TreeMix_output.pdf")
treemix.bootstrap("m3_3m_finalrun_1", out.file = "tmp", #stem of TreeMix files (as above) + number of run with highest likelihood and unique topology phylip.file = "Consensus.newick", #consensus tree in newick format (from the bootstrap procedure generated with PHYLIP)
nboot = 300, fill = TRUE, #nboot is the number of bootstraps used pop.color.file = "col.txt", #specify colours with a tab delimited pop.color.file - first column is pop name, second the colour boot.legend.location = "topright")
treemix.drift(in.file = "m3_3m_finalrun_1", #pairwise matrix for drift estimates with specified order of pops pop.order.color.file = "poporder.txt") + title("Drift")
plot_resid("m3_3m_finalrun_1", #pairwise matrix for residuals poporder = "poporder.txt") + title("Residuals")
dev.off()
I also get 2 warning messages:
1st from treemix.drift function:
"In dcast(data = DiffX, formula = B1 ~ B2, fill = 0, value.var = "DELTA") : The dcast generic in data.table has been passed a data.frame and will attempt to redirect to the reshape2::dcast; please note that reshape2 is deprecated, and this redirection is now deprecated as well. Please do this redirection yourself like reshape2::dcast(DiffX). In the next version, this warning will become an error."
2nd from plot_resid.R: "Warning message: In xy.coords(x, y, xlabel, ylabel, log) : NAs introduced by coercion"
Can I ignore them?
Thanks!