This vignette briefly demonstrates multiple features of babette
, without going into much detail.
First, load the library:
This vignette shows how to:
babette
run ‘BEAST2’In all cases, this is done for a short MCMC chain length of 10K:
Also, in all cases, we use the same BEAST2 options:
babette
run ‘BEAST2’For an alignment, we’ll use a babette
example alignment.
if (is_beast2_installed()) {
out <- bbt_run_from_model(
fasta_filename = get_babette_path("anthus_aco_sub.fas"),
inference_model = inference_model,
beast2_options = beast2_options
)
bbt_delete_temp_files(
inference_model = inference_model,
beast2_options = beast2_options
)
}
#> |parameter |value |
#> |:---------------------|:-----------------------------------------------------------|
#> |input_filename |/home/richel/.cache/beastier/beast2_3dfaf36441195.xml |
#> |output_state_filename |/home/richel/.cache/beastier/beast2_3dfaf211b1690.xml.state |
#> |rng_seed |NA |
#> |n_threads |NA |
#> |use_beagle |FALSE |
#> |overwrite |TRUE |
#> |beast2_path |/home/richel/.local/share/beast/lib/launcher.jar |
#> |verbose |TRUE |
#> Running command: '/usr/lib/jvm/java-11-openjdk-amd64/bin/java -cp /home/richel/.local/share/beast/lib/launcher.jar beast.app.beastapp.BeastLauncher -validate /home/richel/.cache/beastier/beast2_3dfaf36441195.xml'
#> BEAST v2.6.6, 2002-2021 Bayesian Evolutionary Analysis Sampling Trees Designed and developed by Remco Bouckaert, Alexei J. Drummond, Andrew Rambaut & Marc A. Suchard Centre for Computational Evolution University of Auckland r.bouckaert@auckland.ac.nz alexei@cs.auckland.ac.nz Institute of Evolutionary Biology University of Edinburgh a.rambaut@ed.ac.uk David Geffen School of Medicine University of California, Los Angeles msuchard@ucla.edu Downloads, Help & Resources: http://beast2.org/ Source code distributed under the GNU Lesser General Public License: http://github.com/CompEvol/beast2 BEAST developers: Alex Alekseyenko, Trevor Bedford, Erik Bloomquist, Joseph Heled, Sebastian Hoehna, Denise Kuehnert, Philippe Lemey, Wai Lok Sibon Li, Gerton Lunter, Sidney Markowitz, Vladimir Minin, Michael Defoin Platel, Oliver Pybus, Tim Vaughan, Chieh-Hsi Wu, Walter Xie Thanks to: Roald Forsberg, Beth Shapiro and Korbinian StrimmerRandom number seed: 1660292517403File: beast2_3dfaf36441195.xml seed: 1660292517403 threads: 1Loading package NS v1.1.0Loading package BEAST v2.6.6Loading package BEASTLabs v1.9.7Loading package MODEL_SELECTION v1.5.3Loading package BEAST v2.6.661430_aco: 78 4626029_aco: 78 4630116_aco: 78 4630210_aco: 78 4B25702_aco: 78 4Alignment(anthus_aco_sub) 5 taxa 78 sites 9 patternsFailed to load BEAGLE library: no hmsbeagle-jni in java.library.path: [/usr/lib/jvm/java-11-openjdk-amd64/lib/server, /usr/lib/jvm/java-11-openjdk-amd64/lib, /usr/lib/jvm/java-11-openjdk-amd64/../lib, /usr/lib/R/lib, /usr/lib/x86_64-linux-gnu, /usr/lib/jvm/default-java/lib/server, /usr/lib/R/lib, /usr/lib/x86_64-linux-gnu, /usr/lib/jvm/default-java/lib/server, /usr/lib/R/lib, /usr/lib/x86_64-linux-gnu, /usr/lib/jvm/default-java/lib/server, /usr/java/packages/lib, /usr/lib/x86_64-linux-gnu/jni, /lib/x86_64-linux-gnu, /usr/lib/x86_64-linux-gnu, /usr/lib/jni, /lib, /usr/lib]TreeLikelihood(treeLikelihood.anthus_aco_sub0) uses BeerLikelihoodCore4 Alignment(anthus_aco_sub): [taxa, patterns, sites] = [5, 9, 78]===============================================================================Citations for this model:Bouckaert, Remco, Timothy G. Vaughan, Joëlle Barido-Sottani, Sebastián Duchêne, Mathieu Fourment, Alexandra Gavryushkina, Joseph Heled, Graham Jones, Denise Kühnert, Nicola De Maio, Michael Matschiner, Fábio K. Mendes, Nicola F. Müller, Huw A. Ogilvie, Louis du Plessis, Alex Popinga, Andrew Rambaut, David Rasmussen, Igor Siveroni, Marc A. Suchard, Chieh-Hsi Wu, Dong Xie, Chi Zhang, Tanja Stadler, Alexei J. Drummond BEAST 2.5: An advanced software platform for Bayesian evolutionary analysis. PLoS computational biology 15, no. 4 (2019): e1006650.===============================================================================Done!
#> cmd: /usr/lib/jvm/java-11-openjdk-amd64/bin/java -cp /home/richel/.local/share/beast/lib/launcher.jar beast.app.beastapp.BeastLauncher -statefile /home/richel/.cache/beastier/beast2_3dfaf211b1690.xml.state -overwrite /home/richel/.cache/beastier/beast2_3dfaf36441195.xml
#> Creating the beautier temprary folder '/home/richel/.cache/beautier'for BEAST2 tracelog, treelog and screenlog files
#> Creating folder '/home/richel/.cache/beastier'for BEAST2 .xml.state output file
#> error_code: 0
#> stdout:
#>
#> BEAST v2.6.6, 2002-2021
#> Bayesian Evolutionary Analysis Sampling Trees
#> Designed and developed by
#> Remco Bouckaert, Alexei J. Drummond, Andrew Rambaut & Marc A. Suchard
#>
#> Centre for Computational Evolution
#> University of Auckland
#> r.bouckaert@auckland.ac.nz
#> alexei@cs.auckland.ac.nz
#>
#> Institute of Evolutionary Biology
#> University of Edinburgh
#> a.rambaut@ed.ac.uk
#>
#> David Geffen School of Medicine
#> University of California, Los Angeles
#> msuchard@ucla.edu
#>
#> Downloads, Help & Resources:
#> http://beast2.org/
#>
#> Source code distributed under the GNU Lesser General Public License:
#> http://github.com/CompEvol/beast2
#>
#> BEAST developers:
#> Alex Alekseyenko, Trevor Bedford, Erik Bloomquist, Joseph Heled,
#> Sebastian Hoehna, Denise Kuehnert, Philippe Lemey, Wai Lok Sibon Li,
#> Gerton Lunter, Sidney Markowitz, Vladimir Minin, Michael Defoin Platel,
#> Oliver Pybus, Tim Vaughan, Chieh-Hsi Wu, Walter Xie
#>
#> Thanks to:
#> Roald Forsberg, Beth Shapiro and Korbinian Strimmer
#>
#> Writing state to file /home/richel/.cache/beastier/beast2_3dfaf211b1690.xml.state
#> Random number seed: 1660292518660
#>
#> 61430_aco: 78 4
#> 626029_aco: 78 4
#> 630116_aco: 78 4
#> 630210_aco: 78 4
#> B25702_aco: 78 4
#> Alignment(anthus_aco_sub)
#> 5 taxa
#> 78 sites
#> 9 patterns
#>
#> TreeLikelihood(treeLikelihood.anthus_aco_sub0) uses BeerLikelihoodCore4
#> Alignment(anthus_aco_sub): [taxa, patterns, sites] = [5, 9, 78]
#> ===============================================================================
#> Citations for this model:
#>
#> Bouckaert, Remco, Timothy G. Vaughan, Joëlle Barido-Sottani, Sebastián Duchêne, Mathieu Fourment,
#> Alexandra Gavryushkina, Joseph Heled, Graham Jones, Denise Kühnert, Nicola De Maio, Michael Matschiner,
#> Fábio K. Mendes, Nicola F. Müller, Huw A. Ogilvie, Louis du Plessis, Alex Popinga, Andrew Rambaut,
#> David Rasmussen, Igor Siveroni, Marc A. Suchard, Chieh-Hsi Wu, Dong Xie, Chi Zhang, Tanja Stadler,
#> Alexei J. Drummond
#> BEAST 2.5: An advanced software platform for Bayesian evolutionary analysis.
#> PLoS computational biology 15, no. 4 (2019): e1006650.
#>
#> ===============================================================================
#> Start likelihood: -412.0608952631166
#>
#> Operator Tuning #accept #reject Pr(m) Pr(acc|m)
#> ScaleOperator(YuleBirthRateScaler.t:anthus_aco_sub) 0.75000 110 9 0.04000 0.92437 Try setting scaleFactor to about 0.562
#> ScaleOperator(YuleModelTreeScaler.t:anthus_aco_sub) 0.50000 75 56 0.04000 0.57252 Try setting scaleFactor to about 0.25
#> ScaleOperator(YuleModelTreeRootScaler.t:anthus_aco_sub) 0.50000 51 65 0.04000 0.43966 Try setting scaleFactor to about 0.272
#> Uniform(YuleModelUniformOperator.t:anthus_aco_sub) - 643 563 0.40000 0.53317
#> SubtreeSlide(YuleModelSubtreeSlide.t:anthus_aco_sub) 1.00000 8 594 0.20000 0.01329 Try decreasing size to about 0.5
#> Exchange(YuleModelNarrow.t:anthus_aco_sub) - 220 369 0.20000 0.37351
#> Exchange(YuleModelWide.t:anthus_aco_sub) - 10 124 0.04000 0.07463
#> WilsonBalding(YuleModelWilsonBalding.t:anthus_aco_sub) - 8 96 0.04000 0.07692
#>
#> Tuning: The value of the operator's tuning parameter, or '-' if the operator can't be optimized.
#> #accept: The total number of times a proposal by this operator has been accepted.
#> #reject: The total number of times a proposal by this operator has been rejected.
#> Pr(m): The probability this operator is chosen in a step of the MCMC (i.e. the normalized weight).
#> Pr(acc|m): The acceptance probability (#accept as a fraction of the total proposals for this operator).
#>
#>
#> Total calculation time: 0.495 seconds
#> stderr:
#> File: beast2_3dfaf36441195.xml seed: 1660292518660 threads: 1
#> Loading package NS v1.1.0
#> Loading package BEAST v2.6.6
#> Loading package BEASTLabs v1.9.7
#> Loading package MODEL_SELECTION v1.5.3
#> Loading package BEAST v2.6.6
#> Failed to load BEAGLE library: no hmsbeagle-jni in java.library.path: [/usr/lib/jvm/java-11-openjdk-amd64/lib/server, /usr/lib/jvm/java-11-openjdk-amd64/lib, /usr/lib/jvm/java-11-openjdk-amd64/../lib, /usr/lib/R/lib, /usr/lib/x86_64-linux-gnu, /usr/lib/jvm/default-java/lib/server, /usr/lib/R/lib, /usr/lib/x86_64-linux-gnu, /usr/lib/jvm/default-java/lib/server, /usr/lib/R/lib, /usr/lib/x86_64-linux-gnu, /usr/lib/jvm/default-java/lib/server, /usr/java/packages/lib, /usr/lib/x86_64-linux-gnu/jni, /lib/x86_64-linux-gnu, /usr/lib/x86_64-linux-gnu, /usr/lib/jni, /lib, /usr/lib]
#> WARNING: If nothing seems to be happening on screen this is because none of the loggers give feedback to screen.
#> WARNING: This happens when a filename is specified for the 'screenlog' logger.
#> WARNING: To get feedback to screen, leave the filename for screenlog blank.
#> WARNING: Otherwise, the screenlog is saved into the specified file.
#> Writing file /home/richel/.cache/beautier/tracelog_3dfaf72f775a7.log
#> Writing file /home/richel/.cache/beautier/screenlog_3dfafd2a8f52.csv
#> Writing file /home/richel/.cache/beautier/treelog_3dfaf6d352e33.trees
#> End likelihood: -126.77270757772963
#> [1] TRUE
Effective sample sizes, with 20% burn-in removed:
if (is_beast2_installed()) {
traces <- remove_burn_ins(
traces = out$estimates,
burn_in_fraction = 0.2
)
esses <- t(
calc_esses(
traces,
sample_interval = inference_model$mcmc$tracelog$log_every
)
)
colnames(esses) <- "ESS"
knitr::kable(esses)
}
ESS | |
---|---|
posterior | 4 |
likelihood | 4 |
prior | 4 |
treeLikelihood | 4 |
TreeHeight | 4 |
YuleModel | 4 |
birthRate | 4 |
For a reliable inference, use an ESS of at least 200.
if (is_beast2_installed()) {
sum_stats <- t(
calc_summary_stats(
traces$posterior,
sample_interval = inference_model$mcmc$tracelog$log_every
)
)
colnames(sum_stats) <- "Statistic"
knitr::kable(sum_stats)
}
Statistic | |
---|---|
mean | -199.683 |
stderr_mean | 61.32501 |
stdev | 141.624 |
variance | 20057.37 |
median | -130.172 |
mode | n/a |
geom_mean | n/a |
hpd_interval_low | -412.0609 |
hpd_interval_high | -126.3272 |
act | 1000 |
ess | 4 |