Most Probable Dynamics of the Single-Species with Allee Effect under Jump-diffusion Noise

We investigate the most probable phase portrait (MPPP) of a stochastic single-species model with the Allee effect using the non-local Fokker-Planck equation. This stochastic model is driven by non-Gaussian as well as Gaussian noise, and it has three fixed points. One of them is the unstable state which lies between the two stable equilibria. We focus on the transition pathways from the extinction state to the upper fixed stable state for the transcription factor activator in a single-species model. This helps us to study the biological behavior of species. The most probable path is obtained from the solution of the non-local Fokker-Planck equation corresponding to the population system of the single-species model, and the corresponding maximum possible stable equilibrium state is determined. We also obtain the Onsager-Machlup (OM) function for the stochastic model and solve the corresponding most probable paths. The numerical simulation shows that: (i) When non-Gaussian noise is presented in the system, the maximum of the stationary density function is located at the most probable stable equilibrium state; (ii) If the initial value increases from extinction state to the upper stable state, the most probable trajectory goes to the maximal likely equilibrium state, in our case it lies between 9 and 10; (iii) The most probable paths increase to stable state quickly, then maintain a nearly constant level, and approach to the upper stable equilibrium state as time goes on. These numerical experiment findings accelerate growth for further experimental study, in order to achieve good knowledge about dynamical systems in biology.


Introduction
Single-species dynamics is one of the core research areas in theoretical ecology.Research about singlespecies dynamics enables the researcher to find out the conditions of extinction and persistence of the species.The strong motivation for the researchers to develop mathematical models is to understand the cause of cycles, like populations [57].
Population modeling is very important for species management, for example, in developing recovery plans for species threatened by extinction, managing fisheries for the highest possible sustainable yield, and trying to contain or prevent the spread of invasive species [3,11,12].
In literature, one can find several models of the dynamical single-species growth system.Gompertz growth model [29], Verhulst growth model with or without Allee effect [30], power law growth model [15], the interconnections between deterministic and stochastic system [6], Gilpin-Ayala model [56] are only a few to mention.In this study, we compute a single-species model focusing on the Verhulst growth model with the Allee effect developed by Y. Jin [21].
The dynamics of biological phenomena, particularly that of populations of living beings, besides some clear trends, are frequently influenced by unpredictable components due to the complexity and variability of environmental conditions [34].Extensive researchers in modeling and analysis of random fluctuations [5,38,39] in biological dynamical systems have been ongoing for a long time now.The studies of events in the population such as persistence stationary distribution, and extinction in stochastic single-species models become an interesting and important research field.One of the hot issues in population dynamics is developing sufficient conditions for the persistence of biological species as mentioned in [13,26,49] and the references therein.
The population may be affected by sudden environmental noises [33,52].For example, severe acute respiratory syndrome (SARS), human immunodeficiency virus (HIV), the smoking habit [58], and the recent COVID-19 [43], earthquakes [51], temperature [54], and hurricanes [50].These sudden environmental perturbations may bring substantial social and economic losses.Stochastic single-species model perturbed by Brownian motion has been researched extensively by many authors [16,25,36,40,42].However, stochastic extension of population process driven by Gaussian noise cannot explain the aforementioned random and intermittent environmental perturbations.Introducing a Lévy process into the underlying population dynamics would explain the impact of these random jumps [4].There have been a few studies that investigated dynamical systems where the noise source is a Lévy process [32,45,48,56].Implying the Lévy noise into the biological system to simulate the effect caused by the external environment is more effective and nearer to reality than the Gaussian noise.The investigation of the single-species model is still in its infancy even though noisy fluctuations naturally portray random intermittent jumps.Lévy noise is widely applied in studying natural and man-made phenomena in science, among which we mention biology [20], physics [17,47], and economics [27,37].
Under this research, we consider the population dynamics of a single-species growth model with Allee effect perturbed by stable Lévy fluctuations.We also analyze the influence of Lévy noise fluctuation on the system (1).Investigating the impact of noisy fluctuations acts a pivotal part in demonstrating the intricate interactions between the single-species models and their complex surroundings.We study how Allee effects and stochasticity combine to affect population persistence in here.To find the numerical solutions for the Fokker-Planck equation determined by non-local differential equation with symmetric α-stable Lévy motion, we apply a finite difference method probed by Gao et al. [14].
The most probable phase portrait was first proposed by Duan [10,Section 5.3.3].Cheng et al. [8] obtained the analytical results of the MPPP and showed that the MPPP can give useful information about the propagation of stochastic dynamics in the one-dimensional model.Wang et al. [44] studied the stochastic bifurcation by using the qualitative changes of the MPPP to a stochastic system driven by multiplicative stable Lévy noise.In Ref. [19], the scholars investigated the most probable trajectories of the tumor growth system with immune surveillance under correlated Gaussian noises, and derived analytical solution of the most probable steady state by utilizing the extremum theory with the local Fokker-Planck equation (FPE) in the system.A function which summarizes about the behavior of the dynamics of a continuous stochastic process was defined as the Onsager-Machlup function [35].The Onsager-Machlup function for stochastic models driven by both non-Gaussian and Gaussian noises was established in [7].The authors also examined the corresponding MPPP of the stochastic dynamical systems.Cheng et al. [9] focused on the impact of Gaussian noise and jump stable Lévy noise in a genetic regulatory system, and they minimized the OM action functional for the stochastic dynamics driven by Gaussian and obtained the most probable transition pathway.This inspired us to study the MPPP of the single-species model.
Therefore, our goal is about to investigate how the most probable trajectories escape from the singlespecies state to the extinction state more quickly.
Consider the following stochastic single-species growth model with Allee effect: for t ≥ 0 and X 0 = x 0 , where X t− is the left limit of the population size X t .The attack rate γ 4 Represents the handling time of predator M = s/γ 2 The carrying capacity t Time Stochastic force Ñ(dt, dy) = N(dt, dy) − ν α (dy)dt is a compensated Poisson random measure with associated Poisson random measure N(dt, dy) and intensity measure ν α (dy)dt, in which ν α (dy) is Lévy measure on a measurable subset Y of (0, ∞) with ν α (Y) < ∞; see [53].
The following restriction on system (1) is natural for biological meaning: When ǫ(y) > 0, the perturbation stands for the increasing species, e.g.planting, while ǫ(y) < 0 represents that the species is decreasing, e.g.harvesting and epidemics.
The main aim of this study is to investigate stochastic dynamics of single-species biological populations in random environments.We model the evolution of these populations with first order ordinary autonomous differential equations bringing in the coefficients and inputs which are stochastic processes.The two stochastic processes germane to this study are Brownian motion and Lévy process.Brownian motion describes random fluctuations that are continuous in time; see Subsection 2.1.Lévy process, of which Brownian motion is a special case, is used to model random fluctuations which may have discontinuities or jumps; see Subsection 2.2.
Here, we develop the stochastic single-species model with the Allee effect influenced by Gaussian and non-Gaussian noises.Firstly, we review the deterministic model, calculate its equilibrium solutions and describe the behavior of the fixed points.Secondly, we get the highest possible paths, and the corresponding maximum possible stable states attracting the nearby maximum possible paths of the stochastic system (1).We do this by finding the stationary density function which is the solution of the non-local FPE.To solve the non-local partial differential equation, we use the finite difference method proposed in [14].This method helps us explore some dynamical behaviors of the single-species system under the impact of non-Gaussian Lévy noise.
This study is organized as follows.In the second section, we recall the definitions of the one-dimensional Brownian motion B t and symmetric α-stable Lévy motion L α t .In the third section, we discuss the formulation and analysis of the deterministic model ( 2) of the single-species with Allee effect .In the fourth section, we explain the analysis of the stochastic single-species model (1) with Allee effect.We also review the definition of the Onsager-Machlup function and most probable phase portraits in the subsections 4.1 and 4.2, respectively.The numerical results and the biological implication of our experimental findings are presented in the fifth section.We conclude our research by giving a brief summary in the last section.

Preliminaries
Under this section, we define the one-dimensional Brownian motion starting at time t = 0 as a process B t and α-stable Lévy motion L α t , which constitute a class of stochastic processes that have independent and stationary increments as defined below.Throughout this study, we denote R + = (0, ∞), R = (−∞, ∞), and X t ∈ R + , for t ≥ 0.

Brownian motion
Brownian motion B t (also called Wiener process) is a one-dimensional stochastic process defined on complete probability space (Ω, F , F t , P), which has independent and stationary increments [31,41].Brownian motion B t satisfies the following conditions: (i) B t has continuous paths, and its paths are nowhere differentiable almost surely; (ii) B t has stationary increments, i.e., B t − B s is normally distributed with mean 0 and variance t − s for any 0 ≤ s ≤ t; (iii) The process starts at the origin, i.e., B 0 = 0 almost surely; (iv) B t has independent increments, i.e., B t − B s is independent of the past for s < t.

The α-Stable Lévy motion
Lévy motions L t are a class of non-Gaussian stochastic processes.A Lévy motion L t having values in R is determined by a drift coefficient b ∈ R, Q ≥ 0 and a Borel measure ν defined on R \ {0}.The triplet ( b, Q, ν) is the so-called generating triplet of Lévy motion L t .A Lévy motion can be written as linear combination of time t, a Brownian motion and a pure jumping process [28,41], i.e., L t can be expressed as The Lévy-Khinchin formula for Lévy motion has a specific form of its characteristic function. where A stable distribution S α (θ, β, γ) is the distribution for a stable random variable, where the stability index α ∈ (0, 2), the skewness β ∈ (0, ∞), the shift γ ∈ (−∞, ∞), and scale index θ ≥ 0. A α-stable Lévy motion L α t [1,23,24] is a non-Gaussian stochastic process satisfying (i) the random variables t has stochastically continuous sample paths, i.e., for 0 ≤ s ≤ t and δ > 0, the probability P( , 0, 0).In the case of a one-dimensional isotropic α-stable Lévy motion, the Lévy triplet has the drift factor b = 0 and the diffusion coefficient Q = 0.In this study, we focus on jump process with a specific size in generating triplet (0, 0, ν α ) for the random distribution S α which can be defined by , and Γ is the Gamma function.

Dynamical analysis of the deterministic model
The deterministic form of the nonlinear model ( 1) without noise is given as This system can be written as dX dt = − dU(X) dX , where U(X) is the potential function given by The single-species model (2) with Allee effect has equilibrium points X 1 = 0 and where β = 4γ 2 γ 3 γ 4 (sγ 3 γ 4 −γ 2 ) 2 (γ 3 − s).If β < 1, then the equilibrium states of system (2) are X 1 = 0, an extinction equilibrium; , a lower unstable equilibrium; , an upper stable equilibrium.
If β = 1, the single-species deterministic model (2) has only two equilibria states: stable state X 1 = 0, and unstable state For simplicity and convenience of discussion, we choose the parameters , and X 4 = 1−γ 2 2 γ 2 .For β < 1, the extinction state X 1 = 0, and the equilibrium solution X 3 are stable, but X 2 is unstable.Fig. 1(b) shows that when the value of attack rate γ 3 increases, the unstable state X 2 and stable state X 3 get more close to each other, then become one solution and finally disappear indicating the occurrence of the saddle-node bifurcation.The critical value of attack rate γ c = 2.67 of the deterministic single-species system (2) with Allee effect is obtained by solving the equation U(X 1 ) = U(X 3 ).This value is an indication to the transition phenomena between the unstable and stable state for deterministic single-species growth model.The steady state (extinction state) X 1 is stable if γ 3 > γ c , and the steady state X 3 exhibits the stability property for γ 3 < γ c .

Dynamical analysis of the stochastic system
In this section, we discuss the behavior of the solution of the stochastic system (1).Firstly, we recall the definition of the Onsager-Machlup function for the stochastic differential equation driven by jump noise.This helps to measure OM induced by the jump process.Secondly, we examine the corresponding most probable paths.Finally, we present the numerical experiment findings using finite difference method [14].Hence, the numerical solution of the stochastic model provides useful information for understanding the dynamical behavior of the system (1).

Onsager-Machlup functional
The Onsager-Machlup functional defines a probability density for a stochastic process in which the probability density is estimated implicitly.It can be used for purposes of reweighting and sampling trajectories, as well as determining the most probable trajectory based on variational arguments.The most probable transition pathway can be obtained by minimizing the Onsager-Machlup function.The whole procedure enables us to detect the dynamics of the most probable path [18].The stochastic single-species system (1) with Allee effect, as proved by Jin [21], has a unique global and positive solution with the initial condition X 0 = x 0 .The jump-diffusion process X t is adapted and càdlàg; see Fig. 2. Denote the space of càdlàg paths starting at x 0 of a solution process X = {X t , t ≥ 0} of (1) by This space equipped with Skorokhod's J 1 -topology generated by the metric d R + is a Polish space [55].For functions for every t, s ≥ 0 and some λ ∈ Λ R + , We consider the corresponding jump-diffusion process X t (ω) := ω(t), t ∈ [0, T ] defined on the canonical probability space (R [0,T ] , B(R) [0,T ] , P T ).Since the paths of X are càdlàg, we identify X t on the space (D T x 0 , B T x 0 , P) instead of (R [0,T ] , B(R) [0,T ] , P T ), where D T x 0 is defined similarly as the space D x 0 on the time interval [0, T ].The associated Borel σ-algebra is B T x 0 = B(R) [0,T ] ∩ D T x 0 , and then (D T x 0 , B T x 0 ) is a separable metric space [2, Section A.2].The probability measure P is generated by P(A ∩ D T x 0 ) := P T (A) for each A ∈ B(R) [0,T ] .Because every càdlàg function on [0, T ] is bounded, we equip D T x 0 with the uniform norm Hence, D T x 0 is a Banach space.In order to find the most probable tube of X t , we should determine the probability that paths lie within the closed tube It is a subset of the space D T x 0 of càdlàg functions on the interval from 0 to T containing a function z together with its ε-neighborhood.Define the measure µ X on B(R) [0,T ] induced by the solution process X t for the stochastic nonlinear model (1) via For sufficiently small ε > 0, the main contribution of the above probability is given by the measure of the trajectories in the ε-tube of z ∈ D T x 0 : where K(z, ε) ∈ B(R) [0,T ] .As the ε-tube K(z, ε) depends on the reference path z, it is necessary for us to look for the "most probable" trajectory z which maximizes the measure µ X (K(z, ε)) in Eq. ( 4).When we focus on the differentiable functions z ∈ D T x 0 , we have the following meaningful definition.
Definition 1.Let 0 < ε ≪ 1 be given.For a ε-tube surrounding a reference path z(t), the probability of the solution process X t , t ∈ [0, T ] lying in this tube is estimated by where the integrand OM(ż, z) is called Onsager-Machulup function and ∝ denotes the equivalence relation for ε small enough.The intergral T 0 OM(ż, z)dt is the Onsager-Machulup functional.Remark 2. The Onsager-Machulup function is similar to the Lagrangian function of a dynamical system in classical mechanics, and the OM functional would correspond to the action functional.In particular, for an SDE with pure jump Lévy noise, Definition 1 is still applicable, and the minimizer of the OM functional T 0 OM(ż, z)dt gives the most probable path for this non-Gaussian stochastic system.Moreover, the minimizer z may be chosen from a more general function space.
Our main result about the expression of the OM function for a jump-diffusion process is clearly presented in the basic theorem.
Theorem 1.For the stochastic nonlinear system (1) with the jump measure satisfying Y ǫ(y)ν α (dy) < ∞, the Onsager-Machlup function [35] is characterized, up to an additive constant, by: where z ∈ D T x 0 is a differentiable function.The contribution of pure jump Lévy noise to the OM function is the third term.When the jump measure is absent, we cover the OM function for the case of diffusion.In terms of OM function, the measure of tube K(z, ε) defined in (3) can be approximated as follows: where Y c t is defined by The proof of Theorem 1 is given in [7,Theorem 4.1].
In Gaussian noise case (ǫ(y) = 0), the stochastic single-species model (1) becomes We apply Lamperti transforms for solving SDE driven by multiplicative noise [10,Example 6.48].This method allows us to transform the multiplicative noise into additive noise.Because numerically solving an additive-noise SDE is usually easier than solving a multiplicative-noise SDE as in Eq. ( 5).Assume g ∈ C 2 (R) and define Y t = g = ln(X t ).Then the new SDE has the following form: where Since the most probable transition path for a stochastic single-species model is the minimizer of the Onsager-Machlup action functional, denoted by Z m , it can be obtained from the following least action principle where the integrand function (Onsager-Machlup function) [9] is given by Thus Eq. ( 7) satisfies the following Euler-Lagrange equation The most probable transition pathway Z m (t) of system ( 6) is characterized by To solve two-point boundary value problem in Eq. ( 9), we apply the shooting method depicted in Ref. [22].

Most probable phase portraits
As for the solution of the Fokker-Planck equation, the probability density function p(X, t) is a surface in the (X, t, p)-space.For a given time t, the maximizer X m (t) for p(X, t) (i.e., X m (t) = max X∈(0,∞) p(X, t)) shows the most probable (i.e., maximal likely) location of this orbit at time t.The orbit traced out by X m (t) is called a most probable orbit starting at x 0 .Thus, the deterministic orbit X m (t) follows the top ridge of the surface in the (X, t, p)-space as time goes on.

Non-local Fokker-Plank equation
The Fokker-Planck equation describes the time evolution of the probability density function, but it can be solved analytically only in special cases.We are interested in the steady-state probability distribution (equilibrium distribution), and want to express the stationary solution of the non-local Fokker-Planck equation.This makes the estimate of the most probable phase portrait possible in Lévy noise case numerically and algorithmically.
Let f : R → R be a smooth function.Suppose that the solution X t of system (1) has a conditional probability density p(X, t|x 0 , 0).For convenience, we drop the initial condition and simply denote it by p(X, t).On one hand, On the other hand, by virtue of Itô's formula, Taking expectation on both sides of (10), we gain Noting that the infinitesimal generator of the solution X t for system (1) is The equation ( 11) is rewritten as As a result, the Fokker-Planck equation for the stochastic nonlinear system (1) of the solution process To simulate the non-local Fokker-Planck equation ( 13), we apply a numerical finite difference method given in Gao et al. [14].
If the Lévy motion is replaced by Brownian motion, then the local Fokker-Planck equation has the following form: The stationary probability density function p s (X) of Eq. ( 14) can be solved by or equivalently, Due to the complexity of stationary solution, we take the extrema of the stationary probability density function (spdf) located at x s directly.In other words, the spdf satisfies ∂ X ( p s (x s )) = 0. Since p s (x s ) 0, Eq. ( 16) becomes Eq. ( 17) is completely different from the equilibrium state of the deterministic model (2) because of the presence of noise with λ term.The numerical solution of Eq. ( 17) is plotted in Fig. 3

Numerical results and its biological implications
In order to make the readers understand our results much better, we perform some numerical simulations to illustrate our theoretical results.Based on the finite difference method [14], the numerical simulations are very useful in the study of real examples of population.In the present section, we define the bifurcation time as the time between the changes in number of maximal likely equilibrium states.It is a time scale for the birth of a new most probable stable equilibrium state.We also show the intervals where there exist one or two maximal likely stable equilibrium states, the value of the equilibrium states, and the point where the number of metastable states of the stochastic single-species model (1) varies.Since the numerical solutions of a model depend on the values of all its deterministic parameters and noise intensities.Here, we discuss the effect of the parameters in Table 1 on the investigated system (1).For simplicity, we choose four maximum likely pathways together with the initial conditions selected in different intervals.
The potential function denoted by U(X) in Fig. 1(a) has two stable steady states X 1 and X 3 , and an unstable steady state X 2 for β < 1.This function has a maximum value at the unstable equilibrium solution X 2 .At the stable fixed points X 1 and X 3 , the potential function attains its minima.For the value of β > 1, the nonlinear system (2) has only one equilibrium point, which is the trivial point X 1 = 0.
In Fig. 1(b), we sketch the equilibrium states versus attack rate γ 3 .For β < 1, there exist two stable equilibrium states X 1 and X 3 and one unstable equilibrium state X 2 .While β > 1, X 1 = 0 is the unique equilibrium state that is stable.Thus, the parameter β = 1 is the bifurcation parameter value.
The distance between the unstable equilibrium X 2 and the stable fixed point X 3 becomes very small when β approaches to 1.This indicates that the expected time to extinction may be too short, as clarified in Fig. 1(b).
Fig. 2 displays the numerical simulation of the stochastic single-species model (1) with Allee effect when it is persistent or extinct at different value of initial condition x 0 .This figure proves that the solutions of the stochastic nonlinear system (1) are positive, and extinction species occurs when the initial condition is less than the value of X 2 , as demonstrated in Fig. 2(b).While the initial condition is greater than the value of X 2 , there is stochastic persistence.
In Fig. 3(a), we depict the most probable transition pathways Z m (t) of system ( 6) for different values of λ.This figure tells us that as time grows, the most probable paths Z m (t) increase to the stable state X 3 quickly, and remain at a nearly constant level, then approach to the high stable equilibrium state.Fig. 3(b) demonstrates the curves for the most probable steady state x s of the stochastic single-species model (1) with ǫ(y) = 0 driven by Gaussian noise at different values of the noise intensity λ.The steady state curves exhibit a bi-stability in the interval (a 1 , a 2 ).For γ 3 > a 1 , the stable steady state stays at the extinction state.While for γ 3 < a 1 , it is located at the stable equilibrium state.Because of the presence of Gaussian noise with λ term, the numerical result in Fig. 3(b) is completely different from the numerical result in Fig. 1(b).
Fig. 5(b) draws the MPPP for different initial values x 0 .From this figure, we observe that the maximum value of the stationary density function p(X, t) is located at the maximum likely stable state X m (t) = 9.0846 with the initial condition p(X, 0) = 40 π e −40(X−x 0 ) 2 .As the initial condition x 0 increases, it raises the peak point of the stationary density function p(X, t).This shows that the extinction of the species may not happen, and the high peak occurs at the maximum likely stable state X m (t).x 0 =10 x 0 =11 x 0 =14 x 0 =15 Figure 4: Most probable orbits and most probable equilibrium states for stochastic nonlinear system (1).(a) When the initial condition x 0 is less than the unstable equilibrium state X 2 , i.e., x 0 < X 2 .(b) When the initial condition is between X 2 and X 3 , i.e., X 2 < x 0 < X 3 .(c) When the initial condition is greater than X 3 , i.e., x 0 > X 3 .Parameters s = 1, γ 2 = 0.1, γ 3 = 2.67, γ 4 = 1, α = 1.5, ǫ = 0.5, β = 0.27 < 1, λ = 0, and bifurcation time at 1.13 (dot vertical line).
The most probable trajectories of the stochastic single-species model (1) with Allee effect are plotted graphically in figures 4 and 5 (a).Here the values of the noise intensities are set up as λ = 0 and ǫ = 0.5, respectively.We choose the stability index α = 1.5, and the interval D = (0, 15).These figures evolve as the initial value x 0 changes, and they tell us that the maximal likely equilibrium state (maximizer) X m (t) lies between 9 and 10 at the bifurcation time 1.13.In other words, the maximizer in high concentration is between 9 and 10, it is different from the deterministic equilibrium stable solution X 3 = 9.0846 due to the effect of external noises.As seen in Fig. 5(a), the most probable growth state is attracted to the maximal likely equilibrium state in the extinction state, and then it leads to the maximal likely equilibrium state in the high concentration as time moves forward.

Conclusion
In the present work, we have studied the Onsager-Machlup functional and most probable phase portraits for the stochastic growth model (1) for single-species population with strong Allee effects driven by Lévy noise.We have focused on the effect of different values of the initial condition on the MPPP of the nonlinear dynamical system.Small disturbances may cause a transition between the extinction stable state X 1 and the upper equilibrium state X 3 , thus we have used a deterministic quantity, namely the maximal likely trajectory to analyze the transition phenomena in a jump stochastic environment.
In order to find the most likely pathways in transition phenomena, we have calculated the most probable paths of stochastic differential equation in (1) using the stationary density function of the non-local Fokker-Planck equation associated with a non-local partial differential equation.We have investigated the impact of the deterministic parameters, noise intensities and domain size on the FPE.We also have studied the dependence of the probability density on the initial condition x 0 .Our finding has displayed that the maximum of the stationary density function is located at the most probable stable equilibrium state X m .
In conclusion, the most probable path has been used as an indicator that helps the researcher to understand the stochastic dynamics of the single-species model (1) based on the evolution of the probability density function over time.

Figure 3 :
Figure 3: (a) Most probable transition pathways Z m (t) starting at the extinction state X 1 = 0 and ending at the upper equilibrium stable state X 3 = 9.0846 under white noise with respect to the different values of λ.(b) The most probable steady state x s versus the attack rate γ 3 for different values of Gaussian noise intensity with λ term.

Table 1 :
Biological meaning of the parameters and variables in the single-species model.