Bio Phylo Species Trees

Version Compatibility

Reference examples tested with: ASTER 1.15+ (provides astral/ASTRAL-III, wastral, astral-pro), IQ-TREE 2.2+ (concordance factors), PAUP* 4.0a168+ (SVDQuartets), BPP 4+.

Before using code patterns, verify installed versions match. If versions differ:

• CLI: astral --version then astral --help to confirm flags
• CLI: iqtree2 --version then iqtree2 --help to confirm flags
• CLI: PAUP* prints its version at startup; BPP reads a control file

If code throws an error, introspect the installed version and adapt flags rather than retrying.

The modern binary is astral from the ASTER package; input is a file of gene-tree Newick strings, one tree per line. The classic-Java ASTRAL uses -t for ANNOTATION level, but the ASTER astral uses -t for THREADS and -u for annotation -- copying a Java -t 8 (quartet support) into ASTER silently requests 8 threads with no annotation. Confirm which interface is installed before scripting.

Coalescent Species-Tree Estimation -- A Species Tree Is Not a Gene Tree, and Concatenation Is Inconsistent in the Anomaly Zone

"Estimate the species tree from my multi-locus phylogenomic data" -> Treat the dataset as a distribution of genealogies, not one tree, and estimate the species tree most consistent with that distribution under the multispecies coalescent.

• CLI: per-locus gene trees (modern-tree-inference) then astral -i gene_trees.nwk -o species.tre (ASTER)
• CLI: gene and site concordance factors via iqtree2 ... --gcf ... --scfl (shared with modern-tree-inference)

Scope: choosing and running a coalescent species-tree method, contracting noisy gene trees, and reading discordance honestly. Per-locus gene-tree inference and model selection -> modern-tree-inference. Computing and plotting gCF/sCF -> modern-tree-inference, tree-visualization. Single-copy vs multi-copy orthology upstream -> comparative-genomics/ortholog-inference. Dated species trees -> divergence-dating. Rooting -> tree-manipulation.

The Single Most Important Modern Insight

A species tree is not a gene tree. Each locus has its own genealogy, and under the multispecies coalescent (MSC) those genealogies disagree with the species tree -- and with each other -- by incomplete lineage sorting (ILS) alone, with zero estimation error. The disagreement is not noise to average away; it is signal generated by the coalescent running inside the species tree. Three load-bearing facts:

1. The most probable gene tree can be the wrong species tree. When two short speciation events occur in quick succession, the single most common gene-tree topology is not the species tree -- this is the anomaly zone (Degnan and Rosenberg 2006), and such gene trees are anomalous gene trees. A plurality vote across loci then converges on the wrong tree.
2. Concatenation is statistically inconsistent and positively misleading under ILS. It pretends all loci share one tree and behaves like a vote dominated by the most-supported site pattern; in the anomaly zone that pattern reflects the anomalous genealogy, so adding loci raises bootstrap support for the wrong topology toward 100% (Kubatko and Degnan 2007; Roch and Steel 2015). More data makes the error more confident, not less. The support is the trap.
3. The honest question is the distribution, not the tree. Coalescent summary methods (ASTRAL family) stay consistent because they ask "what species tree is most consistent with having generated this gene-tree distribution under the MSC?" -- and they work at the quartet level, where the matching topology is always the plurality and no anomaly zone exists. Concatenation answers a different, wrong question, and both outputs look identically resolved and supported.

When Concatenation Fails

Concatenation is consistent and efficient only when ILS is low. It is actively wrong in three coupled regimes, all driven by short internodes measured in coalescent units (generations / 2Ne):

1. Short internodes plus large ancestral Ne. A short branch gives lineages little time to coalesce before the next deeper population, so many survive and sort randomly. For a rooted triple, P(gene tree matches) = 1 - (2/3)e^(-tau) in coalescent units tau: at tau = 1 about 25% of gene trees are discordant, at tau = 0.3 about 50%, at tau = 0.1 about 60%. Large Ne shortens tau for the same number of years, inflating ILS.
2. The anomaly zone (Degnan and Rosenberg 2006). With two or more short successive internal branches (both well under ~1, especially under ~0.1-0.2 coalescent units), the most probable gene tree differs from the species tree. The zone grows with taxon number and with shorter, more clustered internodes.
3. Rapid radiations. Successive speciation in quick succession is built from exactly those short internodes, so radiations routinely show a majority of gene trees conflicting with the species tree. This is the canonical case for a coalescent method.

Method Taxonomy

Summary methods take estimated gene trees as input and assume they are true -- the central caveat (see failure modes). Site-based (SVDQuartets) and full-likelihood (BPP, StarBEAST2) skip that assumption but pay in scale: roughly tens of taxa for the full-likelihood methods. The speed/scale order (most scalable first) is concatenation < ASTRAL/wASTRAL < SVDQuartets < StarBEAST2 < BPP; model completeness runs the opposite way. Pick the most scalable method whose assumptions the data do not badly violate. Introgression detection (HyDe, Blischak 2018, Syst Biol; QuIBL, Edelman 2019, Science; the D-statistic) lives downstream of the symmetry diagnostic below.

Concatenation vs Coalescent Decision

Do not reflexively coalescent-everything: where ILS is low the two trees agree and concatenation pools weak per-locus signal robustly. The deliverable is concordance, not a single best method -- run both, compute gCF/sCF, and trust the coalescent tree exactly where the two disagree on low-gCF branches.

Estimate the Species Tree from Gene Trees

Goal: Turn per-locus gene trees into a species tree that is consistent under the MSC, after removing gene-tree noise that would bias the quartet counts.

Approach: Infer one gene tree per locus with branch support (modern-tree-inference), contract branches below ~10% support to polytomies (which ASTRAL handles correctly), then run wASTRAL as the primary estimate and astral for the classic localPP / polytomy-test workflow.

# Per-locus gene trees with support live upstream in modern-tree-inference:
#   iqtree2 -S loci_dir -m MFP -B 1000 -T AUTO --prefix loci   # one tree per locus
# Contract gene-tree branches below 10% support (Newick tools / nw_ed) -> polytomies,
# then concatenate into one file, one Newick per line:
#   cat loci/*.contracted.treefile > gene_trees.nwk

# Primary estimate: wASTRAL weights quartets by gene-tree support+length (noisy gene trees)
wastral -i gene_trees.nwk -o species_wastral.tre 2> species_wastral.log

# Classic ASTRAL-III for localPP + polytomy test. In ASTER, -t is THREADS, -u is annotation:
astral -t 8 -i gene_trees.nwk -o species.tre 2> species.log     # -t 8 = 8 threads here
astral -u 2 -i gene_trees.nwk -o species_annot.tre              # full localPP + q1/q2/q3 for all 3 resolutions
astral --root OUTGROUP -i gene_trees.nwk -o species_rooted.tre  # rooting improves branch-length estimation

# Multi-copy gene families (no pre-orthology) -> ASTRAL-Pro:
astral-pro -i family_trees.nwk -o species_pro.tre

ASTRAL returns an UNROOTED tree; root with an outgroup afterward. Branch lengths are in coalescent units (not time, not substitutions), and tip lengths are undefined in that mode. The classic Java astral.5.7.8.jar inverts the flags: there -t is the annotation level (-t 2 full, -t 8 quartet support, -t 10 polytomy test) and there is no -u.

Compute and Read Concordance Factors

Goal: Replace a single bootstrap/localPP per branch with how much of the actual data agrees, and read ILS-vs-introgression off the quartet symmetry.

Approach: Compute gene and site concordance factors against the species tree, then inspect the two minority quartet frequencies (q2, q3) at each contested branch.

# gCF = % of decisive gene trees containing a branch; sCF = % of decisive sites supporting it
iqtree2 -te species.tre --gcf gene_trees.nwk -s concat.fasta --scfl 100 --prefix cf   # -te fixes the topology for likelihood sCF
# cf.cf.tree carries gCF/sCF on node labels; cf.cf.stat has per-branch q1/q2/q3 counts

Then read each contested branch: under pure ILS the two minority topologies are equally probable (q2 ~ q3); introgression breaks that symmetry by inflating the one matching the direction of gene flow. A high-bootstrap branch with gCF ~ 25 is a branch where the data are screaming disagreement that the bootstrap hides.

Interpreting Concordance and Discordance

• gCF / sCF (Minh 2020) are the honest support. Bootstrap and localPP saturate near their ceiling at genome scale even when most loci disagree; gCF/sCF report the actual conflict. Report them on every phylogenomic tree.
• localPP (Sayyari and Mirarab 2016) is NOT a bootstrap. It is a coalescent posterior computed from the quartet frequencies around a branch. localPP = 1.0 is common and expected on resolved branches and does not mean "100 of 100 bootstraps"; localPP ~ 0.33 means the three quartet resolutions are tied (no resolution).
• The polytomy test (Sayyari and Mirarab 2018) is the principled alternative to staring at low support. It tests the null q1 = q2 = q3 (an effectively simultaneous speciation); failing to reject means a resolved bifurcation cannot be claimed -- collapse the branch.
• Quartet symmetry separates ILS from introgression. Symmetric minority quartets (q2 ~ q3) = ILS; handle with a coalescent tree. Asymmetric (q2 != q3, significant D-statistic / ABBA-BABA, HyDe gamma != 0, or QuIBL favoring the introgression component) = gene flow; a bifurcating tree cannot represent it -- model the reticulation. Confirm the asymmetry before invoking introgression.

Per-Method Failure Modes

Gene-Tree Estimation Error Biasing ASTRAL

Trigger: Short, low-information loci (e.g. <500 bp exons) fed to ASTRAL/wASTRAL as if their gene trees were truth. Mechanism: MSC consistency assumes the input gene trees are correct; estimation error is a non-coalescent noise source that biases quartet counts in finite samples and can make the coalescent tree worse than concatenation. Symptom: Coalescent tree noisier or worse than concatenation; unstable across locus subsets; low localPP everywhere; gCF far below sCF. Fix: Contract gene-tree branches below ~10% bootstrap before ASTRAL (they become polytomies that contribute no spurious quartet similarity); better, use wASTRAL (continuous weighting); or drop to SVDQuartets (no gene trees) or BPP (integrates gene-tree uncertainty).