Estimating Variable Effective Population Sizes from Multiple Genomes: A Sequentially Markov Conditional Sampling Distribution Approach

Throughout history, the population size of modern humans has varied considerably due to changes in environment, culture, and technology. More accurate estimates of population size changes, and when they occurred, should provide a clearer picture of human colonization history and help remove confounding effects from natural selection inference. Demography influences the pattern of genetic variation in a population, and thus genomic data of multiple individuals sampled from one or more present-day populations contain valuable information about the past demographic history. Recently, Li and Durbin developed a coalescent-based hidden Markov model, called the pairwise sequentially Markovian coalescent (PSMC), for a pair of chromosomes (or one diploid individual) to estimate past population sizes. This is an efficient, useful approach, but its accuracy in the very recent past is hampered by the fact that, because of the small sample size, only few coalescence events occur in that period. Multiple genomes from the same population contain more information about the recent past, but are also more computationally challenging to study jointly in a coalescent framework. Here, we present a new coalescent-based method that can efficiently infer population size changes from multiple genomes, providing access to a new store of information about the recent past. Our work generalizes the recently developed sequentially Markov conditional sampling distribution framework, which provides an accurate approximation of the probability of observing a newly sampled haplotype given a set of previously sampled haplotypes. Simulation results demonstrate that we can accurately reconstruct the true population histories, with a significant improvement over the PSMC in the recent past. We apply our method, called diCal, to the genomes of multiple human individuals of European and African ancestry to obtain a detailed population size change history during recent times.     

Estimating the number of founder lineages from haplotypes of closely linked SNPs

We consider an isolated population founded by a small number of individuals randomly chosen from a source population of known genetic composition at a known time in the past. We develop a Monte-Carlo maximum-likelihood method for estimating the number of founding individuals from the haplotype frequencies at several SNP (single nucleotide polymorphism) loci in a sample. We assume the isolated population was founded recently enough that that mutation can be ignored and that haplotype frequencies in the source population have not changed. We apply the method to simulated data and show that it is unbiased. With a reasonable number of individuals sampled, it is possible to estimate the number of founders within a factor of 2. We show that the performance of the method is not degraded substantially if the frequencies of the rare haplotypes in the source are not known precisely and if there is some recombination. We illustrate the use of our method by applying it to a previously published data set from a recently founded population of wolves (Canis lupus) in Scandinavia.     

The effect of population history on the lengths of ancestral chromosome segments

An isolated population is a group of individuals who are descended from a founding population who lived some time ago. If the founding individuals are assumed to be noninbred and unrelated, a chromosome sampled from the population can be represented as a mosaic of segments of the original ancestral types. A population in which chromosomes are made up of a few long segments will exhibit linkage disequilibrium due to founder effect over longer distances than a population in which the chromosomes are made up of many short segments. We study the length of intact ancestral segments by obtaining the expected number of junctions (points where DNA of two distinct ancestral types meet) in a chromosome. Assuming random mating, we study analytically the effects of population age, growth patterns, and internal structure on the expected number of junctions in a chromosome. We demonstrate that the type of growth a population has experienced can influence the expected number of junctions, as can population subdivision. These effects are substantial only when population sizes are very small. We also develop an approximation to the variance of the number of junctions and show that the variance is large.     

Using molecular markers with high mutation rates to obtain estimates of relative population size and to distinguish the effects of gene flow and mutation: a demonstration using data from endemic Mauritian skinks

We propose a method of analysing genetic data to obtain separate estimates of the size (N(p)) and migration rate (m(p)) for the sampled populations, without precise prior knowledge of mutation rates at each locus ( micro(L)). The effects of migration and mutation can be distinguished because high migration has the effect of reducing genetic differentiation across all loci, whereas a high mutation rate will only affect the locus in question. The method also takes account of any differences between the spectra of immigrant alleles and of new mutant alleles. If the genetic data come from a range of population sizes, and the loci have a range of mutation rates, it is possible to estimate the relative sizes of the different N(p) values, and likewise the m(p) and the micro(L). Microsatellite loci may also be particularly appropriate because loci with a high mutation rate can reach mutation-drift-migration equilibrium more quickly, and because the spectra of mutants arriving in a population can be particularly distinct from the immigrants. We demonstrate this principle using a microsatellite data set from Mauritian skinks. The method identifies low gene flow between a putative new species and populations of its sister species, whereas the differentiation of two other populations is attributed to small population size. These distinct interpretations were not readily apparent from conventional measures of genetic differentiation and gene diversity. When the method is evaluated using simulated data sets, it correctly distinguishes low gene flow from small population size. Loci that are not at mutation-migration-drift equilibrium can distort the parameter estimates slightly. We discuss strategies for detecting and overcoming this effect.     

The effects of sample size on population genetic diversity estimates in song sparrows Melospiza melodia

To empirically determine the effects of sample size on commonly used measures of average genetic diversity, we genotyped 200 song sparrows Melospiza melodia from two populations, one genetically depauperate (n=100) and the other genetically diverse (n=100), using eight microsatellite loci. These genotypes were used to randomly create 10,000 datasets of differing sizes (5 to 50) for each population to determine what the effects of sample size might be on several estimates of genetic diversity (number of alleles per locus, average observed heterozygosity, and unbiased average expected heterozygosity) in natural populations of conservation concern. We found that at small sample sizes of 5 to 10 individuals, estimates of unbiased heterozygosity outperformed those based on observed heterozygosity or allelic diversity for both low- and high-diversity populations. We also found that when comparing across populations in which different numbers of individuals were sampled, rarefaction provided a useful way to compare estimates of allelic diversity. We recommend that standard errors should be reported for all diversity estimators, especially when sample sizes are small. We also recommend that at least 20 to 30 individuals be sampled in microsatellite studies that assess genetic diversity when working in a population that has an unknown level of diversity. However, research on critically endangered populations (where large sample sizes are impossible or extremely difficult to obtain) should include measures of genetic diversity even if sample sizes are less than ideal. These estimates can be useful in assessing the genetic diversity of the population.     

Sibship within samples of brown trout (<Emphasis Type="Italic">Salmo trutta</Emphasis>) and implications for supportive breeding

We analysed family relationships among brown trout from two small tributary populations that have been suggested as a source of individuals for supportive breeding, using variation at eight microsatellite loci. As a control, we analysed a sample of supposedly unrelated individuals representing a large anadromous population, and we simulated unrelated individuals based on the allelic distributions in all three samples. Two different approaches were used: (1) pairwise estimates of relatedness between individuals and (2) a method for partitioning individuals into half-sib and full-sib families. The anadromous population did not show evidence of a significant number of closely related individuals. In both tributary populations, however, the distributions of pairwise relatedness estimates suggested the presence of several related individuals, and sibship reconstruction suggested fewer families consisting of more individuals than were observed for the simulated individuals. The expected increase of inbreeding coefficient in the two samples due to family structure was 0.026 and 0.030 respectively. Moreover, tests for recent bottlenecks yielded significant outcomes in both populations suggesting a history of low effective population sizes. Depending on the effective population size of captive spawners and past effective population sizes in the populations it could be beneficial to conduct sib-avoidance matings, though this cannot eliminate inbreeding but only delay it. Alternatively, individuals from different populations could be crossed. Sibship reconstruction provided the clearest evidence for family structure, but pairwise relatedness is the best measure for designing mating schemes, as it allows for mating as unrelated individuals as possible rather than just avoiding mating between sibs.     

Joint inference of population assignment and demographic history

A new approach to assigning individuals to populations using genetic data is described. Most existing methods work by maximizing Hardy-Weinberg and linkage equilibrium within populations, neither of which will apply for many demographic histories. By including a demographic model, within a likelihood framework based on coalescent theory, we can jointly study demographic history and population assignment. Genealogies and population assignments are sampled from a posterior distribution using a general isolation-with-migration model for multiple populations. A measure of partition distance between assignments facilitates not only the summary of a posterior sample of assignments, but also the estimation of the posterior density for the demographic history. It is shown that joint estimates of assignment and demographic history are possible, including estimation of population phylogeny for samples from three populations. The new method is compared to results of a widely used assignment method, using simulated and published empirical data sets.     

Detecting Small Amounts of Gene Flow from Phylogenies of Alleles

The method of coalescents is used to find the probability that none of the ancestors of alleles sampled from a population are immigrants. If that is the case for samples from two or more populations, then there would be concordance between the phylogenies of those alleles and the geographic locations from which they are drawn. This type of concordance has been found in several studies of mitochondrial DNA from natural populations. It is shown that if the number of sequences sampled from each population is reasonably large (10 or more), then this type of concordance suggests that the average number of individuals migrating between populations is likely to be relatively small (Nm less than 1) but the possibility of occasional migrants cannot be excluded. The method is applied to the data of E. Bermingham and J. C. Avise on mtDNA from the bowfin, Amia calva.     

Species trees from gene trees: reconstructing Bayesian posterior distributions of a species phylogeny using estimated gene tree distributions

The desire to infer the evolutionary history of a group of species should be more viable now that a considerable amount of multilocus molecular data is available. However, the current molecular phylogenetic paradigm still reconstructs gene trees to represent the species tree. Further, commonly used methods of combining data, such as the concatenation method, are known to be inconsistent in some circumstances. In this paper, we propose a Bayesian hierarchical model to estimate the phylogeny of a group of species using multiple estimated gene tree distributions, such as those that arise in a Bayesian analysis of DNA sequence data. Our model employs substitution models used in traditional phylogenetics but also uses coalescent theory to explain genealogical signals from species trees to gene trees and from gene trees to sequence data, thereby forming a complete stochastic model to estimate gene trees, species trees, ancestral population sizes, and species divergence times simultaneously. Our model is founded on the assumption that gene trees, even of unlinked loci, are correlated due to being derived from a single species tree and therefore should be estimated jointly. We apply the method to two multilocus data sets of DNA sequences. The estimates of the species tree topology and divergence times appear to be robust to the prior of the population size, whereas the estimates of effective population sizes are sensitive to the prior used in the analysis. These analyses also suggest that the model is superior to the concatenation method in fitting these data sets and thus provides a more realistic assessment of the variability in the distribution of the species tree that may have produced the molecular information at hand. Future improvements of our model and algorithm should include consideration of other factors that can cause discordance of gene trees and species trees, such as horizontal transfer or gene duplication.     

Inferring recent migration rates from individual genotypes

We present a novel and straightforward method for estimating recent migration rates between discrete populations using multilocus genotype data. The approach builds upon a two-step sampling design, where individual genotypes are sampled before and after dispersal. We develop a model that estimates all pairwise backwards migration rates ( m_ij , the probability that an individual sampled in population i is a migrant from population j ) between a set of populations. The method is validated with simulated data and compared with the methods of BayesAss and Structure. First, we use data for an island model and then we consider more realistic data simulations for a metapopulation of the greater white-toothed shrew ( Crocidura russula ). We show that the precision and bias of estimates primarily depend upon the proportion of individuals sampled in each population. Weak sampling designs may particularly affect the quality of the coverage provided by 95% highest posterior density intervals. We further show that it is relatively insensitive to the number of loci sampled and the overall strength of genetic structure. The method can easily be extended and makes fewer assumptions about the underlying demographic and genetic processes than currently available methods. It allows backwards migration rates to be estimated across a wide range of realistic conditions.     

Effects of sample size on estimates of population growth rates calculated with matrix models

Background

Matrix models are widely used to study the dynamics and demography of populations. An important but overlooked issue is how the number of individuals sampled influences estimates of the population growth rate (λ) calculated with matrix models. Even unbiased estimates of vital rates do not ensure unbiased estimates of λ–Jensen''s Inequality implies that even when the estimates of the vital rates are accurate, small sample sizes lead to biased estimates of λ due to increased sampling variance. We investigated if sampling variability and the distribution of sampling effort among size classes lead to biases in estimates of λ.

Methodology/Principal Findings

Using data from a long-term field study of plant demography, we simulated the effects of sampling variance by drawing vital rates and calculating λ for increasingly larger populations drawn from a total population of 3842 plants. We then compared these estimates of λ with those based on the entire population and calculated the resulting bias. Finally, we conducted a review of the literature to determine the sample sizes typically used when parameterizing matrix models used to study plant demography.

Conclusions/Significance

We found significant bias at small sample sizes when survival was low (survival = 0.5), and that sampling with a more-realistic inverse J-shaped population structure exacerbated this bias. However our simulations also demonstrate that these biases rapidly become negligible with increasing sample sizes or as survival increases. For many of the sample sizes used in demographic studies, matrix models are probably robust to the biases resulting from sampling variance of vital rates. However, this conclusion may depend on the structure of populations or the distribution of sampling effort in ways that are unexplored. We suggest more intensive sampling of populations when individual survival is low and greater sampling of stages with high elasticities.     

An instantaneous coalescent method insensitive to population structure

Conventional coalescent inferences of population history make the critical assumption that the population under examination is panmictic. However, most populations are structured. This complicates the prevailing coalescent analyses and sometimes leads to inaccurate estimates. To develop a coalescent method unhampered by population structure, we perform two analyses. First, we demonstrate that the coalescent probability of two randomly sampled alleles from the immediate preceding generation(one generation back)is independent of population structure. Second, motivated by this finding, we propose a new coalescent method: i-coalescent analysis. The i-coalescent analysis computes the instantaneous coalescent rate by using a phylogenetic tree of sampled alleles. Using simulated data, we broadly demonstrate the capability of i-coalescent analysis to accurately reconstruct population size dynamics of highly structured populations, although we find this method often requires larger sample sizes for structured populations than for panmictic populations. Overall, our results indicate i-coalescent analysis to be a useful tool, especially for the inference of population histories with intractable structure such as the developmental history of cell populations in the organs of complex organisms.     

Life histories, neighbourhood sizes, and variance structure in some North American conifers

Neighourhood sizes and variance structure were determined for ten North American coniferous forest tree species, using information on seed dispersal, patterns of life history, and demography. A synthetic cohort of reproductive individuals belonging to a critical age group was constructed using normal yield tables. Data on reproductive individuals in a synthetic cohort were used to determine the neighbourhood sizes for each species. The neighbourhood sizes varied from 1244 to 50 118 individuals. Effective population sizes were calculated from data on neighbourhood sizes and replacement patterns of individuals in various species, the interdeme differentiation, F _st ranged from 0.0002 to 0.007. The results suggest that large neighbourhood and effective population sizes, and hence panmixia, may be virtually common among the species examined. These conclusions are consistent with recent findings on population differentiation and estimates of gene flow in confers using allozyme markers.     

A genetic demographic analysis of Lake Malawi rock‐dwelling cichlids using spatio‐temporal sampling

We estimated the effective population sizes (N_e) and tested for short‐term temporal demographic stability of populations of two Lake Malawi cichlids: Maylandia benetos, a micro‐endemic, and Maylandia zebra, a widespread species found across the lake. We sampled a total of 351 individuals, genotyped them at 13 microsatellite loci and sequenced their mitochondrial D‐loop to estimate genetic diversity, population structure, demographic history and effective population sizes. At the microsatellite loci, genetic diversity was high in all populations. Yet, genetic diversity was relatively low for the sequence data. Microsatellites yielded mean N_e estimates of 481 individuals (±99 SD) for M. benetos and between 597 (±106.3 SD) and 1524 (±483.9 SD) individuals for local populations of M. zebra. The microsatellite data indicated no deviations from mutation–drift equilibrium. Maylandia zebra was further found to be in migration–drift equilibrium. Temporal fluctuations in allele frequencies were limited across the sampling period for both species. Bayesian Skyline analyses suggested a recent expansion of M. zebra populations in line with lake‐level fluctuations, whereas the demographic history of M. benetos could only be estimated for the very recent past. Divergence time estimates placed the origin of M. benetos within the last 100 ka after the refilling of the lake and suggested that it split off the sympatric M. zebra population. Overall, our data indicate that micro‐endemics and populations in less favourable habitats have smaller N_e, indicating that drift may play an important role driving their divergence. Yet, despite small population sizes, high genetic variation can be maintained.     

50,000 years of ice and seals: Impacts of the Last Glacial Maximum on Antarctic fur seals

Ice is one of the most important drivers of population dynamics in polar organisms, influencing the locations, sizes, and connectivity of populations. Antarctic fur seals, Arctocephalus gazella, are particularly interesting in this regard, as they are concomitantly reliant on both ice‐associated prey and ice‐free coastal breeding areas. We reconstructed the history of this species through the Last Glacial Maximum (LGM) using genomic sequence data from seals across their range. Population size trends and divergence events were investigated using continuous‐time size estimation analysis and divergence time estimation models. The combined results indicated that a panmictic population present prior to the LGM split into two small refugial populations during peak ice extent. Following ice decline, the western refugial population founded colonies at the South Shetlands, South Georgia, and Bouvetøya, while the eastern refugial population founded the colony on Iles Kerguelen. Postglacial population divergence times closely match geological estimates of when these coastal breeding areas became ice free. Given the predictions regarding continued future warming in polar oceans, these responses of Antarctic fur seals to past climate variation suggest it may be worthwhile giving conservation consideration to potential future breeding locations, such as areas further south along the Antarctic Peninsula, in addition to present colony areas.     

The effect of single nucleotide polymorphism identification strategies on estimates of linkage disequilibrium

At present there is tremendous interest in characterizing the magnitude and distribution of linkage disequilibrium (LD) throughout the human genome, which will provide the necessary foundation for genome-wide LD analyses and facilitate detailed evolutionary studies. To this end, a human high-density single-nucleotide polymorphism (SNP) marker map has been constructed. Many of the SNPs on this map, however, were identified by sampling a small number of chromosomes from a single population, and inferences drawn from studies using such SNPs may be influenced by ascertainment bias (AB). Through extensive simulations, we have found that AB is a potentially significant problem in estimating and comparing LD within and between populations. Specifically, the magnitude of AB is a function of the SNP discovery strategy, number of chromosomes used for SNP discovery, population genetic characteristics of the particular genomic region considered, amount of gene flow between populations, and demographic history of the populations. We demonstrate that a balanced SNP discovery strategy (where equal numbers of chromosomes are sampled from multiple subpopulations) is the optimal study design for generating broadly applicable SNP resources. Finally, we validate our theoretical predictions by comparing our results to publicly available data from ten genes sequenced in 24 African American and 23 European American individuals.     

Recovery of mutations of different sizes from a population sample of DNA sequences under variable mutation rates across sites.

Mutations may be classified according to their positions of occurrence in the genealogy of the sampled DNA sequences from a population. A mutation is said to be of size i if it has i descendants in the sample. Such classifications for mutations may yield detailed insights into the evolutionary history and properties of the population. Statistical methods based on such classification have been developed and shown to be efficient and powerful. However, the utility of these statistical methods critically depends on reliable and robust recovery of mutations of different sizes. We investigated the distributional changes of mutations of different sizes due to genealogy reconstruction using the unweighted pair-group method with arithmetic mean (UPGMA) and the performance of maximum-parsimony method in inferring mutations of different sizes on a given topology. Genealogy reconstruction by UPGMA was found to change the distribution of mutations of different sizes on constructed topologies. Multiple hits at some nucleotide sites made it difficult to infer mutations of different sizes with the maximum-parsimony method, even when the true topology was designated. These results suggest that while the newly developed statistical methods employing information on mutations of different sites are powerful, they also impose significant new challenges for developing methods to accurately recover mutations of different sizes from population DNA sequence data.     

On the number of New World founders: a population genetic portrait of the peopling of the Americas

The founding of New World populations by Asian peoples is the focus of considerable archaeological and genetic research, and there persist important questions on when and how these events occurred. Genetic data offer great potential for the study of human population history, but there are significant challenges in discerning distinct demographic processes. A new method for the study of diverging populations was applied to questions on the founding and history of Amerind-speaking Native American populations. The model permits estimation of founding population sizes, changes in population size, time of population formation, and gene flow. Analyses of data from nine loci are consistent with the general portrait that has emerged from archaeological and other kinds of evidence. The estimated effective size of the founding population for the New World is fewer than 80 individuals, approximately 1% of the effective size of the estimated ancestral Asian population. By adding a splitting parameter to population divergence models it becomes possible to develop detailed portraits of human demographic history. Analyses of Asian and New World data support a model of a recent founding of the New World by a population of quite small effective size.     

Estimation of the ancestral effective population sizes of African great apes under different selection regimes

Reliable estimates of ancestral effective population sizes are necessary to unveil the population-level phenomena that shaped the phylogeny and molecular evolution of the African great apes. Although several methods have previously been applied to infer ancestral effective population sizes, an analysis of the influence of the selective regime on the estimates of ancestral demography has not been thoroughly conducted. In this study, three independent data sets under different selective regimes were used were composed to tackle this issue. The results showed that selection had a significant impact on the estimates of ancestral effective population sizes of the African great apes. The inference of the ancestral demography of African great apes was affected by the selection regime. The effects, however, were not homogeneous along the ancestral populations of great apes. The effective population size of the ancestor of humans and chimpanzees was more impacted by the selection regime when compared to the same parameter in the ancestor of humans, chimpanzees and gorillas. Because the selection regime influenced the estimates of ancestral effective population size, it is reasonable to assume that a portion of the discrepancy found in previous studies that inferred the ancestral effective population size may be attributable to the differential action of selection on the genes sampled.     

Genetic monitoring of wild and repatriated populations of endangered razorback sucker (Xyrauchen texanus, Catostomidae, Teleostei) in Lake Mohave, Arizona-Nevada

The Native Fishes Work Group, formed in 1991, developed and implemented a protocol to enhance the dwindling razorback sucker population in Lake Mohave, Arizona-Nevada. This large, genetically diverse population is severely reduced in size as a result of recruitment failure associated with predation on larvae. To circumvent this problem, wild larvae are captured, reared in protective custody until they are large enough to escape predation, and then released back into the lake. We present results of a monitoring program designed to assess the effectiveness of the sampling design in transmitting the high genetic diversity found in wild adults. Variation in a fragment from the mitochondrial DNA gene cytochrome b was examined by analysis of single-stranded polymorphisms and direct sequencing. Samples were characterized from three life history stages. Characterization of wild adults verified previous results that identified considerable diversity and provided baseline data. Samples of larvae from several temporal collections from throughout the spawning season and four geographical areas were characterized for 7 years (1997-2003) to assess the transmission of genetic variation from wild adults to larvae. Several analyses identified significant differences among temporal collections, resulting from sampling errors associated with finite number of females spawning at a given time and place. Comparisons among areas and years failed to identify significant variation, indicating that pooled collections for each year possess the same levels and patterns of genetic variation. Examination of repatriates representing 11 years (1992-2002) also failed to identify significant differences among cohorts; however, some sample sizes were small and the amova may lack sufficient power to detect differences. Contrasts of wild adults, larvae, and repatriates identified statistically significant differences among collections within these three groups; however, levels of variation are small and not biologically meaningful. More importantly, this analysis failed to detect significant differences among adults, larvae, and repatriates indicating that the program has been achieving its goal of transmitting variation from adults through the larvae and into the repatriate population. The reproductive capability of repatriates has not been examined, so it is unknown if the program will maintain genetic variation found in the original adult population. This will be most easily achieved by periodic monitoring of genetic variation in larval samples. If levels of variation become reduced in repatriates, levels and patterns of diversity in larvae are also expected to become reduced, and deviations in estimates of genetic diversity may become larger and more frequent. If this is the case, intervention may be necessary to ensure that certain individuals are not over-represented in the repatriate population.