Review: Maximum Likelihood Phylogenetic Estimation from DNA Sequences with Variable Rates over Sites: Approximate Methods¶
Citation
- Yang, Z. (1994). Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: approximate methods. Journal of Molecular Evolution, 39(3), 306–314.
- DOI
Abstract¶
Among-site rate variation is modeled using a Gamma distribution. An approximate method based on discretizing the distribution into a small number of rate categories is developed and shown to be accurate and computationally efficient. The approach is incorporated into maximum likelihood estimation of phylogenetic trees from DNA sequences.
The substitution models described in Tavaré (1986) and Felsenstein (1981) assume that every site in the alignment evolves at the same rate. In reality, many sites in a protein-coding gene or ribosomal RNA gene evolve at very different rates: some positions are under strong purifying selection and rarely change; others are nearly neutral and change rapidly. Ignoring this variation leads to biased branch length estimates and can cause well-known artefacts such as long-branch attraction.
Yang (1994) introduced the practical solution now used in virtually all modern phylogenetic software: model among-site rate variation with a Gamma distribution, then approximate it with a small number of discrete rate categories. The Gamma distribution \( \Gamma(\alpha, \beta) \) with shape parameter \( \alpha \) controls the degree of heterogeneity: small \( \alpha \) (near zero) corresponds to extreme rate variation (most sites invariant, a few evolving very fast); large \( \alpha \) corresponds to near-uniform rates. Setting \( \beta = \alpha \) fixes the mean rate to 1, preserving the branch-length scale.
The discrete approximation divides the Gamma distribution into \( K \) equiprobable categories (typically \( K = 4 \)) and assigns each category the mean rate of its portion of the distribution, computed from the incomplete Gamma function. The site likelihood is then the weighted average over the \( K \) rate categories:
where \( r_k \) are the discrete rates and \( L_p(r_k) \) is the likelihood of site \( p \) under rate \( r_k \). The cost is \( K \) likelihood evaluations per site instead of one, typically increasing runtime by a factor of \( K \).
Hifuku supports Gamma rate heterogeneity as an option for substitution model parameterization. The discrete categories and rates are computed following Yang (1994) and are passed to the pruning kernel as an additional rate multiplier on each branch. This allows the same GPU kernel to handle both uniform-rate and Gamma-rate models without code duplication.
Key result¶
For \( K = 4 \) discrete categories, Yang (1994) shows (Table 2) that the approximation is accurate to within 0.1 log-likelihood units compared to the continuous Gamma, while reducing the required numerical integration from continuous to a 4-point sum. This accuracy justifies the universal adoption of \( K = 4 \) as the default in IQ-TREE 2, RAxML, MrBayes, and Hifuku.