Extinction of metastable stochastic populations

Metapopulation dynamics on the brink of extinction

We analyse metapopulation dynamics in terms of an individual-based, stochastic model of a finite metapopulation. We suggest a new approach, using the number of patches in the population as a large parameter. This approach does not require that the number of individuals per patch is large, neither is it necessary to assume a time-scale separation between local population dynamics and migration. Our approach makes it possible to accurately describe the dynamics of metapopulations consisting of many small patches. We focus on metapopulations on the brink of extinction. We estimate the time to extinction and describe the most likely path to extinction. We find that the logarithm of the time to extinction is proportional to the product of two vectors, a vector characterising the distribution of patch population sizes in the quasi-steady state, and a vector-related to Fisher's reproduction vector-that quantifies the sensitivity of the quasi-steady state distribution to demographic fluctuations. We compare our analytical results to stochastic simulations of the model, and discuss the range of validity of the analytical expressions. By identifying fast and slow degrees of freedom in the metapopulation dynamics, we show that the dynamics of large metapopulations close to extinction is approximately described by a deterministic equation originally proposed by Levins (1969). We were able to compute the rates in Levins' equation in terms of the parameters of our stochastic, individual-based model. It turns out, however, that the interpretation of the dynamical variable depends strongly on the intrinsic growth rate and carrying capacity of the patches. Only when the local growth rate and the carrying capacity are large does the slow variable correspond to the number of patches, as envisaged by Levins. Last but not least, we discuss how our findings relate to other, widely used metapopulation models.

Mixing times in evolutionary game dynamics

Without mutation and migration, evolutionary dynamics ultimately leads to the extinction of all but one species. Such fixation processes are well understood and can be characterized analytically with methods from statistical physics. However, many biological arguments focus on stationary distributions in a mutation-selection equilibrium. Here, we address the mixing time required to reach stationarity in the presence of mutation. We show that mixing times in evolutionary games have the opposite behavior from fixation times when the intensity of selection increases: in coordination games with bistabilities, the fixation time decreases, but the mixing time increases. In coexistence games with metastable states, the fixation time increases, but the mixing time decreases. Our results are based on simulations and the Wentzel-Kramers-Brillouin approximation of the master equation.

Dynamics of the chemical master equation, a strip of chains of equations in d-dimensional space

We investigate the multichain version of the chemical master equation, when there are transitions between different states inside the long chains, as well as transitions between (a few) different chains. In the discrete version, such a model can describe the connected diffusion processes with jumps between different types. We apply the Hamilton-Jacobi equation to solve some aspects of the model. We derive exact (in the limit of infinite number of particles) results for the dynamic of the maximum of the distribution and the variance of distribution.

Spontaneous clearance of viral infections by mesoscopic fluctuations

Spontaneous disease extinction can occur due to a rare stochastic fluctuation. We explore this process, both numerically and theoretically, in two minimal models of stochastic viral infection dynamics. We propose a method that reduces the complexity in models of viral infections so that the remaining dynamics can be studied by previously developed techniques for analyzing epidemiological models. Using this technique, we obtain an expression for the infection clearance time as a function of kinetic parameters. We apply our theoretical results to study stochastic infection clearance for specific stages of HIV and HCV dynamics. Our results show that the typical time for stochastic clearance of a viral infection increases exponentially with the size of the population, but infection still can be cleared spontaneously within a reasonable time interval in a certain population of cells. We also show that the clearance time is exponentially sensitive to the viral decay rate and viral infectivity but only linearly dependent on the lifetime of an infected cell. This suggests that if standard drug therapy fails to clear an infection then intensifying therapy by adding a drug that reduces the rate of cell infection rather than immune modulators that hasten infected cell death may be more useful in ultimately clearing remaining pockets of infection.

Stochastic formulation of ecological models and their applications

The increasing use of computer simulation by theoretical ecologists started a move away from models formulated at the population level towards individual-based models. However, many of the models studied at the individual level are not analysed mathematically and remain defined in terms of a computer algorithm. This is not surprising, given that they are intrinsically stochastic and require tools and techniques for their study that may be unfamiliar to ecologists. Here, we argue that the construction of ecological models at the individual level and their subsequent analysis is, in many cases, straightforward and leads to important insights. We discuss recent work that highlights the importance of stochastic effects for parameter ranges and systems where it was previously thought that such effects would be negligible.

Multiple extinction routes in stochastic population models

Isolated populations ultimately go extinct because of the intrinsic noise of elementary processes. In multipopulation systems extinction of a population may occur via more than one route. We investigate this generic situation in a simple predator-prey (or infected-susceptible) model. The predator and prey populations may coexist for a long time, but ultimately both go extinct. In the first extinction route the predators go extinct first, whereas the prey thrive for a long time and then also go extinct. In the second route the prey go extinct first, causing a rapid extinction of the predators. Assuming large subpopulation sizes in the coexistence state, we compare the probabilities of each of the two extinction routes and predict the most likely path of the subpopulations to extinction. We also suggest an effective three-state master equation for the probabilities to observe the coexistence state, the predator-free state, and the empty state.

Switching between phenotypes and population extinction

Many types of bacteria can survive under stress by switching stochastically between two different phenotypes: the "normals" who multiply fast, but are vulnerable to stress, and the "persisters" who hardly multiply, but are resilient to stress. Previous theoretical studies of such bacterial populations have focused on the fitness: the asymptotic rate of unbounded growth of the population. Yet for an isolated population of established (and not very large) size, a more relevant measure may be the population extinction risk due to the interplay of adverse extrinsic variations and intrinsic noise of birth, death and switching processes. Applying a WKB approximation to the pertinent master equation of such a two-population system, we quantify the extinction risk, and find the most likely path to extinction under both favorable and adverse conditions. Analytical results are obtained both in the biologically relevant regime when the switching is rare compared with the birth and death processes, and in the opposite regime of frequent switching. We show that rare switches are most beneficial in reducing the extinction risk.

Exact results in the large system size limit for the dynamics of the chemical master equation, a one dimensional chain of equations

We apply the Hamilton-Jacobi equation (HJE) formalism to solve the dynamics of the chemical master equation (CME). We found exact analytical expressions (in large system-size limit) for the probability distribution, including explicit expression for the dynamics of variance of distribution. We also give the solution for some simple cases of the model with time-dependent rates. We derived the results of the Van Kampen method from the HJE approach using a special ansatz. Using the Van Kampen method, we give a system of ordinary differential equations (ODEs) to define the variance in a two-dimensional case. We performed numerics for the CME with stationary noise. We give analytical criteria for the disappearance of bistability in the case of stationary noise in one-dimensional CMEs.

Noise-assisted interactions of tumor and immune cells

We consider a three-state model comprising tumor cells, effector cells, and tumor-detecting cells under the influence of noises. It is demonstrated that inevitable stochastic forces existing in all three cell species are able to suppress tumor cell growth completely. Whereas the deterministic model does not reveal a stable tumor-free state, the auto-correlated noise combined with cross-correlation functions can either lead to tumor-dormant states, tumor progression, as well as to an elimination of tumor cells. The auto-correlation function exhibits a finite correlation time τ, while the cross-correlation functions shows a white-noise behavior. The evolution of each of the three kinds of cells leads to a multiplicative noise coupling. The model is investigated by means of a multivariate Fokker-Planck equation for small τ. The different behavior of the system is, above all, determined by the variation of the correlation time and the strength of the cross-correlation between tumor and tumor-detecting cells. The theoretical model is based on a biological background discussed in detail, and the results are tested using realistic parameters from experimental observations.

Velocity fluctuations of population fronts propagating into metastable states

The position of propagating population fronts fluctuates because of the discreteness of the individuals and stochastic character of processes of birth, death, and migration. Here we consider a Markov model of a population front propagating into a metastable state, and focus on the weak noise limit. For typical, small fluctuations the front motion is diffusive, and we calculate the front diffusion coefficient. We also determine the probability distribution of rare, large fluctuations of the front position and, for a given average front velocity, find the most likely population density profile of the front. Implications of the theory for population extinction risk are briefly considered.

Determining the stability of genetic switches: explicitly accounting for mRNA noise

Cells use genetic switches to shift between alternate gene-expression states, e.g., to adapt to new environments or to follow a developmental pathway. Here, we study the dynamics of switching in a generic-feedback on-off switch. Unlike protein-only models, we explicitly account for stochastic fluctuations of mRNA, which have a dramatic impact on switch dynamics. Employing the WKB theory to treat the underlying chemical master equations, we obtain accurate results for the quasistationary distributions of mRNA and protein copy numbers and for the mean switching time, starting from either state. Our analytical results agree well with Monte Carlo simulations. Importantly, one can use the approach to study the effect of varying biological parameters on switch stability.

Rare events in population genetics: stochastic tunneling in a two-locus model with recombination

We study the evolution of a population in a two-locus genotype space, in which the negative effects of two single mutations are overcompensated in a high-fitness double mutant. We discuss how the interplay of finite population size N and sexual recombination at rate r affects the escape times t(esc) to the double mutant. For small populations demographic noise generates massive fluctuations in t(esc). The mean escape time varies nonmonotonically with r, and grows exponentially as lnt(esc)∼N(r-r(*))(3/2) beyond a critical value r(*).

Fixation of a deleterious allele under mutation pressure and finite selection intensity

The mean fixation time of a deleterious mutant allele is studied beyond the diffusion approximation. As in Kimura's classical work [M. Kimura, Proc. Natl. Acad. Sci. USA. 77, 522 (1980)], that was motivated by the problem of fixation in the presence of amorphic or hypermorphic mutations, we consider a diallelic model at a single locus comprising a wild-type A and a mutant allele A' produced irreversibly from A at small uniform rate v. The relative fitnesses of the mutant homozygotes A'A', mutant heterozygotes A'A and wild-type homozygotes AA are 1-s, 1-h and 1, respectively, where it is assumed that v<<s. Here, we employ a WKB theory and directly treat the underlying Markov chain (formulated as a birth-death process) obeyed by the allele frequency (whose dynamics is prescribed by the Moran model). Importantly, this approach allows to accurately account for effects of large fluctuations. After a general description of the theory, we focus on the case of a deleterious mutant allele (i.e. s>0) and discuss three situations: when the mutant is (i) completely dominant (s=h); (ii) completely recessive (h=0), and (iii) semi-dominant (h=s/2). Our theoretical predictions for the mean fixation time and the quasi-stationary distribution of the mutant population in the coexistence state, are shown to be in excellent agreement with numerical simulations. Furthermore, when s is finite, we demonstrate that our results are superior to those of the diffusion theory, while the latter is shown to be an accurate approximation only when N(e)s(2)<<1, where N(e) is the effective population size.

Extinction rates of established spatial populations

This paper deals with extinction of an isolated population caused by intrinsic noise. We model the population dynamics in a "refuge" as a Markov process which involves births and deaths on discrete lattice sites and random migrations between neighboring sites. In extinction scenario I, the zero population size is a repelling fixed point of the on-site deterministic dynamics. In extinction scenario II, the zero population size is an attracting fixed point, corresponding to what is known in ecology as the Allee effect. Assuming a large population size, we develop a WKB (Wentzel-Kramers-Brillouin) approximation to the master equation. The resulting Hamilton's equations encode the most probable path of the population toward extinction and the mean time to extinction. In the fast-migration limit these equations coincide, up to a canonical transformation, with those obtained, in a different way, by Elgart and Kamenev [Phys. Rev. E 70, 041106 (2004)]. We classify possible regimes of population extinction with and without an Allee effect and for different types of refuge, and solve several examples analytically and numerically. For a very strong Allee effect, the extinction problem can be mapped into the overdamped limit of the theory of homogeneous nucleation due to Langer [Ann. Phys. (NY) 54, 258 (1969)]. In this regime, and for very long systems, we predict an optimal refuge size that maximizes the mean time to extinction.

Speeding up disease extinction with a limited amount of vaccine

We consider optimal vaccination protocol where the vaccine is in short supply. In this case, the endemic state remains dynamically stable; disease extinction happens at random and requires a large fluctuation, which can come from the intrinsic randomness of the population dynamics. We show that vaccination can exponentially increase the disease extinction rate. For a time-periodic vaccination with fixed average rate, the optimal vaccination protocol is model independent and presents a sequence of short pulses. The effect can be resonantly enhanced if the vaccination pulse period coincides with the characteristic period of the disease dynamics or its multiples. This resonant effect is illustrated using a simple epidemic model. The analysis is based on the theory of fluctuation-induced population extinction in periodically modulated systems that we develop. If the system is strongly modulated (for example, by seasonal variations) and vaccination has the same period, the vaccination pulses must be properly synchronized; a wrong vaccination phase can impede disease extinction.

Time-resolved extinction rates of stochastic populations

Extinction of a long-lived isolated stochastic population can be described as an exponentially slow decay of quasistationary probability distribution of the population size. We address extinction of a population in a two-population system in the case when the population turnover-renewal and removal--is much slower than all other processes. In this case there is a time-scale separation in the system which enables one to introduce a short-time quasistationary extinction rate W1 and a long-time quasistationary extinction rate W2, and to develop a time-dependent theory of the transition between the two rates. It is shown that W1 and W2 coincide with the extinction rates when the population turnover is absent and present, but very slow, respectively. The exponentially large disparity between the two rates reflects fragility of the extinction rate in the population dynamics without turnover.