# Dynamical Analysis of a Parasite-Host Model within Fluctuating Environment

Authors: Kun Xu; Ziang Zhou; Huitao Zhao
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

A parasite-host model within fluctuating environment is proposed. Firstly, the positivity and boundedness of solutions of the model within deterministic environment are discussed, and, then, the asymptotical stability and global stability of equilibria of deterministic model are investigated. Secondly, we show that the stochastic model has a unique global positive solution; furthermore, we show that the stochastic model has a stationary distribution under certain conditions. Finally, we give some numerical simulations to illustrate our analytical results.


---

## Body

## 1. Introduction

Mathematical model plays an ever more important role in the study of epidemiology, which provides understanding of the underlying mechanisms that influence the spread of disease, and, in the process, it suggests control strategies. Classical epidemic model has been extended in many ways to understand mechanisms of disease transmissions; for example, to understand how parasites regulate host populations is an interesting topic in the study of infectious diseases.It is now widely believed that diseases and parasites were responsible for a number of extinctions on islands and on large land masses (see [1] and the references therein). As a result, ecologists acknowledge the importance of disease and parasites in the dynamics of population. Recently, theory on the effects of parasites on host population dynamics has received much attention, such as how the parasite induced reduction of host fecundity and survival rates changes the host population dynamics and how such dynamics are applied to predict threats to biodiversity in general and endangered species in particular [2–4]. In order to understand how six microparasites regulateDaphnia populations and drive the populations to extinction, Ebert et al. [5] proposed the following epidemiological microparasite model:(1)
d
x
d
t
=
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
,
d
y
d
t
=
β
x
y
-
d
+
α
y
,where x and y are the densities of uninfected and infected hosts, respectively; a
>
0 is the maximum per capita birth rate of uninfected hosts; f
>
0 is the relative fecundity of an infected host; c
>
0 measures the per capita density-dependent reduction in birth rate; d
>
0 is the parasite-independent host background mortality; β
>
0 is the infection rate constant; and α
>
0 is the parasite-induced excess death rate.Model (1) uses the classical mass action incidence and has the usual asymptotic behavior. However, this model is deterministic and does consider the effect of environmental noise. As is known to us, real life is full of randomness and stochasticity. So it is important to know whether or not the long-time behavior of the solution for deterministic dynamics system can be changed by stochastic perturbations. In this paper, following the idea of [6–8], we introduce random noise to model (1) as follows:(2)
d
x
t
=
a
x
t
+
f
y
t
1
-
c
x
t
+
y
t
-
d
x
t
-
β
x
t
y
t
d
t
+
σ
1
x
t
d
B
1
t
,
d
y
t
=
β
x
t
y
t
-
d
+
α
y
t
d
t
+
σ
2
y
t
d
B
2
t
.Here we assume that B
i
t (i
=
1,2) are independent Brownian motion and σ
i (i
=
1,2) are the coefficients of the effects of environmental stochastic perturbations on x
(
t
) and y
(
t
).In the following, unless otherwise specified, we assume that(
Ω
,
F
,
{
F
t
}
t
⩾
0
,
P
) is a complete probability space with a filtration {
F
t
}
t
⩾
0 satisfying the usual conditions (i.e., it is increasing and is right continuous while F
0 contains all P-null sets). Let B
i
(
t
) be the Brownian motions defined on this probability space, i
=
1,2. Also let R
+
n
=
{
x
∈
R
n
,
x
i
>
0
f
o
r
a
l
l
1
⩽
i
⩽
n
} and x
(
t
)
=
(
x
t
,
y
t
)
T.The rest of the paper is organized as follows. In Section2, we show the positivity and the boundedness of solutions of deterministic model (1) with positive initial condition; the local and global stability of the equilibria of system (1) is also investigated in this section. In Section 3, we show that the stochastic model (2) has a unique global positive solution; furthermore, we show that the stochastic model (2) has a stationary distribution under certain conditions. In Section 4, we will give some numerical simulations to support the theoretical prediction. In Section 5, a brief discussion is given.


## 2. Deterministic Model

In this section, we first discuss some basic dynamical properties of the deterministic model (1), which is subjected to positive initial conditions: (3)
x
0
≡
x
0
>
0
,
y
0
≡
y
0
>
0
.


### 2.1. Positivity and Boundedness

In this subsection, we study the positivity and boundedness of solutions of system (1) with initial condition (3). Using the analysis methods in [1], we can have the following results.Lemma 1.
Let(
x
(
t
)
,
y
(
t
)
) be the solution of system (1) with initial condition (3). Then (
x
(
t
)
,
y
(
t
)
) is positive and ultimately bounded for all t
⩾
0.Proof.
First, let us considery
(
t
) for t
⩾
0. From the second equation of system (1), we have (4)
y
t
=
y
0
e
∫
0
t
β
x
s
-
d
+
α
d
s
.From the initial condition (3) and (4), we have y
(
t
)
>
0 for all t
⩾
0. To prove x
(
t
)
>
0, we get from (1) that (5)
d
d
t
y
x
=
x
d
y
/
d
t
-
y
d
x
/
d
t
x
2
=
β
y
-
d
+
α
y
x
+
d
y
x
+
β
y
y
x
-
a
y
x
1
+
f
y
x
1
-
c
x
+
y
⩾
α
y
x
+
β
y
y
x
-
a
y
x
1
+
f
y
x
1
-
c
x
+
y
.A comparison argument shows that y
/
x
>
0 for t
⩾
0. It follows that x
(
t
)
>
0 for t
⩾
0.
Then from (1), we have (6)
d
x
+
y
d
t
=
a
x
t
+
f
y
t
1
-
c
x
t
+
y
t
-
d
x
+
y
-
α
y
⩽
a
x
t
+
f
y
t
1
-
c
x
t
+
y
t
.Therefore, (
x
+
y
)
<
1
/
c
+
ϵ for all large t, where ϵ is an arbitrarily small positive constant. Thus x
(
t
), y
(
t
) are ultimately bounded.Lemma 2.
Let(
x
(
t
)
,
y
(
t
)
) be the solution of system (1) with initial condition (3). Then there is t
0
>
0 such that, for t
>
t
0, we have (
x
(
t
)
,
y
(
t
)
)
∈
Δ
=
{
(
x
(
t
)
,
y
(
t
)
)
:
x
>
0
,
y
>
0
,
x
+
y
<
1
/
c
}.Proof.
It is easy to see that if0
<
x
0
+
y
0
<
1
/
c, then 0
<
x
(
t
)
+
y
(
t
)
<
1
/
c for all t
>
0. Moreover, if there exists t
0
>
0, such that 0
<
x
(
t
0
)
+
y
(
t
0
)
<
1
/
c, then 0
<
x
(
t
)
+
y
(
t
)
<
1
/
c for all t
>
t
0. Hence, we assume below that x
0
+
y
0
>
1
/
c.
If the conclusion of the proposition is false, thenx
(
t
)
+
y
(
t
)
⩾
1
/
c for all t
>
0. It follows from the first equation of system (1) that (7)
d
x
d
t
=
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
⩽
-
d
x
.Hence lim
t
→
∞
⁡
x
(
t
)
=
0. Thus, for any sufficient small 0
<
ϵ
<
d
+
α
/
β, there is T
>
0, such that, for t
>
T, we have x
(
t
)
<
ϵ. Then, from the second equation of system (1), we have (8)
d
y
d
t
=
β
x
y
-
d
+
α
y
⩽
ϵ
β
-
d
-
α
y
.Because ϵ
β
-
d
-
α
<
0 for t
>
T, we have lim
t
→
∞
⁡
y
(
t
)
=
0. Therefor lim
t
→
∞
⁡
(
x
(
t
)
+
y
(
t
)
)
=
0
≪
1. This is a contradiction to x
(
t
)
+
y
(
t
)
⩾
1
/
c for all t
>
0. This ends the proof.Lemma 3.
In the interior of the setΔ, system (1) has no closed orbit.Proof.
Let(9)
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
≕
P
x
,
y
,
β
x
y
-
d
+
α
y
≕
Q
x
,
y
.Take B
(
x
,
y
)
=
1
/
x
y as a Dulac function, if follows from (9) that, for (
x
,
y
)
∈
int
⁡
Δ, (10)
∂
B
P
∂
x
+
∂
B
Q
∂
y
=
-
a
c
y
-
a
f
1
-
c
y
x
2
<
0
.From Bendixson-Dulac Theorem [9], we know that system (1) does not have a limit cycle in int
⁡
Δ.


### 2.2. Equilibria and Their Stability

It is easy to see thatO
(
0,0
) is always an equilibrium of system (1). Direct calculation shows that when d
<
a, system has a semisingular equilibrium E
0
(
a
-
d
/
a
c
,
0
).SupposeE
∗
(
x
∗
,
y
∗
) a positive equilibrium of system (1). Then x
∗
=
d
+
α
/
β and y
∗ satisfies the following equation: (11)
A
y
2
+
B
y
+
C
=
0
,where A
=
a
c
f, B
=
(
a
c
+
a
c
f
+
β
)
x
∗
-
a
f, and C
=
x
∗
(
a
c
x
∗
+
d
-
a
). Obviously, if d
<
β
/
c
-
α and d
<
a
[
1
-
c
(
d
+
α
)
/
β
], then C
<
0; thus system (1) has unique positive equilibrium E
∗
(
x
∗
,
y
∗
). If a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a and f
⩽
1, then C
>
0 and B
>
0; thus system (1) has no positive equilibrium.Next, we always assume thatf
⩽
1. The stability of equilibria E
0 and E
∗ is determined by the eigenvalues of the Jacobian matrix (12)
J
x
,
y
=
a
11
x
,
y
a
12
x
,
y
a
21
x
,
y
a
22
x
,
y
,where (13)
a
11
x
,
y
=
a
-
d
-
2
a
c
x
-
a
c
+
a
c
f
+
β
y
,
a
12
x
,
y
=
a
f
-
a
c
+
a
c
f
+
β
x
-
2
a
c
f
y
,
a
21
x
,
y
=
β
y
,
a
22
x
,
y
=
β
x
-
d
-
α
.For equilibriumE
0, the Jacobian matrix in it is (14)
J
E
0
=
-
a
-
d
a
12
x
0
,
0
0
β
a
-
d
a
c
-
d
-
αwith a
12
(
x
0
,
0
)
=
a
f
-
(
a
c
+
a
c
f
+
β
)
a
-
d
/
a
c. It is easy to see that if a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a, J
(
E
0
) have two negative eigenvalues, hence E
0 is asymptotically stable. Notice that the equilibrium O
(
0,0
) is always unstable when d
<
a, and, combined with Lemma 3, we have E
0 is globally stable with a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a. However, if d
<
a
[
1
-
c
(
d
+
α
)
/
β
], J
(
E
0
) has one positive eigenvalue and one negative eigenvalue; hence E
0 is a saddle and is unstable.For equilibriumE
∗, the Jacobian matrix in it is (15)
J
E
∗
=
a
11
x
∗
,
y
∗
a
12
x
∗
,
y
∗
β
y
∗
0
,where (16)
a
11
x
∗
,
y
∗
=
a
-
d
-
2
a
c
x
∗
-
a
c
+
a
c
f
+
β
y
∗
,
a
12
x
∗
,
y
∗
=
a
f
-
a
c
+
a
c
f
+
β
x
∗
-
2
a
c
f
y
∗
.The characteristic equation of J
(
E
∗
) is as follows: (17)
λ
2
-
a
11
x
∗
,
y
∗
λ
-
β
y
∗
a
12
x
∗
,
y
∗
=
0
.On the other hand, from the first equation of system (1), we get (
x
∗
,
y
∗
) satisfies (18)
a
c
+
a
c
f
+
β
x
∗
y
∗
=
a
-
d
-
a
c
x
∗
x
∗
+
a
f
y
∗
1
-
c
y
∗
.From the first equation of (16) and (18), we have (19)
a
11
x
∗
,
y
∗
x
∗
=
a
-
d
-
2
a
c
x
∗
x
∗
-
a
c
+
a
c
f
+
β
x
∗
y
∗
=
a
-
d
-
a
c
x
∗
x
∗
-
a
c
x
∗
2
-
a
-
d
-
a
c
x
∗
x
∗
-
a
f
y
∗
1
-
c
y
∗
=
-
a
c
x
∗
2
-
a
f
y
∗
1
-
c
y
∗
.Noting that x
∗
+
y
∗
⩽
1
/
c, thus we have a
11
(
x
∗
,
y
∗
)
x
∗
<
0; it follows that a
11
(
x
∗
,
y
∗
)
<
0.From the second equation of (16) and (18), we have (20)
a
12
x
∗
,
y
∗
y
∗
=
a
f
y
∗
-
a
c
+
a
c
f
+
β
x
∗
y
∗
-
2
a
c
f
y
∗
2
=
a
f
y
∗
-
2
a
c
f
y
∗
2
-
a
f
y
∗
1
-
c
y
∗
-
a
-
d
-
a
c
x
∗
x
∗
=
-
a
c
f
y
∗
2
-
a
-
d
-
a
c
x
∗
x
∗
.Note that E
∗ coexists with d
<
a
[
1
-
c
(
d
+
α
)
/
β
]. It follows that a
12
(
x
∗
,
y
∗
)
y
∗
<
0, which implies -
β
y
∗
a
12
(
x
∗
,
y
∗
)
>
0. Hence, we conclude that all the roots of (17) have negative real part. It means that the equilibrium E
∗ is asymptotically stable.Notice that the equilibriumO
(
0,0
) is always unstable when d
<
a. For d
<
a
[
1
-
c
(
d
+
α
)
/
β
], E
0 is unstable; combined with Lemma 3, we have E
∗ is globally stable in R
+
2 when d
<
a
[
1
-
c
(
d
+
α
)
/
β
].Based on the above analysis, we have the following.Theorem 4.
For system (1), we have the following results:(i)
Ifd
<
a
[
1
-
c
(
d
+
α
)
/
β
], E
0 is a saddle and unstable. If a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a, E
0 is globally stable in R
+
2.
(ii)
Iff
⩽
1 and d
<
a
[
1
-
c
(
d
+
α
)
/
β
], E
∗ is globally stable in R
+
2.


## 2.1. Positivity and Boundedness

In this subsection, we study the positivity and boundedness of solutions of system (1) with initial condition (3). Using the analysis methods in [1], we can have the following results.Lemma 1.
Let(
x
(
t
)
,
y
(
t
)
) be the solution of system (1) with initial condition (3). Then (
x
(
t
)
,
y
(
t
)
) is positive and ultimately bounded for all t
⩾
0.Proof.
First, let us considery
(
t
) for t
⩾
0. From the second equation of system (1), we have (4)
y
t
=
y
0
e
∫
0
t
β
x
s
-
d
+
α
d
s
.From the initial condition (3) and (4), we have y
(
t
)
>
0 for all t
⩾
0. To prove x
(
t
)
>
0, we get from (1) that (5)
d
d
t
y
x
=
x
d
y
/
d
t
-
y
d
x
/
d
t
x
2
=
β
y
-
d
+
α
y
x
+
d
y
x
+
β
y
y
x
-
a
y
x
1
+
f
y
x
1
-
c
x
+
y
⩾
α
y
x
+
β
y
y
x
-
a
y
x
1
+
f
y
x
1
-
c
x
+
y
.A comparison argument shows that y
/
x
>
0 for t
⩾
0. It follows that x
(
t
)
>
0 for t
⩾
0.
Then from (1), we have (6)
d
x
+
y
d
t
=
a
x
t
+
f
y
t
1
-
c
x
t
+
y
t
-
d
x
+
y
-
α
y
⩽
a
x
t
+
f
y
t
1
-
c
x
t
+
y
t
.Therefore, (
x
+
y
)
<
1
/
c
+
ϵ for all large t, where ϵ is an arbitrarily small positive constant. Thus x
(
t
), y
(
t
) are ultimately bounded.Lemma 2.
Let(
x
(
t
)
,
y
(
t
)
) be the solution of system (1) with initial condition (3). Then there is t
0
>
0 such that, for t
>
t
0, we have (
x
(
t
)
,
y
(
t
)
)
∈
Δ
=
{
(
x
(
t
)
,
y
(
t
)
)
:
x
>
0
,
y
>
0
,
x
+
y
<
1
/
c
}.Proof.
It is easy to see that if0
<
x
0
+
y
0
<
1
/
c, then 0
<
x
(
t
)
+
y
(
t
)
<
1
/
c for all t
>
0. Moreover, if there exists t
0
>
0, such that 0
<
x
(
t
0
)
+
y
(
t
0
)
<
1
/
c, then 0
<
x
(
t
)
+
y
(
t
)
<
1
/
c for all t
>
t
0. Hence, we assume below that x
0
+
y
0
>
1
/
c.
If the conclusion of the proposition is false, thenx
(
t
)
+
y
(
t
)
⩾
1
/
c for all t
>
0. It follows from the first equation of system (1) that (7)
d
x
d
t
=
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
⩽
-
d
x
.Hence lim
t
→
∞
⁡
x
(
t
)
=
0. Thus, for any sufficient small 0
<
ϵ
<
d
+
α
/
β, there is T
>
0, such that, for t
>
T, we have x
(
t
)
<
ϵ. Then, from the second equation of system (1), we have (8)
d
y
d
t
=
β
x
y
-
d
+
α
y
⩽
ϵ
β
-
d
-
α
y
.Because ϵ
β
-
d
-
α
<
0 for t
>
T, we have lim
t
→
∞
⁡
y
(
t
)
=
0. Therefor lim
t
→
∞
⁡
(
x
(
t
)
+
y
(
t
)
)
=
0
≪
1. This is a contradiction to x
(
t
)
+
y
(
t
)
⩾
1
/
c for all t
>
0. This ends the proof.Lemma 3.
In the interior of the setΔ, system (1) has no closed orbit.Proof.
Let(9)
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
≕
P
x
,
y
,
β
x
y
-
d
+
α
y
≕
Q
x
,
y
.Take B
(
x
,
y
)
=
1
/
x
y as a Dulac function, if follows from (9) that, for (
x
,
y
)
∈
int
⁡
Δ, (10)
∂
B
P
∂
x
+
∂
B
Q
∂
y
=
-
a
c
y
-
a
f
1
-
c
y
x
2
<
0
.From Bendixson-Dulac Theorem [9], we know that system (1) does not have a limit cycle in int
⁡
Δ.


## 2.2. Equilibria and Their Stability

It is easy to see thatO
(
0,0
) is always an equilibrium of system (1). Direct calculation shows that when d
<
a, system has a semisingular equilibrium E
0
(
a
-
d
/
a
c
,
0
).SupposeE
∗
(
x
∗
,
y
∗
) a positive equilibrium of system (1). Then x
∗
=
d
+
α
/
β and y
∗ satisfies the following equation: (11)
A
y
2
+
B
y
+
C
=
0
,where A
=
a
c
f, B
=
(
a
c
+
a
c
f
+
β
)
x
∗
-
a
f, and C
=
x
∗
(
a
c
x
∗
+
d
-
a
). Obviously, if d
<
β
/
c
-
α and d
<
a
[
1
-
c
(
d
+
α
)
/
β
], then C
<
0; thus system (1) has unique positive equilibrium E
∗
(
x
∗
,
y
∗
). If a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a and f
⩽
1, then C
>
0 and B
>
0; thus system (1) has no positive equilibrium.Next, we always assume thatf
⩽
1. The stability of equilibria E
0 and E
∗ is determined by the eigenvalues of the Jacobian matrix (12)
J
x
,
y
=
a
11
x
,
y
a
12
x
,
y
a
21
x
,
y
a
22
x
,
y
,where (13)
a
11
x
,
y
=
a
-
d
-
2
a
c
x
-
a
c
+
a
c
f
+
β
y
,
a
12
x
,
y
=
a
f
-
a
c
+
a
c
f
+
β
x
-
2
a
c
f
y
,
a
21
x
,
y
=
β
y
,
a
22
x
,
y
=
β
x
-
d
-
α
.For equilibriumE
0, the Jacobian matrix in it is (14)
J
E
0
=
-
a
-
d
a
12
x
0
,
0
0
β
a
-
d
a
c
-
d
-
αwith a
12
(
x
0
,
0
)
=
a
f
-
(
a
c
+
a
c
f
+
β
)
a
-
d
/
a
c. It is easy to see that if a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a, J
(
E
0
) have two negative eigenvalues, hence E
0 is asymptotically stable. Notice that the equilibrium O
(
0,0
) is always unstable when d
<
a, and, combined with Lemma 3, we have E
0 is globally stable with a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a. However, if d
<
a
[
1
-
c
(
d
+
α
)
/
β
], J
(
E
0
) has one positive eigenvalue and one negative eigenvalue; hence E
0 is a saddle and is unstable.For equilibriumE
∗, the Jacobian matrix in it is (15)
J
E
∗
=
a
11
x
∗
,
y
∗
a
12
x
∗
,
y
∗
β
y
∗
0
,where (16)
a
11
x
∗
,
y
∗
=
a
-
d
-
2
a
c
x
∗
-
a
c
+
a
c
f
+
β
y
∗
,
a
12
x
∗
,
y
∗
=
a
f
-
a
c
+
a
c
f
+
β
x
∗
-
2
a
c
f
y
∗
.The characteristic equation of J
(
E
∗
) is as follows: (17)
λ
2
-
a
11
x
∗
,
y
∗
λ
-
β
y
∗
a
12
x
∗
,
y
∗
=
0
.On the other hand, from the first equation of system (1), we get (
x
∗
,
y
∗
) satisfies (18)
a
c
+
a
c
f
+
β
x
∗
y
∗
=
a
-
d
-
a
c
x
∗
x
∗
+
a
f
y
∗
1
-
c
y
∗
.From the first equation of (16) and (18), we have (19)
a
11
x
∗
,
y
∗
x
∗
=
a
-
d
-
2
a
c
x
∗
x
∗
-
a
c
+
a
c
f
+
β
x
∗
y
∗
=
a
-
d
-
a
c
x
∗
x
∗
-
a
c
x
∗
2
-
a
-
d
-
a
c
x
∗
x
∗
-
a
f
y
∗
1
-
c
y
∗
=
-
a
c
x
∗
2
-
a
f
y
∗
1
-
c
y
∗
.Noting that x
∗
+
y
∗
⩽
1
/
c, thus we have a
11
(
x
∗
,
y
∗
)
x
∗
<
0; it follows that a
11
(
x
∗
,
y
∗
)
<
0.From the second equation of (16) and (18), we have (20)
a
12
x
∗
,
y
∗
y
∗
=
a
f
y
∗
-
a
c
+
a
c
f
+
β
x
∗
y
∗
-
2
a
c
f
y
∗
2
=
a
f
y
∗
-
2
a
c
f
y
∗
2
-
a
f
y
∗
1
-
c
y
∗
-
a
-
d
-
a
c
x
∗
x
∗
=
-
a
c
f
y
∗
2
-
a
-
d
-
a
c
x
∗
x
∗
.Note that E
∗ coexists with d
<
a
[
1
-
c
(
d
+
α
)
/
β
]. It follows that a
12
(
x
∗
,
y
∗
)
y
∗
<
0, which implies -
β
y
∗
a
12
(
x
∗
,
y
∗
)
>
0. Hence, we conclude that all the roots of (17) have negative real part. It means that the equilibrium E
∗ is asymptotically stable.Notice that the equilibriumO
(
0,0
) is always unstable when d
<
a. For d
<
a
[
1
-
c
(
d
+
α
)
/
β
], E
0 is unstable; combined with Lemma 3, we have E
∗ is globally stable in R
+
2 when d
<
a
[
1
-
c
(
d
+
α
)
/
β
].Based on the above analysis, we have the following.Theorem 4.
For system (1), we have the following results:(i)
Ifd
<
a
[
1
-
c
(
d
+
α
)
/
β
], E
0 is a saddle and unstable. If a
[
1
-
c
(
d
+
α
)
/
β
]
<
d
<
a, E
0 is globally stable in R
+
2.
(ii)
Iff
⩽
1 and d
<
a
[
1
-
c
(
d
+
α
)
/
β
], E
∗ is globally stable in R
+
2.


## 3. Stochastic Model

### 3.1. Global Positive Solution

In this subsection we show the solution of system (2) is global and nonnegative. As we have known, in order for a stochastic differential equation to have a unique global (i.e., no explosion in a finite time) solution for any given initial value, the coefficients of the equation are generally required to satisfy the linear growth condition and local Lipschitz condition [10]. However, the coefficients of (2) do not satisfy the linear growth condition, though they are locally Lipschitz continuous, so the solution of (2) may explode at a finite time [10]. In this subsection, using the Lyapunov analysis method (mentioned in [10]), we will show the solution of (2) is positive and global.Theorem 5.
For any given initial value(
x
(
0
)
,
y
(
0
)
)
∈
R
+
2, there is unique positive solution (
x
(
t
)
,
y
(
t
)
) of model (2) on t
⩾
0 and the solution will remain in R
+
2 with probability 1; namely, (
x
(
t
)
,
y
(
t
)
)
∈
R
+
2 for all t
⩾
0 almost surely.Proof.
Since the coefficients of (2) are locally Lipschitz continuous, for any given initial value (
x
(
0
)
,
y
(
0
)
)
∈
R
+
2, there is a unique local solution (
x
(
t
)
,
y
(
t
)
) on t
∈
[
0
,
τ
e
), where τ
e is the explosion time. To show this solution is global, we need to show that τ
e
=
∞ a.s.
First, we provex
(
t
) and y
(
t
) do not explode to infinity in a finite time. Let k
0
>
0 be sufficiently large for x
(
0
)
∈
1
/
k
0
,
k
0 and y
(
0
)
∈
1
/
k
0
,
k
0. For each integer k
⩾
k
0, define the stopping time (21)
τ
k
=
inf
⁡
t
∈
0
,
τ
e
:
x
t
∉
1
k
,
k
o
r
y
t
∉
1
k
,
k
,where throughout this paper we set inf
⁡
∅
=
∞ (∅ denotes the empty set). Clearly, τ
k is increasing as k
→
∞. Set τ
∞
=
lim
k
→
∞
⁡
τ
k; therefore τ
∞
⩽
τ
e a.s. If we can show that τ
∞
=
∞ a.s., then τ
e
=
∞ and (
x
(
t
)
,
y
(
t
)
)
∈
R
+
2 a.s. for all t
⩾
0. In other words, to complete the proof it is required to show that τ
∞
=
∞ a.s. If this statement is false, then there are a pair of constants T
>
0 and ϵ
∈
(
0,1
) such that (22)
P
τ
∞
⩽
T
>
ϵ
.Hence, there is an integer k
1
>
k
0 such that (23)
P
τ
k
⩽
T
⩾
ϵ
∀
k
⩾
k
1
.
Define aC
2-function V
:
R
+
2
→
R
-
+ by (24)
V
x
,
y
=
x
-
1
-
ln
⁡
x
+
y
-
1
-
ln
⁡
y
,and the nonnegativity of this function can be seen from (25)
u
-
1
-
ln
⁡
u
⩾
0
f
o
r
u
>
0
.Let k
⩾
k
0 and T
>
0 be arbitrary. Applying Itô’s formula, we obtain (26)
d
V
x
,
y
=
L
V
x
,
y
d
t
+
σ
1
x
-
1
d
B
1
t
+
σ
2
y
-
1
d
B
2
t
,where L
V
:
R
+
2
→
R
-
+ is defined by (27)
L
V
x
,
y
=
1
-
1
x
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
+
1
-
1
y
β
x
y
-
d
+
α
y
+
1
2
σ
1
2
+
1
2
σ
2
2
⩽
d
+
α
-
a
+
1
2
σ
1
2
+
1
2
σ
2
2
+
a
-
d
+
a
c
-
β
x
-
a
c
x
2
+
a
f
+
a
c
+
a
c
f
+
β
-
d
-
α
y
-
1
2
a
c
f
y
2
-
a
f
y
x
+
1
2
a
c
f
y
2
x
2
≕
K
.Thus, (28)
d
V
x
,
y
⩽
K
d
t
+
σ
1
x
-
1
d
B
1
t
+
σ
2
y
-
1
d
B
2
t
.We can now integrate both sides of (28) from 0 to τ
k
∧
T and then take the expectations (29)
E
V
x
τ
k
∧
T
,
y
τ
k
∧
T
⩽
x
S
0
,
I
0
+
K
E
τ
k
∧
Tso we have (30)
E
V
x
τ
k
∧
T
,
y
τ
k
∧
T
⩽
V
x
0
,
y
0
+
K
T
.Let Ω
k
=
{
τ
k
⩽
T
} for k
⩾
k
1 and, by (23), P
(
Ω
k
)
⩾
ϵ. Note that, for every ω
∈
Ω
k, there is x
(
τ
k
,
ω
) or y
(
τ
k
,
ω
) equaling either k or 1
/
k and therefore V
(
x
(
τ
k
,
ω
)
,
y
(
τ
k
,
ω
)
) is not less than either k
-
1
-
ln
⁡
k or 1
/
k
-
1
-
ln
⁡
1
/
k
=
1
/
k
-
1
+
ln
⁡
K. Hence (31)
V
x
τ
k
,
ω
,
y
τ
k
,
ω
⩾
k
-
1
-
ln
⁡
k
∧
1
k
-
1
+
ln
⁡
K
.It follows from (30) that (32)
V
x
0
,
y
0
+
K
T
⩾
E
1
Ω
k
ω
V
x
τ
k
,
ω
,
y
τ
k
,
ω
⩾
ϵ
k
-
1
-
ln
⁡
k
∧
1
k
-
1
+
ln
⁡
K
,where 1
Ω
k is the indicator function of Ω
k. Letting k
→
∞, we have (33)
∞
>
V
x
0
,
y
0
+
K
T
=
∞
,and it is a contradiction; then we must have τ
∞
=
∞. Therefore, it implies x
(
t
) and y
(
t
) will not explode in a finite time with probability one.


### 3.2. Stochastic Persistence

In this subsection, we are intended to prove the stochastic persistence of the model system (2) under certain parametric restrictions. Stochastic persistence means that if we start from a positive initial condition, that is, from an interior point of the first quadrant, then solution trajectories of the stochastic model will always remain within the interior of the first quadrant and remain bounded at all future time. There are several concepts of stochastic persistence [11]; here we use the notion of stochastic persistence in mean. Before proving the main results of this section we define stochastic persistence in mean.Definition 6.
The population ofx
(
t
) is said to be strongly persistent in the mean if 〈
x
(
t
)
〉
∗
>
0, where (34)
x
t
≔
1
t
∫
0
t
x
s
d
s
,
x
t
∗
≔
l
i
m
i
n
f
t
→
+
∞
1
t
∫
0
t
x
s
d
s
.
The proof of strong persistence result for the stochastic model (2) is based upon the following lemma (see Lemma 4 in [11]).Lemma 7.
Supposex
(
t
)
∈
C
[
Ω
×
R
+
,
R
+
0
], where R
+
0
=
{
a
:
a
>
0
,
a
∈
R
}.(i)
If there are positive constantsμ, T, and λ
⩾
0 such that (35)
ln
⁡
x
t
⩽
λ
t
-
μ
∫
0
t
x
s
d
s
+
∑
i
=
1
n
β
i
B
i
tfor t
⩾
T, where β
i (i
=
1,2
,
…
,
n) are constants, then 〈
x
(
t
)
〉
∗
⩽
λ
/
μ, a
.
s
.
(ii)
If there are positive constantsμ, T, and λ
⩾
0 such that (36)
ln
⁡
x
t
⩾
λ
t
-
μ
∫
0
t
x
s
d
s
+
∑
i
=
1
n
β
i
B
i
tfor t
⩾
T, where β
i (i
=
1,2
,
…
,
n) are constants, then 〈
x
(
t
)
〉
∗
⩾
λ
/
μ, a
.
s
.Remark 8.
〈
x
(
t
)
〉
∗ is defined by (37)
x
t
∗
≔
l
i
m
s
u
p
t
→
+
∞
1
t
∫
0
t
x
s
d
s
.Theorem 9.
Assume thatf
⩽
1. If d
<
a
f
-
α
+
σ
1
2
/
2
+
σ
2
2
/
2, then the total population of x
(
t
) and y
(
t
) of stochastic model system (2) starting from any interior point of first quadrant is strongly persistent in mean. If d
⩾
a
+
σ
1
2
/
2
+
σ
2
2
/
2, then the total population of x
(
t
) and y
(
t
) of stochastic model system (2) starting from any interior point of first quadrant will tend to extinction.Proof.
LetV
(
x
(
t
)
+
y
(
t
)
)
=
ln
⁡
(
x
(
t
)
+
y
(
t
)
) for (
x
,
y
)
∈
R
+
2. Then using Itô’s formula, we get from the first equation of system (2) (38)
d
ln
⁡
x
t
+
y
t
=
a
x
t
+
a
f
y
t
x
t
+
y
t
-
a
c
x
t
+
f
y
t
-
d
-
α
y
t
x
t
+
y
t
+
σ
1
2
2
+
σ
2
2
2
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
⩾
a
f
-
d
-
α
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
.Integrating both sides from 0 to t and dividing by t, we have (39)
l
n
x
t
+
y
t
/
x
0
+
y
0
t
⩾
a
f
-
d
-
α
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
+
σ
1
B
1
t
+
σ
2
B
2
t
t
.From Lemma 7, we have the following inequality: (40)
x
t
+
y
t
∗
⩾
a
f
-
d
-
α
+
σ
1
2
/
2
+
σ
2
2
/
2
a
c
.Hence x
t
+
y
t
∗
>
0, whenever a
f
-
d
-
α
+
σ
1
2
/
2
+
σ
2
2
/
2
>
0.
On the other hand(41)
d
ln
⁡
x
t
+
y
t
=
a
x
t
+
a
f
y
t
x
t
+
y
t
-
a
c
x
t
+
f
y
t
-
d
-
α
y
t
x
t
+
y
t
+
σ
1
2
2
+
σ
2
2
2
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
⩽
a
-
d
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
.Integrating both sides from 0 to t and dividing by t, we have (42)
ln
⁡
x
t
+
y
t
/
x
0
+
y
0
t
⩽
a
-
d
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
+
σ
1
B
1
t
+
σ
2
B
2
t
t
.From Lemma 7, we have the following inequality: (43)
x
t
+
y
t
∗
⩽
a
-
d
+
σ
1
2
/
2
+
σ
2
2
/
2
a
c
.It follows that 〈
(
x
(
t
)
+
y
(
t
)
)
〉
∗
=
0, whenever a
-
d
+
σ
1
2
/
2
+
σ
2
2
/
2
⩽
0.


### 3.3. Existence of Stationary Distribution

In this section, we prove the existence of stationary distribution of prey and predator populations. For this purpose we find the stationary distribution for solutions of system (2), which in turn implies the stability in stochastic sense. Before proving the main theorem related to the stationary distribution we state a useful lemma from [12] which will be useful to prove the theorem.LetX
(
t
) be a homogeneous Markov process defined in the l-dimensional Euclidean space, denoted by E
l, and described by the following system of stochastic differential equation: (44)
d
X
t
=
b
X
d
t
+
∑
r
=
1
k
f
r
X
d
B
r
t
.We assume there exists a bounded domainU
∈
E
l with regular boundary Γ, having the following properties:(P1)
In the domainU and some neighborhood thereof, the smallest eigenvalue of the diffusion matrix A
3
(
x
) is bounded away from zero.
(P2)
Ifx
∈
E
l
∖
U, the mean time τ at which a path emerging from x reaches the set U is finite, and sup
x
∈
S
⁡
E
x
τ
<
∞ for every compact subset S
⊂
E
l.Lemma 10.
If assumptions (P1) and (P2) hold, then the Markov processX
(
t
) has a stationary distribution μ
·. Let g
· be a function integrable with respect to the measure μ. Then (45)
P
X
l
i
m
T
→
∞
∫
0
T
g
X
t
d
t
=
∫
E
l
g
x
μ
d
x
=
1for all x
∈
E
l.Theorem 11.
Assume thatd
<
β
/
c
-
α and d
<
a
[
1
-
c
(
d
+
α
)
/
β
]. Then there exists a stationary distribution μ
(
·
) for the stochastic system (2).Proof.
Ifd
<
β
/
c
-
α and d
<
a
[
1
-
c
(
d
+
α
)
/
β
], then system (1) has unique positive equilibrium E
∗
(
x
∗
,
y
∗
).
Define a positive definite functionV
:
E
2
→
R
+, where E
2
=
I
n
t
(
R
+
2
) as follows: (46)
V
x
,
y
=
x
-
x
∗
-
x
∗
ln
⁡
x
x
∗
+
y
-
y
∗
-
y
∗
ln
⁡
y
y
∗
≔
V
1
+
V
2
.Applying Itô’s formula, we can calculate (47)
d
V
1
=
1
-
x
∗
x
d
x
+
1
2
x
∗
x
2
d
x
2
=
-
a
c
x
-
x
∗
2
+
a
f
x
-
x
∗
y
x
-
y
∗
x
∗
-
a
c
+
a
c
f
+
β
x
-
x
∗
y
-
y
∗
+
x
∗
2
σ
1
2
-
a
c
f
x
-
x
∗
y
2
x
-
y
∗
2
x
∗
d
t
+
σ
1
x
-
x
∗
d
B
1
=
-
a
c
+
a
f
y
x
x
∗
+
a
c
f
y
2
x
x
∗
x
-
x
∗
2
-
a
c
+
a
c
f
+
β
+
a
f
x
∗
+
a
c
f
y
∗
x
∗
x
-
x
∗
y
-
y
∗
+
x
∗
2
σ
1
2
d
t
+
σ
1
x
-
x
∗
d
B
1
,
d
V
2
=
1
-
y
∗
y
d
y
+
1
2
y
∗
y
2
d
y
2
=
β
x
-
x
∗
y
-
y
∗
+
y
∗
2
σ
2
2
d
t
+
σ
2
y
-
y
∗
d
B
2
.
Therefore(48)
d
V
=
d
V
1
+
d
V
2
=
L
V
d
t
+
σ
1
x
-
x
∗
d
B
1
+
σ
2
y
-
y
∗
d
B
2
,where (49)
L
V
=
-
a
c
+
a
f
y
x
x
∗
+
a
c
f
y
2
x
x
∗
x
-
x
∗
2
-
a
c
+
a
c
f
+
a
f
x
∗
+
a
c
f
y
∗
x
∗
x
-
x
∗
y
-
y
∗
+
M
1
.By Young’s inequality, we have (50)
x
-
x
∗
y
-
y
∗
⩽
x
-
x
∗
2
2
ϵ
+
ϵ
y
-
y
∗
2
2
.Then, we have (51)
L
V
⩽
-
a
c
+
a
f
y
x
x
∗
+
a
c
f
y
2
x
x
∗
x
-
x
∗
2
+
a
c
+
a
c
f
+
a
f
x
∗
+
a
c
f
y
∗
x
∗
x
-
x
∗
y
-
y
∗
+
M
1
⩽
-
a
c
-
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
ϵ
x
-
x
∗
2
+
ϵ
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
y
-
y
∗
2
+
M
1
.Then the hyperbola (52)
a
c
-
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
ϵ
x
-
x
∗
2
-
ϵ
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
y
-
y
∗
2
=
M
1lies in I
n
t
(
R
+
2
). We can take U to be a neighborhood of the hyperbola with U
-
⊆
E
2
=
I
n
t
(
R
+
2
), where U
- is the compact closure of U. So for x
∈
U
∖
E
2, L
V
<
0, which implies condition (P2) in Lemma 10 is satisfied.
Furthermore, there is(53)
M
=
min
⁡
σ
1
2
x
2
,
σ
2
2
y
2
,
x
,
y
∈
U
-
>
0such that the diffusion matrix of system (2) has the following:(54)
∑
i
,
j
=
1
2
a
i
j
ξ
i
ξ
j
=
σ
1
2
x
2
ξ
1
2
+
σ
2
2
y
2
ξ
2
2
⩾
M
ξ
2for all (
x
,
y
)
∈
U
-, ξ
∈
R
2, which indicates that condition (P1) in Lemma 10 is also satisfied.
Hence, we conclude that stochastic system (2) has a stationary distribution μ
(
·
).


## 3.1. Global Positive Solution

In this subsection we show the solution of system (2) is global and nonnegative. As we have known, in order for a stochastic differential equation to have a unique global (i.e., no explosion in a finite time) solution for any given initial value, the coefficients of the equation are generally required to satisfy the linear growth condition and local Lipschitz condition [10]. However, the coefficients of (2) do not satisfy the linear growth condition, though they are locally Lipschitz continuous, so the solution of (2) may explode at a finite time [10]. In this subsection, using the Lyapunov analysis method (mentioned in [10]), we will show the solution of (2) is positive and global.Theorem 5.
For any given initial value(
x
(
0
)
,
y
(
0
)
)
∈
R
+
2, there is unique positive solution (
x
(
t
)
,
y
(
t
)
) of model (2) on t
⩾
0 and the solution will remain in R
+
2 with probability 1; namely, (
x
(
t
)
,
y
(
t
)
)
∈
R
+
2 for all t
⩾
0 almost surely.Proof.
Since the coefficients of (2) are locally Lipschitz continuous, for any given initial value (
x
(
0
)
,
y
(
0
)
)
∈
R
+
2, there is a unique local solution (
x
(
t
)
,
y
(
t
)
) on t
∈
[
0
,
τ
e
), where τ
e is the explosion time. To show this solution is global, we need to show that τ
e
=
∞ a.s.
First, we provex
(
t
) and y
(
t
) do not explode to infinity in a finite time. Let k
0
>
0 be sufficiently large for x
(
0
)
∈
1
/
k
0
,
k
0 and y
(
0
)
∈
1
/
k
0
,
k
0. For each integer k
⩾
k
0, define the stopping time (21)
τ
k
=
inf
⁡
t
∈
0
,
τ
e
:
x
t
∉
1
k
,
k
o
r
y
t
∉
1
k
,
k
,where throughout this paper we set inf
⁡
∅
=
∞ (∅ denotes the empty set). Clearly, τ
k is increasing as k
→
∞. Set τ
∞
=
lim
k
→
∞
⁡
τ
k; therefore τ
∞
⩽
τ
e a.s. If we can show that τ
∞
=
∞ a.s., then τ
e
=
∞ and (
x
(
t
)
,
y
(
t
)
)
∈
R
+
2 a.s. for all t
⩾
0. In other words, to complete the proof it is required to show that τ
∞
=
∞ a.s. If this statement is false, then there are a pair of constants T
>
0 and ϵ
∈
(
0,1
) such that (22)
P
τ
∞
⩽
T
>
ϵ
.Hence, there is an integer k
1
>
k
0 such that (23)
P
τ
k
⩽
T
⩾
ϵ
∀
k
⩾
k
1
.
Define aC
2-function V
:
R
+
2
→
R
-
+ by (24)
V
x
,
y
=
x
-
1
-
ln
⁡
x
+
y
-
1
-
ln
⁡
y
,and the nonnegativity of this function can be seen from (25)
u
-
1
-
ln
⁡
u
⩾
0
f
o
r
u
>
0
.Let k
⩾
k
0 and T
>
0 be arbitrary. Applying Itô’s formula, we obtain (26)
d
V
x
,
y
=
L
V
x
,
y
d
t
+
σ
1
x
-
1
d
B
1
t
+
σ
2
y
-
1
d
B
2
t
,where L
V
:
R
+
2
→
R
-
+ is defined by (27)
L
V
x
,
y
=
1
-
1
x
a
x
+
f
y
1
-
c
x
+
y
-
d
x
-
β
x
y
+
1
-
1
y
β
x
y
-
d
+
α
y
+
1
2
σ
1
2
+
1
2
σ
2
2
⩽
d
+
α
-
a
+
1
2
σ
1
2
+
1
2
σ
2
2
+
a
-
d
+
a
c
-
β
x
-
a
c
x
2
+
a
f
+
a
c
+
a
c
f
+
β
-
d
-
α
y
-
1
2
a
c
f
y
2
-
a
f
y
x
+
1
2
a
c
f
y
2
x
2
≕
K
.Thus, (28)
d
V
x
,
y
⩽
K
d
t
+
σ
1
x
-
1
d
B
1
t
+
σ
2
y
-
1
d
B
2
t
.We can now integrate both sides of (28) from 0 to τ
k
∧
T and then take the expectations (29)
E
V
x
τ
k
∧
T
,
y
τ
k
∧
T
⩽
x
S
0
,
I
0
+
K
E
τ
k
∧
Tso we have (30)
E
V
x
τ
k
∧
T
,
y
τ
k
∧
T
⩽
V
x
0
,
y
0
+
K
T
.Let Ω
k
=
{
τ
k
⩽
T
} for k
⩾
k
1 and, by (23), P
(
Ω
k
)
⩾
ϵ. Note that, for every ω
∈
Ω
k, there is x
(
τ
k
,
ω
) or y
(
τ
k
,
ω
) equaling either k or 1
/
k and therefore V
(
x
(
τ
k
,
ω
)
,
y
(
τ
k
,
ω
)
) is not less than either k
-
1
-
ln
⁡
k or 1
/
k
-
1
-
ln
⁡
1
/
k
=
1
/
k
-
1
+
ln
⁡
K. Hence (31)
V
x
τ
k
,
ω
,
y
τ
k
,
ω
⩾
k
-
1
-
ln
⁡
k
∧
1
k
-
1
+
ln
⁡
K
.It follows from (30) that (32)
V
x
0
,
y
0
+
K
T
⩾
E
1
Ω
k
ω
V
x
τ
k
,
ω
,
y
τ
k
,
ω
⩾
ϵ
k
-
1
-
ln
⁡
k
∧
1
k
-
1
+
ln
⁡
K
,where 1
Ω
k is the indicator function of Ω
k. Letting k
→
∞, we have (33)
∞
>
V
x
0
,
y
0
+
K
T
=
∞
,and it is a contradiction; then we must have τ
∞
=
∞. Therefore, it implies x
(
t
) and y
(
t
) will not explode in a finite time with probability one.


## 3.2. Stochastic Persistence

In this subsection, we are intended to prove the stochastic persistence of the model system (2) under certain parametric restrictions. Stochastic persistence means that if we start from a positive initial condition, that is, from an interior point of the first quadrant, then solution trajectories of the stochastic model will always remain within the interior of the first quadrant and remain bounded at all future time. There are several concepts of stochastic persistence [11]; here we use the notion of stochastic persistence in mean. Before proving the main results of this section we define stochastic persistence in mean.Definition 6.
The population ofx
(
t
) is said to be strongly persistent in the mean if 〈
x
(
t
)
〉
∗
>
0, where (34)
x
t
≔
1
t
∫
0
t
x
s
d
s
,
x
t
∗
≔
l
i
m
i
n
f
t
→
+
∞
1
t
∫
0
t
x
s
d
s
.
The proof of strong persistence result for the stochastic model (2) is based upon the following lemma (see Lemma 4 in [11]).Lemma 7.
Supposex
(
t
)
∈
C
[
Ω
×
R
+
,
R
+
0
], where R
+
0
=
{
a
:
a
>
0
,
a
∈
R
}.(i)
If there are positive constantsμ, T, and λ
⩾
0 such that (35)
ln
⁡
x
t
⩽
λ
t
-
μ
∫
0
t
x
s
d
s
+
∑
i
=
1
n
β
i
B
i
tfor t
⩾
T, where β
i (i
=
1,2
,
…
,
n) are constants, then 〈
x
(
t
)
〉
∗
⩽
λ
/
μ, a
.
s
.
(ii)
If there are positive constantsμ, T, and λ
⩾
0 such that (36)
ln
⁡
x
t
⩾
λ
t
-
μ
∫
0
t
x
s
d
s
+
∑
i
=
1
n
β
i
B
i
tfor t
⩾
T, where β
i (i
=
1,2
,
…
,
n) are constants, then 〈
x
(
t
)
〉
∗
⩾
λ
/
μ, a
.
s
.Remark 8.
〈
x
(
t
)
〉
∗ is defined by (37)
x
t
∗
≔
l
i
m
s
u
p
t
→
+
∞
1
t
∫
0
t
x
s
d
s
.Theorem 9.
Assume thatf
⩽
1. If d
<
a
f
-
α
+
σ
1
2
/
2
+
σ
2
2
/
2, then the total population of x
(
t
) and y
(
t
) of stochastic model system (2) starting from any interior point of first quadrant is strongly persistent in mean. If d
⩾
a
+
σ
1
2
/
2
+
σ
2
2
/
2, then the total population of x
(
t
) and y
(
t
) of stochastic model system (2) starting from any interior point of first quadrant will tend to extinction.Proof.
LetV
(
x
(
t
)
+
y
(
t
)
)
=
ln
⁡
(
x
(
t
)
+
y
(
t
)
) for (
x
,
y
)
∈
R
+
2. Then using Itô’s formula, we get from the first equation of system (2) (38)
d
ln
⁡
x
t
+
y
t
=
a
x
t
+
a
f
y
t
x
t
+
y
t
-
a
c
x
t
+
f
y
t
-
d
-
α
y
t
x
t
+
y
t
+
σ
1
2
2
+
σ
2
2
2
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
⩾
a
f
-
d
-
α
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
.Integrating both sides from 0 to t and dividing by t, we have (39)
l
n
x
t
+
y
t
/
x
0
+
y
0
t
⩾
a
f
-
d
-
α
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
+
σ
1
B
1
t
+
σ
2
B
2
t
t
.From Lemma 7, we have the following inequality: (40)
x
t
+
y
t
∗
⩾
a
f
-
d
-
α
+
σ
1
2
/
2
+
σ
2
2
/
2
a
c
.Hence x
t
+
y
t
∗
>
0, whenever a
f
-
d
-
α
+
σ
1
2
/
2
+
σ
2
2
/
2
>
0.
On the other hand(41)
d
ln
⁡
x
t
+
y
t
=
a
x
t
+
a
f
y
t
x
t
+
y
t
-
a
c
x
t
+
f
y
t
-
d
-
α
y
t
x
t
+
y
t
+
σ
1
2
2
+
σ
2
2
2
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
⩽
a
-
d
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
d
t
+
σ
1
d
B
1
t
+
σ
2
d
B
2
t
.Integrating both sides from 0 to t and dividing by t, we have (42)
ln
⁡
x
t
+
y
t
/
x
0
+
y
0
t
⩽
a
-
d
+
σ
1
2
2
+
σ
2
2
2
-
a
c
x
t
+
y
t
+
σ
1
B
1
t
+
σ
2
B
2
t
t
.From Lemma 7, we have the following inequality: (43)
x
t
+
y
t
∗
⩽
a
-
d
+
σ
1
2
/
2
+
σ
2
2
/
2
a
c
.It follows that 〈
(
x
(
t
)
+
y
(
t
)
)
〉
∗
=
0, whenever a
-
d
+
σ
1
2
/
2
+
σ
2
2
/
2
⩽
0.


## 3.3. Existence of Stationary Distribution

In this section, we prove the existence of stationary distribution of prey and predator populations. For this purpose we find the stationary distribution for solutions of system (2), which in turn implies the stability in stochastic sense. Before proving the main theorem related to the stationary distribution we state a useful lemma from [12] which will be useful to prove the theorem.LetX
(
t
) be a homogeneous Markov process defined in the l-dimensional Euclidean space, denoted by E
l, and described by the following system of stochastic differential equation: (44)
d
X
t
=
b
X
d
t
+
∑
r
=
1
k
f
r
X
d
B
r
t
.We assume there exists a bounded domainU
∈
E
l with regular boundary Γ, having the following properties:(P1)
In the domainU and some neighborhood thereof, the smallest eigenvalue of the diffusion matrix A
3
(
x
) is bounded away from zero.
(P2)
Ifx
∈
E
l
∖
U, the mean time τ at which a path emerging from x reaches the set U is finite, and sup
x
∈
S
⁡
E
x
τ
<
∞ for every compact subset S
⊂
E
l.Lemma 10.
If assumptions (P1) and (P2) hold, then the Markov processX
(
t
) has a stationary distribution μ
·. Let g
· be a function integrable with respect to the measure μ. Then (45)
P
X
l
i
m
T
→
∞
∫
0
T
g
X
t
d
t
=
∫
E
l
g
x
μ
d
x
=
1for all x
∈
E
l.Theorem 11.
Assume thatd
<
β
/
c
-
α and d
<
a
[
1
-
c
(
d
+
α
)
/
β
]. Then there exists a stationary distribution μ
(
·
) for the stochastic system (2).Proof.
Ifd
<
β
/
c
-
α and d
<
a
[
1
-
c
(
d
+
α
)
/
β
], then system (1) has unique positive equilibrium E
∗
(
x
∗
,
y
∗
).
Define a positive definite functionV
:
E
2
→
R
+, where E
2
=
I
n
t
(
R
+
2
) as follows: (46)
V
x
,
y
=
x
-
x
∗
-
x
∗
ln
⁡
x
x
∗
+
y
-
y
∗
-
y
∗
ln
⁡
y
y
∗
≔
V
1
+
V
2
.Applying Itô’s formula, we can calculate (47)
d
V
1
=
1
-
x
∗
x
d
x
+
1
2
x
∗
x
2
d
x
2
=
-
a
c
x
-
x
∗
2
+
a
f
x
-
x
∗
y
x
-
y
∗
x
∗
-
a
c
+
a
c
f
+
β
x
-
x
∗
y
-
y
∗
+
x
∗
2
σ
1
2
-
a
c
f
x
-
x
∗
y
2
x
-
y
∗
2
x
∗
d
t
+
σ
1
x
-
x
∗
d
B
1
=
-
a
c
+
a
f
y
x
x
∗
+
a
c
f
y
2
x
x
∗
x
-
x
∗
2
-
a
c
+
a
c
f
+
β
+
a
f
x
∗
+
a
c
f
y
∗
x
∗
x
-
x
∗
y
-
y
∗
+
x
∗
2
σ
1
2
d
t
+
σ
1
x
-
x
∗
d
B
1
,
d
V
2
=
1
-
y
∗
y
d
y
+
1
2
y
∗
y
2
d
y
2
=
β
x
-
x
∗
y
-
y
∗
+
y
∗
2
σ
2
2
d
t
+
σ
2
y
-
y
∗
d
B
2
.
Therefore(48)
d
V
=
d
V
1
+
d
V
2
=
L
V
d
t
+
σ
1
x
-
x
∗
d
B
1
+
σ
2
y
-
y
∗
d
B
2
,where (49)
L
V
=
-
a
c
+
a
f
y
x
x
∗
+
a
c
f
y
2
x
x
∗
x
-
x
∗
2
-
a
c
+
a
c
f
+
a
f
x
∗
+
a
c
f
y
∗
x
∗
x
-
x
∗
y
-
y
∗
+
M
1
.By Young’s inequality, we have (50)
x
-
x
∗
y
-
y
∗
⩽
x
-
x
∗
2
2
ϵ
+
ϵ
y
-
y
∗
2
2
.Then, we have (51)
L
V
⩽
-
a
c
+
a
f
y
x
x
∗
+
a
c
f
y
2
x
x
∗
x
-
x
∗
2
+
a
c
+
a
c
f
+
a
f
x
∗
+
a
c
f
y
∗
x
∗
x
-
x
∗
y
-
y
∗
+
M
1
⩽
-
a
c
-
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
ϵ
x
-
x
∗
2
+
ϵ
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
y
-
y
∗
2
+
M
1
.Then the hyperbola (52)
a
c
-
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
ϵ
x
-
x
∗
2
-
ϵ
a
c
+
a
c
f
+
a
f
/
x
∗
+
a
c
f
y
∗
/
x
∗
2
y
-
y
∗
2
=
M
1lies in I
n
t
(
R
+
2
). We can take U to be a neighborhood of the hyperbola with U
-
⊆
E
2
=
I
n
t
(
R
+
2
), where U
- is the compact closure of U. So for x
∈
U
∖
E
2, L
V
<
0, which implies condition (P2) in Lemma 10 is satisfied.
Furthermore, there is(53)
M
=
min
⁡
σ
1
2
x
2
,
σ
2
2
y
2
,
x
,
y
∈
U
-
>
0such that the diffusion matrix of system (2) has the following:(54)
∑
i
,
j
=
1
2
a
i
j
ξ
i
ξ
j
=
σ
1
2
x
2
ξ
1
2
+
σ
2
2
y
2
ξ
2
2
⩾
M
ξ
2for all (
x
,
y
)
∈
U
-, ξ
∈
R
2, which indicates that condition (P1) in Lemma 10 is also satisfied.
Hence, we conclude that stochastic system (2) has a stationary distribution μ
(
·
).


## 4. Numerical Simulation Results

In this section we provide numerical simulation results to substantiate the analytical findings for the deterministic and stochastic model reported in the previous sections. To find the approximate strong solution of system (2) in Itô sense with positive initial condition we use the Milstein’s method having strong order of convergence γ
=
1 [13]. This scheme is obtained from Euler-Maruyama scheme by incorporating a correction term for stochastic increment.For numerical simulations of the stochastic model (2), by the method mentioned in [13], we consider the discretized equation (55)
x
k
+
1
=
x
k
+
a
x
k
+
f
y
k
1
-
c
x
k
-
c
y
k
-
d
x
k
-
β
x
k
y
k
Δ
t
+
σ
1
x
k
ξ
1
j
Δ
t
+
1
2
x
k
σ
1
2
Δ
t
ξ
1
j
2
-
1
,
y
k
+
1
=
y
k
+
β
x
k
y
k
-
d
+
α
y
k
Δ
t
+
σ
2
y
k
ξ
2
j
Δ
t
+
1
2
y
k
σ
2
2
Δ
t
ξ
2
j
2
-
1
,where ξ
1
j and ξ
2
j are two independent Gaussian random variables N
(
0,1
) for j
=
1,2
,
…, and we choose Δ
t
=
0.01.Example 12.
We first choosea
=
0.2, c
=
0.2, d
=
0.14, f
=
0.2, β
=
0.2, and α
=
0.2. It is easy to compute that the deterministic model (1) has semisingular equilibrium E
0
(
1.5,0
). From Theorem 4, E
0
(
1.5,0
) is globally stable in this case. But, the white noise may make system (2) appear to be different phenomena. In detail, Figure 1 shows that the trajectory of x
(
t
) is fluctuating in a neighborhood of x
0
=
1.5; however, the trajectory of y
(
t
) tends to 0 when σ
1
=
0.01, σ
2
=
0.02. The stationary distribution of population of x
(
t
) is also provided in Figure 1(b). From stationary distribution of the population of x
(
t
), it is clear that the population of x
(
t
) is distributed normally around the mean value 1.5.
Figure2 also shows that the trajectory of x
(
t
) is fluctuating in a neighborhood of x
0
=
1.5, and the trajectory of y
(
t
) also tends to 0 when σ
1
=
0.1, σ
2
=
0.05, but we can see from Figures 1 and 2 that when the intensity of noise is amplified, the amplitude of fluctuation is also amplified. Furthermore, if we take σ
1
=
0.4, σ
2
=
0.2, then, from Theorem 9, 〈
(
x
(
t
)
+
y
(
t
)
)
〉
∗
=
0, it implies that the total population of x
(
t
) and y
(
t
) will tend to extinction (see Figure 3).Figure 1
(a) Results of one-simulation run of the stochastic system (2) (red curves) and deterministic model (1) (blue curves) with parameter values a
=
0.2, c
=
0.2, d
=
0.14, f
=
0.2, β
=
0.2, α
=
0.2, and σ
1
=
0.01, σ
2
=
0.02. (b) Stationary distribution of uninfected and infected hosts from 10000-simulation run.
(a)
(b)Figure 2
(a) Results of one-simulation run of the stochastic system (2) (red curves) and deterministic model (1) (blue curves) with parameter values a
=
0.2, c
=
0.2, d
=
0.14, f
=
0.2, β
=
0.2, α
=
0.2, and σ
1
=
0.1, σ
2
=
0.05. (b) Stationary distribution of uninfected and infected hosts from 10000-simulation run.
(a)
(b)Figure 3
Results of one-simulation run of the stochastic system (2) (red curves) and deterministic model (1) (blue curves) with parameter values (1) (blue curves) with parameter values a
=
0.2, c
=
0.2, d
=
0.14, f
=
0.2, β
=
0.2, α
=
0.2, and σ
1
=
0.4, σ
2
=
0.2.Example 13.
We choosea
=
0.2, c
=
0.2, d
=
0.04, f
=
0.2, β
=
0.2, and α
=
0.2. We compute that the deterministic model (1) has positive equilibrium E
∗
(
1.2,0.5135
). From Theorem 4, E
∗
(
1.2,0.5135
) is globally stable in this case. But the white noise may make system (2) appear to be different phenomena. In detail, Figure 4(a) shows that the trajectory of x
(
t
) and y
(
t
) is fluctuating in a neighborhood of (
1.2,0.5135
) when σ
1
=
0.01, σ
2
=
0.02. The stationary distributions of population of x
(
t
) and y
(
t
) are also provided in Figure 4(b). From stationary distributions of two populations it is clear that they are distributed normally around the mean values 1.2 and 0.5153, respectively.
Figure5 also shows that the trajectory of x
(
t
) and y
(
t
) is fluctuating in a neighborhood of (
1.2,0.5135
) when σ
1
=
0.1, σ
2
=
0.05, but we can see from Figures 4 and 5 that when the intensity of noise is amplified, the amplitude of fluctuation is also amplified. Furthermore, if we take σ
1
=
0.6, σ
2
=
0.3, then, from Theorem 9, x
t
+
y
t
∗
=
0, it implies that the total population of x
(
t
) and y
(
t
) will tend to extinction (see Figure 6).Figure 4
(a) Results of one-simulation run of the stochastic system (2) (red curves) and deterministic model (1) (blue curves) with parameter values a
=
0.2, c
=
0.2, d
=
0.04, f
=
0.2, β
=
0.2, α
=
0.2, and σ
1
=
0.01, σ
2
=
0.02. (b) Stationary distribution of uninfected and infected hosts from 10000-simulation run.
(a)
(b)Figure 5
(a) Results of one-simulation run of the stochastic system (2) (red curves) and deterministic model (1) (blue curves) with parameter values a
=
0.2, c
=
0.2, d
=
0.04, f
=
0.2, β
=
0.2, α
=
0.2, and σ
1
=
0.1, σ
2
=
0.05. (b) Stationary distribution of uninfected and infected hosts from 10000-simulation run.
(a)
(b)Figure 6
Results of one-simulation run of the stochastic system (2) (red curves) and deterministic model (1) (blue curves) with parameter values a
=
0.2, c
=
0.2, d
=
0.04, f
=
0.2, β
=
0.2, α
=
0.2, and σ
1
=
0.6, σ
2
=
0.3.


## 5. Discussion

In this paper, we investigated a parasite-host model within fluctuating environment. We first investigated the positivity and boundedness of the solution of deterministic model (1); we show that the solution of deterministic model (1) with the initial condition from the first quadrant is positive and bounded. Our results also show that when the parasite-independent host background mortality is less than the maximum per capita birth rate of uninfected hosts, then the population of infected hosts maybe tend to 0. However, when the parasite-independent host background mortality is bigger than the maximum per capita birth rate of uninfected hosts, then the uninfected and infected hosts will coexist.For the stochastic model (2), by using suitable Lyapunov functions, we show that the solution of stochastic model is positive and global, and the total population of the solutions will be persistent or tend to extinction under certain conditions. Further, we can conclude that the stochastic model (2) has a stationary distribution when the deterministic model (1) has a positive equilibrium. Our simulation results show that the solutions of stochastic model oscillate around the equilibria of deterministic model under certain conditions. That is to say, if the effects of environmental stochastic perturbations are smaller enough than the natural death rate, the solution of stochastic model will oscillate around the equilibria of deterministic model, and we observe that the amplitude of fluctuation of population distribution is amplified with the intensity of noise amplified. However, the numerical results suggest that the populations of uninfected and infected hosts become extinct after some initial large amplitude oscillation.



---

*Source: 2972956-2016-03-01.xml*