Bayesian Phylogenetics with BEAST for Molecular Clock Dating

Introduction to Bayesian Phylogenetic Inference

Bayesian phylogenetic inference provides a probabilistic framework for reconstructing evolutionary relationships and estimating divergence times from molecular sequence data. This approach integrates prior knowledge about evolutionary processes with observed sequence data to generate posterior probability distributions over tree topologies, branch lengths, and substitution model parameters [1, 2]. The software package BEAST (Bayesian Evolutionary Analysis Sampling Trees) implements Markov chain Monte Carlo (MCMC) sampling algorithms that enable simultaneous estimation of phylogeny and evolutionary timescales from nucleotide or amino acid sequence alignments [1, 3].

The fundamental principle underlying molecular clock dating is the observation that genetic substitutions accumulate at a relatively constant rate over time in many lineages [4, 2]. However, empirical evidence demonstrates that substitution rates can vary substantially among lineages, necessitating relaxed clock models that accommodate rate heterogeneity [3, 5, 6]. Bayesian approaches are particularly well suited to this challenge because they can incorporate flexible prior distributions on rate variation while jointly estimating tree topology and divergence times [7, 8].

Molecular Clock Models in BEAST

Strict Clock Model

The strict molecular clock assumes a single constant substitution rate across all branches in the phylogeny [2, 6]. This model is appropriate for datasets where rate variation among lineages is minimal, such as closely related populations or rapidly evolving viruses sampled over short time periods [9, 10]. The strict clock is parameterized by a single rate parameter mu, which represents the expected number of substitutions per site per unit time [2]. While computationally efficient, the strict clock is frequently rejected by likelihood ratio tests or Bayes factors when applied to datasets spanning diverse taxonomic groups or long evolutionary timescales [2, 6].

Relaxed Clock Models

Relaxed clock models relax the assumption of rate constancy by allowing substitution rates to vary across branches according to a specified probability distribution [3, 5, 6]. BEAST implements two primary classes of relaxed clocks: uncorrelated and autocorrelated models.

The uncorrelated relaxed clock assumes that rates on different branches are drawn independently from a common distribution, typically a lognormal or exponential distribution [3, 6]. This model is parameterized by a mean rate and a standard deviation that controls the degree of rate variation among branches [6]. The uncorrelated lognormal (UCLN) model is widely used because it accommodates substantial rate heterogeneity without assuming any correlation between adjacent branches [3, 6].

The autocorrelated relaxed clock assumes that substitution rates evolve gradually along the tree, such that rates on adjacent branches are more similar than rates on distant branches [5, 6]. This model is biologically plausible for datasets where rate variation is driven by gradual changes in generation time or metabolic rate [5]. The random local clock model represents an intermediate approach that allows discrete rate shifts at specific nodes in the phylogeny, with constant rates maintained within each rate regime [5].

Model Selection and Comparison

Selection among competing clock models is performed using marginal likelihood estimation and Bayes factor comparison [1, 3]. Path sampling and stepping-stone sampling are recommended for accurate marginal likelihood estimation in BEAST, as these methods are less sensitive to the choice of prior distributions than harmonic mean estimators [1]. Simulation studies have demonstrated that model selection accuracy depends on the degree of rate variation present in the data and the sampling density of taxa [2, 3].

Prior Calibration and Fossil Constraints

Calibration Prior Specification

Molecular clock dating requires external calibration information to convert relative substitution rates into absolute time estimates [7, 11, 8]. Calibration priors are typically derived from fossil evidence, biogeographic events, or known sampling dates for heterochronous sequences [7, 12]. In veterinary virology, calibration may be based on the earliest documented isolation dates of a pathogen, historical records of disease emergence, or radiocarbon dating of ancient DNA samples from archaeological specimens [11, 12].

Calibration priors are specified as probability distributions on the ages of specific nodes in the phylogeny [7, 8]. The choice of prior distribution shape (e.g., uniform, lognormal, gamma) reflects the degree of uncertainty in the calibration estimate [7]. Lognormal priors are commonly used for fossil calibrations because they accommodate asymmetric uncertainty, with a hard minimum bound and a soft maximum bound [7, 8].

Temporal Sampling Schemes

For rapidly evolving pathogens such as RNA viruses, temporal calibration can be achieved through heterochronous sampling, where sequences are collected at different time points [9, 12, 10]. The tip-dating approach directly incorporates sampling dates as calibrations, allowing estimation of substitution rates and divergence times without external fossil constraints [9, 12]. The precision of tip-dating estimates depends on the temporal span of sampling, the number of sequences per time point, and the evolutionary rate of the pathogen [9, 12].

Studies have shown that increasing the temporal range of sampling improves the accuracy of rate estimates, while dense sampling within a narrow time window may lead to biased estimates if rate variation is present [9, 12]. The empirical calibrated radiocarbon sampler provides a method for incorporating radiocarbon date uncertainty into Bayesian phylogenetic analyses of ancient DNA, accounting for both calibration curve error and laboratory dating error [11].

MCMC Convergence and Diagnostic Assessment

Chain Mixing and Effective Sample Size

Bayesian phylogenetic inference with BEAST relies on MCMC sampling to approximate the posterior distribution [1, 3]. Adequate chain mixing is essential for obtaining reliable parameter estimates. The effective sample size (ESS) is a diagnostic metric that quantifies the number of independent samples in the MCMC chain after accounting for autocorrelation [1]. ESS values above 200 are generally considered acceptable for continuous parameters, while tree topology parameters may require higher ESS thresholds [1].

Convergence Diagnostics

Multiple independent MCMC chains should be run to assess convergence to the stationary distribution [1]. The potential scale reduction factor (PSRF) compares within-chain and between-chain variance to determine whether chains have converged to the same distribution [1]. Trace plots of log-likelihood values and key model parameters should be visually inspected for stationarity and the absence of trends [1].

Adaptive Proposals and Computational Efficiency

Recent developments in BEAST have introduced adaptive proposal mechanisms that automatically tune MCMC proposal distributions during the burn-in phase [3]. These adaptive proposals improve mixing efficiency for relaxed clock models by adjusting step sizes based on acceptance rates [3]. The adaptive dating approach specifically addresses the challenge of sampling branch rates in large phylogenies by proposing rate changes that are correlated with branch lengths [3].

Applications in Veterinary Virology

Estimating Emergence Times of Pathogens

Bayesian molecular clock dating has been applied extensively to estimate the emergence times of viral pathogens in animal populations [9, 12]. For example, analyses of avian influenza virus sequences have been used to infer the timing of introduction events into poultry populations and the subsequent diversification of lineages [9]. The accuracy of these estimates depends on the availability of temporally structured sequence data and appropriate model specification [9, 12].

Phylodynamic Inference

Phylodynamic approaches combine phylogenetic inference with epidemiological modeling to reconstruct population dynamics from genetic data [1]. BEAST implements coalescent-based population size estimation, allowing inference of effective population size trajectories through time [1]. These methods have been applied to study the population dynamics of rabies virus in canine populations, foot-and-mouth disease virus in livestock, and African swine fever virus in wild boar populations.

Host Range and Spillover Events

Molecular clock dating can be used to investigate the timing of cross-species transmission events and the evolutionary history of host adaptation [9]. By estimating the divergence times of viral lineages sampled from different host species, researchers can infer the temporal sequence of host range expansion and identify potential reservoir populations [9]. These analyses require careful consideration of sampling bias and rate variation among host lineages [2, 9].

Workflow for Bayesian Phylogenetic Analysis

The following Mermaid diagram illustrates a typical workflow for Bayesian phylogenetic analysis with BEAST:

flowchart TD
    A[Sequence Alignment], > B[Model Selection]
    B, > C[Prior Specification]
    C, > D[Calibration Priors]
    C, > E[Clock Model Priors]
    C, > F[Tree Prior]
    D, > G[BEAST XML Generation]
    E, > G
    F, > G
    G, > H[MCMC Sampling]
    H, > I[Convergence Assessment]
    I, > J{ESS > 200?}
    J, >|No| H
    J, >|Yes| K[Posterior Summary]
    K, > L[Maximum Clade Credibility Tree]
    K, > M[Parameter Estimates]
    K, > N[HPD Intervals]

Model Parameterization and Interpretation

Substitution Models

The choice of nucleotide substitution model affects both topological inference and divergence time estimates [1, 2]. BEAST supports a wide range of substitution models, including the general time-reversible (GTR) model with gamma-distributed rate heterogeneity among sites [1]. Model selection is typically performed using information criteria such as AIC or BIC, or through Bayesian model averaging [1].

Rate Priors and Their Influence

The prior distribution on the mean substitution rate can substantially influence posterior estimates, particularly when the data contain limited temporal information [2, 7]. For tip-dating analyses, a diffuse prior on the substitution rate is often appropriate, allowing the data to inform the rate estimate [9]. For fossil-calibrated analyses, the rate prior should be specified based on prior knowledge from related taxa or experimental measurements [7, 8].

Node Age Estimates and Credible Intervals

Posterior distributions of node ages are summarized using highest posterior density (HPD) intervals, which represent the shortest interval containing a specified probability mass (typically 95%) [1, 7]. The width of HPD intervals reflects the combined uncertainty from phylogenetic reconstruction, substitution rate estimation, and calibration prior specification [2, 7]. Factors that increase HPD width include sparse taxon sampling, limited temporal range, and conflict between the molecular data and calibration priors [2, 12].

Frequently Asked Questions

What is the difference between strict and relaxed molecular clocks?

The strict molecular clock assumes a single constant substitution rate across all branches in the phylogeny, while relaxed clock models allow rates to vary among lineages according to a specified probability distribution [2, 6]. Strict clocks are appropriate for closely related sequences with minimal rate variation, whereas relaxed clocks are necessary for datasets spanning diverse taxonomic groups or long evolutionary timescales [3, 5].

How are calibration priors specified in BEAST?

Calibration priors are specified as probability distributions on the ages of specific nodes in the phylogeny, typically derived from fossil evidence, biogeographic events, or known sampling dates [7, 8]. The shape of the prior distribution reflects the degree of uncertainty in the calibration estimate, with lognormal distributions commonly used for fossil calibrations [7].

What is the minimum temporal sampling range required for tip dating?

The minimum temporal sampling range required for accurate tip dating depends on the substitution rate of the organism and the number of sequences analyzed [9, 12]. Simulation studies suggest that a temporal range spanning at least 10% of the total evolutionary timescale is generally sufficient for reliable rate estimation, although longer ranges improve precision [9, 12].

How can I assess MCMC convergence in BEAST?

MCMC convergence is assessed using effective sample size (ESS) values, trace plots, and multiple independent chains [1]. ESS values above 200 for continuous parameters indicate adequate sampling, while trace plots should show stationarity without trends [1]. The potential scale reduction factor (PSRF) compares within-chain and between-chain variance to confirm convergence [1].

What factors influence the accuracy of divergence time estimates?

Accuracy of divergence time estimates is influenced by the quality and quantity of sequence data, the appropriateness of the clock model, the informativeness of calibration priors, and the degree of rate variation among lineages [2, 7, 12]. Sparse taxon sampling, limited temporal range, and model misspecification can all lead to biased estimates [2, 9].

References

[1] Baele G, Ji X, Hassler GW, et al. BEAST X for Bayesian phylogenetic, phylogeographic and phylodynamic inference. Nat Methods. 2025. https://pubmed.ncbi.nlm.nih.gov/40624354/

[2] Ritchie AM, Hua X, Bromham L. Investigating the reliability of molecular estimates of evolutionary time when substitution rates and speciation rates vary. BMC Ecol Evol. 2022. https://pubmed.ncbi.nlm.nih.gov/35538412/

[3] Douglas J, Zhang R, Bouckaert R. Adaptive dating and fast proposals: Revisiting the phylogenetic relaxed clock model. PLoS Comput Biol. 2021. https://pubmed.ncbi.nlm.nih.gov/33529184/

[4] Osozawa S. Spermatophyta Molecular Clock: Time Drift and Recent Acceleration. Plant Environ Interact. 2025. https://pubmed.ncbi.nlm.nih.gov/40978090/

[5] Drummond AJ, Suchard MA. Bayesian random local clocks, or one rate to rule them all. BMC Biol. 2010. https://pubmed.ncbi.nlm.nih.gov/20807414/

[6] Brown RP, Yang Z. Bayesian dating of shallow phylogenies with a relaxed clock. Syst Biol. 2010. https://pubmed.ncbi.nlm.nih.gov/20525625/

[7] Grimm GW, Kapli P, Bomfleur B, et al. Using more than the oldest fossils: dating osmundaceae with three Bayesian clock approaches. Syst Biol. 2015. https://pubmed.ncbi.nlm.nih.gov/25503771/

[8] Gustafsson AL, Verola CF, Antonelli A. Reassessing the temporal evolution of orchids with new fossils and a Bayesian relaxed clock, with implications for the diversification of the rare South American genus Hoffmannseggella (Orchidaceae: Epidendroideae). BMC Evol Biol. 2010. https://pubmed.ncbi.nlm.nih.gov/20546585/

[9] Dearlove B, Tovanabutra S, Owen CL, et al. Factors influencing estimates of HIV-1 infection timing using BEAST. PLoS Comput Biol. 2021. https://pubmed.ncbi.nlm.nih.gov/33524022/

[10] Yang Z, O'Brien JD, Zheng X, et al. Tree and rate estimation by local evaluation of heterochronous nucleotide data. Bioinformatics. 2007. https://pubmed.ncbi.nlm.nih.gov/17110369/ *** Disclaimer: This article is for educational and informational purposes only. It is not intended to substitute for professional veterinary advice, diagnosis, treatment, or regulatory guidance. Always consult a licensed veterinarian or qualified specialist regarding animal health, disease diagnosis, and therapeutic decisions.

[11] Molak M, Suchard MA, Ho SY, et al. Empirical calibrated radiocarbon sampler: a tool for incorporating radiocarbon-date and calibration error into Bayesian phylogenetic analyses of ancient DNA. Mol Ecol Resour. 2015. https://pubmed.ncbi.nlm.nih.gov/24964386/

[12] Molak M, Lorenzen ED, Shapiro B, et al. Phylogenetic estimation of timescales using ancient DNA: the effects of temporal sampling scheme and uncertainty in sample ages. Mol Biol Evol. 2013. https://pubmed.ncbi.nlm.nih.gov/23024187/