# Most Probable Dynamics of the Single-Species with Allee Effect under Jump-
diffusion Noise 11footnotemark: 1

Almaz Tesfay<EMAIL_ADDRESS>Shenglan Yuan<EMAIL_ADDRESS>Daniel Tesfay<EMAIL_ADDRESS>James Brannan<EMAIL_ADDRESS>School
of Mathematics and Statistics & Center for Mathematical Sciences, Huazhong
University of Science and Technology, Wuhan 430074, China  Department of
Mathematics, Mekelle University, P.O.Box 231, Mekelle, Ethiopia  Department of
Mathematical Sciences, Clemson University, Clemson, South Carolina 29634, USA

###### Abstract

For a stochastic single-species model with Allee effect, we compute the most
probable phase portrait (MPPP) using the non-local Fokker-Planck equation.
This bistable stochastic model is driven by multiplicative Lévy noise as well
as white noise. The two fixed points are stable equilibria and one is unstable
state between them. In order to study biological behavior of species, 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. By
calculating the solution of the non-local Fokker-Planck equation corresponding
to the population system of single-species model, the maximum possible path is
obtained and 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. Our
numerical simulation shows that: (i) The maximum of the stationary density
function is located at the most probable stable equilibrium state when non-
Gaussian noise is presented in the system; (ii) When the initial value
increases from extinction state to the upper stable state, the most probable
trajectories converge to the maximal likely equilibrium state (maximizer), in
our case lies between 9 and 10; (iii) As time goes on, the most probable paths
increase to stable state quickly, and remain a nearly constant level, then
approach to the upper stable equilibrium state. These numerical experiment
findings help researchers for further experimental study, in order to achieve
good knowledge about dynamic systems in biology.

###### keywords:

Single-species model, most probable phase portrait, jump-diffusion processes,
Onsager–Machlup function, Fokker-Planck equation, extinction probability.

###### MSC:

[2020] 39A50, 45K05, 65N22.

††journal: …

##  1 Introduction

Single-species dynamics have always been a core research area in theoretical
ecology. Studying single-species dynamics enables us to identify the
conditions of extinction and the persistence of the species. One of the
original and still strong motivations for developing mathematical models is to
understand the cause of cycles in particular populations [1].

Population modeling is also important for species management: managing
fisheries for the highest possible sustainable yield, developing recovery
plans for species threatened by extinction, and trying to contain or prevent
the spread of invasive species [2, 3, 4].

In literature, one can find several models of the dynamical single-species
growth system. Gompertz growth model [5], Verhulst growth model with or
without Allee effect [6], power law growth model [7], the interconnections
between deterministic and stochastic system [8], Gilpin–Ayala model [9] are
only a few to mention. In this study, we investigate a single-species system
focusing on the Verhulst growth model with Allee effect which was developed by
Y. Jin [10].

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 [11]. Extensive investigations in modeling and
analysis of random fluctuations [12, 14] in biological dynamical systems have
been on going for long time now. The studies of events in population such as
persistence, extinction, and stationary distribution in stochastic single-
species models constitute an interesting and important research field. One of
the hot issues in population dynamics is developing sufficient condition for
the persistence of biological species as mentioned in [15, 16, 17] and the
references therein.

The population may be affected by sudden environmental noises [18]. For
example, such as severe acute respiratory syndrome (SARS), human
immunodeficiency virus (HIV), the smoking habit [19], and the recent COVID-19
[20], earthquakes, temperature, and hurricanes. 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 [21, 22, 23, 24, 25]. 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 [26]. To the best of our knowledge,
there have been only few studies conducted so far on dynamical systems in
which the noise source is a Lévy process [27, 28, 9]. Introducing the non-
Gaussian noise into the biological system to simulate impulses caused by
external disturbances (mentioned above) is more appropriate and more near to
reality than the white noise. The investigation of the single-species model is
still in its infancy even though noisy fluctuations that naturally portray
random intermittent jumps. Lévy noise is widely applied in studying natural
and man made phenomenon in science, among which we mention, biology [29, 30],
physics [31], and economics [32, 33].

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 impact of Lévy noise fluctuation on the system (1). Investigating
the impact of noisy fluctuations serves in demonstrating the intricate
interactions between the single-species models and their complex surroundings.
In here, we study how Allee effects and stochasticity combine to affect
population persistence. To find the numerical solutions for the Fokker-Planck
equation determined by non-local differential equation with symmetric
$\alpha$-stable Lévy motion, we apply a finite difference method that is
studied by [34].

The most probable phase portrait was first proposed by J. Duan [35, Section
5.3.3]. Cheng [36] obtained the analytical results of the MPPP and showed that
the MPPP can give some information about the propagation of stochastic
dynamics in 1-dimensional model. Wang [37] studied the stochastic bifurcation
by using the qualitative changes of the MPPP to a stochastic system driven by
a stable Lévy noise. In [38], the scholars investigated the most probable
trajectories of the tumor system growth with immune surveillance under
correlated Gaussian noises, and derived analytical solution of the most
probable steady state by using the extremum theory with the local
Fokker–Planck equation in the system. The Onsager–Machlup function for a
stochastic model driven by non-Gaussian noise and Gaussian noise was
established in [39]. The authors also examined the corresponding MPPP of the
stochastic dynamical systems. Cheng et al. [40] focused on the impact of white
noise and jump stable Lévy noise in a genetic regulatory system. Furthermore,
they obtained the most probable transition pathway by minimizing the
Onsager–Machlup action functional for the stochastic model under Gaussian
noise.

Motivated by the above works, we are interested in studying the most probable
trajectories of the single-species growth system under the Gaussian and non-
Gaussian noises. Consequently, the most probable paths can be introduced to
explore the dynamics of the stochastic system. We are able to take the MPPP as
a deterministic geometric tool and give the evolution trajectory of high
probability events (common events). Therefore in here, our main focus will be
investigating how the most probable trajectories escape to the extinction
state from the single-species state more quickly.

Consider the following stochastic single-species growth model with Allee
effect:

$dX_{t}=X_{t-}\left[\left(s-\gamma_{2}\,X_{t-}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t-}+1}\right)dt+\lambda
dB_{t}+\int_{\mathbb{Y}}\epsilon(y)\tilde{N}(dt,dy)\right],\quad t\geq 0,\quad
X_{0}=x_{0},$ (1)

where $X_{t-}$ is the left limit of the population size $X_{t}$.

Table 1: Biological meaning of the parameters and variables in the single-species model. Parameter | Definition
---|---
$s$ | The growth rate
$\gamma_{2}$ | Intraspecific competition rate
$\gamma_{3}$ | The attack rate
$\gamma_{4}$ | Represents the handling time of predator
$M=s/\gamma_{2}$ | The carrying capacity
$t$ | Time

Stochastic force $\tilde{N}(dt,dy)=N(dt,dy)-\nu(dy)dt$ is a compensated
Poisson random measure with associated Poisson random measure $N(dt,dy)$ and
intensity measure $\nu(dy)dt$, in which $\nu(dy)$ is Lévy measure on a
measurable subset $\mathbb{Y}$ of $(0,\infty)$ with $\nu(\mathbb{Y})<\infty.$

The following restriction on system (1) is natural for biological meaning:

$1+\epsilon(y)>0,\quad y\in\mathbb{Y}.$

When $\epsilon(y)>0$, the perturbation stands for increasing of the species
(e.g., planting), while $\epsilon(y)<0$ stands for decreasing (e.g.,
harvesting and epidemics).

The main aim of this study is to investigation 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 which the coefficients and inputs are stochastic processes. The
two stochastic processes germane to this study are Brownian motion and the
Lévy process. Brownian motion describes random fluctuations that are
continuous in time; see Subsection 2.1. The 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 Allee effect influenced by Gaussian and non-Gaussian
noises. Firstly, we review the deterministic model, calculate its equilibrium
solutions (fixed points) and describe the behavior of the fixed points.
Secondly, we obtain the maximum possible paths, and the corresponding maximum
possible stable equilibrium states, that is, 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 Fokker-
Planck equation. To solve the non-local partial differential equation, we use
finite difference method that was proposed in [34]. This method help 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
$\alpha$-stable Lévy motion $L_{t}^{\alpha}$. In the third section, we discuss
the formulation and analysis of the deterministic model of the single-species
model (3) with Allee effect . In the fourth section, we explain the analysis
of the stochastic single-species model with Allee effect. We also review the
definition of the Onsager-Machlup function and most probable phase portraits
in Subsection 4.1 and 4.2, respectively. The numerical results and the
biological implication of our experimental findings are given in the fifth
section. We conclude our research by giving a brief summary in the last
section.

##  2 Preliminaries

Under this section, we define the one-dimensional Brownian motion starting at
time $t=0$ to be a process $B_{t}$ and $\alpha$-stable Lévy motion
$L_{t}^{\alpha}$, which constitute a class of stochastic processes that have
independent and stationary increments as defined below. Throughout this study,
we denote $\mathbb{R}_{+}=(0,\infty)$, $\mathbb{R}=(-\infty,\infty)$, and
$X_{t}\in\mathbb{R}_{+}$, for $t\geq 0$.

###  2.1 Brownian motion

Brownian motion $B_{t}$ (also called Wiener process), [41, P. 15],[42, 43] is
a one-dimensional stochastic process, defined on complete probability space
$(\Omega,\mathfrak{F},\mathfrak{F}_{t},\mathbb{P})$ that also has independent
and stationary increments. 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\leq s\leq t$;
(iii) The process starts at the origin, i.e., $B_{0}=0$ almost surely;
(iv) $B_{t}$ has independent increments, $B_{t}-B_{s}$, for $s<t$, is
independent of the past.

###  2.2 The $\alpha$-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 $\mathbb{R}$ is defined by a drift coefficient
$\hat{b}\in\mathbb{R}$, $\hat{Q}\geq 0$ and a Borel measure $\nu$ defined on
${\mathbb{R}}\setminus{\\{0\\}}$. The triplet $(\hat{b},\hat{Q},\nu)$ 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 $B_{t}$ and a
pure jumping process [35, 44, 45], i.e., $L_{t}$ can be expressed as

$L_{t}=\hat{b}t+B_{\hat{Q}}(t)+\int_{|y|<1}y\tilde{N}(t,dy)+\int_{|y|\geq
1}yN(t,dy),$

where $N(t,dy)$ is the independent Poisson random measure on
$\mathbb{R}^{+}\times{{\mathbb{R}}\setminus{\\{0\\}}}$,
$\tilde{N}(t,dy)=N(t,dy)-\nu(dy)dt$ is the compensated Poisson random measure,
$\nu(S)=\mathbb{E}(N(1,S))$ is the jump measure, and $B_{t}$ is the
independent standard Brownian motion.

The Lévy-Khinchin formula for Lévy motion has a specific form for its
characteristic function.

$\mathbb{E}[e^{(i\,\xi{L_{t}})}]=e^{(t\phi(\xi))},\quad 0\leq t<\infty,$

where

$\phi(\xi)=i\,\xi\hat{b}-\frac{\hat{Q}}{2}\,\xi^{2}+\int_{{\mathbb{R}}\setminus{\\{0\\}}}(e^{i\,\xi
z}-1-i\,\xi zI{{}_{|z|<1}})\nu_{\alpha}(dz),\quad\xi\in\mathbb{R}.$

A stable distribution $S_{\alpha}(\theta,\beta,\gamma)$ is the distribution
for a stable random variable, where the stability index $\alpha\in(0,2)$, the
skewness $\beta\in(0,\infty)$, the shift $\gamma\in(-\infty,\infty)$, and
scale index $\theta\geq 0$. A $\alpha$-stable Lévy motion $L_{t}^{\alpha}$
[13, 46, 47] is a non-Gaussian stochastic process satisfying
(i) the random variables $L_{t_{i+1}}^{\alpha}-L_{t_{i}}^{\alpha}$ are
independent for $0\leq t_{1}<t_{2}<...<t_{i-1}<t_{i}<\infty,$ and for each
$i=1,2,...$;
(ii) $L_{t}^{\alpha}$ has stochastically continuous sample paths, i.e., for
$0\leq s\leq t$ and $\delta>0$, the probability
$\mathbb{P}(|L_{t}^{\alpha}-L_{s}^{\alpha}|>\delta)\rightarrow 0$ as
$t\rightarrow s$;
(iii) $L_{0}^{\alpha}=0$, almost surely;
(iv) $L_{t}^{\alpha}-L_{s}^{\alpha}$ and $L_{t-s}^{\alpha}$ have the same
distribution $S_{\alpha}((t-s)^{1/{\alpha}},0,0)$.

In the case of a one-dimensional isotropic $\alpha$-stable Lévy motion, the
Lévy triplet has the drift vector $\hat{b}=0$ and the diffusion $\hat{Q}=0$.
In this study, we focus on jump process with a specific size in generating
triplet $(0,0,\nu_{\alpha})$ for the random variable $S_{\alpha}$ which can be
defined by $\Delta
L_{t}^{\alpha}=L_{t}^{\alpha}-L_{t-}^{\alpha}<\infty,\,\,t\geq 0,$ where
$L_{t-}^{\alpha}$ is the left limit of the $\alpha$-stable Lévy motion in
$\mathbb{R}$ at any time $t$. Here
$\nu_{\alpha}(dz)=c({\alpha})\frac{1}{|z|^{1+\alpha}}dz$ is Lévy measure with
$c_{\alpha}=\alpha\frac{\Gamma(\frac{1+\alpha}{2})}{{2^{1-\alpha}\pi^{\frac{1}{2}}}\Gamma{(1-\frac{\alpha}{2})}}$,
and $\Gamma$ is the Gamma function.

###### Remark 1

A special case of $\alpha$-stable Lévy motion is Brownian motion when $\alpha$
= 2. Poisson process, $\alpha$-stable process, compound Poisson process, etc.
are also examples of Lévy processes [44].

##  3 Dynamical analysis of the deterministic model

The deterministic form of the model (1) without noise is given as

$\frac{dX_{t}}{dt}=X_{t}\left(s-\gamma_{2}\,X_{t}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t}+1}\right)=F(X_{t}),\quad
t\geq 0,\quad X_{0}=x_{0}.$ (2)

This system can be written as $\frac{dX_{t}}{dt}=-\frac{dU(X)}{dX}$, where
$U(X)$ is the potential function given by

$\displaystyle
U(X)=-\frac{s\,X^{2}}{2}+\frac{\gamma_{2}\,X^{3}}{3}+\frac{\gamma_{3}}{(\gamma_{3}\,\gamma_{4})^{2}}[\gamma_{3}\,\gamma_{4}\,X+1-\ln(\gamma_{3}\,\gamma_{4}\,X+1)].$

The single species model (2) with Allee effect has equilibrium points
$X_{1}=0$ and

$\displaystyle X_{2,3}=$
$\displaystyle\frac{(s\gamma_{3}\gamma_{4}-\gamma_{2})\pm\sqrt{(s\gamma_{3}\gamma_{4}-\gamma_{2})^{2}-4\gamma_{2}\gamma_{3}\gamma_{4}(\gamma_{3}-s)}}{2\gamma_{2}\gamma_{3}\gamma_{4}}$
$\displaystyle=\frac{(s\gamma_{3}\gamma_{4}-\gamma_{2})\pm(s\gamma_{3}\gamma_{4}-\gamma_{2})\,\sqrt{1-\beta}}{2\gamma_{2}\gamma_{3}\gamma_{4}}$
$\displaystyle=\frac{(s\gamma_{3}\gamma_{4}-\gamma_{2})\,\Big{(}1\pm\,\sqrt{1-\beta}\,\,\Big{)}}{2\gamma_{2}\gamma_{3}\gamma_{4}},$

where
$\beta=\frac{4\gamma_{2}\gamma_{3}\gamma_{4}}{(s\gamma_{3}\gamma_{4}-\gamma_{2})^{2}}\,(\gamma_{3}-s).$

If $\beta<1$, then the equilibrium states of system (2) are

$\displaystyle X_{1}=0,\quad\text{an extinction equilibrium;}$ $\displaystyle
X_{2}=\frac{(s\gamma_{3}\gamma_{4}-\gamma_{2})\,\Big{(}1-\sqrt{1-\beta}\Big{)}}{2\gamma_{2}\gamma_{3}\gamma_{4}},\quad\text{a
lower unstable equilibrium;}$ $\displaystyle
X_{3}=\frac{(s\gamma_{3}\gamma_{4}-\gamma_{2})\,\Big{(}1+\sqrt{1-\beta}\Big{)}}{2\gamma_{2}\gamma_{3}\gamma_{4}},\quad\text{an
upper stable equilibrium}.$

If $\beta=1$ the single-species deterministic model has only two equilibria
states:

(a)

(b)

Figure 1:  (a) Numerical simulation of the potential function of the
deterministic model (2). Dashed black lines: local unstable and stable
equilibria at $X_{2}=2.6159$ and $X_{3}=6.3841$, respectively. (b) The
phaselines of the single-species model (2). Parameters
$s=1,\,\gamma_{2}=0.1,\,\gamma_{3}=2.67,\,\gamma_{4}=\frac{1}{\gamma_{3}}.$

$\displaystyle\text{stable
state}~{}~{}X_{1}=0,\quad\text{and}\quad\text{unstable
state}~{}~{}X_{4}=\frac{s\gamma_{3}\gamma_{4}-\gamma_{2}}{2\gamma_{2}\gamma_{3}\gamma_{4}}.$

The derivative of $F(X_{t})$ is

$\displaystyle
s-2\gamma_{2}X-\gamma_{3}\frac{1}{(\gamma_{3}\,\gamma_{4}\,X+1)^{2}}.$

For simplicity and convenience of discussion, we choose the parameters
$\gamma_{3}\,\gamma_{4}=1$, $s=1$, $0<\gamma_{2}<1$,
$0<\gamma_{3}<\frac{(1+\gamma_{2})^{2}}{4\,\gamma_{2}}$, then
$\beta=\frac{4\gamma_{2}(\gamma_{3}-1)}{(1-\gamma_{2})^{2}},$ and
$X_{4}=\frac{1-\gamma_{2}}{2\,\gamma_{2}}$. For $\beta<1$, the extinction
state, $X_{1}=0$, and the equilibrium solution $X_{3}$ are stable but $X_{2}$
is unstable. Figure 1(b) shows that when the value of attack rate $\gamma_{3}$
increases, the unstable state $X_{2}$ and stable state $X_{3}$ get more close
to each other, and become one solution and finally disappear indicating the
occurrence of the saddle-node bifurcation.

The critical value of attack rate $\gamma_{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 for deterministic single-
species growth model. The steady state (extinction state) $X_{1}$ is stable if
$\gamma_{3}>\gamma_{c}$, and the steady state $X_{3}$ exhibits the stability
property for $\gamma_{3}<\gamma_{c}$.

##  4 Dynamical analysis of the stochastic system

In this section, we discuss the behavior of the solution of the system
presented in Eq. (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 put the numerical experiment
findings using finite difference method [34]. Hence, the numerical solution of
the stochastic model provides useful information for understanding the
dynamical behavior of the system.

###  4.1 Onsager-Machlup functional

The Onsager–Machlup functional defines a probability density for a stochastic
process in which the probability density is estimated implicitly, and it is
similar to the Lagrangian of a dynamical system. It can be used for purposes
of reweighting and sampling trajectories, as well as for 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.

The stochastic single-species system (1) with Allee effect, as proved by Jin
[10], 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$\grave{\textup{a}}$dl$\grave{\textup{a}}$g; see Figure 2. Denote the space
of càdlàg paths starting at $x_{0}$ of a solution process $X=\\{X_{t},t\geq
0\\}$ of (1) by

$\mathcal{D}_{x_{0}}=\\{X:\ \mbox{for any}~{}t\geq 0,\ \lim_{s\uparrow
t}X_{s}=X_{t-},\ \lim_{s\downarrow t}X_{s}=X_{t}\ \mbox{exist}\ \mbox{and}\
X_{0}=x_{0}\\}.$

This space equipped with Skorokhod’s $\mathcal{J}_{1}$-topology generated by
the metric $\rm{d}_{\mathbb{R}^{+}}$ is a Polish space [55]. For functions
$x_{1},x_{2}\in\mathcal{D}_{x_{0}}$, define

$\displaystyle{\rm
d}_{\mathbb{R}^{+}}(x_{1},x_{2})=\inf\Big{\\{}\varepsilon>0:|x_{1}(t)-x_{2}(\bar{\lambda}\,t)|\leq\varepsilon,\
\left|\ln\frac{\arctan(\bar{\lambda}\,t)-\arctan(\bar{\lambda}\,s)}{\arctan(t)-\arctan(s)}\right|\leq\varepsilon,\Big{.}$
$\displaystyle\left.\mbox{for every}~{}t,s\geq 0\ \mbox{and some}\
\bar{\lambda}\in\Lambda^{\mathbb{R}^{+}}\right\\},$

where

$\Lambda^{\mathbb{R}^{+}}=\\{\bar{\lambda}:\mathbb{R}^{+}\rightarrow\mathbb{R};\bar{\lambda}\
\mbox{is injective increasing},\ \lim\limits_{t\rightarrow
0}\bar{\lambda}(t)=0,\
\lim\limits_{t\rightarrow\infty}\bar{\lambda}(t)=\infty\\}.$

We consider the corresponding jump-diffusion process
$X_{t}(\omega):=\omega(t),t\in[0,T]$ defined on the canonical probability
space $({\mathbb{R}}^{[0,T]},\mathcal{B}(\mathbb{R})^{[0,T]},\mathbb{P}_{T})$.
Since the paths of $X$ are càdlàg, we identify $X_{t}$ on the space
$(\mathcal{D}_{x_{0}}^{T},\mathcal{B}_{x_{0}}^{T},\mathbb{P})$ instead of
$({\mathbb{R}}^{[0,T]},\mathcal{B}(\mathbb{R})^{[0,T]},\mathbb{P}_{T})$, where
$\mathcal{D}_{x_{0}}^{T}$ is defined similarly as the space
$\mathcal{D}_{x_{0}}$ on the time interval $[0,T]$. The associated Borel
$\sigma$-algebra is
$\mathcal{B}_{x_{0}}^{T}=\mathcal{B}(\mathbb{R})^{[0,T]}\cap\mathcal{D}_{x_{0}}^{T}$,
and then $(\mathcal{D}_{x_{0}}^{T},\mathcal{B}_{x_{0}}^{T})$ is a separable
metric space [48, Section A.2]. The probability measure $\mathbb{P}$ is
generated by $\mathbb{P}(A\cap\mathcal{D}_{x_{0}}^{T}):=\mathbb{P}_{T}(A)$ for
each $A\in\mathcal{B}(\mathbb{R})^{[0,T]}$. Because every càdlàg function on
$[0,T]$ is bounded, we equip $\mathcal{D}_{x_{0}}^{T}$ with the uniform norm

$\|x\|=\sup_{t\in[0,T]}|x(t)|,\quad x(t)\in\mathcal{D}_{x_{0}}^{T}.$

Hence, $\mathcal{D}_{x_{0}}^{T}$ is a Banach space. We are interested in
finding the most probable tube of $X_{t}$, so it is meaningful to determine
the probability that paths lie within the closed tube

$K(z,\varepsilon)=\\{x\in\mathcal{D}_{x_{0}}^{T}:\,\parallel
x-z\parallel\leq\varepsilon,\,z\in\mathcal{D}_{x_{0}}^{T},\,\varepsilon>0\\}.$
(3)

It is worthwhile to note that this probability can be calculated or estimated
by using the induced measure in function space [39]. Define the measure
$\mu_{X}$ on $\mathcal{B}(\mathbb{R})^{[0,T]}$ induced by the solution process
$X_{t}$ for the model (1) via

$\mu_{X}(B)=\mathbb{P}(\\{w:X_{t}(\omega)\in
B\\}),\quad\text{for}\,\,B\in\mathcal{B}(\mathbb{R})^{[0,T]}.$

Using the induced measure, once an $\varepsilon>0$ is given, we are able to
investigate the probabilities of tubes of the same ‘thickness’ for all
$z\in\mathcal{D}_{x_{0}}^{T}$ by

$\mu_{X}(K(z,\varepsilon))=\mathbb{P}(\\{w:X_{t}(\omega)\in
K(z,\varepsilon)\\}),$ (4)

when $K(z,\varepsilon)\in\mathcal{B}(\mathbb{R})^{[0,T]}$. As the tube
$K(z,\varepsilon)$ depends on the choice of a function $z(t)$, we have to look
for that function $z(t)$ which maximizes Eq. (4). If we restrict ourselves on
differentiable functions $z(t)$, then the following definition makes sense
[49].

###### Definition 1

Let $\varepsilon>0$ be given. Consider a tube surrounding a reference path
$z(t)$. If for $\varepsilon$ sufficiently small we estimate the probability of
the solution process $X_{t},t\in[0,T]$ lying in this tube in the form :

$\mathbb{P}(\parallel X-z\parallel\leq\varepsilon)\varpropto
C(\varepsilon)\exp\left\\{-\frac{1}{2}\int_{0}^{T}\textup{OM}(\dot{z},z)dt\right\\},$

then the integrand $\textup{OM}(\dot{z},z)$ is called Onsager-Machulup
function. Where $\varpropto$ denotes the equivalence relation for
$\varepsilon$ small enough. We call $\int_{0}^{T}\textup{OM}(\dot{z},z)dt$ the
Onsager-Machulup functional. In analogy to classical mechanics, we also call
the OM function the Lagrangian function, and the OM functional the action
functional.

###### Remark 2

In particular, for an SDE with pure jump Lévy noise, Definition 1 would be
still applicable, and the minimizer of the Onsager-Machulup functional
$\int_{0}^{T}\textup{OM}(\dot{z},z)dt$ gives a notion of most probable path
for this system. In addition, the minimizer $z(t)$ 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 present in the following theorem.

###### Theorem 1

For the stochastic system (1) with the jump measure satisfying
$\int_{\mathbb{Y}}\epsilon(y)\nu(dy)<\infty$, the Onsager-Machlup function
[39] is characterized, up to an additive constant, by:

$\textup{OM}(\dot{z},z)=\left[\frac{\dot{z}-z\left(s-\gamma_{2}\,z-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,z+1}\right)}{\lambda
z}\right]^{2}+s-2\gamma_{2}z-\frac{\gamma_{3}}{(\gamma_{3}\,\gamma_{4}\,z+1)^{2}}+2\frac{\dot{z}-z\left(s-\gamma_{2}\,z-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,z+1}\right)}{\lambda^{2}z}\int_{\mathbb{Y}}\epsilon(y)\nu_{\alpha}(dy),$

where $z(t)\in\mathcal{D}_{x_{0}}^{T}$ 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,\varepsilon)$
defined as (3) can be approximated as follows:

$\mu_{X}(K(z,\varepsilon))\varpropto\mu_{Y^{c}}(K(0,\varepsilon))\exp\left\\{-\frac{1}{2}\int_{0}^{T}\textup{OM}(\dot{z},z)dt\right\\},$

(a) Persistent sample paths of model (1).

(b) Extinct sample paths of model (1).

Figure 2: The numerical simulation of the system (1) when it is persistent or
extinct at different value initial condition $x_{0}$. (a) The initial
condition is 5. (b) The initial condition is 0.3. Parameters $s=1,$
$\gamma_{2}=0.1$, $\gamma_{3}=2.67$, $\gamma_{4}=1$, $\alpha=1.5$,
$\beta=0.2700<1$.

where $Y_{t}^{c}$ is defined by

$dY_{t}^{c}=Y_{t}^{c}\left\\{\lambda\,Y_{t}dB_{t}+\int_{\mathbb{Y}}\epsilon(y)\tilde{N}(dt,dy)\right\\},\quad
t\in[0,T].$

The proof of this theorem is given in [39, Theorem 4.1].

In Gaussian noise case ($\epsilon=0)$, the stochastic single-species model (1)
becomes

$dX_{t}=X_{t}\left[\left(s-\gamma_{2}\,X_{t}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t}+1}\right)dt+\lambda
dB_{t}\right],\quad t\geq 0,\quad X_{0}=x_{0}.$ (5)

We apply Lamperti transforms for solving multiplicative noise SDE [35, Example
6.48]. This method helps us to transform the multiplicative noise into
additive noise. Because numerically solving an additive-noise SDE is usually
easier than solving a SDE in Eq. (5).

Assume $g\in C^{2}(\mathbb{R})$ and define $Y_{t}=g=\ln(X_{t})$. Then the new
SDE has the following form:

$dY_{t}=G(Y_{t})dt+\lambda\,dB_{t},$ (6)

where

$G(Y_{t})=\left(h(X_{t})-\frac{\lambda^{2}}{2}\right)\Big{|}_{\,X_{t}\,=\,\exp(Y_{t})}$

and

$h(X_{t})=s-\gamma_{2}\,X_{t}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t}+1}$

.

(a)

(b)

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 at different values of $\lambda$. (b) The
most probable steady state $x_{m}$ versus the attack rate $\gamma_{3}$ for
different values of Gaussian noise intensity with $\lambda$ term.

The most probable transition pathway $Z_{m}(t)$ can be obtained the following
least action principle

$\displaystyle\delta\int_{0}^{T}\textup{OM}(\dot{z},z)dt=0,$

where the integrand function (Onsager-Machlup function) [40] is given by

$\textup{OM}(\dot{z},z)=\left(\frac{G(z)-\dot{z}}{\lambda}\right)^{2}+\dot{G}(z).$
(7)

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}(t)$, thus Eq. (7) satisfies the following Euler–Lagrange equation

$\displaystyle\frac{d}{dt}\frac{\partial\textup{OM}(\dot{z},z)}{\partial\dot{z}}=\frac{\partial\textup{OM}(\dot{z},z)}{\partial
z}.$ (8)

The most probable transition pathway $Z_{m}(t)$ of system (6) is given by

$\displaystyle\ddot{Z}_{m}(t)=\frac{\lambda^{2}}{2}\ddot{G}(Z_{m})+\dot{G}(Z_{m})\,G(Z_{m}),\qquad
0<t<T,$ $\displaystyle Z_{m}(0)=X_{1},\,\,\,Z_{m}(T)=X_{3}.$ (9)

To solve two-point boundary value problem in Eq. (4.1), we apply the shooting
method proposed in [50].

###  4.2 Most probable phase portraits

As 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)$
($X_{m}(t)=\mathrm{max}_{X\in(0,\infty)}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.

####  4.2.1 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 express the stationary solution of the non-local
Fokker–Planck equation. This makes the estimate of the most probable phase
portrait numerically and algorithmically possible in Lévy noise case.

Let $f:\mathbb{R}\rightarrow\mathbb{R}$ be a smooth function. Assume that the
solution $X_{t}$ of system (1) has a conditional probability density
$p(X,t|x_{0},0)$ [51, 52, 53, 54]. For convenience, we drop the initial
condition and simply denote it by $p(X,t)$. On one hand,

(a) $x_{0}<X_{2}$

(b) $X_{2}<x_{0}<X_{3}$

(c) $x_{0}>X_{m}$.

Figure 4: Most probable orbits and ‘most probable equilibrium states’ for
system (1). (a) When the initial condition $x_{0}$ is less than the unstable
equilibrium state $X_{2}$. (b) When initial condition is between the unstable
equilibrium state and $X_{3}$. (c) When initial condition is greater than
$X_{3}$. Parameters $s=1,$ $\gamma_{2}=0.1$, $\gamma_{3}=2.67$,
$\gamma_{4}=1$, $\alpha=1.5$, $\epsilon=0.5,$ $\beta=0.2700<1$, $\lambda=0$,
and bifurcation time at 1.13 (dot vertical line).

$\mathbb{E}f(X_{t})=\int_{\mathbb{R}}f(X)p(X,t)dX,$

and thus

$\frac{d}{dt}\mathbb{E}f(X_{t})=\int_{\mathbb{R}}f(X)\frac{\partial}{\partial
t}p(X,t)dX.$

On the other hand, by virtue of Itô’s formula,

$\displaystyle df(X_{t})$
$\displaystyle=X_{t-}\left(s-\gamma_{2}\,X_{t-}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t-}+1}\right)f^{\prime}(X_{t})dt+\int_{\mathbb{Y}}\big{(}f(X_{t}+\epsilon(y)X_{t-})-f(X_{t})-\epsilon(y)X_{t-}f^{\prime}(X_{t})\big{)}\nu_{\alpha}(dy)dt.$
(10)

Taking expectation on both sides of (10), we gain

$\displaystyle d\mathbb{E}f(X_{t})$
$\displaystyle=\mathbb{E}\Big{[}X_{t-}\left(s-\gamma_{2}\,X_{t-}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t-}+1}\right)f^{\prime}(X_{t})dt+\int_{\mathbb{Y}}\big{(}f(X_{t}+\epsilon(y)X_{t-})-f(X_{t})-\epsilon(y)X_{t-}f^{\prime}(X_{t})\big{)}\nu_{\alpha}(dy)dt\Big{]}.$
(11)

Noting that the infinitesimal generator of the solution $X_{t}$ for system (1)
is

$\displaystyle Ap(X,t)$
$\displaystyle=X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)\partial_{X}p(X,t)+\int_{\mathbb{Y}}\big{(}f(X+\epsilon(y)X)-f(X)-\epsilon(y)X\partial_{X}p(X,t)\big{)}\nu_{\alpha}(dy).$

The equation (11) is rewritten as

$\displaystyle\frac{d}{dt}\mathbb{E}f(X_{t})$
$\displaystyle=\mathbb{E}\Big{[}X_{t}\left(s-\gamma_{2}\,X_{t}-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X_{t}+1}\right)f^{\prime}(X_{t})+\int_{\mathbb{Y}}\big{(}f(X_{t}+\epsilon(y)X_{t})-f(X_{t})-\epsilon(y)X_{t}f^{\prime}(X_{t})\big{)}\nu_{\alpha}(dy)\Big{]}$
$\displaystyle=\int_{\mathbb{R}}\Big{[}X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)f^{\prime}(X)+\int_{\mathbb{Y}}\big{(}f(X+\epsilon(y)X)-f(X)-\epsilon(y)X\,f^{\prime}(X)\big{)}\nu_{\alpha}(dy)\Big{]}p(X,t)dX.$
(12)

Figure 5:  (a) Most probable orbits and ‘most probable equilibrium states’ for
system (1) with equilibrium state $X_{m}$ between 9 and 10. (b) The solution
of Fokker-Planck equation of the mode (1). The stationary density function of
the FPE has maximum value at the equilibrium state $X_{m}$. The other
parameters are fixed: $s=1,$ $\gamma_{2}=0.1$, $\gamma_{3}=2.67$,
$\gamma_{4}=1$, $\alpha=1.5$, $\epsilon=0.5,$ $\beta=0.2700<1$, $\lambda=0$
and $x_{0}\in(0,10]$, bifurcation time is 1.13.

As a result, the Fokker-Planck equation for the system (1) for the solution
process $X=\\{X_{t},t\geq 0\\}$ with initial condition
$p(X,0)=\sqrt{\frac{40}{\pi}}e^{-40\left(X-x_{0}\right)^{2}}$ is

$\displaystyle\partial_{t}p(X,t)$
$\displaystyle=-\left(s-2\gamma_{2}X-\frac{\gamma_{3}}{(\gamma_{3}\,\gamma_{4}\,X+1)^{2}}\right)p(X,t)-X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)\partial_{X}p(X,t)$
$\displaystyle+\int_{\mathbb{Y}}\big{(}f(X+\epsilon(y)X)-f(X)-\epsilon(y)X\,f^{\prime}(X)\big{)}\nu_{\alpha}(dy)p(X,t).$
(13)

To simulate the non-local Fokker-Planck equation (4.2.1), we apply a numerical
finite difference method given in Gao et al. [34].

If the Lévy motion is replaced by Brownian motion, then local Fokker-Planck
equation has the following form:

$\displaystyle\partial_{t}p(X,t)=-\partial_{X}\left[X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)p(X,t)\right]+\frac{\lambda^{2}}{2}\partial_{XX}[X^{2}\,p(X,t)].$
(14)

The stationary probability density function $p_{s}(X)$ of Eq. (14) can be
solved by

$\displaystyle
0=-\partial_{X}\left[X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)p_{s}(X)\right]+\frac{\lambda^{2}}{2}\partial_{XX}[X^{2}\,p_{s}(X)],$
(15)

or equivalently,

$\displaystyle 0$
$\displaystyle=-\left[X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)p_{s}(X)\right]+\frac{\lambda^{2}}{2}\partial_{X}[X^{2}\,p_{s}(X)],$
$\displaystyle=\left[X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)-\lambda^{2}\,X\right]p_{s}(X)-\frac{\lambda^{2}\,X^{2}}{2}\partial_{X}\,p_{s}(X).$
(16)

Due to the complexity of stationary probability density solution, we take the
extrema of the stationary probability density function (spdf) located at
$X=x_{s}$ directly. In other words, the spdf satisfies
$\partial_{X}(\,p_{s}(x_{s}))=0$. Since, $p_{s}(x_{s})\neq 0,$ Eq. (4.2.1)
becomes

$\displaystyle
X\left(s-\gamma_{2}\,X-\frac{\gamma_{3}}{\gamma_{3}\,\gamma_{4}\,X+1}\right)-\lambda^{2}\,X=0.$
(17)

Eq. (17) is completely different from the equilibrium state of deterministic
model (2) because of the presence of noise with $\lambda$ term. The numerical
solution of Eq. (17) is plotted in Figure 3(b).

##  5 Numerical results and its biological implications

In order to make the readers understand our results more better, we perform
some numerical simulations to illustrate our theoretical results. The
numerical simulations which are very useful in the study of real examples of
population, are given by the following finite difference method [34]. In the
present section, we define 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. 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 these parameters on the
investigated system. For simplicity, we choose four maximum likely pathways
together with the initial conditions. The initial conditions are selected in
different intervals.

While we plot the above Figures, we fix the deterministic parameters $s=1$,
$\gamma_{2}=0.1$, $\gamma_{3}=2.67,$ $\gamma_{4}=1/\gamma_{3},$ the noise
intensity $\epsilon=0.5,$ and the stability index $\alpha=1.5$.

The potential function denoted by $U(X)$ in Figure 1(a) has two stable steady
states $X_{1}$ and $X_{3}$, and an unstable steady state $X_{2}$ for
$\beta<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 $\beta>1,$ the system
(2) has only one equilibrium point, which is the trivial point $X_{1}=0$.

In Figure 1(b), we plot the equilibrium states versus attack rate
$\gamma_{3}$. For $\beta<1$, there exist two stable equilibrium states $X_{1}$
and $X_{3}$ and one unstable equilibrium state $X_{2}.$ While $\beta>1$,
$X_{1}=0$ is the only equilibrium state that is stable. Thus, the parameter
$\beta=1$ is bifurcation parameter value.

The distance between the unstable $X_{2}$ and the stable $X_{3}$ equilibrium
points becomes very small when $\beta$ approaches to 1. This indicates that
the expected time to extinction may be too short; see Figure 1(b).

Figure 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 system are positive, and extinction species occurs when the
initial condition is less than $X_{2}$; see Figure 2(b). While the initial
condition is greater than the value of $X_{2}$, there is stochastic
persistence.

In Figure 3(a), we plot the most probable transition pathways, denoted by
$Z_{m}(t)$ of system (2.2) for different values of $\lambda$. This Figure
tells us that as time increases, the most probable $Z_{m}(t)$ increases to
stable state $X_{3}$ quickly, and remains a nearly constant level, then decays
(approaches) to the high stable equilibrium state. Figure 3(b) demonstrates
the curves for the most probable steady state $x_{s}$ of the stochastic
single-species model with Allee effect driven by Gaussian noise at different
values of the Gaussian noise intensity $\lambda.$ The steady state curves
exhibit a bi-stability in the interval $(a_{1},,a_{2})$. For
$\gamma_{3}>a_{1}$, the stable steady state stays at the extinction state.
While for $\gamma_{3}<a_{1}$, it remains at the stable equilibrium state.
Because of the presence of noise with $\lambda$ term, the numerical result in
Figure 3(b) is completely different from the numerical result in Figure 1(b).

The most probable trajectories of the stochastic single-species (1) with Allee
effect are plotted in Figures 4 and 5 (a). Here the values of the noise
intensities $\lambda=0$ and $\epsilon=0.5,$ respectively, and stability index
$\alpha=1.5,$ and we choose the interval $D=(0,15)$. These Figures evolve as
initial value $x_{0}$ changes, and they tell us that the maximal likely
equilibrium state (maximizer) $X_{m}(t)$ which 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 impact of the noises. Figure 5(a) tells us
that 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 goes on.

Figure 5(b) plots 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)=\sqrt{\frac{40}{\pi}}e^{-40\left(X-x_{0}\right)^{2}}$. As the initial
condition $x_{0}$ increases, the peak point of the stationary density function
$p(X,t)$ increases. This exhibits extinction of the species may not happen,
and the high peak occurs at the maximum likely stable state $X_{m}(t)$.

##  6 Conclusion

A function which summarizes about the behavior of the dynamics of a continuous
stochastic process is known as the Onsager–Machlup function.

In the present work, we have studied the Onsager–Machlup functional and most
probable phase portraits for the stochastic model for single-species
population growth 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 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 our stochastic differential equation in
(1) using the stationary density function of the non-local Fokker-Planck
equation which is associated with a non-local partial differential equation.
We have investigated on the effect of deterministic parameters, noise
intensities and the domain size on the FPE, i.e., we have studied the
dependence of the probability density on the initial condition $x_{0}$. Our
finding displays 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 dynamics of the single-species model based on
the evolution of the probability density function over time.

## Acknowledgments

This work was partly supported by the NSFC grant 12001213.

## References

  * [1] Y. Zhang, J. Lv, X. Zou. ”Dynamics of stochastic single‐species models.” _Math. Methods Appl. Sci._ 43.15, 8728-8735, 2020.
  * [2] L. J. Allen, _An introduction to stochastic processes with applications to biology_. CRC Press, 2010.
  * [3] B. Dennis and L. Assas, S. Elaydi, E. Kwessi, and G. Livadiotis, “Allee effects and resilience in stochastic populations,” _Theor. Ecol._ , vol. 9, no. 3, pp. 323–335, 2016.
  * [4] B. Dennis, “Allee effects: population growth, critical density, and the chance of extinction,” _Nat. Resour. Model._ , vol. 3, no. 4, pp. 481–538, 1989.
  * [5] C. Lo, Stochastic Gompertz model of tumour cell growth, _J. Theoret. Biol_. 248 (2), 317–321, 2007.
  * [6] M. Liu and C. Bai, “A remark on a stochastic logistic model with Lévy jumps,” _Appl. Math. Comput._ , vol. 251, pp. 521–526, 2015.
  * [7] C. Guiot, P.P. Delsanto, A. Carpinteri, N. Pugno, Y. Mansury, and T.S. Deisboeck, ”The dynamic evolution of the power exponent in a universal growth model of tumors”, _J. Theoret. Biol_. 240 (3) (2006) 459–463.
  * [8] R. Baratti, J. Alvarez, S. Tronci, M. Grosso, and A. Schaum, ”Characterization with Fokker-Planck theory of the nonlinear stochastic dynamics of a class of two-state continuous bioreactors”, _J. Process Control_ , 102(20):66-84, 2021
  * [9] X. Zhang, K. Wang, “Stability analysis of a stochastic Gilpin–Ayala model driven by Lévy noise,” _Commun. Nonlinear Sci. Numer. Simul._ , vol. 19, no. 5, pp. 1391–1399, 2014.
  * [10] Y. Jin, “Analysis of a stochastic single species model with Allee effect and jump-diffusion,” _Adv. Differ. Equations_ ,vol. 2020, no. 1, pp. 1–11, 2020.
  * [11] A. Noor, A. Barnawi, and R. Nour, A. Assiri, and M. El-Beltagy, “Analysis of the stochastic population modelwith random parameters,” _Entropy_ , vol. 22, no. 5, p. 562, 2020.
  * [12] C. Braumann, _Introduction to Stochastic Differential Equations with Applications to Modelling in Biology and Finance_ , John Wiley & Sons, 2019.
  * [13] F. Klebaner, _Introduction to stochastic calculus with applications_. World Scientific Publishing Company, 2005.
  * [14] Z. Sun, J. Lv, and X. Zou, “Dynamical analysis on two stochastic single-species models,” _Appl. Math. Lett._ , vol. 99, p. 105982, 2020.
  * [15] S. Elaydi and R. J. Sacker, “Population models with Allee effect: A new model,” _J. Biol. Dyn._ , vol. 4, no. 4, pp. 397–408, 2010.
  * [16] Q. Yang and D. Jiang, “A note on asymptotic behaviors of stochastic population model with Allee effect,” _Appl. Math. Model._ , vol. 35, no. 9, pp. 4611–4619, 2011.
  * [17] M. Krstic and M.Jovanovi, “On stochastic population model with the Allee effect,” _Math. Comput. Model._ , vol. 52, no. 1-2, pp. 370–379, 2010.
  * [18] X. Mao, G. Marion and E. Renshaw, Environmental Brownian noise suppresses explosions in population dynamics, _Stoch. Process. Appl_. 97, 95–110, 2002.
  * [19] A. Zeb, S. Kumar, A. Tesfay, and A. Kumar, ”A stability analysis on a smoking model with stochastic perturbation”, _Int. J. Numerical Methods for Heat & Fluid Flow_, Vol. ahead-of-print No. ahead-of-print. https://doi.org/10.1108/HFF-02-2021-0140, 2021.
  * [20] A. Tesfay, T. Saeed, A. Zeb, D. Tesfay, A. Khalaf, J Brannan, Dynamics of a stochastic COVID-19 epidemic model with jump-diffusion, _Adv. Differ. Equations_ , https://doi.org/10.1186/s13662-021-03396-8, 2021:228, 2021.
  * [21] A. Tesfay, D. Tesfay, J. Brannan, J. Duan. ”A logistic-harvest model with Allee effect under multiplicative noise. _Stochastics Dyn._ , https://doi.org/10.1142/S0219493721500441, 2021..
  * [22] D. Jiang, N. Shi, and Y. Zhao, Mathematical modeling of fish resources harvesting with predator at maximum sustainable yield, _Math. Theory Model_. 6, 22245804, 2016.
  * [23] D. Rahmani and M. Saraj, The logistic modeling population: Having harvesting factor, Yugosl. _J. Oper. Res_. 25, 107–115, 2015.
  * [24] M. Kot, _Elements of Mathematical Ecology_ , Cambridge Univ. Press, 2001.
  * [25] A. Tesfay, D. Tesfay, S. Yuan, J. Brannan, J. Duan, Stochastic bifurcation in single-species model induced by $\alpha$-stable Lévy noise, J. Stat.Mech. Theory and Experiment, http://dx.doi.org/10.1088/1742-5468/ac2255, arxiv.org/abs/2010.06410, 2021.
  * [26] D. Applebaum and M. Siakalli, Asymetric stability properties of stochastic differential equation driven by Lévy noise, _J. Appl. Probab_ 46, 16–29, 2009.
  * [27] L. Meng and Z. Baichuan, A remark on stochastic logistic model with Lévy jumps, _Appl. Math. Comp_. 25, 521–526, 2015.
  * [28] R. Wu, X. Zou and K. Wang, Dynamics of logistic system driven by Lévy noise under regime switching, Electron. _J. Diff. Equ_. 76, 1072–6691, 2014.
  * [29] N. Humphries, N. Queiroz, J. R. Dyer, N. G. Pade, M. K. Musyl, K. M. Schaefer, & D. W. Sims, Environmental context explains Lévy and Brownian movement patterns of marine predators. _Nat. Lett._ , 465(7301), 1066-1069, 2010.
  * [30] A. Cartea and D. del Castillo-Negrete. Fractional diffusion models of option prices in markets with jumps. Phys. A, _Stat. Mech. Appl._ , 374(2):749–763, 2007.
  * [31] M. Hao, J. Duan, R. Song, and W. Xu. Asymmetric non-Gaussian effects in a tumor growth model with immunization. _Appl. Math. Model_ , 38:4428–4444, 2014.
  * [32] T. Koren, A. Chechkin, and J. Klafter. On the first passage time and leapover properties of Lévy mmotion. _Physica A_ , 379:10–22, 2007.
  * [33] T. Srokowski. Asymmetric Lévy flights in nonhomogeneous environments. _J. Stat. Mech. Theory E._ , 2014(5):P05024, 2014.
  * [34] T. Gao, J. Duan, and X. Li, “Fokker–planck equations for stochastic dynamical systems with symmetric Lévy motions,” _Appl. Math. Comput._ , vol. 278, pp. 1–20, 2016.
  * [35] J. Duan, _An introduction to stochastic dynamics_ , vol. 51. Cambridge University Press, 2015.
  * [36] Z. Cheng, J. Duan, and L. Wang, “Most probable dynamics of some nonlinear systems under noisy fluctuations,” _Commun. Nonlinear ence Numer. Simul._ , vol. 30, no. 1-3, pp. 108–114, 2016.
  * [37] H. Wang, X. Chen, and J. Duan, “A stochastic pitchfork bifurcation in most probable phase portraits,” _Int. J. Bifurc. Chaos_ , vol. 28, no. 1, 2018.
  * [38] P. Han, W. Xu, L. Wang, H. Zhang, and S. Ma, “Most probable dynamics of the tumor growth model with immune surveillance under cross-correlated noises,” _Phys. A Statal Mech. its Appl._ , vol. 547, 2020.
  * [39] Y. Chao and J. Duan, “The onsager–machlup function as lagrangian for the most probable path of a jump-diffusion process,” _Nonlinearity_ , vol. 32, no. 10, p. 3715, 2019.
  * [40] X. Cheng, H. Wang, X. Wang, J. Duan, and X. Li, “Most probable transition pathways and maximal likely trajectories in agenetic regulatory system,” _Phys. A Statal Mech. its Appl._ , vol. 531, p. 121779, 2019.
  * [41] X. Mao, _Stochastic differential equations and applications_. Elsevier, 2007.
  * [42] J. Meiss, _Differential dynamical systems_. SIAM, 2007.
  * [43] B. Oksendal, _Stochastic differential equations: an introduction with applications_. Springer Science & Business Media, 2013.
  * [44] A. Tesfay, D. Tesfay, A. Khalaf, and J. Brannan, “Mean exit time and escape probability for the stochastic logistic growth model with multiplicative $\alpha$-stable Lévy noise,” _Stochastics Dyn._ ,https://doi.org/10.1142/S0219493721500167, 2020..
  * [45] A. Khalaf, A. Tesfay, and X. Wang, “Impulsive stochastic volterra integral equations driven by Lévyy noise,” _Bull. Iran. Math. Soc._ , https://doi.org/10.1007/s41980-020-00465-8, 2020..
  * [46] S. Ken-Iti, _Lévy processes and infinitely divisible distributions_. Cambridge university press, 1999.
  * [47] D. Applebaum, _Lévy processes and stochastic calculus_. Cambridge university press, 2009.
  * [48] L. Arnold, “ _Random dynamical systems_ ,” in Dynamical systems, pp. 1–43, Springer, 1995.
  * [49] L. Onsager and S. Machlup, “Fluctuations and irreversible processes,” _Physical Review_ , vol. 91, no. 6, p. 1505, 1953.
  * [50] H. Keller,_Numerical Solution of Two Point Boundary Value Problems, Society for Industrial and Applied Mathematics_ , England, 1976.
  * [51] Y. Zhang, X. Wang, Q. Huang, J. Duan, and T. Li, “Numerical analysis and applications of Fokker-Plancke quations for stochastic dynamical systems with multiplicative $\alpha$-stable noises,” _Appl. Math. Model._ , 2020.
  * [52] X. Sun and J. Duan, “Fokker-planck equations for nonlinear dynamical systems driven by non-Gaussian Lévy processes,” _J. Math. Phys._ , vol. 53, no. 7, p. 072701, 2012.
  * [53] M. Wang and J. Duan, “Existence and regularity of linear nonlocal Fokker–Planck equation with growing drift,” _J. Math. Anal. Appl._ , vol. 449, no. 1, 2016.
  * [54] J. He, J. Duan, and H. Gao, “A nonlocal Fokker–Planck equation for non-Gaussian stochastic dynamical systems,” _Appl. Math. Lett._ , vol. 49, pp. 1–6, 2015.
  * [55] S. Yuan, J. Hu, X. Liu, and J. Duan, “Slow manifolds for dynamical systems with non-Gaussian stable Lévy noise,” _Anal. Appl._ , vol. 17, no. 03, pp. 477–511, 2019.