The free energy method and the Wright–Fisher model with 2 alleles

We systematically investigate the Wright–Fisher model of population genetics with the free energy functional formalism of statistical mechanics and in the light of recent mathematical work on the connection between Fokker–Planck equations and free energy functionals. In statistical physics, entropy increases, or equivalently, free energy decreases, and the asymptotic state is given by a Gibbs-type distribution. This also works for the Wright–Fisher model when rewritten in divergence to identify the correct free energy functional. We not only recover the known results about the stationary distribution, that is, the asymptotic equilibrium state of the model, in the presence of positive mutation rates and possibly also selection, but can also provide detailed formulae for the rate of convergence towards that stationary distribution. In the present paper, the method is illustrated for the simplest case only, that of two alleles.     

Evolution of linkage disequilibrium of the founders in exponentially growing populations

Evolution of linkage disequilibrium of the founders in exponentially growing populations was studied using a time-inhomogeneous It? process model. The model is an extension of the diffusion approximation of the Wright-Fisher model. As a measure of linkage disequilibrium, the squared standard linkage deviation, which is defined by a ratio of the moments, was considered. A system of ordinary differential equations that these moments obey was obtained. This system can be solved numerically. By simulations, it was shown that the squared standard linkage deviation gives a good approximation of the expectation of the squared correlation coefficient of gamete frequencies. In addition, a perturbative solution was obtained when the growth rate is not large. By using the perturbation, an asymptotic formula for the squared standard linkage deviation after a large number of generations was obtained. According to the formula, the squared standard linkage deviation tends to be 1/(4Nc), where N is the current size of the population and c is the recombination fraction between two loci. It is dependent on neither the initial effective size, the growth rate, nor the mutation rate. In exponentially growing populations, linkage disequilibrium will be asymptotically the same as that in a constant size population, the effective size of which is the current effective size.     

Statistical Mechanics and the Evolution of Polygenic Quantitative Traits

The evolution of quantitative characters depends on the frequencies of the alleles involved, yet these frequencies cannot usually be measured. Previous groups have proposed an approximation to the dynamics of quantitative traits, based on an analogy with statistical mechanics. We present a modified version of that approach, which makes the analogy more precise and applies quite generally to describe the evolution of allele frequencies. We calculate explicitly how the macroscopic quantities (i.e., quantities that depend on the quantitative trait) depend on evolutionary forces, in a way that is independent of the microscopic details. We first show that the stationary distribution of allele frequencies under drift, selection, and mutation maximizes a certain measure of entropy, subject to constraints on the expectation of observable quantities. We then approximate the dynamical changes in these expectations, assuming that the distribution of allele frequencies always maximizes entropy, conditional on the expected values. When applied to directional selection on an additive trait, this gives a very good approximation to the evolution of the trait mean and the genetic variance, when the number of mutations per generation is sufficiently high (4Nμ > 1). We show how the method can be modified for small mutation rates (4Nμ → 0). We outline how this method describes epistatic interactions as, for example, with stabilizing selection.     

Exact moment calculations for genetic models with migration,mutation, and drift

Using properties of moment stationarity we develop exact expressions for the mean and covariance of allele frequencies at a single locus for a set of populations subject to drift, mutation, and migration. Some general results can be obtained even for arbitrary mutation and migration matrices, for example: (1) Under quite general conditions, the mean vector depends only on mutation rates, not on migration rates or the number of populations. (2) Allele frequencies covary among all pairs of populations connected by migration. As a result, the drift, mutation, migration process is not ergodic when any finite number of populations is exchanging genes. In addition, we provide closed-form expressions for the mean and covariance of allele frequencies in Wright's finite-island model of migration under several simple models of mutation, and we show that the correlation in allele frequencies among populations can be very large for realistic rates of mutation unless an enormous number of populations are exchanging genes. As a result, the traditional diffusion approximation provides a poor approximation of the stationary distribution of allele frequencies among populations. Finally, we discuss some implications of our results for measures of population structure based on Wright's F-statistics.     

The frequency of shifts between alternative equilibria

We derive a formula giving the frequency with which random drift shifts a population between alternative equilibria. This formula is valid when such shifts are rare (Ns much greater than 1), and applies over a wide range of mutation rates. When the number of mutations entering the population is low (4 N mu much less than 1), the rate of stochastic shifts reduces to the product of the mutation rate and the probability of fixation of a single mutation. However, when many mutations enter the population in each generation (4 N mu much greater than 1), the rate is higher than would be expected if mutations were established independently, and converges to that given by a gaussian approximation. We apply recent results on bistable systems to extend this formula to the general multidimensional case. This gives an explicit expression for the frequency of stochastic shifts, which depends only on the equilibrium probability distribution near the saddle point separating the alternative stable states. The plausibility of theories of speciation through random drift are discussed in the light of these results.     

Cumulant dynamics of a population under multiplicative selection, mutation, and drift

We revisit the classical population genetics model of a population evolving under multiplicative selection, mutation, and drift. The number of beneficial alleles in a multilocus system can be considered a trait under exponential selection. Equations of motion are derived for the cumulants of the trait distribution in the diffusion limit and under the assumption of linkage equilibrium. Because of the additive nature of cumulants, this reduces to the problem of determining equations of motion for the expected allele distribution cumulants at each locus. The cumulant equations form an infinite dimensional linear system and in an authored appendix Adam Prügel-Bennett provides a closed form expression for these equations. We derive approximate solutions which are shown to describe the dynamics well for a broad range of parameters. In particular, we introduce two approximate analytical solutions: (1) Perturbation theory is used to solve the dynamics for weak selection and arbitrary mutation rate. The resulting expansion for the system's eigenvalues reduces to the known diffusion theory results for the limiting cases with either mutation or selection absent. (2) For low mutation rates we observe a separation of time-scales between the slowest mode and the rest which allows us to develop an approximate analytical solution for the dominant slow mode. The solution is consistent with the perturbation theory result and provides a good approximation for much stronger selection intensities.     

Inferring complex DNA substitution processes on phylogenies using uniformization and data augmentation

A new method is developed for calculating sequence substitution probabilities using Markov chain Monte Carlo (MCMC) methods. The basic strategy is to use uniformization to transform the original continuous time Markov process into a Poisson substitution process and a discrete Markov chain of state transitions. An efficient MCMC algorithm for evaluating substitution probabilities by this approach using a continuous gamma distribution to model site-specific rates is outlined. The method is applied to the problem of inferring branch lengths and site-specific rates from nucleotide sequences under a general time-reversible (GTR) model and a computer program BYPASSR is developed. Simulations are used to examine the performance of the new program relative to an existing program BASEML that uses a discrete approximation for the gamma distributed prior on site-specific rates. It is found that BASEML and BYPASSR are in close agreement when inferring branch lengths, regardless of the number of rate categories used, but that BASEML tends to underestimate high site-specific substitution rates, and to overestimate intermediate rates, when fewer than 50 rate categories are used. Rate estimates obtained using BASEML agree more closely with those of BYPASSR as the number of rate categories increases. Analyses of the posterior distributions of site-specific rates from BYPASSR suggest that a large number of taxa are needed to obtain precise estimates of site-specific rates, especially when rates are very high or very low. The method is applied to analyze 45 sequences of the alpha 2B adrenergic receptor gene (A2AB) from a sample of eutherian taxa. In general, the pattern expected for regions under negative selection is observed with third codon positions having the highest inferred rates, followed by first codon positions and with second codon positions having the lowest inferred rates. Several sites show exceptionally high substitution rates at second codon positions that may represent the effects of positive selection.     

On the asymptotic trajectory of the roots of Lotka's equation.

A mathematical formula is given for the asymptotic trajectory of the complex roots of Lotka's equation. This formula is obtained by use of a Taylor expansion of the net fertility function in the neighborhood of the age beyond which fertility is zero. The approximate trajectory is compared with an exact trajectory obtained by use of a computational algorithm for finding complex roots suggested by Turner. For two examples, the agreement is surprisingly good.     

Optimal nonsingular control of fed-batch fermentation

Presented is a new simple method for multidimensional optimization of fed-batch fermentations based on the use of the orthogonal collocation technique. Considered is the problem of determination of optimal programs for fermentor temperature, substrate concentration in feed, feeding profile, and process duration. By reformulation of the state and control variables is obtained a nonsingular form of the optimization problem which has considerable advantage over the singular case since a complicated procedure for determination of switching times for feeding is avoided. The approximation of the state variables by Lagrange polynomials enables simple incorporation of split boundary conditions in the approximation, and the use of orthogonal collocations provides stability for integration of state and costate variables. The interpolation points are selected to obtain highest accuracy for approximation of the objective functional by the Radau-Lobatto formula. The control variables are determined by optimization of the Hamiltonian at the collocation points with the DFP method. Constraints are imposed on state and control variables.The method is applied for a homogeneous model of fermentation with volume, substrate, biomass, and product concentrations as the state variables. Computer study shows considerable simplicity of the method, its high accuracy for low order of approximation, and efficient convergence.     

Neutral evolution in spatially continuous populations

We introduce a general recursion for the probability of identity in state of two individuals sampled from a population subject to mutation, migration, and random drift in a two-dimensional continuum. The recursion allows for the interactions induced by density-dependent regulation of the population, which are inevitable in a continuous population. We give explicit series expansions for large neighbourhood size and for low mutation rates respectively and investigate the accuracy of the classical Malécot formula for these general models. When neighbourhood size is small, this formula does not give the identity even over large scales. However, for large neighbourhood size, it is an accurate approximation which summarises the local population structure in terms of three quantities: the effective dispersal rate, sigma(e); the effective population density, rho(e); and a local scale, kappa, at which local interactions become significant. The results are illustrated by simulations.     

Allele frequencies in a multidimensional Wright-Fisher model with general mutation

A formula is obtained for the probability that two genes at a single locus, sampled at random from a population at time t, are of particular types. The model assumed is a diffusion approximation to a neutral Wright-Fisher model in which mutation is general and not necessarily symmetric. An example is given of a population in which one allele has a high mutation rate, and the others have an equal, low mutation rate. The matrix Q, with elements given by the probability of sampling two alleles of particular types, is calculated exactly and approximately for this case. A formula is given for the distribution of the number of segregating sites occurring in two randomly sampled finite sequences of completely linked sites, with general mutation at a site and identical mutation structure between sites.     

A note on the accuracy of PAC-likelihood inference with microsatellite data

Stephens and Donnelly have introduced a simple yet powerful importance sampling scheme for computing the likelihood in population genetic models. Fundamental to the method is an approximation to the conditional probability of the allelic type of an additional gene, given those currently in the sample. As noted by Li and Stephens, the product of these conditional probabilities for a sequence of draws that gives the frequency of allelic types in a sample is an approximation to the likelihood, and can be used directly in inference. The aim of this note is to demonstrate the high level of accuracy of "product of approximate conditionals" (PAC) likelihood when used with microsatellite data. Results obtained on simulated microsatellite data show that this strategy leads to a negligible bias over a wide range of the scaled mutation parameter theta. Furthermore, the sampling variance of likelihood estimates as well as the computation time are lower than that obtained with importance sampling on the whole range of theta. It follows that this approach represents an efficient substitute to IS algorithms in computer intensive (e.g. MCMC) inference methods in population genetics.     

EVALUATING A SIMPLE APPROXIMATION TO MODELING THE JOINT EVOLUTION OF SELF‐FERTILIZATION AND INBREEDING DEPRESSION

A comprehensive understanding of plant mating system evolution requires detailed genetic models for both the mating system and inbreeding depression, which are often intractable. A simple approximation assuming that the mating system evolves by small infrequent mutational steps has been proposed. We examine its accuracy by comparing the evolutionarily stable selfing rates it predicts to those obtained from an explicit genetic model of the selfing rate, when inbreeding depression is caused by partly recessive deleterious mutations at many loci. Both models also include pollen limitation and pollen discounting. The approximation produces reasonably accurate predictions with a low or moderate genomic mutation rate to deleterious alleles, on the order of U = 0.02–0.2. However, for high mutation rates, the predictions of the full genetic model differ substantially from those of the approximation, especially with nearly recessive lethal alleles. This occurs because when a modifier allele affecting the selfing rate is rare, homozygous modifiers are produced mainly by selfing, which enhances the opportunity for purging nearly recessive lethals and increases the marginal fitness of the allele modifying the selfing rate. Our results confirm that explicit genetic models of selfing rate and inbreeding depression are required to understand mating system evolution.     

A stochastic theory for rate genes in large populations with overlapping generations

In this paper there is developed a stochastic theory for rare and nonrecessive genes in large populations that may have individuals of several age groups present at one time. The analysis is based on an age-dependent branching process due to Goodman. An approximate formula for the probability of extinction of a line of mutant genes, originating in an ancestral heterozygote in age group 0, is calculated. Expressions are also given for the asymptotic rates of approach of the probabilities of extinction of lines at finite times to their limiting values. These expressions apply regardless of the age of the ancestral heterozygote or whether the line has a positive probability of surviving indefinitely. Mean frequencies at equilibrium are calculated when there is recurrent mutation to an unfavorable gene.     

Statistical Methods for Analyzing Drosophila Germline Mutation Rates

Most studies of mutation rates implicitly assume that they remain constant throughout development of the germline. However, researchers recently used a novel statistical framework to reveal that mutation rates differ dramatically during sperm development in Drosophila melanogaster. Here a general framework is described for the inference of germline mutation patterns, generated from either mutation screening experiments or DNA sequence polymorphism data, that enables analysis of more than two mutations per family. The inference is made more rigorous and flexible by providing a better approximation of the probabilities of patterns of mutations and an improved coalescent algorithm within a single host with realistic assumptions. The properties of the inference framework, both the estimation and the hypothesis testing, were investigated by simulation. The refined inference framework is shown to provide (1) nearly unbiased maximum-likelihood estimates of mutation rates and (2) robust hypothesis testing using the standard asymptotic distribution of the likelihood-ratio tests. It is readily applicable to data sets in which multiple mutations in the same family are common.     

Confidence intervals for spectral mean and ratio statistics

We propose a new method, to construct confidence intervals forspectral mean and related ratio statistics of a stationary process,that avoids direct estimation of their asymptotic variances.By introducing a bandwidth, a self-normalization procedure isadopted and the distribution of the new statistic is asymptoticallynuisance-parameter free. The bandwidth is chosen using informationcriteria and a moving average sieve approximation. Through asimulation study, we demonstrate good finite sample performanceof our method when the sample size is moderate, while a comparisonwith an empirical likelihood-based method for ratio statisticsis made, confirming a wider applicability of our method.     

An analysis of neutral-alleles and variable-environment diffusion models

Three diffusion models are formulated for the evolution of a diploid population with K alleles at one locus with completely symmetric mutation and random genetic drift, a variable-environment, and all the above mechanisms. For the diallelic case, the transient behavior is studied by solving the corresponding diffusion equations by an asymptotic method valid for short time intervals. The transient behavior of the three models is compared for the case when their stationary distributions are identical. The expected amount of heterozygosity is computed using the asymptotic solution and is compared to an exact result. The asymptotic results are extended to the general case with K alleles at the locus for the symmetric mutation and variable-environment models.Research supported by the National Science Foundation under Grant MCS 79-01718     

How small are small mutation rates?

We consider evolutionary game dynamics in a finite population of size N. When mutations are rare, the population is monomorphic most of the time. Occasionally a mutation arises. It can either reach fixation or go extinct. The evolutionary dynamics of the process under small mutation rates can be approximated by an embedded Markov chain on the pure states. Here we analyze how small the mutation rate should be to make the embedded Markov chain a good approximation by calculating the difference between the real stationary distribution and the approximated one. While for a coexistence game, where the best reply to any strategy is the opposite strategy, it is necessary that the mutation rate μ is less than N ^−1/2exp[−N] to ensure that the approximation is good, for all other games, it is sufficient if the mutation rate is smaller than (N ln N)⁻¹. Our results also hold for a wide class of imitation processes under arbitrary selection intensity.