Source-linked AI summary
Fast coalescent-based computation of local branch support from quartet frequencies
Erfan Sayyari, Siavash Mirarab
TL;DR
Quantifying support for species-tree branches remains underdeveloped, despite risks of high support for incorrect branches. This paper derives fast quartet-based local posterior probabilities and branch lengths, showing reliable support with few incorrect highly supported branches and high support for most correct branches given enough genes.
Problem
Despite progress in species-tree reconstruction, methods for calculating branch support remain limited, and supermatrix analyses can highly support incorrect branches.
Method
The method uses four-taxon multispecies-coalescent properties and quartet frequencies to compute branch support, local posterior probabilities, and coalescent-unit lengths.
Results
Very few highly supported branches are incorrect, and with sufficient genes most correct branches receive high support.
Takeaways & Limitations
The approach provides a fast, reliable measure for assessing species-tree branch accuracy without requiring extensive computation.
Takeaways & Limitations
Experiments violated some assumptions, and unmodeled discordance from paralogy or horizontal gene transfer remains for future examination.
Abstract
from arXiv · showhide
Species tree reconstruction is complicated by effects of Incomplete Lineage Sorting (ILS), commonly modeled by the multi-species coalescent model. While there has been substantial progress in developing methods that estimate a species tree given a collection of gene trees, less attention has been paid to fast and accurate methods of quantifying support. In this paper, we propose a fast algorithm to compute quartet-based support for each branch of a given species tree with regard to a given set of gene trees. We then show how the quartet support can be used in the context of the multi-species coalescent model to compute i) the local posterior probability that the branch is in the species tree and ii) the length of the branch in coalescent units. We evaluate the precision and recall of the local posterior probability on a wide set of simulated and biological data, and show that it has very high precision and improved recall compared to multi-locus bootstrapping. The estimated branch lengths are highly accurate when gene trees have little error, but are underestimated when gene tree estimation error increases. Computation of both branch length and local posterior probability is implemented as a new feature in ASTRAL.
Introduction
The paper addresses limited attention to support estimation for MSC-based species trees by introducing fast quartet-frequency methods for local posterior probabilities and coalescent branch lengths. The approach is evaluated as a reliable support measure and implemented in ASTRAL.
- Motivation: MSC-based species-tree methods have advanced, but calculating support remains comparatively underdeveloped.Concatenation can produce high support for incorrect branches, motivating improved support quantification.
- Contribution: Quartet-based support values are more precise, more powerful, and much faster to compute than MLBS under the MSC.Quartet gene trees lack anomaly zones, so the most probable gene-tree topology matches the species-tree topology for any quartet.
- Method: The method computes quartet-based support for each species-tree branch from unrooted gene trees and derives Bayesian local posterior probabilities for the three alternatives.For each internal branch, the four sides are assumed correct, leaving three possible topologies around that branch.
- Method: The fast algorithm computes quartet support for all three alternatives in Θ(nl) time, where l is the number of species and n is the number of genes.Quartet frequencies also provide estimates of internal branch lengths in coalescent units.
- Results: Simulated and empirical datasets show that local posterior probability is reliable: very few highly supported branches are incorrect, and sufficient genes give most correct branches high support.The method remains reliable when model assumptions are violated.
- Implementation: ASTRAL implements the method and now estimates species-tree topologies, branch lengths, and local posterior probabilities.The implementation is available through the ASTRAL project.
New Approaches · Calculation of local posterior probability
The paper defines quartet-based local support around a scored branch and computes its posterior probability under MSC-based assumptions. The method uses averaged quartet frequencies, a conservative dependence model, and an efficient Θ(nl) computation, with posterior values increasing with gene-tree count but also depending on alternative-topology frequencies.
- New Approaches: Each internal branch induces a quadripartition, and quartets around it are formed by selecting one leaf from each of four clusters.These quartets provide the local observations used to score the branch.
- New Approaches: The method scores a quadripartition by estimating both the probability that it belongs to the true species tree and its coalescent-unit branch length.The locality assumption treats all four surrounding clusters as correct, making computation tractable for large datasets.
- New Approaches: The theoretical setup assumes true gene trees differ from the species tree only through ILS under the MSC and come from an unbiased sample.Real data instead involve estimated gene trees, potentially biased samples, and possible violations of the MSC or locality assumptions.
- Calculation of local posterior probability: For each quartet around a branch, the MSC predicts one dominant topology and two equally probable alternatives, with quartet frequencies modeled using multinomial distributions.The dominant-topology probability is θ=1−2/3 for a quartet whose intervening branch length is d in coalescent units.
- Calculation of local posterior probability: Because quartets around a branch are dependent, treating them as independent would inflate confidence by replacing n observations with m×n observations.The method instead models a shared hidden vector of true quartet frequencies and treats each quartet’s frequencies as noisy estimates.
- Calculation of local posterior probability: An efficient Θ(nl) algorithm computes the average quartet-frequency vector used to estimate local posterior probability.The posterior calculation assumes the scored quadripartition’s four clusters are each present in the species tree.
- Calculation of local posterior probability: A branch appearing in 40% of gene trees has 66.1% pp with 50 genes, 93.0% pp with 200 genes, and 99.7% pp with 500 genes.Thus, a branch with 60% gene-tree discordance can still receive high confidence given enough genes.
- Calculation of local posterior probability: When a branch appears in 40% of gene trees but one alternative appears in 45%, its local posterior probability is only 1.90%.Posterior probability therefore depends on the frequencies of both alternative topologies, not only the scored topology.
Calculation of branch length
The section frames branch-length calculation in coalescent units through a true parameter θ and a theorem conditioned on the branch being present in the species tree.
- Branch length is expressed in coalescent units using the true parameter θ for a correct branch.
- Theorem 2 analyzes the calculation under Theorem 1’s conditions, assuming the branch represented by Q is in the species tree.
Calculation of quartet support
The method avoids enumerating all quartet counts around each branch by mapping quadripartitions to gene-tree tripartitions and evaluating 12 cluster matchings. A postorder traversal reduces per-branch computation to Θ(nl), yielding Θ(l2n) time for all branch lengths and local posterior probabilities.
- Computational challenge: Naively enumerating quartet counts can be very slow because a single branch may have Θ(l4) surrounding quartets.The worst-case number of quartets around one branch is Θ(l4).
- Cluster-matching computation: The algorithm maps pairs of sister clusters in a species-tree quadripartition to gene-tree tripartition clusters and evaluates all 12 possible matchings.For each matching, it multiplies sizes of corresponding cluster intersections and sums the results to count shared quartet topologies.
- Efficient implementation: Although direct cluster-intersection computation requires O(l2n) time, a postorder traversal computes the same sum in Θ(nl).The traversal avoids separately scoring each gene-tree tripartition.
- Overall complexity: Θ(l2n) time is required to compute branch lengths and local posterior probabilities for an entire tree.This follows from scoring l−3 branches after computing the per-branch quantity efficiently.
Other considerations · Materials and Methods · Datasets
The methods are implemented in ASTRAL with special handling for missing data and unresolved gene trees. Evaluation uses simulated and biological datasets spanning varied species counts, ILS levels, gene-tree error, and published empirical collections.
- Other considerations: ASTRAL implements the methods using the Colt package for numerical computations.Colt is cited as the numerical-computation package used in the implementation.
- Other considerations: Missing data require setting m to the average number of quartets present in gene trees rather than the number of quartets around Q.This adjustment is used when computing (7).
- Other considerations: For each branch, n equals the number of genes containing at least one quartet around Q, while unresolved polytomies are scored across all surrounding tripartitions.The branch-specific n accounts for genes missing every quartet around Q; polytomies of degree d require scoring all associated tripartitions.
- Datasets: The evaluation uses both simulated and biological datasets.The datasets support testing across controlled simulations and published empirical data.
- Simulated data: 201 taxa and three ILS levels characterize the A-200 simulations, with tree lengths of 500k, 2M, and 10M generations.The simulations test accuracy under heterogeneous conditions with many species.
- Simulated data: ASTRAL species-tree error ranges from 3% to 15% across the simulated datasets, depending on dataset and gene-tree error.A-200 errors range between average 3% and 15% across gene-count settings, while avian errors range from 5% to 15%.
- Simulated data: 48 taxa, 20 replicates, and 1000 genes define each avian simulation condition, varying gene length from 250bp to 1500bp.The four site-count conditions create varying levels of gene-tree estimation error, with 200 bootstrap replicates.
- Biological dataset:: Four published biological datasets contain 103, 46, 48, and 201 taxa, with 424, 310, 2022, and 256 genes, respectively.They cover plant, angiosperm, and avian datasets, including the 48-taxon binned supergene tree collection.
Evaluation procedure
The evaluation tests accuracy under model assumptions, robustness to violations, and comparison with site-only MLBS. It scores true and estimated gene and species trees using branch-level support and length metrics.
- Evaluation questions: The evaluation asks how accurate branch lengths and support values are under the model, how violations affect results, and how local posteriors compare with site-only MLBS.These questions are studied using true and estimated gene trees and species trees.
- Scoring design: For every internal species-tree branch, the procedure scores the focal topology and its two alternatives.For AB|CD, the alternatives are AD|BC and AC|BD.
- Support metrics: Support evaluation varies thresholds and counts true and false branches using precision, recall, false positive rate, and ROC curves.Precision is the percentage of supported branches that are correct; recall is the percentage of true branches meeting the threshold; FPR is the percentage of false branches meeting it.
- MLBS comparison: Because MLBS and posterior values are not directly comparable, their performance is compared with ROC curves rather than matched thresholds.MLBS comparison is feasible only for the avian dataset because bootstrapping was computationally infeasible for A-200.
- Branch-length evaluation: Branch-length accuracy compares true and estimated lengths on correct branches using logarithmic error and root mean squared error.On Avian datasets, the procedure compares MP-EST and ASTRAL errors.
Results · A-200 dataset · Avian
On the A-200 dataset, local posterior probability achieves very high precision, while recall decreases with greater ILS, fewer genes, and gene-tree error. Branch-length accuracy is strong with true gene trees but worsens and becomes biased toward underestimation with estimated gene trees; avian analyses compare this support with site-only MLBS under varying gene-tree estimation error.
- Results · A-200 dataset: 100% precision is achieved for branches with local posterior probability ≥0.99 across all model conditions, and precision is at least 99.8% at the 0.95 threshold.Lowering the threshold to 0.9 preserves high precision.
- A-200 dataset: At the 0.95 threshold, recall ranges from 98.7% with low ILS and 1000 genes to 55.0% with high ILS and 50 genes.Increasing ILS and reducing gene number lowers recall while maintaining high precision.
- A-200 dataset: Gene-tree estimation error reduces precision to 98.4–99.8% and can reduce recall by as much as 10%.Thus, estimation error has a small effect on precision but a substantial effect on recall.
- A-200 dataset: Changing the threshold from 0.99 to 0.95 improves recall, for example from 39% to 48% under high ILS with 50 genes, while reducing precision by at most 0.3%.ROC curves examine the tradeoff between recall and false-positive branches.
- Branch Length: With 1000 true gene trees, logarithmic branch-length error ranges from 0.03 to 0.10, corresponding to branches 7% or 25% shorter or longer than true branches.With true gene trees, error remains below 0.25 as gene number decreases and is largely unbiased except for extreme branch lengths.
- Branch Length: With 1000 estimated gene trees and low ILS, logarithmic branch-length error is 0.42, corresponding to branches estimated on average 2.6 times too short or long.Estimated gene trees substantially increase error, with a bias toward underestimation, especially under low ILS.
- Avian: The avian analysis compares local posterior probabilities with site-only MLBS using estimated gene trees and ASTRAL species trees while reducing gene length from 1500bp to assess estimation error.This experiment studies how increasing gene-tree estimation error affects support comparisons.
Posterior and MLBS · Biological datasets
Local posterior probabilities provide highly precise branch support with better recall than MLBS, while branch-length estimates are accurate with true gene trees but increasingly underestimated as gene-tree error grows. Biological datasets show strong support for some relationships, substantial uncertainty and short branches in difficult radiations, and effects of collapsing poorly supported gene-tree branches.
- Branch Length: 0.06 log error with true gene trees corresponded to branches about 14% shorter or longer than the true branch.Branch-length accuracy declined as gene-tree estimation error increased: 1500bp genes had 0.17 log error, whereas 250bp genes had 0.59 error.
- Branch Length: Estimated gene trees produced branch lengths biased toward underestimation, increasingly so for shorter alignments; ASTRAL gained an advantage over MP-EST as gene-tree error increased.ASTRAL and MP-EST had similar log error for highly accurate gene trees, while RMSE was comparable and slightly favored MP-EST with accurate gene trees.
- Biological datasets: In the angiosperm dataset, Amborella as sister to the remaining angiosperms had 75% BS but 100% posterior and collapsed posterior, with estimated branch length 0.160.The estimated length closely matched the 1KP estimate of 0.156.
- Biological datasets: In the avian genomic dataset, 12 branches at the base of Neoaves were below 0.1 coalescent units, and collapsing low-support gene-tree branches increased lengths by a median of 0.23 units.Eleven branches remained below 0.1 after collapsing, and uncertainty persisted despite 2022 input trees: Hoatzin-adjacent branches had collapsed posterior below 50% and posterior below 75%.
Discussions
The method computes local posterior probabilities rapidly without extensive MCMC sampling or bootstrapping, but relies on assumptions that can be violated by real-data gene-tree and species-tree error. Despite conservative support estimates, it maintains very high precision and generally high recall, with commonly used thresholds providing strong precision.
- Computational efficiency: Posterior computation takes only 10 seconds for 200 taxa and 1000 genes, and about three minutes for 1000 taxa and 1000 genes on a laptop.The speed depends on avoiding extensive MCMC sampling or bootstrapping.
- Assumptions and limitations: The method assumes true MSC-generated gene trees and locality around each internal branch, but both assumptions can be violated on real data.Simulations therefore included gene-tree estimation error from average RF distances of 25% to 67% and species-tree error.
- Support reliability: The method produces very few false positives with high support, even under high gene-tree estimation error, although it may underestimate support.This conservatism is presented as less serious than abundant false branches receiving high support.
- Support reliability: Local posterior probabilities are less conservative than MLBS and generally retain high recall except with very few genes for a given amount of ILS.Increasing gene-tree estimation error negatively affected recall but retained very high precision.
- Support thresholds: The commonly used 0.95 threshold yields very high precision with moderately high recall, while thresholds of 0.9 or even 0.7 also give high precision.The appropriate threshold depends on the desired false-positive rate and the amount of gene-tree error.
Gene tree estimation error · High support despite high discordance
Gene tree estimation error affects quartet-support precision and especially branch-length accuracy, while the method remains precise even with high estimation error. High confidence can nevertheless support correctly inferred relationships under substantial discordance, provided enough genes are available, with requirements increasing under higher ILS.
- Gene tree estimation error: Estimated gene trees yield more precise support under high ILS than low ILS at a given threshold, unlike true gene trees.The authors relate this pattern to gene-tree estimation error contributing more to observed discordance under low ILS.
- Gene tree estimation error: Branch-length accuracy degrades more from true to estimated gene trees under low ILS than high ILS, and estimation error causes underestimation.Branch lengths showed encouraging consistency across biological datasets but should be interpreted cautiously when gene trees are estimated.
- Gene tree estimation error: Bootstrapped gene trees show higher discordance than true or ML gene trees, giving MLBS a biased, noisier input that may explain underperformance.The bias arises because site-resampled bootstrap trees increase discordance with the species tree and among themselves, especially for genes with limited phylogenetic signal.
- Gene tree estimation error: The method avoids bootstrapping and uses available gene trees, remaining precise despite high gene-tree estimation error, although it does not directly model uncertainty.Collapsing poorly resolved gene-tree branches may help, but aggressive filtering can introduce bias; the authors collapse branches below 33% support.
- High support despite high discordance: A branch appearing in only 40% of gene trees can still be resolved with high confidence when enough genes are available, such as around 500.This illustrates that high support can arise despite high discordance for a correctly inferred relationship.
- High support despite high discordance: In the 1KP tree, branches appearing in 49% and 50% of collapsed gene-tree quartets both received posterior probabilities of 1.0.The branches placed Zygnematales sister to land plants and made Bryophytes monophyletic, respectively.
- High support despite high discordance: The number of genes needed for high resolution and low false-positive rates depends on model conditions; higher ILS requires more genes.The method can be extended to estimate how many genes may be required to resolve a tree.
Limitations and future work … 2 Supplementary methods
The paper identifies gene-tree uncertainty, prior specification, model violations, and unmodeled biological discordance as limitations, while outlining future methodological comparisons and extensions. Supplementary material documents dataset properties, support and branch-length analyses, and algorithms for quartet-frequency calculation.
- Limitations and future work: Gene-tree uncertainty could be incorporated by weighting quartets using SH-like support, BS, or concordance measures.These approaches are proposed as promising directions for improving local posterior probabilities.
- Limitations and future work: Future work should compare the method with Bayesian co-estimation methods on small datasets where those methods can run.
- Limitations and future work: The study did not examine changing λ prior parameters or alternative priors such as Dirichlet distributions and birth-death processes.These prior choices are left for future investigation.
- Limitations and future work: Experiments violated some assumptions because real-data model violations can bias gene-tree estimates, and simulations included discordance only from ILS and estimation error.Undetected paralogy and horizontal gene transfer were not simulated.
- Limitations and future work: Future work should assess how other biological sources of discordance affect local posterior probability reliability.The paper specifically identifies undetected paralogy and horizontal gene transfer as examples.
- 1 Supplementary Figures and Tables: The supplementary figures and tables report support accuracy, dataset discordance, local-posterior precision and recall, and branch-length accuracy across simulated and biological datasets.They include analyses using true, NJST, and concatenation species trees, as well as MLBS and local posterior support.
- 1 Supplementary Figures and Tables: The supplementary material presents support comparisons for angiosperm and avian datasets and branch-length analyses involving ASTRAL and MPEST.Additional avian analyses vary sites per gene to control gene-tree estimation error.
- 2 Supplementary methods: The supplementary methods provide a frequency-calculation algorithm that takes gene trees and a quadripartition and returns frequencies for three quartet topologies.The implementation uses bitset-indexed species clusters and includes an explicit formula for the calculation.
3 Proofs
The proofs derive a closed-form local posterior probability for a quartet topology under the Yule-process multispecies coalescent model and characterize maximum-likelihood and maximum-a-posteriori branch-length estimates. The branch-length estimates are zero below a quartet-frequency threshold, with the MAP derivation using the Yule-induced exponential prior.
- Local posterior probability: Theorem 1 expresses the local posterior probability of quartet topology Q as h(z1) divided by h(z1) + 2z2−z1h(z2) + 2z3−z1h(z3).The result assumes gene trees are generated by the MSC on a Yule-process species tree and z1 corresponds to Q's topology.
- Local posterior probability: The posterior derivation uses multinomial quartet counts and the absence of anomaly zones, so Q belongs to the species tree exactly when θ1 > 1 3.The Yule-process branch-length assumption enables the closed-form calculation from the lemma’s equations.
- Branch-length estimation: The ML branch-length estimate is −ln 3 2(1− z1 n+2λ) when 3z1 ≥n; otherwise, both ML and MAP estimates are zero.The ML objective is shown to be concave, and the estimate is obtained by maximizing its log likelihood over nonnegative branch lengths.
- Branch-length estimation: For MAP estimation, the Yule process supplies an exponential branch-length prior with rate 2λ, and the resulting posterior objective is also concave.The waiting time between consecutive speciation events has mean 1 2λ, yielding the exponential prior.
4 Commands and version numbers
The section specifies ASTRAL versions and commands for scoring species trees, computing branch lengths, and estimating local posterior probabilities. It also documents MP-EST settings for branch-length estimation on a fixed topology.
- ASTRAL: ASTRAL 4.9.1 was used for scoring with the command `java -Xmx2000M -jar astral.4.9.1.jar -i [GENE TREES] -q [SPECIES TREE] -t 4`.The scoring version was available from the posteval branch.
- ASTRAL: ASTRAL 4.9.8 computed main-species-tree branch lengths using the MAP estimate with `java -Xmx2000M -jar astral.4.9.8.jar -i [GENE TREES] -q [SPECIES TREE] -t 2`.This command was specifically used for branch-length computation.
- ASTRAL: ASTRAL 4.10.0 provides a command to score inferred species trees and compute branch lengths and local posterior probabilities: `java -Xmx2000M -jar astral.4.10.0 -i [GENE TREES] -q [SPECIES TREE] -t 3`.The section directs users to the updated code repository for this functionality.
- MP-EST: MP-EST 1.5 estimated branch lengths on a fixed topology by running twice with different random seeds and selecting the highest-likelihood tree.A custom shell script, available in the cited repository, performed the repeated runs and selection.