Extinctions of coupled populations, and rare event dynamics under non-Gaussian noise

The survival of natural populations may be greatly affected by environmental conditions that vary in space and time. We look at a population residing in two locations (patches) coupled by migration, in which the local conditions fluctuate in time. We report on two findings. First, we find that, unlike rare events in many other systems, here the histories leading to a rare extinction event are not dominated by a single path. We develop the appropriate framework, which turns out to be a hybrid of the standard saddle-point method, and the Donsker-Varadhan formalism which treats rare events of atypical averages over a long time. It provides a detailed description of the statistics of histories leading to the rare event and the mean time to extinction. The framework applies to rare events in a broad class of systems driven by non-Gaussian noise. Second, applying this framework to the population-dynamics model, we find a phase transition in its extinction behavior. Strikingly, a patch which is a sink (where individuals die more than are born) can nonetheless reduce the probability of extinction, even if it lowers the average population's size and growth rate.

Figure 1

(a) Extinction trajectories and probabilities in a single patch are correctly predicted by the IM. A realization of a system history is shown that ends with extinction for an isolated single patch. The dashed line is the IM prediction for the decline rate $r_d=r_1=0.1$ where $r_1$ is the growth rate. (b) When two patches are coupled (here $r_2=r_1$), extinction proceeds at a faster decline rate given by our theoretical prediction (17) $r_d=0.226$ (black dashed line). The IM prediction for the decline rate which is given by $r_1$ (orange doted line) is incorrect here. In the regime where the population is stabilized by noise-induced stabilization alone ($r_{1,2}<0$) the IM is altogether meaningless, see Sec. 4. Other parameters: $T_1=lnK/r_1$, $lnK=10$, and $D=0.2$. (c) Sink habitats can help to protect against extinction: When a source is coupled to a sink, the extinction probability per unit time is lower than when it is isolated, see Sec. 6b. Here extinction probability for a source patch alone (upper circles) is shown along with a source coupled to a sink (lower circles, indicating lower probabilities). Both are calculated by numerical solutions of the Fokker-Planck equation corresponding to the dynamics (1) with a logistic regulating term. Black lines are the analytical predictions (14) for large $lnK$, in perfect agreement with the numerics. The orange doted dashed line is the (incorrect) IM evaluation for the two patch case given by the slope $-2r_1=-3.8$. Here $D=1$.

Figure 2

(a) The stationary joint probability distribution $P_s(x_{+},x_{-})$ corresponding to the dynamics (1), with logistic regulating term. Arrows show the direction of the noiseless dynamics. The orange dotted arrow marks the incorrect IM prediction for the extinction trajectory, which simply follows the noiseless dynamics. Even though $r_{1,2}<0$, the average growth rate is positive, $r_{+}+\langle g\rangle >0$, and extinction is a rare process, which proceeds via an ensemble of trajectories. Parameters: $r_1=r_2=-0.1$, $lnK=10$, and $D=0.3$. (b) The $x_{-}$ distribution during typical growth (7) (solid line) is narrower during an extinction event (dashed line). This means that the abundances in the two patches tend to match more closely during an extinction event. Both obtained analytically at $r_1=r_2=2.402$.

Figure 3

(a) The extinction exponent $W$ (14) defined as $\text{MTE}\sim K^W$, as a function of the coupling strength $D$ for $r_1=2.2$ and $r_2=0.2$, together with the high (J5) and low (H5) asymptotic expressions in blue and magenta dashed lines respectively. The lower dotted orange line is the incorrect IM prediction obtained by Eqs. (D1) and (D2). The green dotted line is the upper bound in (18) given by the expression (8). (b) The rate function (12) in black solid line for $D={10}^2$. The blue dashed line is the high $D$ asymptotics (J3), showing a perfect match.

Figure 4

(a) the stationary distribution $P\left(x_{-}^c\right)$ of $x_{-}^c$, the $x_{-}$ coordinate conditioned on extinction in black solid lines. The dashed green and dashed red are the analytical predictions at small $D$, Eqs. (28) and (H2) respectively, in the extinction sink regime, with $r_1=4,r_2=-4.41$ (dashed red), and extinction source regime, $r_1=4,r_2=-3.47$ (dashed green). (a) The rate function (12) in black solid line with $r_{-}=1$. The blue dashed line is the low $D$ asymptotics (21). In both panels $D={10}^{-2}$.

Figure 5

Extinction at small coupling $D$, showing three regimes. (a) Migration into patch 2 during extinction, normalized by population size, for $r_1=2$. (b) The extinction probability exponent $W$, decline rate $r_d$, and growth rate. Note that in the extinction sink regime $r_2<-r_1$, patch 2 acts as a sink also during extinction, in that migration is non-negligible, and $W=4$ as for patch 1 alone; while in the extinction source regime $-r_1<r_2$, migration into patch 2 vanishes and the extinction probability is suppressed $W>4$. This happens both when patch 2 is a sink, $r_2<0$ (light red), or a source (light orange). The thin blue lines are analytical predictions for $D\to 0^{+}$, and the thicker gray lines for $D={10}^{-3}$.

Figure 6

A realization of extinction trajectories for $x_1^c$ in blue and $x_2^c$ in magenta generated from the process (19) and (20). (a) The extinction sink regime ($r_1=4,r_2=-4.41$). (b) extinction source regime, for a sink patch 2, $r_2=-3.47<0$. Dashed lines show the predicted decline rate (17). In both $D={10}^{-3}$ and $lnK=20$.

Figure 7

The confining potential $V_k$ (C3) [(a) and (c)] and its associated ground state ${\psi }_k$ [(b) and (d)] in thick black lines for $k=-15$. The dashed blue line correspond to $k=0$, that is, the unconditional dynamics. In the left panels (a) and (b) $r_{-}=10$, where the unconditional distribution is peaked away from the origin. For $k<0$, the ground-state eigenfunction's peak is advanced toward the origin, thus suppressing the value of $G$ (10). In the right panels (c) and (d) $r_{-}=0$, where the unconditional distribution is symmetric around the origin. The effect of $k<0$ here is to decrease the width of the eigenfunction which again suppresses the value of $G$. In all panels $D=0.1$.