450 likes | 627 Views
Monte Carlo methods for estimating population genetic parameters. Rasmus Nielsen University of Copenhagen. Outline. Idiosyncratic history and background on ML estimation of demographic parameters based on DNA sequence data. A new computational approach/modification.
E N D
Monte Carlo methods for estimating population genetic parameters Rasmus Nielsen University of Copenhagen
Outline • Idiosyncratic history and background on ML estimation of demographic parameters based on DNA sequence data. • A new computational approach/modification. • Idiosyncratic history and background on ML estimation of demographic parameters based on SNP data. • Ascertainment and large scale SNP data sets.
Felsenstein’s Equation So where Gi,i=1,2,…k, has been simulated from p(G|Q).
Importance Sampling So where Gi,i=1,2,…k, has been simulated from h(G).
Griffiths and Tavare Recursion Simulate mutation (coalescent) from and correct using importance sampling.
Example (Nielsen 1998) • Infinite sites model • Estimation of T • Estimation of population phylogenies
Integro-recursion Ugliest equation ever published in a biological journal…
MCMC Set up a Markov chain on state space on all supported values of Q and G and with stationary distribution p(Q, G | X). Now since this can easily be done using Metropolis-Hastings sampling, i.e. updates to Q and G are proposed from a proposal distribution q(Q , G→ Q’ , G’) and accepted with probability
Some problems… • Histogram estimator or other smoothing must be used. • Likelihood ratios hard to estimate (e.g. M=0).
A new method • It is possible to calculate the marginal prior probability of a genealogy • It turns out that this math is doable, for most components of Θ such as q and M. • The we can sample from the marginal posterior of G using the previously discussed MCMC procedures. Slide inspired by Jody Hey
We then recover the posterior for Q using Approximated by Slide inspired by Jody Hey
Advantages • Eliminates problems with covariance between parameters leading to mixing problems. • Provides a smooth posterior/likelihood function useful for optimization and likelihood ratio estimation. Disadvantages • Requires more calculation in each MCMC iteration
Likelihood ratio estimation 6 loci, 15 gene copies, H0: m1=m2
Other approaches • Kuhner and Felsenstein use a combination of MCMC and importance sampling to estimate surfaces (no prior for the parameters). • PAC methods suggested by Stephens and Donnelly samples from a close approximation to to estimate an approximate likelihood. • ABC (Beaumont, Pritchard, Tavare and others) methods are a very popular and promising class of methods based on (1) reducing the data to summary statistics, (2) simulate new data from the prior, (3) accepting the parameter value under which the data was simulated if the difference between simulated and true statistics is less thand.
A more efficient method..Griffiths and Tavare (1998), Nielsen (2000)
A more efficient method..Griffiths and Tavare (1998), Nielsen (2000)
Ascertainment Sample vs. Typed Sample Ascertainment sample Typed sample
1.0 0.9 0.8 E[D'] 0.7 0.6 0.5 0 1 2 3 4 5 6 7 8 9 10 r =2Nc no ascertainment bias ascertainment bias
Correcting for ascertainment biases Now, for simplicity, consider the case without a sweep, then where (in the simplest possible case) and
In this simple case, the maximum likelihood estimate of P is simply given by , k = 1, 2, …, n – 1, where nk is the number of SNPs with allele frequency k. Selective sweeps: Similarly define
Hudson’s (2001) Estimator when n = 100, m = 5, r = 5, and #SNP pairs = 200. Uncorrected Corrected
Complications • Double-hit ascertainment (HapMap) • Ascertainment based on chimpanzee (HapMap) • Panel depth may vary among SNPs and/or among regions (HapMap). • Ascertainment method may vary among SNPs (HapMap). • Population structure (HapMap). • Loss of information regarding asc. scheme (HapMap??).
HapMap ascertainment depth distrb. (ignores many important components)
Perlegen HapMap
Data Directly sequenced polymorphism data from 20 European-Americans, 19 African-Americans and one chimpanzee from 9,316 protein coding genes Data set previously described in Bustamante, C.D. et al. 2005. Natural selection on protein-coding genes in the human genome. Nature 437, 1153-7.
Demographic model European-Americans African-Americans migration Admixture Population growth Bottleneck
Estimation Sampling probabilities from the 2D frequency spectrum , Number of SNPs with pattern j in the 2D frequency spectrum SNPs within a gene are correlated. But estimator is consistent. The estimate has the same properties as a real likelihood estimator except that it converges slightly slower because of the correlation (Wiuf 2006).
European-Americans Godness-of-fit: p = 0.6
Acknowledgements Jody Hey, John Wakeley, Melissa Hubisz, Andy Clark, Carlos Bustamante, Scott Williamson, Aida Andres, Amit Andip, Adam Boyko, Anders Albrechtsen,Mark Adams, Michelle Cargill and other staff at Celera Genomics and Applied Biosystems.