Inferring the joint demographic history of multiple populations from multidimensional SNP frequency data

Demographic models built from genetic data play important roles in illuminating prehistorical events and serving as null models in genome scans for selection. We introduce an inference method based on the joint frequency spectrum of genetic variants within and between populations. For candidate models we numerically compute the expected spectrum using a diffusion approximation to the one-locus two-allele Wright-Fisher process, involving up to three simultaneous populations. Our approach is a composite likelihood scheme, since linkage between neutral loci alters the variance but not the expectation of the frequency spectrum. We thus use bootstraps incorporating linkage to estimate uncertainties for parameters and significance values for hypothesis tests. Our method can also incorporate selection on single sites, predicting the joint distribution of selected alleles among populations experiencing a bevy of evolutionary forces, including expansions, contractions, migrations, and admixture. As applications, we model human expansion out of Africa and the settlement of the New World, using 5 Mb of noncoding DNA resequenced in 68 individuals from 4 populations (YRI, CHB, CEU, and MXL) by the Environmental Genome Project. We also combine our demographic model with a previously estimated distribution of selective effects among newly arising amino acid mutations to accurately predict the frequency spectrum of nonsynonymous variants across three continental populations (YRI, CHB, CEU).

💡 Research Summary

The paper presents a comprehensive framework for inferring complex demographic histories of multiple populations using the joint allele frequency spectrum (JFS) derived from genome‑wide single‑nucleotide polymorphism (SNP) data. The authors build on the classic Wright‑Fisher model for a single locus with two alleles, approximating its discrete‑generation dynamics with a continuous‑time diffusion equation. By imposing appropriate boundary conditions, they incorporate a wide range of evolutionary forces—population size changes, bidirectional migration, admixture events, and even selection on individual sites—into the diffusion framework. The diffusion equation is solved numerically using finite‑difference methods, allowing the expected JFS to be computed for up to three populations simultaneously.

Parameter estimation proceeds via a composite likelihood approach: the likelihood of the observed JFS is approximated as the product of marginal likelihoods for each SNP, assuming independence. Because real SNPs are linked, this assumption underestimates variance. To correct for linkage disequilibrium, the authors employ a block bootstrap strategy that preserves the correlation structure by resampling contiguous genomic blocks. Re‑estimating parameters across many bootstrap replicates yields robust standard errors, confidence intervals, and hypothesis‑test p‑values.

A major methodological advance is the integration of selection. The authors treat the selection coefficient as a time‑dependent parameter that can be added to the diffusion operator. By coupling this with an externally estimated distribution of fitness effects for new nonsynonymous mutations, they predict the joint frequency distribution of selected alleles across populations undergoing arbitrary demographic scenarios. This extension successfully reproduces the excess of high‑frequency nonsynonymous variants that neutral models cannot explain.

The framework is applied to real data from the Environmental Genome Project: 68 individuals from four populations (YRI, CHB, CEU, and MXL) sequenced across 5 Mb of non‑coding DNA. The inferred model includes an early African effective size, a rapid out‑of‑Africa expansion, distinct growth rates for East Asian and European lineages, asymmetric migration between Asia and Europe, and a recent admixture component in the Mexican (MXL) sample. Bootstrap‑derived confidence intervals demonstrate that the data reject a simple single‑step out‑of‑Africa model in favor of a more gradual, multi‑phase expansion.

Finally, the authors combine their demographic estimates with a previously derived distribution of selective effects for amino‑acid‑changing mutations. Using this joint model, they accurately predict the observed frequency spectra of nonsynonymous variants in YRI, CHB, and CEU, achieving high concordance across low‑ and intermediate‑frequency bins. This result underscores the importance of jointly modeling demography and selection when interpreting patterns of functional variation.

In summary, the study delivers a versatile, numerically tractable method for simultaneous inference of demographic parameters and selection pressures from multi‑population SNP data. Its strengths lie in (1) the ability to handle up to three populations in a unified diffusion‑based JFS calculation, (2) rigorous uncertainty quantification via block bootstrapping, and (3) seamless incorporation of selection to explain functional variant spectra. Limitations include the current three‑population ceiling and the computational intensity of solving high‑dimensional diffusion equations, which will require further algorithmic optimization for whole‑genome applications. Nonetheless, the approach sets a new standard for reconstructing human evolutionary history and provides a robust null model for genome‑wide scans of adaptation.