Search

In this article, discrete and stochastic changes in (effective) population size are incorporated into the spectral representation of a biallelic diffusion process for drift and small mutation rates. A forward algorithm inspired by Hidden-Markov-Model (HMM) literature is used to compute exact sample allele frequency spectra for three demographic scenarios: single changes in (effective) population size, boom-bust dynamics, and stochastic fluctuations in (effective) population size. An approach for fully agnostic demographic inference from these sample allele spectra is explored, and sufficient statistics for step-wise changes in population size are found. Further, convergence behaviours of the polymorphic sample spectra for population size changes on different time scales are examined and discussed within the context of inference of the effective population size. Joint visual assessment of the sample spectra and the temporal coefficients of the spectral decomposition of the forward diffusion process is found to be important in determining departure from equilibrium. Stochastic changes in (effective) population size are shown to shape sample spectra particularly strongly.

Inference of the marginal likelihood of sample allele configurations using backward algorithms yields identical results with the Kingman coalescent, the Moran model, and the diffusion model (up to a scaling of time). For inference of probabilities of ancestral population allele frequencies at any given point in the past - either of discrete ancestral allele configurations as in the coalescent, or of ancestral allele proportions as in the backward diffusion - backward approaches need to be combined with corresponding forward ones. This is done in so-called forward-backward algorithms. In this article, we utilize orthogonal polynomials in forward-backward algorithms. They enable efficient calculation of past allele configurations of an extant sample and probabilities of ancestral population allele frequencies in equilibrium and in non-equilibrium. We show that the genealogy of a sample is fully described by the backward polynomial expansion of the marginal likelihood of its allele configuration.

In population genetics, extant samples are usually used for inference of past population genetic forces. With the Kingman coalescent and the backward diffusion equation, inference of the marginal likelihood proceeds from an extant sample backward in time. Conditional on an extant sample, the Moran model can also be used backward in time with identical results, up to a scaling of time. In particular, all three approaches -- the coalescent, the backward diffusion, and the Moran model -- lead to the identical marginal likelihood of the sample. If probabilities of ancestral states are also inferred, either of discrete ancestral allele particle configurations, as in the coalescent, or of ancestral population allele proportions, as in the backward diffusion, the backward algorithm needs to be combined with the corresponding forward algorithm to the forward-backward algorithm. Generally orthogonal polynomials, solving the diffusion equation, are numerically simpler than the other approaches: they implicitly sum over many intermediate ancestral particle configurations; furthermore, while the Moran model requires iterative matrix multiplication with a transition matrix of a dimension of the population size squared, expansion of the polynomials is only necessary up to the sample size. For discrete samples, forward-in-time moving pure birth processes similar to the Polya- or Hoppe-urn models complement the backward-looking coalescent. Because, the sample size is a random variable forward in time, pure-birth processes are unsuited to model population demography given extant samples. With orthogonal polynomials, however, not only ancestral allele proportions but also probabilities of ancestral particle configurations can be calculated easily. Assuming only mutation and drift, the use of orthogonal polynomials is numerically advantageous over alternative strategies., Comment: 26 pages, 2 figures

In this article, a biallelic reversible mutation model with linear and quadratic selection is analyzed. The approach reconnects to one proposed by Kimura ( Possibility of extensive neutral evolution under stabilizing selection with special reference to nonrandom use of codons (PNAS,1981)), who starts from a diffusion model and derives its equilibrium distribution up to a constant. We use a boundary-mutation Moran model, which approximates a general mutation model for small effective mutation rates, and derive its equilibrium distribution for polymorphic and monomorphic variants in small to moderately sized populations. Using this model, we show that biased mutation rates and linear selection alone can cause patterns of polymorphism rates within and substitution rates between populations that are usually ascribed to balancing or overdominant selection. We illustrate this using a data set of short introns and fourfold degenerate sites from Drosophila simulans and Drosophila melanogaster.

The stationary sampling distribution of a neutral decoupled Moran or Wright-Fisher diffusion with neutral mutations is known to first order for a general rate matrix with small but otherwise unconstrained mutation rates. Using this distribution as a starting point we derive results for maximum likelihood estimates of scaled mutation rates from site frequency data under three model assumptions: a twelve-parameter general rate matrix, a nine-parameter reversible rate matrix, and a six-parameter strand-symmetric rate matrix. The site frequency spectrum is assumed to be sampled from a fixed size population in equilibrium, and to consist of allele frequency data at a large number of unlinked sites evolving with a common mutation rate matrix without selective bias. We correct an error in a previous treatment of the same problem (Burden and Tang, 2017) affecting the estimators for the general and strand-symmetric rate matrices. The method is applied to a biological dataset consisting of a site frequency spectrum extracted from short autosomal introns in a sample of Drosophila melanogaster individuals., Comment: 39 pages, 4 figures, simulation to test accuracy of the model added

With the multiallelic parent-independent mutation-drift model, the equilibrium proportions of alleles are known to be Dirichlet distributed. A special case is the biallelic model, in which the proportions are beta distributed. A sample taken from these models is then Dirichlet-multinomially or beta-binomially distributed, respectively. Maximum likelihood (ML) estimators for the mutation parameters of the biallelic parent-independent mutation model are available via an expectation maximization algorithm. Assuming small scaled mutation rates, the distribution of a sample of size $M$ can be expanded in a Taylor series of first order. Then the ML estimators for the two parameters in the biallelic model can be expressed using the site frequency spectrum. In this article, we go beyond parent-independent mutation and analyse a strand-symmetric mutation model with six scaled mutation parameters that deviates from parent independent mutation and, generally, from detailed balance. We derive ML estimators for these six parameters assuming mutation-drift equilibrium and small scaled mutation rates. This is the first time that ML estimators are provided for a mutation model more complex than parent-independent mutation.

In a review on the contribution of J.B.S. Haldane to the development of the Bayes factor hypothesis test (arXiv:1511.08180), Etz and Wagenmakers focus on Haldane's proposition of a mixture prior in a genetic example (Haldane 1932, A note on inverse probability. Mathematical Proceedings of the Cambridge Philosophical Society, 28, 55-61.). As Haldane never followed up on these ideas, it is difficult to gauge his motivation and intentions. I argue that contrary to Haldane's stated intention of replacing flat priors with more reasonable assumptions, he actually chose in this example an unreasonable flat prior. Considering the information available to Haldane, I derive a superior prior and compare to Haldane's flat prior. Haldane's main intent with his article seems to have been to explore the different parameter regions of the binomial and the conjugate beta. Furthermore, I agree with Etz and Wagenmakers that Haldane serendipitously adopted a mixture prior comprising a point mass and smooth distribution in his genetic example.

The evolution of the allelic proportion $x$ of a biallelic locus subject to the forces of mutation and drift is investigated in a diffusion model, assuming small scaled mutation rates. The overall scaled mutation rate is parametrized with $\theta=(\mu_1+\mu_0)N$ and the ratio of mutation rates with $\alpha=\mu_1/(\mu_1+\mu_0)=1-\beta$. The equilibrium density of this process is beta with parameters $\alpha\theta$ and $\beta\theta$. Away from equilibrium, the transition density can be expanded into a series of modified Jacobi polynomials. If the scaled mutation rates are small, i.e., $\theta \ll 1$, it may be assumed that polymorphism derives from mutations at the boundaries. A model, where the interior dynamics conform to the pure drift diffusion model and the mutations are entering from the boundaries is derived. In equilibrium, the density of the proportion of polymorphic alleles, \ie\ $x$ within the polymorphic region $[1/N,1-1/N]$, is $\alpha\beta\theta(\tfrac1x+\tfrac1{1-x})=\tfrac{\alpha\beta\theta}{x(1-x)}$, while the mutation bias $\alpha$ influences the proportion of monomorphic alleles at 0 and 1. Analogous to the expansion with modified Jacobi polynomials, a series expansion of the transition density is derived, which is connected to Kimura's well known solution of the pure drift model using Gegenbauer polynomials. Two temporal and two spatial regions are separated. The eigenvectors representing the spatial component within the polymorphic region depend neither on the on the scaled mutation rate $\theta$ nor on the mutation bias $\alpha$. Therefore parameter changes, e.g., growing or shrinking populations or changes in the mutation bias, can be modeled relatively easily, without the change of the eigenfunctions necessary for the series expansion with Jacobi polynomials.