Predictability of escape for a stochastic saddle-node bifurcation: when rare events are typical

TL;DR

This paper explores how stochastic fluctuations influence the predictability of escape events near a saddle-node bifurcation, revealing a transition from rare to typical escapes governed by a control parameter.

Contribution

It introduces a unified framework describing the interplay between deterministic bifurcations and stochastic noise, identifying regimes where escapes are rare or typical, with distinct trajectories.

Findings

Small noise regime: escapes are rare, described by Freidlin-Wentzell theory.

Large noise regime: escapes become typical, with peaked escape time distribution.

Different escape trajectories: algebraic divergence in rare events, exponential in typical events.

Abstract

Transitions between multiple stable states of nonlinear systems are ubiquitous in physics, chemistry, and beyond. Two types of behaviors are usually seen as mutually exclusive: unpredictable noise-induced transitions and predictable bifurcations of the underlying vector field. Here, we report a new situation, corresponding to a fluctuating system approaching a bifurcation, where both effects collaborate. We show that the problem can be reduced to a single control parameter governing the competition between deterministic and stochastic effects. Two asymptotic regimes are identified: when the control parameter is small (e.g. small noise), deviations from the deterministic case are well described by the Freidlin-Wentzell theory. In particular, escapes over the potential barrier are very rare events. When the parameter is large (e.g. large noise), such events become typical. Unlike pure noise-induced transitions, the distribution of the escape time is peaked around a value which is asymptotically predicted by an adiabatic approximation. We show that the two regimes are characterized by qualitatively different reacting trajectories, with algebraic and exponential divergence, respectively.

Figures

Equations

$$d{x}_t=-{\partial }_{x}V(x_t,t)dt+\sqrt{2\sigma }d{W}_t,$$

$${\tau }_M(x_0,t_0)=\inf \{t\ge t_0|x_t\ge M\}.$$

$$\dot{x}=x^2+t+2p,\dot{p}=-2xp.$$

$${\mathbb{E}}_{EK}[\delta ({\tau }_M-t)]=\frac{\sqrt{-t}}{\pi }{e}^{-\frac{4(-t{)}^{3/2}}{(3ϵ)}}\exp \left[-\frac{ϵ}{2\pi }{e}^{-\frac{4(-t{)}^{3/2}}{(3ϵ)}}\right].$$

$${\mathbb{E}}_{EK}[{\tau }_M^n]\sim \frac{(-1{)}^n}{2\pi }{\left(\frac{3ϵ\ln ϵ}{4}\right)}^{2n/3}.$$

Full text

Motivated by the possibility to predict the approach of a tipping point, many earlier studies have focused on early-warnings, i.e. features of a time-series which change before the transition occurs. In that framework, deterministic bifurcations are announced by phenomena such as increasing autocorrelation or variance, which are absent in noise-induced transitions Scheffer et al. (2009); Ditlevsen and Johnsen (2010); Thompson and Sieber (2011). Here, we adopt a different point of view, and study the universal statistical and dynamical features of transitions occurring under the joint effect of loss of stability and stochastic forcing. In the classical noise-induced case, transitions are completely random (they follow a Poisson distribution), the reaction rate satisfies the Arrhenius law Hänggi et al. (1990), and typical reacting trajectories follow the optimal path minimizing the action, as predicted by Freidlin-Wentzell theory Freidlin and Wentzell (1998). On the other hand, the transition is completely governed by the deterministic behavior in the bifurcation scenario. In this paper, we study when and how the transition occurs under a sweeping of the bifurcation parameter in time, in the presence of noise. An important motivation is to understand how much can be learned about the transition from observations of trajectories such as those represented in Fig. 1. What is the parameter governing the competition between deterministic and stochastic behavior? Can we distinguish trajectories corresponding to these two regimes?

We show that the transition time in a system approaching loss of stability is controlled by a single parameter and is always more predictable than the static noise-induced case. When the noise is small, the behavior is close to deterministic, the particle escapes slightly after the bifurcation and follows a universal trajectory with algebraic divergence. In the large noise regime, the escape time is determined by a balance between deterministic (lowering of the potential barrier) and stochastic effects, similarly to stochastic resonance Benzi et al. (1981); Gammaitoni et al. (1998). We show that the probability distribution of the escape time reaches a peak well before the bifurcation time, and can be predicted by an adiabatic approximation, corresponding to an Eyring-Kramers regime. Typical reacting trajectories leave the attractor in an exponential manner, and they show no imprint of the saddle-node, unlike the standard time-independent case.

The model.— Let us consider an overdamped Langevin particle in a time-dependent potential $V(x,t)$, undergoing a saddle-node bifurcation at $t=0$. The system is described by the stochastic differential equation:

[TABLE]

where $W_{t}$ is the standard Brownian motion. The most simple such potential has the form $V(x,t)=-a^{3}x^{3}/3-a\omega tx$, where the spatial scale $a$ and the time-dependent bifurcation parameter $\omega t$ determine the height of the potential barrier $\Delta V=4(-\omega t)^{3/2}/3$ and the width of the potential well $\sqrt{-3\omega t/a^{1/3}}$. By a proper choice of units, all the relevant parameters (geometry of the potential, speed of approach of the bifurcation, and noise amplitude) can be absorbed into a single, non-dimensional parameter $\epsilon=\sigma/(a^{2/3}\sqrt{\omega})$. With the rescaled variables, the potential now reads $V(x,t)=-x^{3}/3-tx$, and the stochastic differential equation is the normal form for the saddle-node bifurcation (with time-dependent bifurcation parameter) perturbed by noise: $dx_{t}=(x_{t}^{2}+t)dt+\sqrt{2\epsilon}dW_{t}$. We should keep in mind that the universality of the saddle-node normal form is only valid close to the bifurcation. Therefore, we expect our results to be universal for slow enough bifurcation parameter drift for arbitrary potentials. We shall denote by $x_{\pm}(t)=\pm\sqrt{-t}$ the fixed points for the stationary problem, which exist for $t<0$ only. The particle, initially lying in the stable state ($x_{0}=x_{-}(t_{0})$), may escape over the potential barrier under the influence of noise, or simply follow the deterministic dynamics and escape after the potential barrier has been removed by the bifurcation (see Fig. 2) Berglund and Gentz (2002, 2006); Kuehn (2011); Miller and Shaw (2012).

The single control parameter $\epsilon$ governs the competition between stochastic and deterministic effects. It decreases with the speed of the bifurcation and the potential stiffness, and increases with noise.

To give a precise meaning to the notion of escape, we shall compute the probability distribution of the first passage time, defined by

[TABLE]

Given the shape of the potential, the results do not depend on $M$ for $M$ large enough. For homogeneous Markov processes, a closed set of equations for the moments $\mathbb{E}[\tau_{M}^{n}]$ may be obtained, which leads to an explicit quadrature formula for the mean first-passage time for a 1D system Gardiner (2009); Risken (1989). Since the stochastic process defined by Eq. 1 is not time-homogeneous, these theoretical results do not apply here. We will discuss the behavior of the random variable $\tau_{M}$ using numerical results obtained with standard Monte-Carlo simulations and numerical solutions of the Fokker-Planck equation associated to Eq. 1, as well as theoretical arguments in the two limiting regimes $\epsilon\ll 1$ and $\epsilon\gg 1$.

Deterministic and small-noise behavior.— In the deterministic case ($\epsilon=0$), we have the dynamical saddle-node bifurcation, for which an analytical solution for the trajectory $x(t;x_{0},t_{0})$ with initial conditions $x_{0},t_{0}$ can be found in terms of Airy functions. In particular, the attractor $\bar{x}(t)=\lim_{t_{0}\to-\infty}x(t;x_{-}(t_{0}),t_{0})$ simply reads $\bar{x}(t)=\text{Ai}^{\prime}(-t)/\text{Ai}(-t)$. When $-t$ is large, it follows the stationary solution $x_{-}(t)=-\sqrt{-t}$. At a time of order one before the bifurcation ($t=0$), the trajectory detaches, and diverges to infinity after the bifurcation (see Fig. 2). The singularity occurs at a time $t_{\star}\approx 2.33811$, which is the opposite of the largest root of the Airy function Ai. The divergence is algebraic: $\bar{x}(t)\sim(t_{\star}-t)^{-1}$, and the deterministic first-passage time is easily related to the singularity: $\bar{\tau}_{M}=t_{\star}-1/M+o(1/M)$.

When the noise amplitude $\epsilon$ is small, escapes over the potential barrier have so low (albeit non-vanishing) probability that the behavior of the system is dominated by escapes after the bifurcation occurs Berglund and Gentz (2002, 2006); Kuehn (2011); Miller and Shaw (2012). This regime is close to the deterministic behavior.

The \acpdf of the first-passage time $\tau_{M}$, computed numerically, is shown in Fig. 3. When $\epsilon$ is small, the \acpdf is close to Gaussian. As $\epsilon$ increases, the \acpdf becomes more and more skewed, and heavy tails develop on the left, due to the presence of early exits. This can be interpreted in the framework of large deviation theory Freidlin and Wentzell (1998). Let us assume that the \acpdf of $\tau_{M}$ satisfies a large deviation principle: $p_{M}(\tau)\asymp e^{-I_{M}(\tau)/\epsilon}$ for $\epsilon\to 0$. If the rate function $I_{M}$ possesses a single zero $\tau_{\star}$ (also a global minimum), then the Gaussian behavior of $\tau_{M}$ corresponds to a quadratic approximation around $\tau_{\star}$. Besides, the integral defining the mean first-passage time $\mathbb{E}[\tau_{M}]$ can be evaluated with a saddle-point approximation. At first order, we expect the deterministic first-passage time $\bar{\tau}_{M}$, and the error associated with the approximation is linear in $\epsilon$ Erdélyi (1956): $\mathbb{E}[\tau_{M}]=\bar{\tau}_{M}(1+O(\epsilon))$. Based on similar considerations, the asymptotic behavior of the standard deviation is expected to be $\sqrt{\mathbb{E}[[\tau_{M}^{2}]]}\sim\sqrt{\epsilon}$. These provide a good fit of numerical results, as shown in Fig. 4.

In the limit $M\to\infty$, Ref. Miller and Shaw (2012) provides an exact result which shows that the random variable $\tau_{\infty}$ (the time at which the trajectory becomes unbounded) satisfies a large deviation principle, with rate function $I_{\infty}(\tau)=\text{Ai}^{2}(-\tau)/[4\pi^{2}\text{Bi}^{2}(-\tau)\int_{-\infty}^{T}\text{Ai}^{4}(-t)dt]$, when $\tau$ lies in the range between the first two zeros of Bi. In this interval, the only root of $I_{\infty}$ is $t_{\star}$, and the above reasoning applies, with $\mathbb{E}[\tau_{\infty}]\to t_{\star}$.

More generally, the large deviation property for $p_{M}$ can be understood at a formal level by writing the probability of fluctuations within the path integral formalism, following the pioneering work of Onsager and Machlup Onsager and Machlup (1953); *Machlup1953. In the $\epsilon\ll 1$ regime, the probability of a path $x(t)$ satisfies $P[x]\asymp e^{-\mathcal{A}[x]/\epsilon}$, introducing the Freidlin-Wentzell action functional $\mathcal{A}[x]=\frac{1}{4}\int dt(\dot{x}+V^{\prime}(x,t))^{2}$ Freidlin and Wentzell (1998). This probability distribution is dominated by the deterministic attractor $\bar{x}(t)$, for which $\mathcal{A}[\bar{x}]=0$. By contraction, the probability to reach $M$ at time $\tau$ is dominated by another action minimizer $x^{\star}_{M,\tau}(t)$ such that $\mathcal{A}[x^{\star}_{M,\tau}]=A_{M}(\tau)\equiv\inf_{x}\{\mathcal{A}[x]|x(t_{0})=x_{0},x(\tau)=M\}$. As we shall see below, such optimal paths follow the deterministic attractor as long as possible, at no cost in the action functional, then detach from the deterministic trajectory in a monotonously increasing manner. The monotonicity of the action minimizers can be proved directly by writing the Hamilton equations associated to the variational problem:

[TABLE]

As a consequence, there is no optimal path which reaches $M$ more than once, because of the cost in the action functional associated with climbing the gradient of the potential. Hence, the trajectories $x^{\star}_{M,\tau}(t)$ dominate the \acpdf of the first-passage time $\tau_{M}$, which satisfies a large deviation principle, with rate function $I_{M}(\tau)=A_{M}(\tau)$. Then, the Gaussian behavior of $\tau_{M}$ corresponds to a quadratic approximation of $A_{M}(\tau)$ around its minimum $\bar{\tau}_{M}$.

In fact, although the \acpdf of $\tau_{M}$ already exhibits substantial deviation from Gaussianity for $\epsilon>0.1$ (see Fig. 3), the above approach describes accurately the first two moments up to order one values of $\epsilon$. For these low-order statistics, a sharp transition between the asymptotic regimes $\epsilon\ll 1$ and $\epsilon\gg 1$ occurs near $\epsilon=20$ (see Fig. 4). On average, the transition always happens before the deterministic case. It happens before the bifurcation for $\epsilon>\epsilon_{c}\approx 4$, and after for $\epsilon<\epsilon_{c}$. The discrepancies seen in Fig. 4 for the mean first-passage time at very small $\epsilon$ are due to numerical errors.

Adiabatic approximation in the large-noise regime.— When the noise amplitude is large, the range of times for which the escape rate is not too small is long enough for those events to dominate the distribution of the first-passage time. In this regime, escapes over the potential barrier, which are usually rare events, become the typical events. However, that range is also short enough for the distribution of the first-passage time to be peaked around a given value (Fig. 5), determined by the competition between stochastic and deterministic effects. This is very different from the classical Kramers problem, for which the first-passage time is distributed according to an exponential law Hänggi et al. (1990).

Because the relaxation time scale is much smaller than the scale at which the potential evolves, this case can be treated with an adiabatic approximation. We introduce the transition probability $P(x,t|x_{0},t_{0})$, which satisfies the (forward) Fokker-Planck equation with initial condition $P(x,t_{0}|x_{0},t_{0})=\delta(x-x_{0})$. With reflecting boundary condition on the left and absorbing boundary condition at a fixed value $M>x_{+}(t_{0})$, $G(x_{0},t_{0};M,t)=\int_{-\infty}^{M}dxP(x,t|x_{0},t_{0})$ is the probability that a particle initially at $x_{0}$ has not reached $M$ at time $t$. In other words, $\text{Prob}(\tau_{M}>t)=G(x_{0},t_{0};M,t)$. $G$ always satisfies a backward Fokker-Planck equation. For homogeneous Markov processes, because $\partial_{t}G=-\partial_{t_{0}}G$, this partial differential equation allows to compute explicitly the moments of the first-passage time Gardiner (2009). Besides, when the potential barrier height $\Delta V$ is large ($\Delta V\gg\epsilon$), transition times form a Poisson process with transition rate given by the Eyring-Kramers formula: $\lambda=\sqrt{V^{\prime\prime}(x_{a})V^{\prime\prime}(x_{s})}/(2\pi)e^{-\Delta V/\epsilon}$, where $x_{a}$ is the position of the attractor and $x_{s}$ that of the saddle point Berglund (2013). Here, since the potential variations are adiabatic, the transition rate at each time $t$ is well approximated by the Eyring-Kramers formula for the “frozen” potential at fixed $t$: $\partial_{t}G(x_{0},t_{0};M,t)=-\lambda_{EK}(t)G(x_{0},t_{0};M,t)$, with $\lambda_{EK}(t)=\frac{\sqrt{-t}}{\pi}\exp(-\Delta V(t)/\epsilon)$, where $\Delta V(t)=4(-t)^{3/2}/3$. This formula is expected to be valid for times $t\ll-\epsilon^{2/3}$, and initial conditions close to the attractor.

An explicit formula is obtained for the \acpdf of the first-passage time in this regime:

[TABLE]

We show in Fig. 5 that this approximation indeed provides a very good fit of the numerically computed \acpdf of the first-passage time $\tau_{M}$ when $\epsilon$ is large enough (here $\epsilon=100$). From Eq. 4, we deduce the asymptotic behavior of the moments of the first-passage time: when $\epsilon\to+\infty$,

[TABLE]

The mean first-passage time and its standard deviation are shown in Fig. 4; again, the theoretical result fits very well the numerical simulations above a critical $\epsilon$ approximately equal to 20. Besides, for the adiabatic approximation to be self-consistent, we need $\mathbb{E}_{EK}[\tau_{M}]\ll-\epsilon^{2/3}$. This condition is asymptotically verified, but this is only due to the logarithmic corrections in Eq. 5. Hence, the adiabatic approximation converges relatively slowly. This explains why the theoretical result is not very accurate for $\epsilon=10$ for instance (see Fig. 5). For such moderate values of the control parameter, the approximation slightly overestimates early escapes, and makes a dramatic error on escapes occurring later than the average time. Indeed, for such regimes, escapes after the bifurcation occurs, which make no sense in the Eyring-Kramers approximation, are already relatively probable events (about one or two standard deviations away from the mean). Although such events are still unaccounted for at larger $\epsilon$, they are then so improbable that it does not hamper the accuracy of the approximation for low-order moments.

Predictability of the reacting trajectory.— Now, we consider the statistics of the escape dynamics.

For small $\epsilon$, the particle typically escapes after the bifurcation, and the dynamics is then essentially deterministic. Hence, all the reacting trajectories have the same shape as the deterministic attractor, even though escape typically occurs slightly before. This can be illustrated by conditioning trajectories on the first-passage time. Fig. 6 shows that, when conditioning on a typical value for the first-passage time (less than one standard deviation away from the mean), apart from fluctuations of order $\sqrt{\epsilon}$, the reacting trajectories remain close to a trajectory with the same shape as the deterministic attractor $\bar{x}(t)$. That trajectory can be predicted as an instanton: it is a minimizer of the action $\mathcal{A}[x]$ with fixed initial and final points. In particular, it has an algebraic divergence of the form $x(t)\sim(t_{\star}^{\prime}-t)^{-1}$, where $t_{\star}^{\prime}=\tau_{M}+1/M$ for trajectories conditioned on $\tau_{M}$. Even rare transitions look similar to the deterministic attractor and are hardly distinguishable.

The situation is more complex in the large $\epsilon$ case. In the classical case of time-independent potential barrier activation (e.g. the Kramers problem), the instanton is degenerated: it takes an infinite time to leave the attractor and an infinite time to reach the saddle-point Caroli et al. (1981). The time spent by reactive trajectories in the vicinity of the saddle-point is distributed according to a Gumbel law, and scales with $\ln\epsilon$. This is why pure noise-induced transition times are unpredictable. Here, it is very different (Fig. 7): the degeneracy of the instanton (starting from the attractor at $t\to-\infty$) is lifted because of the time dependence. Nevertheless, compared to stochastic trajectories, the instanton triggers slightly before (shifting it in time describes correctly the dynamics away from the saddle-point) and still takes a longer time to pass the saddle-point. This is because the amplitude of typical fluctuations is large enough to smooth the instanton slow-down at the saddle-point. As a consequence, some trajectories remain in the vicinity of the saddle-point for a short while, but the majority of them swing directly to the other side of the potential. This can also be seen as a consequence of the cost, in terms of the action functional, associated to staying at the position of the saddle-node, which is not a deterministic solution of the time-dependent problem. Finally, let us note that because $\epsilon\gg 1$, we are no longer in the large deviation regime for $P[x]$, and there is a priori no reason for the statistics of observables such as reacting trajectories to be dominated by an action-minimizing path. This also holds when considering the appropriate action for finite $\epsilon$, $\mathcal{A}^{\prime}=\mathcal{A}-\epsilon/2\int dt\partial_{x}^{2}V$ Graham (1973). Hence, in the $\epsilon\gg 1$ regime, typical reacting trajectories are more predictable than the instanton in the sense that there is no imprint from the saddle-node, unlike the standard stationary case but similarly to the glacial-interglacial transitions (Fig. 1). We have chosen a typical value for $\tau_{M}$ on Fig. 7, but the conclusions remain true for values which deviate significantly from the mean first-passage time.

Conclusion.— We have given a global picture of the possible scenarios for transitions in a noisy system undergoing loss of stability, and the associated predictability. We have shown that there exist two regimes characterized by a single control parameter $\epsilon$. When $\epsilon$ is small, the escape time only deviates from the deterministic value in a Gaussian manner, and the reacting trajectories have a universal shape with an algebraic divergence. On the contrary, when $\epsilon$ is large, escapes over the potential barrier become typical, but they are different from the standard Kramers problem: their \acpdf is peaked, and can be predicted by an adiabatic approach consistent with large deviation theory. Reacting trajectories leave the attractor exponentially fast and do not stick to the saddle-point. Such trajectories are not described by large deviation theory. These results open new prospects for the analysis of time series exhibiting abrupt transitions such as those encountered in climate dynamics.