Using maximum likelihood to estimate population size from temporal changes in allele frequencies.

We develop a maximum-likelihood framework for using temporal changes in allele frequencies to estimate the number of breeding individuals in a population. We use simulations to compare the performance of this estimator to an F-statistic estimator of variance effective population size. The maximum-likelihood estimator had a lower variance and smaller bias. Taking advantage of the likelihood framework, we extend the model to include exponential growth and show that temporal allele frequency data from three or more sampling events can be used to test for population growth.     

A Coalescent Estimator of the Population Recombination Rate

Population genetic models often use a population recombination parameter 4Nc, where N is the effective population size and c is the recombination rate per generation. In many ways 4Nc is comparable to 4Nu, the population mutation rate. Both combine genome level and population level processes, and together they describe the rate of production of genetic variation in a population. However, 4Nc is more difficult to estimate. For a population sample of DNA sequences, historical recombination can only be detected if polymorphisms exist, and even then most recombination events are not detectable. This paper describes an estimator of 4Nc, hereafter designated γ (gamma), that was developed using a coalescent model for a sample of four DNA sequences with recombination. The reliability of γ was assessed using multiple coalescent simulations. In general γ has low to moderate bias, and the reliability of γ is comparable, though less, than that for a widely used estimator of 4Nu. If there exists an independent estimate of the recombination rate (per generation, per base pair), γ can be used to estimate the effective population size or the neutral mutation rate.     

The inference of stepwise changes in substitution rates using serial sequence samples.

It is frequently true that molecular sequences do not evolve in a strictly clocklike manner. Instead, substitution rate may vary for a number of reasons, including changes in selection pressure and effective population size, as well as changes in mean generation time. Here we present two new methods for estimating stepwise changes in substitution rates when serially sampled molecular sequences are available. These methods are based on multiple rates with dated tips (MRDT) models and allow different rates to be estimated for different intervals of time. These intervals may correspond to the sampling intervals or to a priori--defined intervals that are not coincident with the times the serial samples are obtained. Two methods for obtaining estimates of multiple rates are described. The first is an extension of the phylogeny-based maximum-likelihood estimation procedure introduced by Rambaut. The second is a new parameterization of the pairwise distance least-squares procedure used by Drummond and Rodrigo. The utility of these methods is demonstrated on a genealogy of HIV sequences obtained at five different sampling times from a single patient over a period of 34 months.     

Estimating mutation rate and generation time from longitudinal samples of DNA sequences

We present in this paper a simple method for estimating the mutation rate per site per year which also yields an estimate of the length of a generation when mutation rate per site per generation is known. The estimator, which takes advantage of DNA polymorphisms in longitudinal samples, is unbiased under a number of population models, including population structure and variable population size over time. We apply the new method to a longitudinal sample of DNA sequences of the env gene of human immunodeficiency virus type 1 (HIV-1) from a single patient and obtain 1.62 x 10(-2) as the mutation rate per site per year for HIV-1. Using an independent data set to estimate the mutation rate per generation, we obtain 1.8 days as the length of a generation of HIV-1, which agrees well with recent estimates based on viral load data. Our estimate of generation time differs considerably from a recent estimate by Rodrigo et al. when the same mutation rate per site per generation is used. Some factors that may contribute to the difference among different estimators are discussed.     

Evolution of intrahost HIV-1 genetic diversity during chronic infection

HIV-1 is one of the fastest evolving entities known. Given that census population sizes of HIV-1 within individuals are much greater than the inverse mutation rate, every possible single point mutation in the viral genome occurs each generation. This enormous capability to generate genetic variation allows for escape from immune surveillance and antiviral therapy. However, compared to this potential, populations of HIV-1 within individuals exhibit little genetic variation. This discrepancy between the known mutation rate of HIV-1 and the average level of genetic variation in the env gene observed in vivo is reflected in comparisons of the actual numbers of productively infected cells, estimated as 10(7), and the effective population size, estimated as 10(3). Using approximate Bayesian computation, we evaluated several hypotheses based on a variety of selective and demographic processes to explain the low effective population size of HIV-1. Of the models we examined, the metapopulation model, in which HIV-1 evolves within an individual as a large collection of small subpopulations subject to frequent migration, extinction, and recolonization, was most consistent with the observed levels of genetic variation and the average frequencies of those variants. The metapopulation model links previous studies of viral dynamics and population genetics.     

Coalescent Theory for a Partially Selfing Population

A coalescent theory for a sample of DNA sequences from a partially selfing diploid population and an algorithm for simulating such samples are developed in this article. Approximate formulas are given for the expectation and the variance of the number of segregating sites in a sample of k sequences from n individuals. Several new estimators of the important parameters θ = 4Nμ and the selfing rate s, where N and μ are, respectively, the effective population size and the mutation rate per sequence per generation, are proposed and their sampling properties are studied.     

A robust measure of HIV-1 population turnover within chronically infected individuals

A simple nonparameteric test for population structure was applied to temporally spaced samples of HIV-1 sequences from the gag-pol region within two chronically infected individuals. The results show that temporal structure can be detected for samples separated by about 22 months or more. The performance of the method, which was originally proposed to detect geographic structure, was tested for temporally spaced samples using neutral coalescent simulations. Simulations showed that the method is robust to variation in samples sizes and mutation rates, to the presence/absence of recombination, and that the power to detect temporal structure is high. By comparing levels of temporal structure in simulations to the levels observed in real data, we estimate the effective intra-individual population size of HIV-1 to be between 10(3) and 10(4) viruses, which is in agreement with some previous estimates. Using this estimate and a simple measure of sequence diversity, we estimate an effective neutral mutation rate of about 5 x 10(-6) per site per generation in the gag-pol region. The definition and interpretation of estimates of such "effective" population parameters are discussed.     

A method for accurate inference of population size from serially sampled genealogies distorted by selection

The serial coalescent extends traditional coalescent theory to include genealogies in which not all individuals were sampled at the same time. Inference in this framework is powerful because population size and evolutionary rate may be estimated independently. However, when the sequences in question are affected by selection acting at many sites, the genealogies may differ significantly from their neutral expectation, and inference of demographic parameters may become inaccurate. I demonstrate that this inaccuracy is severe when the mutation rate and strength of selection are jointly large, and I develop a new likelihood calculation that, while approximate, improves the accuracy of population size estimates. When used in a Bayesian parameter estimation context, the new calculation allows for estimation of the shape of the pairwise coalescent rate function and can be used to detect the presence of selection acting at many sites in a sequence. Using the new method, I investigate two sets of dengue virus sequences from Puerto Rico and Thailand, and show that both genealogies are likely to have been distorted by selection.     

Estimating the Timing of Mother-to-Child Transmission of the Human Immunodeficiency Virus Type 1 Using a Viral Molecular Evolution Model

Background

Mother-to-child transmission (MTCT) is responsible for most pediatric HIV-1 infections worldwide. It can occur during pregnancy, labor, or breastfeeding. Numerous studies have used coalescent and molecular clock methods to understand the epidemic history of HIV-1, but the timing of vertical transmission has not been studied using these methods. Taking advantage of the constant accumulation of HIV genetic variation over time and using longitudinally sampled viral sequences, we used a coalescent approach to investigate the timing of MTCT.

Materials and Methods

Six-hundred and twenty-two clonal env sequences from the RNA and DNA viral population were longitudinally sampled from nine HIV-1 infected mother-and-child pairs [range: 277–1034 days]. For each transmission pair, timing of MTCT was determined using a coalescent-based model within a Bayesian statistical framework. Results were compared with available estimates of MTCT timing obtained with the classic biomedical approach based on serial HIV DNA detection by PCR assays.

Results

Four children were infected during pregnancy, whereas the remaining five children were infected at time of delivery. For eight out of nine pairs, results were consistent with the transmission periods assessed by standard PCR-based assay. The discordance in the remaining case was likely confused by co-infection, with simultaneous introduction of multiple maternal viral variants at the time of delivery.

Conclusions

The study provided the opportunity to validate the Bayesian coalescent approach that determines the timing of MTCT of HIV-1. It illustrates the power of population genetics approaches to reliably estimate the timing of transmission events and deepens our knowledge about the dynamics of viral evolution in HIV-infected children, accounting for the complexity of multiple transmission events.     

Estimating Effective Population Size or Mutation Rate Using the Frequencies of Mutations of Various Classes in a Sample of DNA Sequences

Mutations resulting in segregating sites of a sample of DNA sequences can be classified by size and type and the frequencies of mutations of different sizes and types can be inferred from the sample. A framework for estimating the essential parameter θ = 4Nu utilizing the frequencies of mutations of various sizes and types is developed in this paper, where N is the effective size of a population and μ is mutation rate per sequence per generation. The framework is a combination of coalescent theory, general linear model and Monte-Carlo integration, which leads to two new estimators θ(ξ) and θ(η) as well as a general Watterson''s estimator θ(K) and a general Tajima''s estimator θ(π). The greatest strength of the framework is that it can be used under a variety of population models. The properties of the framework and the four estimators θ(K), θ(π), θ(ξ) and θ(η) are investigated under three important population models: the neutral Wright-Fisher model, the neutral model with recombination and the neutral Wright''s finite-islands model. Under all these models, it is shown that θ(ξ) is the best estimator among the four even when recombination rate or migration rate has to be estimated. Under the neutral Wright-Fisher model, it is shown that the new estimator θ(ξ) has a variance close to a lower bound of variances of all unbiased estimators of θ which suggests that θ(ξ) is a very efficient estimator.     

The effective number of a population that varies cyclically in size. I. Discrete generations

We consider a dioecious population having numbers of males and females that vary over time in cycles of length k. It is shown that if k is small in comparison with the numbers of males and females in any generation of the cycle, the effective population number (or size), N(e), is approximately equal to the harmonic mean of the effective population sizes during any given cycle. This result holds whether the locus under consideration is autosomal or sex-linked and whether inbreeding effective population numbers or variance effective population numbers are involved in the calculation of N(e). If, however, only two successive generations in the cycle are considered and the population changes in size between these generations, the inbreeding effective population number, N(eI), differs from the variance effective population number, N(eV). The mutation effective population number turns out to be the same as the number derived using calculations involving probabilities of identity by descent. It is also shown that, at least in one special case, the eigenvalue effective population number is the same as N(eV).     

Evolution of the human immunodeficiency virus envelope gene is dominated by purifying selection

The evolution of the human immunodeficiency virus (HIV-1) during chronic infection involves the rapid, continuous turnover of genetic diversity. However, the role of natural selection, relative to random genetic drift, in governing this process is unclear. We tested a stochastic model of genetic drift using partial envelope sequences sampled longitudinally in 28 infected children. In each case the Bayesian posterior (empirical) distribution of coalescent genealogies was estimated using Markov chain Monte Carlo methods. Posterior predictive simulation was then used to generate a null distribution of genealogies assuming neutrality, with the null and empirical distributions compared using four genealogy-based summary statistics sensitive to nonneutral evolution. Because both null and empirical distributions were generated within a coalescent framework, we were able to explicitly account for the confounding influence of demography. From the distribution of corrected P-values across patients, we conclude that empirical genealogies are more asymmetric than expected if evolution is driven by mutation and genetic drift only, with an excess of low-frequency polymorphisms in the population. This indicates that although drift may still play an important role, natural selection has a strong influence on the evolution of HIV-1 envelope. A negative relationship between effective population size and substitution rate indicates that as the efficacy of selection increases, a smaller proportion of mutations approach fixation in the population. This suggests the presence of deleterious mutations. We therefore conclude that intrahost HIV-1 evolution in envelope is dominated by purifying selection against low-frequency deleterious mutations that do not reach fixation.     

Influence of random genetic drift on human immunodeficiency virus type 1 env evolution during chronic infection

Human immunodeficiency virus type 1 (HIV-1) has high replication and mutation rates that generate large census populations and high levels of genetic variation. We examined the roles of natural selection, population growth, random genetic drift, and recombination in shaping the variation in 1509 C2-V5 env sequences derived from nine men with chronic HIV-1 infection. These sequences were obtained from clinical visits that reflect the first 6-13.7 years of infection. Pairwise comparisons of nonsynonymous and synonymous distances, Tajima's D test, Fu and Li's D* test, and a test of recurrent mutation revealed evidence for episodes of nonneutral evolution in a total of 22 out of 145 blood samples, representing six of the nine individuals. Using three coalescent-based maximum-likelihood estimators, we found viral effective population sizes in all nine individuals to be approximately 10(3). We also show that a previous estimate of the effective population size of approximately 10(5) based on rare haplotype frequencies decreases to approximately 10(3) upon correcting a biased sampling procedure. We conclude that the genetic variation in these data sets can be explained by a predominance of random genetic drift of neutral mutations with brief episodes of natural selection that were frequently masked by recombination.     

Maximum-likelihood estimation of migration rates and effective population numbers in two populations using a coalescent approach.

A new method for the estimation of migration rates and effective population sizes is described. It uses a maximum-likelihood framework based on coalescence theory. The parameters are estimated by Metropolis-Hastings importance sampling. In a two-population model this method estimates four parameters: the effective population size and the immigration rate for each population relative to the mutation rate. Summarizing over loci can be done by assuming either that the mutation rate is the same for all loci or that the mutation rates are gamma distributed among loci but the same for all sites of a locus. The estimates are as good as or better than those from an optimized FST-based measure. The program is available on the World Wide Web at http://evolution.genetics. washington.edu/lamarc.html/.     

Estimating effective population size from samples of sequences: a bootstrap Monte Carlo integration method.

We would like to use maximum likelihood to estimate parameters such as the effective population size N(e) or, if we do not know mutation rates, the product 4N(e) mu of mutation rate per site and effective population size. To compute the likelihood for a sample of unrecombined nucleotide sequences taken from a random-mating population it is necessary to sum over all genealogies that could have led to the sequences, computing for each one the probability that it would have yielded the sequences, and weighting each one by its prior probability. The genealogies vary in tree topology and in branch lengths. Although the likelihood and the prior are straightforward to compute, the summation over all genealogies seems at first sight hopelessly difficult. This paper reports that it is possible to carry out a Monte Carlo integration to evaluate the likelihoods approximately. The method uses bootstrap sampling of sites to create data sets for each of which a maximum likelihood tree is estimated. The resulting trees are assumed to be sampled from a distribution whose height is proportional to the likelihood surface for the full data. That it will be so is dependent on a theorem which is not proven, but seems likely to be true if the sequences are not short. One can use the resulting estimated likelihood curve to make a maximum likelihood estimate of the parameter of interest, N(e) or of 4N(e) mu. The method requires at least 100 times the computational effort required for estimation of a phylogeny by maximum likelihood, but is practical on today's work stations. The method does not at present have any way of dealing with recombination.     

A separation-of-timescales approach to the coalescent in a continuous population

This article presents an analysis of a model of isolation by distance in a continuous, two-dimensional habitat. An approximate expression is derived for the distribution of coalescence times for a pair of sequences sampled from specific locations in a rectangular habitat. Results are qualitatively similar to previous analyses of isolation by distance, but account explicitly for the location of samples relative to the habitat boundaries. A separation-of-timescales approach takes advantage of the fact that the sampling locations affect only the recent coalescent behavior. When the population size is larger than the number of generations required for a lineage to cross the habitat range, the long-term genealogical process is reasonably well described by Kingman's coalescent with time rescaled by the effective population size. This long-term effective population size is affected by the local dispersal behavior as well as the geometry of the habitat. When the population size is smaller than the time required to cross the habitat, deep branches in the genealogy are longer than would be expected under the standard neutral coalescent, similar to the pattern expected for a panmictic population whose population size was larger in the past.     

Stochastic Simulations Suggest that HIV-1 Survives Close to Its Error Threshold

The use of mutagenic drugs to drive HIV-1 past its error threshold presents a novel intervention strategy, as suggested by the quasispecies theory, that may be less susceptible to failure via viral mutation-induced emergence of drug resistance than current strategies. The error threshold of HIV-1, , however, is not known. Application of the quasispecies theory to determine poses significant challenges: Whereas the quasispecies theory considers the asexual reproduction of an infinitely large population of haploid individuals, HIV-1 is diploid, undergoes recombination, and is estimated to have a small effective population size in vivo. We performed population genetics-based stochastic simulations of the within-host evolution of HIV-1 and estimated the structure of the HIV-1 quasispecies and . We found that with small mutation rates, the quasispecies was dominated by genomes with few mutations. Upon increasing the mutation rate, a sharp error catastrophe occurred where the quasispecies became delocalized in sequence space. Using parameter values that quantitatively captured data of viral diversification in HIV-1 patients, we estimated to be substitutions/site/replication, ∼2–6 fold higher than the natural mutation rate of HIV-1, suggesting that HIV-1 survives close to its error threshold and may be readily susceptible to mutagenic drugs. The latter estimate was weakly dependent on the within-host effective population size of HIV-1. With large population sizes and in the absence of recombination, our simulations converged to the quasispecies theory, bridging the gap between quasispecies theory and population genetics-based approaches to describing HIV-1 evolution. Further, increased with the recombination rate, rendering HIV-1 less susceptible to error catastrophe, thus elucidating an added benefit of recombination to HIV-1. Our estimate of may serve as a quantitative guideline for the use of mutagenic drugs against HIV-1.     

Very large long-term effective population size in the virulent human malaria parasite Plasmodium falciparum

It has been proposed that the virulent human malaria parasite Plasmodium falciparum underwent a recent severe population bottleneck. In order to test this hypothesis, we estimated the effective population size of this species from the patterns of nucleotide substitution at 23 nuclear protein-coding loci, using a variety of methods based on coalescent theory. Both simple methods and phylogenetically based maximum-likelihood methods yielded the conclusion that the effective population size of this species has been of the order of at least 10(5) for the past 300,000-400,000 years.     

MAXIMUM-LIKELIHOOD ESTIMATION OF POPULATION DIVERGENCE TIMES AND POPULATION PHYLOGENY IN MODELS WITHOUT MUTATION

In this paper we present a method for estimating population divergence times by maximum likelihood in models without mutation. The maximum-likelihood estimator is compared to a commonly applied estimator based on Wright's F_ST statistic. Simulations suggest that the maximum-likelihood estimator is less biased and has a lower variance than the F_ST-based estimator. The maximum-likelihood estimator provides a statistical framework for the analysis of population history given genetic data. We demonstrate how maximum-likelihood estimates of the branching pattern of divergence of multiple populations may be obtained. We also describe how the method may be applied to test hypotheses such as whether populations have maintained equal population sizes. We illustrate the method by applying it to two previously published sets of human restriction fragment length polymorphism (RFLP) data.     

On tree intensity estimation for forest inventories: Some statistical issues

The effect of the plot shape, number of subplots and their spatial arrangement on the sample variance for spatially explicit point populations is analysed for a simple intensity estimator. We derive the sample variance and covariance for sampling designs involving more than one subplot. Some numerical approximations are also presented. If a clustered point pattern has to be sampled, the best strategy to reduce the sample variance is to consider as many rectangular subplots as possible, for a prescribed total sample area, distributed over a grid. In contrast, if a regular point pattern is to be sampled, then a single circular subplot should be considered. If we assume that the point configuration is Poisson, then we can consider any subplot shape and spatial distribution ensuring no overlapping between the subplots. A case study in forestry is considered to assess the validity of our results.