Ancestral population genomics: the coalescent hidden Markov model approach - PubMed (original) (raw)

Ancestral population genomics: the coalescent hidden Markov model approach

Julien Y Dutheil et al. Genetics. 2009 Sep.

Abstract

With incomplete lineage sorting (ILS), the genealogy of closely related species differs along their genomes. The amount of ILS depends on population parameters such as the ancestral effective population sizes and the recombination rate, but also on the number of generations between speciation events. We use a hidden Markov model parameterized according to coalescent theory to infer the genealogy along a four-species genome alignment of closely related species and estimate population parameters. We analyze a basic, panmictic demographic model and study its properties using an extensive set of coalescent simulations. We assess the effect of the model assumptions and demonstrate that the Markov property provides a good approximation to the ancestral recombination graph. Using a too restricted set of possible genealogies, necessary to reduce the computational load, can bias parameter estimates. We propose a simple correction for this bias and suggest directions for future extensions of the model. We show that the patterns of ILS along a sequence alignment can be recovered efficiently together with the ancestral recombination rate. Finally, we introduce an extension of the basic model that allows for mutation rate heterogeneity and reanalyze human-chimpanzee-gorilla-orangutan alignments, using the new models. We expect that this framework will prove useful for population genomics and provide exciting insights into genome evolution.

PubMed Disclaimer

Figures

F<sc>igure</sc> 1.—

Figure 1.—

The four different types that the genealogy of four sequences can take. The four species are H, human; C, chimpanzee; G, gorilla; and O, orangutan (outgroup). Following H

obolth

et al. (2007), these genealogies are labeled HC1, HC2, HG, and CG.

F<sc>igure</sc> 2.—

Figure 2.—

The genealogy of sequences varies along the genome. This graph shows the changes in topologies by colors (see Figure 1) together with the divergence between human and chimpanzee along 20,000 positions simulated using a coalescent with recombination (see

data and methods

).

F<sc>igure</sc> 3.—

Figure 3.—

The CoalHMM framework. The basic demographic model for human (H), chimpanzee (C), and gorilla (G) with the orangutan (O) as outgroup is used to design a set of candidate genealogies (topology plus coalescence times) and a corresponding transition matrix. These are computed according to the parameters of the model, namely speciation times (τ_X_), ancestral population sizes (θ_X_), and recombination rate ρ. The genealogies and transition probabilities are then used within a standard hidden Markov model to compute the likelihood and perform the posterior decoding.

F<sc>igure</sc> 4.—

Figure 4.—

Population and genealogy parameter notations.

F<sc>igure</sc> 5.—

Figure 5.—

Population parameters recovery: (a) 07 implementation; (b–d) 09 implementation. (a and b) Ancestral recombination graph (ARG) simulated using a coalescent with recombination; (c) ARG simulated from the HMM; (d) ARG with averaged coalescent times. The dashed lines show the true values.

F<sc>igure</sc> 6.—

Figure 6.—

(a) Estimation of recombination rates. Ten data sets were simulated for each value of the recombination rate. The recombination rate was fixed to the same value in all lineages, and the 09 implementation with ρH = ρC = ρG was used for the estimation. Other parameters are set to the same values as in Figure 5. The dashed line shows the 1:1 ratio, and the solid line is fitted to the data. The slope is estimated to be 0.29. (b) Global vs. lineage-specific recombination rates. One hundred simulations were performed in each case, with parameter values identical to those in Figure 5. Dashed lines represent the true value of the recombination rate. Overlapping notches show nonsignificant differences in the medians.

F<sc>igure</sc> 7.—

Figure 7.—

Population parameters recovery when the bias correction is used. The parameter values used in the simulation are the same as in Figure 5.

F<sc>igure</sc> 8.—

Figure 8.—

Distribution of posterior fragment lengths. Shaded bars show the posterior distribution, hatched bars show the distribution from the ancestral recombination graph, and the solid lines show the geometric distribution (expected under the Markov assumption).

F<sc>igure</sc> 9.—

Figure 9.—

Comparison of the posterior decoding with the known ancestral recombination graph in 100 simulated data sets. (a) Joint distribution of reconstructed states (by column) and real states (by row). The numbers indicate the average percentage of each configuration. The proportions in all combinations sum to 100%, and the frequencies on the top-left/bottom-right diagonal correspond to the correctly inferred cases. The sums over each row give the true frequencies of each state, and the sums over columns provide the frequencies inferred by the CoalHMM method. (b) Precision and recall measures. The precision measure corresponds, for each column of panel 9a, to the diagonal proportion divided by the sum of the column (number of true positives over number of positives). The recall measure corresponds, for each row of panel 9a, to the diagonal proportion divided by the sum of the row (number of true positives over actual number of sites in a given genealogy). (c) Comparison of the reconstructed and real frequencies of hidden states. The proportion of incomplete lineage sorting is taken as the frequency of sites in alternative genealogies.

F<sc>igure</sc> 10.—

Figure 10.—

Model comparison on target 1 of the human–chimpanzee–gorilla–orangutan alignment. LRT: likelihood-ratio test. ΔBIC: difference in Bayesian information criterion. The arrows show the direction of the model comparison, pointing toward the alternative model. BIC is the difference in the BIC indexes between the two models. A negative value hence favors the alternative model, whereas a positive value favors the null model. NS: nonsignificant. ***P_-value <1_e − 16. Δ ln L: logarithm of the likelihood compared to the value of the 09 model. Γ: Gamma-distributed mutation rates.

F<sc>igure</sc> 11.—

Figure 11.—

Parameter estimates obtained with different models applied to target 1. The parameter estimate for the gorilla recombination rate with the 09 model is 24 cM/Mb.

Similar articles

Cited by

References

    1. Burgess, R., and Z. Yang, 2008. Estimation of hominoid ancestral population sizes under Bayesian coalescent models incorporating mutation rate variation and sequencing errors. Mol. Biol. Evol. 25: 1979–1994. - PubMed
    1. Chen, F. C., and W. H. Li, 2001. Genomic divergences between humans and other hominoids and the effective population size of the common ancestor of humans and chimpanzees. Am. J. Hum. Genet. 68: 444–456. - PMC - PubMed
    1. Durbin, R., S. Eddy, A. Krogh and G. Mitchison, 1998. Biological Sequence Analysis: Probabilistic Models of Proteins and Nucleic Acids. Cambridge University Press, Cambridge/London/New York.
    1. Dutheil, J., and B. Boussau, 2008. Non-homogeneous models of sequence evolution in the Bio++ suite of libraries and programs. BMC Evol. Biol. 8: 255. - PMC - PubMed
    1. Dutheil, J., S. Gaillard, E. Bazin, S. Glémin, V. Ranwez et al., 2006. Bio++: a set of C++ libraries for sequence analysis, phylogenetics, molecular evolution and population genetics. BMC Bioinformatics 7: 188. - PMC - PubMed

Publication types

MeSH terms

LinkOut - more resources