# Asymptotic Analysis of Impulsive Dispersal Predator-Prey Systems with Markov Switching on Finite-State Space

Authors: Guodong Liu; Zhengbo Chang; Xinzhu Meng
Journal: Journal of Function Spaces
 (2019)

---

## Abstract

In this paper, we investigate the stochastic dynamics of two dispersal predator-prey systems perturbed by white noise, impulsive effect, and regime switching. For the system just interrupted by white noise, we first prove that the stochastic impulsive system has a nontrivial positive periodic solution. Then the sufficient conditions for persistence in mean and extinction of the system are obtained. For the system with Markov regime switching, we verify it is ergodic and has a stationary distribution. And conditions for extinction of the prey species are established. Finally, we provide a series of numerical simulations to illustrate the theoretical analysis.


---

## Body

## 1. Introduction

Due to the ecological effects of human activities and climate changes, e.g., the pollution of the environment and the construction of manufacturing industries and highways, many complete habitats have been increasingly broken into many smaller patches. In some patches, there is less environmental pollution but more predation, and, in other patches, the situation is quite the opposite. As a consequence, the biological species disperse among different patches to search for a better habitat, which has profound effects on both the permanence and evolution of species. Therefore, the effect of dispersion on the possible persistence of single or multiple species has become an imperative subject in population biology and various dispersal models have been investigated by many authors in recent years.Levin [1] firstly proposed a single species dispersal system in 1974. Takeuchi in [2] established the sufficient conditions for persistence of the single dispersal system and verified the existence and uniqueness of the globally positive equilibrium. Then two species or multiple species dispersal system have been investigated by many authors [3–7]. Kuang et al. investigated the relationship between the stability of the population and the dispersal rate of the prey species in [3]. When two species interact, one as a predator and the other as prey, then the Lotka-Volterra predator-prey model is frequently used to describe the dynamics of the biological system [8–12]. Cui [4] explored a nonautonomous dispersal predator-prey system and the sufficient conditions for persistence were obtained. In [5], a threshold for the persistence and extinction of a dispersal predator-prey model with time delay was derived.In view of the environmental pollution, we consider a dispersal predator-prey system with impulsive toxicant input [13–17], which can be described as follows:(1)x1˙t=x1tr1-d1c0t-b1x1t-e1ytx1t+myt+D1x2t-x1t,x2˙t=x2tr2-d2c0t-b2x2t-e2ytx2t+myt+D2x1t-x2t,y˙t=ytr3-b3yt+ε1x1tx1t+myt+ε2x2tx2t+myt,c0˙t=kcet-d+gc0t,ce˙t=-hcet,t≠nT,Δx1t=0,Δx2t=0,Δyt=0,Δc0t=0,Δcet=u,t=nT,n∈Z+,where xi(t)(i=1,2) are the prey density in patch i at time t and y(t) denotes the predator density at time t. The prey species can disperse between two patches and Di(i=1,2) represent the dispersal rates in different patches. Δx1(t)=x1(t+)-x1(t),  Δx2(t)=x2(t+)-x2(t), and Δy(t)=y(t+)-y(t). The environments of two patches are polluted, ce(t) is the toxicant concentration in environment at time t and c0(t) is the toxicant concentration in the organism at time t. Δc0(t)=c0(t+)-c0(t) and Δce(t)=ce(t+)-ce(t). ri(i=1,2,3) are the intrinsic growth rates of each species. bi(i=1,2,3) stand for the intraspecific competition rates of each species. The ratio-dependent functional responses are defined by eiy(t)/xi(t)+my(t)(i=1,2), where m represents the predator interference. ei and εi  (i=1,2) are predation rates and ingestion rates, respectively. di  (i=1,2) are dose-response parameters of the prey in different patches to the toxicant. k represents the toxicant uptake rate of the prey species from environment, and d and g denote the organismal net depuration rate and net ingestion rate, respectively. h stands for the loss rate of toxicant in environment and u is the pulsed input concentration of the toxicant at period nT, n∈Z+. All the parameters above are positive.It is acknowledged that deterministic model cannot describe the species activities completely since the species activities are disturbed by environmental noise frequently. Now let us take white noise [18–34] into account. Suppose that parameters ri(i=1,2,3) are stochastically perturbed by the white noise like (2)ri→ri+σiBi˙t,i=1,2,3,where Bi(t)(i=1,2,3) are mutually independent Brownian motions defined on a complete probability space (Ω,F,{Ft}t≥0,P) with a filtration {Ft}t≥0 satisfying the usual conditions (i.e., it is increasing and right continuous while F0 contains all P-null sets). σi(i=1,2,3) represent the intensities of the white noise. Li et al. in [34] considered a stochastic dispersal population system; then, we perturb system (1) by white noise to the following impulsive stochastic dispersal predator-prey system(3)dx1t=x1tr1-d1c0t-b1x1t-e1ytx1t+myt+D1x2t-x1tdt+σ1x1tdB1t,dx2t=x2tr2-d2c0t-b2x2t-e2ytx2t+myt+D2x1t-x2tdt+σ2x2tdB2t,dyt=ytr3-b3yt+ε1x1tx1t+myt+ε2x2tx2t+mytdt+σ3ytdB3t,dc0t=kcet-d+gc0tdt,dcet=-hcetdt,t≠nT,Δx1t=0,Δx2t=0,Δyt=0,Δc0t=0,Δcet=u,t=nT,n∈Z+.Generally, there are two types of noise in environment. One is white noise, and the other is colored noise, namely, telegraph noise. The telegraph noise is always described as a switching in two or more environmental regimes distinguished by factors such as temperature and humidity. According to the experiment done by Aquaculture Institute of Guizhou Province in 2008 [35], the reproductive effects of African catfish in different season were compared by using local African catfish fed by native fisherman with artificial propagation of oxytocin, in which induced spawning groups were 14 in rainy season and 10 in dry season with male and female ratio 1:2. The results show that the average egg number of per kilogram is 29593 in rainy season, which is 89.06 percent more than that in dry season. It is obvious that different seasons can affect the intrinsic growth rates of the species significantly. Recently, the stochastic population models under regime switching have attracted many authors’ attention [36–46]. Usually, the regime switching between environmental regimes is memoryless and the waiting time for the next switching follows exponential distribution. In this paper, we model the telegraph noise by a continuous-time Markov chain γ(t),t≥0 with finite-state space S={1,2,…,n}, then system (3) becomes(4)dx1t=x1tr1γt-d1γtc0t-b1γtx1t-e1γtytx1t+mγtyt+D1γtx2t-x1tdt+σ1γtx1tdB1t,dx2t=x2tr2γt-d2γtc0t-b2γtx2t-e2γtytx2t+mγtyt+D2γtx1t-x2tdt+σ2γtx2tdB2t,dyt=ytr3γt-b3γtyt+ε1γtx1tx1t+mγtyt+ε2γtx2tx2t+mγtytdt+σ3γtytdB3t,dc0t=kcet-d+gc0tdt,dcet=-hcetdt,t≠nT,Δx1t=0,Δx2t=0,Δyt=0,Δc0t=0,Δcet=u,t=nT,n∈Z+.For each i∈S, all parameters above that related to regime i are positive.This paper is organized as follows. In Section2, we introduce some notations and lemmas, which are necessary for later discussion. In Section 3, we verify the existence of a positive T-periodic solution of system (3) and establish the sufficient conditions for persistence in mean and extinction of the prey species. In Section 4, we use a class of stochastic Lyapunov functions with regime switching to obtain the ergodic property and positive recurrence of system (4), which account for some recurring events of a population system. We also analyze the extinction of system (4). Furthermore, numerical simulations and examples are introduced and carried out to illustrate the theoretical results in Section 5. Finally, we give some discussions in Section 6.


## 2. Preliminaries

Throughout this paper, defineR+3={x=(x1,x2,y)∈R3:xi>0,i=1,2,y>0}.  {γ(t),t≥0} is a right-continuous Markov chain taking values in a finite-state space S with generator Γ=(qij)n×n given by (5)Pγt+Δt=j∣γt=i=qijΔt+oΔt,ifi≠j,1+qiiΔt+oΔt,ifi=j,where Δt>0. qij≥0 are the transition rates from i to j if i≠j while ∑i=1nqij=0. In this paper, we assume that γ(·) and B(·) are independent and the Markov chain γ(t) is irreducible, which means that the system can switch from any regime to any other one. This is equivalent to the condition that for any i,j∈S, there are finite numbers i1,i2,…,il∈S such that qi,i1qi1,i2⋯qil,j>0 and so implies the ergodicity property according to Markov theory for finite states. Note that Γ always has a trivial eigenvalue. The algebraic interpretation of irreducibility is that rank(Γ)=n-1. Under these conditions, the Markov chain has a unique stationary distribution π=(π1,π2,…,πn)∈R1×n which can be determined by solving the linear equation πΓ=0 subject to ∑k=1nπk=1andπk>0,∀k∈S.For convenience, we defineχ^=mini∈S⁡{χ(i)}, χˇ=maxi∈S⁡{χ(i)}, and fu=supt∈0,+∞⁡f(t).Firstly, we consider a nonautonomous stochastic differential equation(6)dφt=ft,φtdt+gt,φtdBt.Lemma 1 (see [42]).
Suppose that the coefficients of (6) are T-periodic in t and there exists a function V(t,φ)∈C2 which is T-periodic in t, and V(t,φ) satisfies the following conditions:
(i)infφ>K⁡V(t,φ)→∞asK→∞;
(ii)LV(t,φ)≤-1 outside some compact set.
Then there exists a solution for system (6) which is a T-periodic Markov process.Then we introduce a lemma which gives a criterion for positive recurrence in terms of Lyapunov functions. Consider the diffusion process(ϕ(t),γ(t)) described by the equation:(7)dϕt=fϕt,γtdt+gϕt,γtdBt,ϕ0=ϕ0,γ0=γ0,where f(·,·):Rn×S→Rn, g(·,·):Rn×S→Rn×d satisfying g(ϕ,k)gT(ϕ,k)=G(ϕ,k). For each k∈S and for any twice continuously differentiable function V(·,k), (ϕ(t),γ(t)) has a generator L given as follows: (8)LV=∑i=1nfiϕ,k∂Vϕ,k∂ϕi+12∑i,j=1nGijϕ,k∂2Vϕ,k∂ϕi∂ϕj+ΞVϕ,·,where (9)ΞVϕ,·k=∑i=1nqkiVϕ,i=∑i≠k∈SqkiVϕ,i-Vϕ,k,k∈S.Lemma 2 (see [43]).
If the following conditions are satisfied:
(i) fori≠j,qij>0;
(ii) for eachk∈S, G(ϕ,k) is symmetric and satisfies (10)ϱζ2≤ζTGϕ,kζ≤ϱ-1ζ2forallζ∈Rn,with some constants ϱ∈(0,1] for all ϕ∈Rn;
(iii) there exists a bounded open subsetD of Rn with a regular (i.e., smooth) boundary satisfying that for each k∈S there exists a nonnegative function V(·,k):Dc→R such that V(·,k) is twice continuously differentiable and that for some ς>0, LV(ϕ,k)≤-ς, for any (ϕ,k)∈Dc×S, the model (7) is ergodic and positive recurrent. That is to say, there exists a unique stationary density μ(·,·) and for any Borel measurable function F(·,·):Rn×S→R such that ∑k∈S∫Rn|F(ϕ,k)|μ(ϕ,k)dϕ<∞, we have P(limt→∞⁡1/t∫0tF(ϕ(s),γ(s))ds=∑k∈S∫RnF(ϕ,k)μ(ϕ,k)dϕ)=1.Now, we consider the subsystem of systems (3) and (4),(11)c0˙t=kcet-d+gc0t,ce˙t=-hcet,t≠nT,Δc0t=0,Δcet=u,t=nT,n∈Z+.Lemma 3 (see [13]).
System (11) has a unique positive solution T-periodic solution (c0∗(t),ce∗(t))T and for each solution (c0(t),ce(t))T of (11), c0(t)→c0∗(t), ce(t)→ce∗(t) as t→+∞, where (12)c0∗t=c0∗0e-d+gt-nT+kue-d+gt-nT-e-ht-nTh-d-g1-e-hT,ce∗t=ue-ht-nT1-e-hT,c0∗0=kue-d+gT-e-hTh-d-g1-e-d+gT1-e-hT,ce∗0=u1-e-hT,for t∈(nT,(n+1)T] and n∈Z+.By Lemma 2.2 of [14], we obtain that(13)limt→∞⁡1t∫0tc0tdt=kuhd+gT≜c¯0.Following the same argument as in [13, 22, 36, 47], we can prove the existence and uniqueness of the global positive solution for systems (3) and (4). The proof is omitted here.


## 3. Analysis of System (3)

By Lemma3, we know that the limit system of (3) is(14)dx1t=x1tr1-d1c0∗t-b1x1t-e1ytx1t+myt+D1x2t-x1tdt+σ1x1tdB1t,dx2t=x2tr2-d2c0∗t-b2x2t-e2ytx2t+myt+D2x1t-x2tdt+σ2x2tdB2t,dyt=ytr3-b3yt+ε1x1tx1t+myt+ε2x2tx2t+mytdt+σ3ytdB3t,with initial value (x1(0),x2(0),y(0)), where c0∗(t) is positive and continuous function of period T.


### 3.1. Existence of Periodic Solution of System (3)

Now, we discuss the existence of periodic solution of system (14). Define τi=1/T∫0T(ri-dic0∗(s)-1/2σi2)ds,  i=1,2.Theorem 4.
Suppose thatτi>0  (i=1,2); then, there exists a positive T-periodic solution for system (14).Proof.
Obviously,c0∗(t) is continuously bounded positive periodic functions in t. cmax=c0∗(0) and cmin=c0∗e-(d+g)T+ku(e-(d+g)T-e-hT)/(h-d-g)(1-e-hT). We need to show that conditions (i) and (ii) in Lemma 1 hold. Define a nonnegative C2-function (15)Vt,x1,x2,y=x1-r1-d1cmaxb1ln⁡x1+r1-d1cmaxb1ϖ1t+x2-r2-d2cmaxb2ln⁡x2+r2-d2cmaxb2ϖ2t+y≜V1+V2+V3,where V1=x1-(r1-d1cmax/b1)ln⁡x1+(r1-d1cmax/b1)ϖ1(t), V2=x2-(r2-d2cmax/b2)ln⁡x2+(r2-d2cmax/b2)ϖ2(t) and V3=y. And ϖi(t)(i=1,2) are functions defined on [0,∞) satisfying ϖi′(t)=ri-dic0∗(t)-1/2σi2-τi and ϖi(0)=0. Obviously, ϖi(t)(i=1,2) are two T-periodic functions on [0,∞). Therefore, V(t,x1,x2,y) is T-periodic in t and satisfies (16)lim⁡infk→∞,x1,x2,y∈R+3\Uk⁡Vt,x1,x2,y=∞,where Uk=(1/k,k)×(1/k,k)×(1/k,k). Therefore, condition (i) of Lemma 1 holds. Then we prove that condition (ii) of Lemma 1 holds.
Applying Itô’s formula to system (14), we have(17)LV1=x1r1-d1c0∗t-b1x1-e1yx1+my+D1x2-x1-r1-d1cmaxb1r1-d1c0∗t-b1x1-e1yx1+my-12σ12-r1-d1cmaxD1b1x2x1-1+r1-d1cmaxb1ϖ1′t≤-b1x12+r1-d1cmin+r1-d1cmaxx1+D1x2+r1-d1cmaxD1+e1/mb1-r1-d1cmaxb1τ1.Similarly,(18)LV2=x2r2-d2c0∗t-b2x2-e2yx2+my+D2x1-x2-r2-d2cmaxb2r2-d2c0∗t-b2x2-e2yx2+my-12σ22-r2-d2cmaxD2b2x1x2-1+r2-d2cmaxb2ϖ2′t≤-b2x22+r2-d2cmin+r2-d2cmaxx2+D2x1+r2-d2cmaxD2+e2/mb2-r2-d2cmaxb2τ2,and(19)LV3=yr3-b3y+ε1x1x1+my+ε2x2x2+my≤-b3y2+r3+ε1+ε2y.Therefore,(20)LV=LV1+LV2+LV3≤-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2-b3y2+r3+ε1+ε2y+r1-d1cmaxD1+e1/mb1+r2-d2cmaxD2+e2/mb2.Define a bounded closed set (21)D=x1,x2,y∈R+3:ϵ<x1<1ϵ,ϵ<x2<1ϵ,ϵ<y<1ϵ,where 0<ϵ<1 is a sufficiently small number such that(22)r1-d1cmin+r1-d1cmax+D2ϵ≤12r1-d1cmaxτ1b1+r2-d2cmaxτ2b2-r2-d2cmin+r2-d2cmax+D124b2-r3+ε1+ε224b3-r1-d1cmaxD1+e1/mb1-r2-d2cmaxD2+e2/mb2≜-K1,where K1<0;(23)r2-d2cmin+r2-d2cmax+D1ϵ≤12r1-d1cmaxτ1b1+r2-d2cmaxτ2b2-r1-d1cmin+r1-d1cmax+D224b1-r3+ε1+ε224b3-r1-d1cmaxD1+e1/mb1-r2-d2cmaxD2+e2/mb2≜-K2,where K2<0;(24)r3+ε1+ε2ϵ≤12r1-d1cmaxτ1b1+r2-d2cmaxτ2b2-r1-d1cmin+r1-d1cmax+D224b1-r2-d2cmin+r2-d2cmax+D124b2-r1-d1cmaxD1+e1/mb1-r2-d2cmaxD2+e2/mb2≜-K3,where K3<0;(25)-b12ϵ2+K4≤-1,where(26)K4=supx1,x2,y∈R+3⁡-b12x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y;-b22ϵ2+K5≤-1,where(27)K5=supx1,x2,y∈R+3-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b222x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y;-b32ϵ2+K6≤-1,where (28)K6=supx1,x2,y∈R+3⁡-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b322y2+r3+ε1+ε2y.
Denote(29)Dϵ1=x1,x2,y∈R+3:0<x1<ϵ,Dϵ2=x1,x2,y∈R+3:0<x2<ϵ,Dϵ3=x1,x2,y∈R+3:0<y<ϵ,Dϵ4=x1,x2,y∈R+3:x1>1ϵ,Dϵ5=x1,x2,y∈R+3:x2>1ϵ,Dϵ6=x1,x2,y∈R+3:y>1ϵ.Note that R+3\D=Dϵ1∪Dϵ2∪Dϵ3∪Dϵ4∪Dϵ5∪Dϵ6. Now, we prove that LV(t,x1,x2,y)≤-1, (x1,x2,y)∈R+3\D.
Case  1. If (x1,x2,y)∈Dϵ1, from (20), we obtain that(30)LV≤r1-d1cmin+r1-d1cmax+D2ϵ-b2x2-r2-d2cmin+r2-d2cmax+D12b22+r2-d2cmin+r2-d2cmax+D124b2-b3y-r3+ε1+ε22b32+r3+ε1+ε224b3-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2≤K1<0.
Case  2. If (x1,x2,y)∈Dϵ2, from (20), one has(31)LV≤r2-d2cmin+r2-d2cmax+D1ϵ-b1x1-r1-d1cmin+r1-d1cmax+D22b12+r1-d1cmin+r1-d1cmax+D224b1-b3y-r3+ε1+ε22b32+r324b3-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2≤K2<0.
Case  3. If (x1,x2,y)∈Dϵ3, from (20), we have(32)LV≤-b1x1-r1-d1cmin+r1-d1cmax+D22b12+r1-d1cmin+r1-d1cmax+D224b1-b2x2-r2-d2cmin+r2-d2cmax+D12b22+r2-d2cmin+r2-d2cmax+D124b2+r3+ε1+ε2ϵ-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2≤K3<0.
Case  4. If (x1,x2,y)∈Dϵ4, from (20), we have(33)LV≤-b12x12-b122x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y≤-b12ϵ2+K4.
Case  5. If (x1,x2,y)∈Dϵ5, from (20), we have(34)LV≤-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b22x22-b222x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y≤-b22ϵ2+K5.
Case  6. If (x1,x2,y)∈Dϵ6, from (20), we have(35)LV≤-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b322y2-b322y2+r3+ε1+ε2y≤-b22ϵ2+K6.Hence, (36)LV≤min⁡K1,K2,K3,-1,x1,x2,y∈R+3\D.
Therefore, the proof of Theorem4 is completed.


### 3.2. Extinction of System (3)

Theorem 5.
For any initial value(x1(0),x2(0),y(0))∈R3 of system (3), if lmax-σmin2/4<0, then the prey is extinct, where lmax=max⁡{r1-d1cmin+D2-D1,r2-d2cmin-D2+D1}, σmin2=min⁡{σ12,σ22}.Proof.
DefineV(x)=x1+x2. Applying Itô’s formula to the limit system (14) of (3) yields(37)dVx=x1r1-d1c0∗t-b1x1-e1yx1+my+D1x2-x1+x2r2-d2c0∗t-b2x2-e2yx2+my+D2x1-x2dt+σ1x1dB1t+σ2x2dB2t=-b1x12+r1-d1c0∗t-D1+D2x1-e1x1yx1+my-b2x22+r2-d2c0∗t-D2+D1x2-e2x2yx2+mydt+σ1x1dB1t+σ1x2dB2t.Applying Itô’s formula again, we obtain that(38)dlnVx=1VxdVx-12V2xdVx2≤1Vx-b1x12+r1-d1cmin+D2-D1x1-e1x1yx1+my-b2x22+r2-d2cmin-D2+D1x2-e2x2yx2+mydt+1Vxσ1x1dB1t+σ2x2dB2t-12V2xσ12x12+σ22x22dt≤1Vxlmaxx1+x2-b1x12-b2x22-e1x1yx1+my-e2x2yx2+mydt-σmin24dt+1Vxσ1x1dB1t+σ2x2dB2t≤lmax-σmin24dt+1Vxσ1x1dB1t+σ2x2dB2t.Integrating both sides from 0 to t, one has(39)ln⁡Vx≤ln⁡Vx0+∫0tlmax-σmin24dt+M1t+M2t,where Mi(t),i=1,2 are local martingales defined respectively by (40)Mit=∫0tσixix1+x2dBiswithMi0=0,i=1,2.By strong law of large number, we have (41)limt→∞Mitt=0,almostsurelyforeachi=1,2.Dividing (39) by t on both sides and then letting t→∞, we have(42)lim⁡supt→∞⁡1tln⁡Vx≤lim⁡supt→∞⁡1t∫0tlmax-σmin24ds<0,a.s.
Therefore,limt→+∞⁡x1+x2=0, which means the prey species is extinct.


### 3.3. Persistence in Mean of System (3)

Theorem 6.
For any initial value(x1(0),x2(0),y(0))∈R3 of system (3), if lmin-emax/m-σmax2/2>0, then the prey species is persistent in mean, where lmin=min⁡{r1-d1cmax+D2-D1,r2-d2cmax-D2+D1}, emax=max⁡{e1,e2}, σmax2=max⁡{σ12,σ22}.Proof.
DefineU(x)=x1+x2. Applying Itô’s formula yields(43)dln⁡Ux=1UxdUx-12U2xdUx2≥1Ux-b1x12+r1-d1cmax+D2-D1x1-e1x1yx1+my-b2x22+r2-d2cmax-D2+D1x2-e2x2yx2+mydt+1Uxσ1x1dB1t+σ2x2dB2t-12U2xσ12x12+σ22x22dt≥1Uxlminx1+x2-bmaxx12+x22-e1x1yx1+my-e2x2yx2+mydt-σmax22dt+1Uxσ1x1dB1t+σ2x2dB2t≥lmin-bmaxUx-emaxm-σmax22dt+1Uxσ1x1dB1t+σ2x2dB2t,where bmax=max{b1,b2}.
Integrating (43) from 0 to t, we get(44)ln⁡Ux≥ln⁡Ux0+∫0tlmin-bmaxUx-emaxm-σmax22dt+N1t+N2t,where Ni(t)=∫0tσixi/U(x)dBi(t),i=1,2.
Dividing (44) by t on both sides and then letting t→∞, we have(45)lim⁡supt→∞⁡1tln⁡Ux≥lmin-emaxm-σmax22-lim⁡inft→∞⁡bmax1t∫0tUxds,a.s.,which implies that(46)lim⁡inft→∞⁡1t∫0tUxds≥-lim⁡supt→∞⁡1bmaxtln⁡Ux+1bmaxlmin-emaxm-σmax22,a.s.According to Corollary 3.1 in [13], we have (47)lim⁡inft→∞⁡1t∫0tUxds≥ε0,where ε0 is a constant, which means the prey species is persistent in mean.


## 3.1. Existence of Periodic Solution of System (3)

Now, we discuss the existence of periodic solution of system (14). Define τi=1/T∫0T(ri-dic0∗(s)-1/2σi2)ds,  i=1,2.Theorem 4.
Suppose thatτi>0  (i=1,2); then, there exists a positive T-periodic solution for system (14).Proof.
Obviously,c0∗(t) is continuously bounded positive periodic functions in t. cmax=c0∗(0) and cmin=c0∗e-(d+g)T+ku(e-(d+g)T-e-hT)/(h-d-g)(1-e-hT). We need to show that conditions (i) and (ii) in Lemma 1 hold. Define a nonnegative C2-function (15)Vt,x1,x2,y=x1-r1-d1cmaxb1ln⁡x1+r1-d1cmaxb1ϖ1t+x2-r2-d2cmaxb2ln⁡x2+r2-d2cmaxb2ϖ2t+y≜V1+V2+V3,where V1=x1-(r1-d1cmax/b1)ln⁡x1+(r1-d1cmax/b1)ϖ1(t), V2=x2-(r2-d2cmax/b2)ln⁡x2+(r2-d2cmax/b2)ϖ2(t) and V3=y. And ϖi(t)(i=1,2) are functions defined on [0,∞) satisfying ϖi′(t)=ri-dic0∗(t)-1/2σi2-τi and ϖi(0)=0. Obviously, ϖi(t)(i=1,2) are two T-periodic functions on [0,∞). Therefore, V(t,x1,x2,y) is T-periodic in t and satisfies (16)lim⁡infk→∞,x1,x2,y∈R+3\Uk⁡Vt,x1,x2,y=∞,where Uk=(1/k,k)×(1/k,k)×(1/k,k). Therefore, condition (i) of Lemma 1 holds. Then we prove that condition (ii) of Lemma 1 holds.
Applying Itô’s formula to system (14), we have(17)LV1=x1r1-d1c0∗t-b1x1-e1yx1+my+D1x2-x1-r1-d1cmaxb1r1-d1c0∗t-b1x1-e1yx1+my-12σ12-r1-d1cmaxD1b1x2x1-1+r1-d1cmaxb1ϖ1′t≤-b1x12+r1-d1cmin+r1-d1cmaxx1+D1x2+r1-d1cmaxD1+e1/mb1-r1-d1cmaxb1τ1.Similarly,(18)LV2=x2r2-d2c0∗t-b2x2-e2yx2+my+D2x1-x2-r2-d2cmaxb2r2-d2c0∗t-b2x2-e2yx2+my-12σ22-r2-d2cmaxD2b2x1x2-1+r2-d2cmaxb2ϖ2′t≤-b2x22+r2-d2cmin+r2-d2cmaxx2+D2x1+r2-d2cmaxD2+e2/mb2-r2-d2cmaxb2τ2,and(19)LV3=yr3-b3y+ε1x1x1+my+ε2x2x2+my≤-b3y2+r3+ε1+ε2y.Therefore,(20)LV=LV1+LV2+LV3≤-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2-b3y2+r3+ε1+ε2y+r1-d1cmaxD1+e1/mb1+r2-d2cmaxD2+e2/mb2.Define a bounded closed set (21)D=x1,x2,y∈R+3:ϵ<x1<1ϵ,ϵ<x2<1ϵ,ϵ<y<1ϵ,where 0<ϵ<1 is a sufficiently small number such that(22)r1-d1cmin+r1-d1cmax+D2ϵ≤12r1-d1cmaxτ1b1+r2-d2cmaxτ2b2-r2-d2cmin+r2-d2cmax+D124b2-r3+ε1+ε224b3-r1-d1cmaxD1+e1/mb1-r2-d2cmaxD2+e2/mb2≜-K1,where K1<0;(23)r2-d2cmin+r2-d2cmax+D1ϵ≤12r1-d1cmaxτ1b1+r2-d2cmaxτ2b2-r1-d1cmin+r1-d1cmax+D224b1-r3+ε1+ε224b3-r1-d1cmaxD1+e1/mb1-r2-d2cmaxD2+e2/mb2≜-K2,where K2<0;(24)r3+ε1+ε2ϵ≤12r1-d1cmaxτ1b1+r2-d2cmaxτ2b2-r1-d1cmin+r1-d1cmax+D224b1-r2-d2cmin+r2-d2cmax+D124b2-r1-d1cmaxD1+e1/mb1-r2-d2cmaxD2+e2/mb2≜-K3,where K3<0;(25)-b12ϵ2+K4≤-1,where(26)K4=supx1,x2,y∈R+3⁡-b12x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y;-b22ϵ2+K5≤-1,where(27)K5=supx1,x2,y∈R+3-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b222x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y;-b32ϵ2+K6≤-1,where (28)K6=supx1,x2,y∈R+3⁡-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b322y2+r3+ε1+ε2y.
Denote(29)Dϵ1=x1,x2,y∈R+3:0<x1<ϵ,Dϵ2=x1,x2,y∈R+3:0<x2<ϵ,Dϵ3=x1,x2,y∈R+3:0<y<ϵ,Dϵ4=x1,x2,y∈R+3:x1>1ϵ,Dϵ5=x1,x2,y∈R+3:x2>1ϵ,Dϵ6=x1,x2,y∈R+3:y>1ϵ.Note that R+3\D=Dϵ1∪Dϵ2∪Dϵ3∪Dϵ4∪Dϵ5∪Dϵ6. Now, we prove that LV(t,x1,x2,y)≤-1, (x1,x2,y)∈R+3\D.
Case  1. If (x1,x2,y)∈Dϵ1, from (20), we obtain that(30)LV≤r1-d1cmin+r1-d1cmax+D2ϵ-b2x2-r2-d2cmin+r2-d2cmax+D12b22+r2-d2cmin+r2-d2cmax+D124b2-b3y-r3+ε1+ε22b32+r3+ε1+ε224b3-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2≤K1<0.
Case  2. If (x1,x2,y)∈Dϵ2, from (20), one has(31)LV≤r2-d2cmin+r2-d2cmax+D1ϵ-b1x1-r1-d1cmin+r1-d1cmax+D22b12+r1-d1cmin+r1-d1cmax+D224b1-b3y-r3+ε1+ε22b32+r324b3-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2≤K2<0.
Case  3. If (x1,x2,y)∈Dϵ3, from (20), we have(32)LV≤-b1x1-r1-d1cmin+r1-d1cmax+D22b12+r1-d1cmin+r1-d1cmax+D224b1-b2x2-r2-d2cmin+r2-d2cmax+D12b22+r2-d2cmin+r2-d2cmax+D124b2+r3+ε1+ε2ϵ-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2≤K3<0.
Case  4. If (x1,x2,y)∈Dϵ4, from (20), we have(33)LV≤-b12x12-b122x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y≤-b12ϵ2+K4.
Case  5. If (x1,x2,y)∈Dϵ5, from (20), we have(34)LV≤-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b22x22-b222x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b3y2+r3+ε1+ε2y≤-b22ϵ2+K5.
Case  6. If (x1,x2,y)∈Dϵ6, from (20), we have(35)LV≤-b1x12+r1-d1cmin+r1-d1cmax+D2x1-r1-d1cmaxb1τ1+r1-d1cmaxD1+e1/mb1-b2x22+r2-d2cmin+r2-d2cmax+D1x2-r2-d2cmaxb2τ2+r2-d2cmaxD2+e2/mb2-b322y2-b322y2+r3+ε1+ε2y≤-b22ϵ2+K6.Hence, (36)LV≤min⁡K1,K2,K3,-1,x1,x2,y∈R+3\D.
Therefore, the proof of Theorem4 is completed.


## 3.2. Extinction of System (3)

Theorem 5.
For any initial value(x1(0),x2(0),y(0))∈R3 of system (3), if lmax-σmin2/4<0, then the prey is extinct, where lmax=max⁡{r1-d1cmin+D2-D1,r2-d2cmin-D2+D1}, σmin2=min⁡{σ12,σ22}.Proof.
DefineV(x)=x1+x2. Applying Itô’s formula to the limit system (14) of (3) yields(37)dVx=x1r1-d1c0∗t-b1x1-e1yx1+my+D1x2-x1+x2r2-d2c0∗t-b2x2-e2yx2+my+D2x1-x2dt+σ1x1dB1t+σ2x2dB2t=-b1x12+r1-d1c0∗t-D1+D2x1-e1x1yx1+my-b2x22+r2-d2c0∗t-D2+D1x2-e2x2yx2+mydt+σ1x1dB1t+σ1x2dB2t.Applying Itô’s formula again, we obtain that(38)dlnVx=1VxdVx-12V2xdVx2≤1Vx-b1x12+r1-d1cmin+D2-D1x1-e1x1yx1+my-b2x22+r2-d2cmin-D2+D1x2-e2x2yx2+mydt+1Vxσ1x1dB1t+σ2x2dB2t-12V2xσ12x12+σ22x22dt≤1Vxlmaxx1+x2-b1x12-b2x22-e1x1yx1+my-e2x2yx2+mydt-σmin24dt+1Vxσ1x1dB1t+σ2x2dB2t≤lmax-σmin24dt+1Vxσ1x1dB1t+σ2x2dB2t.Integrating both sides from 0 to t, one has(39)ln⁡Vx≤ln⁡Vx0+∫0tlmax-σmin24dt+M1t+M2t,where Mi(t),i=1,2 are local martingales defined respectively by (40)Mit=∫0tσixix1+x2dBiswithMi0=0,i=1,2.By strong law of large number, we have (41)limt→∞Mitt=0,almostsurelyforeachi=1,2.Dividing (39) by t on both sides and then letting t→∞, we have(42)lim⁡supt→∞⁡1tln⁡Vx≤lim⁡supt→∞⁡1t∫0tlmax-σmin24ds<0,a.s.
Therefore,limt→+∞⁡x1+x2=0, which means the prey species is extinct.


## 3.3. Persistence in Mean of System (3)

Theorem 6.
For any initial value(x1(0),x2(0),y(0))∈R3 of system (3), if lmin-emax/m-σmax2/2>0, then the prey species is persistent in mean, where lmin=min⁡{r1-d1cmax+D2-D1,r2-d2cmax-D2+D1}, emax=max⁡{e1,e2}, σmax2=max⁡{σ12,σ22}.Proof.
DefineU(x)=x1+x2. Applying Itô’s formula yields(43)dln⁡Ux=1UxdUx-12U2xdUx2≥1Ux-b1x12+r1-d1cmax+D2-D1x1-e1x1yx1+my-b2x22+r2-d2cmax-D2+D1x2-e2x2yx2+mydt+1Uxσ1x1dB1t+σ2x2dB2t-12U2xσ12x12+σ22x22dt≥1Uxlminx1+x2-bmaxx12+x22-e1x1yx1+my-e2x2yx2+mydt-σmax22dt+1Uxσ1x1dB1t+σ2x2dB2t≥lmin-bmaxUx-emaxm-σmax22dt+1Uxσ1x1dB1t+σ2x2dB2t,where bmax=max{b1,b2}.
Integrating (43) from 0 to t, we get(44)ln⁡Ux≥ln⁡Ux0+∫0tlmin-bmaxUx-emaxm-σmax22dt+N1t+N2t,where Ni(t)=∫0tσixi/U(x)dBi(t),i=1,2.
Dividing (44) by t on both sides and then letting t→∞, we have(45)lim⁡supt→∞⁡1tln⁡Ux≥lmin-emaxm-σmax22-lim⁡inft→∞⁡bmax1t∫0tUxds,a.s.,which implies that(46)lim⁡inft→∞⁡1t∫0tUxds≥-lim⁡supt→∞⁡1bmaxtln⁡Ux+1bmaxlmin-emaxm-σmax22,a.s.According to Corollary 3.1 in [13], we have (47)lim⁡inft→∞⁡1t∫0tUxds≥ε0,where ε0 is a constant, which means the prey species is persistent in mean.


## 4. Analysis of System (4)

According to (13), we consider the following system:(48)dx1t=x1tr1γt-d1γtc¯0-b1γtx1t-e1γtytx1t+mγtyt+D1γtx2t-x1tdt+σ1γtx1tdB1t,dx2t=x2tr2γt-d2γtc¯0-b2γtx2t-e2γtytx2t+mγtyt+D2γtx1t-x2tdt+σ2γtx2tdB2t,dyt=ytr3γt-b3γtyt+ε1γtx1tx1t+mγtyt+ε2γtx2tx2t+mγtytdt+σ3γtytdB3t,with initial value (x1(0),x2(0),y(0),γ(0)).


### 4.1. Extinction of System (4)

Theorem 7.
For any initial value(x1(0),x2(0),y(0),γ(0))∈R+3×S of system (4), if  Υ=Σk∈SLmax(k)-Θmin2(k)/4<0, then the prey species will be extinct, where Lmax(k)=max⁡{r1(k)-d1(k)c¯0-D1(k)+D2(k),r2(k)-d2(k)c¯0+D1(k)-D2(k)} and Θmin(k)=min⁡{σ1(k),σ2(k)}.Proof.
Similar to the proof of Theorem4, we define V(x)=x1+x2. Applying Itô’s formula to (48), one has(49)dVx=x1r1γt-d1γtc¯0-b1γtx1-e1γtyx1+mγty+D1γtx2-x1+x2r2γt-d2γtc¯0-b2γtx2-e2γtyx2+mγty+D2γtx1-x2dt+σ1γtx1dB1t+σ2γtx2dB2t=-b1γtx12+r1γt-d1γtc¯0-D1γt+D2γtx1-b1γtx22+r2γt-d2γtc¯0-D2γt+D1γtx2-e1γtx1yx1+mγty-e2γtx2yx2+mγtydt+σ1γtx1dB1t+σ2γtx2dB2t.Applying Itô’s formula again, we can obtain that(50)dln⁡Vx=1VxdVx-12V2xdVx2≤1Vx-b1γtx12+r1γt-d1γtc¯0+D2γt-D1γtx1-b2γtx22+r2γt-d2γtc¯0-D2γt+D1γtx2dt+1Vxσ1γtx1dB1t+σ2γtx2dB2t-12V2xσ12γtx12+σ22γtx22dt≤Lmaxγt-Θmin2γt4dt+1Vxσ1γtx1dB1t+σ2γtx2dB2t.Then,(51)ln⁡Vx-ln⁡Vx0t≤1t∫0tLmaxγt-Θminγt4dt+1tVx∫0tσ1γtx1dB1t+∫0tσ2γtx2dB2t.By the ergodic theory of the Markov chain and the strong law of large number, we have(52)lim⁡supt→∞⁡ln⁡Vxt≤∑k=1mπkLmaxk-Θmin2k4=Υ<0.
Therefore,(53)limt→∞⁡x1+x2=0,a.s.


### 4.2. Unique Stationary Distribution of System (4)

Theorem 8.
Assume that fori≠j, qij>0, ∑k∈Sπkr1k-D1k-d1kc¯0-1/2σ12k-e1(k)/m(k)>0, ∑k∈Sπk(r2(k)-D2(k)-d2(k)c¯0-1/2σ22(k)-e2(k)/m(k))>0 and ∑k∈Sπk(r3(k)-1/2σ32(k))>0, then the stochastic process (x1(t),x2(t),y(t),γ(t)) of model (4) is ergodic and has a unique stationary distribution in R+3×S.Proof.
By the assumptionqij>0 for i≠j, condition (i) in Lemma 2 is satisfied. Let g(x1,x2,y,k)=diag(σ1(k),σ2(k),σ3(k)), then (54)Gx1,x2,y,k=gx1,x2,y,kgTx1,x2,y,k=diagσ12kx12,σ22kx22,σ32ky2.Define a bounded open subset as follows: (55)E=ϵ,1ϵ×ϵ,1ϵ×ϵ,1ϵ⊂R+3,where ϵ∈(0,1) is a constant. Let ϱ=min⁡{1,K7,1/K8}, where K7=minx1,x2,y,k∈E¯×S⁡{σ12(k)x12+σ22(k)x22+σ32(k)y2} and K8=maxx1,x2,y,k∈E¯×S⁡{σ12(k)x12+σ22(k)x22+σ32(k)y2}. For (x1,x2,y,k)∈E×S, we have (56)ϱζ2≤ζTGx1,x2,y,kζ=σ12kx12ζ12+σ22kx22ζ22+σ32ky2ζ32≤ϱ-1ζ2for all ζ∈R3. Thus condition (ii) in Lemma 2 holds. Therefore, it remains for us to verify condition (iii) in Lemma 2.
Define aC2-function on R+3×S, (57)V~x1,x2,y,k=1-νηkx1-ν+1-νδkx2-ν+1-νγky-ν+ε1ke1kx1+ε2ke2kx2+y≜V~1+V~2,where V~1=(1-νηk)x1-ν+(1-νδk)x2-ν+(1-νξk)y-ν and V~2=ε1(k)/e1(k)x1+ε2(k)/e2(k)x2+y. ν is a sufficiently small positive constant satisfying ν<min⁡{1/ηk,1/δk,1/ξk}.
An application of the operatorL to V~1 yields(58)LV1x1,x2,y,k=-ν1-νηkx1-ν-1x1r1k-d1kc¯0-b1kx1-e1kyx1+mky+D1kx2-x1+12ν1+ν1-νηkσ12kx1-ν-ν∑i≠kqkiηi-ηkx1-ν-ν1-νδkx2-ν-1x2r2k-d2kc¯0-b2kx2-e2kyx2+mky+D2kx2-x2+12ν1+ν1-νδkσ22kx2-ν-ν∑i≠kqkiδi-δkx2-ν-ν1-νξky-ν-1yr3k-b3ky+ε1kx1x1+mky+ε1kx1x2+mky+12ν1+ν1-νξkσ32ky-ν-ν∑i≠kqkiξi-ξky-ν≤ν1-νηk-r1k-D1k-d1kc¯0-12σ12k-e1kmk-∑i≠kqkiηi-ηk+12νσ12k-νηk1-νηk∑i≠kqkiηi-ηkx1-ν+ν1-νδk-r2k-D2k-d2kc¯0-12σ22k-e2kmk-∑i≠kqkiδi-δk+12νσ22k-νδk1-νδk∑i≠kqkiδi-δkx2-ν+ν1-νξk-r3k-12σ32k-∑i≠kqkiξi-ξk+12νσ32k-νξk1-νξk∑i≠kqkiξi-ξky-ν+ν1-νηkb1kx11-ν+ν1-νδkb2kx21-ν+ν1-νξkb3ky1-ν.Define the vectors β=(β1,β2,…,βm)T, θ=(θ1,θ2,…,θm)T and ϑ=(ϑ1,ϑ2,…,ϑm)T with βk=-(r1(k)-D1(k)-d1(k)c¯0-1/2σ12(k)-e1(k)/m(k)), θk=-(r2(k)-D2(k)-d2(k)c¯0-1/2σ22(k)) and ϑk=-(r3(k)-1/2σ32(k)-e2(k)/m(k)). As the generator matrix Γ is irreducible, for each βk, θk, and ϑk, there exists η=(η1,η2,…,ηm)T, δ=(δ1,δ2,…,δm)T, and ξ=(ξ1,ξ2,…,ξm)T, respectively, which is a solution of the Poisson system (59)Γηk-βk=-∑j=1mπjβj,Γδk-θk=-∑j=1mπjθj,Γϑk-ϑk=-∑j=1mπjϑj.Therefore we have(60)-r1k-D1k-d1kc¯0-12σ12k-e1kmk-∑i≠kqkiηi-ηk=-∑k∈Sπkr1k-D1k-d1kc¯0-12σ12k-e1kmk,(61)-r2k-D2k-d2kc¯0-12σ22k-e2kmk-∑i≠kqkiδi-δk=-∑k∈Sπkr2k-D2k-d2kc¯0-12σ22k-e2kmk,(62)-r3k-12σ32k-∑i≠kqkiξi-ξk=-∑k∈Sπkr3k-12σ32k.Combining (60), (61) and (62) with (58), we obtain(63)LV~1x1,x2,y,k≤ν1-νηk-∑k∈Sπkr1k-D1k-d1kc¯0-12σ12k-e1kmk+12νσ12k-νηk1-νηk∑i≠kqkiηi-ηkx1-ν+ν1-νδk-∑k∈Sπkr2k-D2k-d2kc¯0-12σ22k-e2kmk+12νσ22k-νδk1-νδk∑i≠kqkiδi-δkx2-ν+ν1-νγk-∑k∈Sπkr3k-12σ32k+12νσ32k-νξk1-νξk∑i≠kqkiξi-ξky-ν+ν1-νηkb1kx11-ν+ν1-νδkb2kx21-ν+ν1-νξkb3ky1-ν.Applying Itô’s formula to V~2(x1,x2,y,k) yields(64)LV~2x1,x2,y,k=ε1ke1kx1r1k-d1kc¯0-b1kx1-e1kyx1+mky+ε1ke1kD1kx2-x1+ε2ke2kx2r2k-d2kc¯0-b2kx2-e1kyx2+mky+ε2ke2kD2kx1-x2+yr3k-b3ky+ε1kx1x1+mky+ε2kx2x2+mky=-b1kx12+ε1ke1kr1k-d1kc¯0-ε1ke1kD1k+ε2ke2kD2kx1-b2kx22+ε2ke2kr2k-d2kc¯0+ε1ke1kD1k-ε2ke2kD2kx2-b3ky2-r3ky≤-b^1x12+εˇ1e^1rˇ1+εˇ2e^2Dˇ2x1-b^2x22+εˇ2e^2rˇ2+εˇ2e^2Dˇ2x2-b^3y2-r^3y.Thus LV~=LV~1+LV~2 can be estimated as follows:(65)LV~≤ν1-νηk-∑k∈Sπkr1k-D1k-d1kc¯0-12σ12k-e1kmk+12νσ12k-νηk1-νηk∑i≠kqkiηi-ηkx1-ν+ν1-νηkb1kx11-ν-b^1x12+εˇ1e^1rˇ1+εˇ2e^2Dˇ2x1+ν1-νδk-∑k∈Sπkr2k-D2k-d2kc¯0-12σ22k-e2kmk+12νσ22k-νδk1-νδk∑i≠kqkiδi-δkx2-ν+ν1-νδkb2kx21-ν-b^2x22+εˇ2e^2rˇ2+εˇ2e^2Dˇ2x2+ν1-νξk-∑k∈Sπkr3k-12σ32k+12νσ32k-νξk1-νξk∑i≠kqkiξi-ξky-ν+ν1-νξkb3ky1-ν-b^3y2-r^3y. It is easy to see that (66)LV~≤f1x1+f2x2+f3x3≤f1x1+f2ux2+f3ux3,ifx1→0orx1→+∞,f1ux1+f2x2+f3ux3,ifx2→0orx2→+∞,f1ux1+f2ux2+f3x3,ifx3→0orx3→+∞.Consequently, we know that, for a sufficiently small ϵ, (67)LV~x1,x2,y,k≤-1,∀x1,x2,y,k∈E¯×S.
By Lemma2, we get that (x1(t),x2(t),y(t),γ(t)) is ergodic and positive recurrent; that is, system (4) has a unique stationary distribution.
This completes the proof of Theorem8.


## 4.1. Extinction of System (4)

Theorem 7.
For any initial value(x1(0),x2(0),y(0),γ(0))∈R+3×S of system (4), if  Υ=Σk∈SLmax(k)-Θmin2(k)/4<0, then the prey species will be extinct, where Lmax(k)=max⁡{r1(k)-d1(k)c¯0-D1(k)+D2(k),r2(k)-d2(k)c¯0+D1(k)-D2(k)} and Θmin(k)=min⁡{σ1(k),σ2(k)}.Proof.
Similar to the proof of Theorem4, we define V(x)=x1+x2. Applying Itô’s formula to (48), one has(49)dVx=x1r1γt-d1γtc¯0-b1γtx1-e1γtyx1+mγty+D1γtx2-x1+x2r2γt-d2γtc¯0-b2γtx2-e2γtyx2+mγty+D2γtx1-x2dt+σ1γtx1dB1t+σ2γtx2dB2t=-b1γtx12+r1γt-d1γtc¯0-D1γt+D2γtx1-b1γtx22+r2γt-d2γtc¯0-D2γt+D1γtx2-e1γtx1yx1+mγty-e2γtx2yx2+mγtydt+σ1γtx1dB1t+σ2γtx2dB2t.Applying Itô’s formula again, we can obtain that(50)dln⁡Vx=1VxdVx-12V2xdVx2≤1Vx-b1γtx12+r1γt-d1γtc¯0+D2γt-D1γtx1-b2γtx22+r2γt-d2γtc¯0-D2γt+D1γtx2dt+1Vxσ1γtx1dB1t+σ2γtx2dB2t-12V2xσ12γtx12+σ22γtx22dt≤Lmaxγt-Θmin2γt4dt+1Vxσ1γtx1dB1t+σ2γtx2dB2t.Then,(51)ln⁡Vx-ln⁡Vx0t≤1t∫0tLmaxγt-Θminγt4dt+1tVx∫0tσ1γtx1dB1t+∫0tσ2γtx2dB2t.By the ergodic theory of the Markov chain and the strong law of large number, we have(52)lim⁡supt→∞⁡ln⁡Vxt≤∑k=1mπkLmaxk-Θmin2k4=Υ<0.
Therefore,(53)limt→∞⁡x1+x2=0,a.s.


## 4.2. Unique Stationary Distribution of System (4)

Theorem 8.
Assume that fori≠j, qij>0, ∑k∈Sπkr1k-D1k-d1kc¯0-1/2σ12k-e1(k)/m(k)>0, ∑k∈Sπk(r2(k)-D2(k)-d2(k)c¯0-1/2σ22(k)-e2(k)/m(k))>0 and ∑k∈Sπk(r3(k)-1/2σ32(k))>0, then the stochastic process (x1(t),x2(t),y(t),γ(t)) of model (4) is ergodic and has a unique stationary distribution in R+3×S.Proof.
By the assumptionqij>0 for i≠j, condition (i) in Lemma 2 is satisfied. Let g(x1,x2,y,k)=diag(σ1(k),σ2(k),σ3(k)), then (54)Gx1,x2,y,k=gx1,x2,y,kgTx1,x2,y,k=diagσ12kx12,σ22kx22,σ32ky2.Define a bounded open subset as follows: (55)E=ϵ,1ϵ×ϵ,1ϵ×ϵ,1ϵ⊂R+3,where ϵ∈(0,1) is a constant. Let ϱ=min⁡{1,K7,1/K8}, where K7=minx1,x2,y,k∈E¯×S⁡{σ12(k)x12+σ22(k)x22+σ32(k)y2} and K8=maxx1,x2,y,k∈E¯×S⁡{σ12(k)x12+σ22(k)x22+σ32(k)y2}. For (x1,x2,y,k)∈E×S, we have (56)ϱζ2≤ζTGx1,x2,y,kζ=σ12kx12ζ12+σ22kx22ζ22+σ32ky2ζ32≤ϱ-1ζ2for all ζ∈R3. Thus condition (ii) in Lemma 2 holds. Therefore, it remains for us to verify condition (iii) in Lemma 2.
Define aC2-function on R+3×S, (57)V~x1,x2,y,k=1-νηkx1-ν+1-νδkx2-ν+1-νγky-ν+ε1ke1kx1+ε2ke2kx2+y≜V~1+V~2,where V~1=(1-νηk)x1-ν+(1-νδk)x2-ν+(1-νξk)y-ν and V~2=ε1(k)/e1(k)x1+ε2(k)/e2(k)x2+y. ν is a sufficiently small positive constant satisfying ν<min⁡{1/ηk,1/δk,1/ξk}.
An application of the operatorL to V~1 yields(58)LV1x1,x2,y,k=-ν1-νηkx1-ν-1x1r1k-d1kc¯0-b1kx1-e1kyx1+mky+D1kx2-x1+12ν1+ν1-νηkσ12kx1-ν-ν∑i≠kqkiηi-ηkx1-ν-ν1-νδkx2-ν-1x2r2k-d2kc¯0-b2kx2-e2kyx2+mky+D2kx2-x2+12ν1+ν1-νδkσ22kx2-ν-ν∑i≠kqkiδi-δkx2-ν-ν1-νξky-ν-1yr3k-b3ky+ε1kx1x1+mky+ε1kx1x2+mky+12ν1+ν1-νξkσ32ky-ν-ν∑i≠kqkiξi-ξky-ν≤ν1-νηk-r1k-D1k-d1kc¯0-12σ12k-e1kmk-∑i≠kqkiηi-ηk+12νσ12k-νηk1-νηk∑i≠kqkiηi-ηkx1-ν+ν1-νδk-r2k-D2k-d2kc¯0-12σ22k-e2kmk-∑i≠kqkiδi-δk+12νσ22k-νδk1-νδk∑i≠kqkiδi-δkx2-ν+ν1-νξk-r3k-12σ32k-∑i≠kqkiξi-ξk+12νσ32k-νξk1-νξk∑i≠kqkiξi-ξky-ν+ν1-νηkb1kx11-ν+ν1-νδkb2kx21-ν+ν1-νξkb3ky1-ν.Define the vectors β=(β1,β2,…,βm)T, θ=(θ1,θ2,…,θm)T and ϑ=(ϑ1,ϑ2,…,ϑm)T with βk=-(r1(k)-D1(k)-d1(k)c¯0-1/2σ12(k)-e1(k)/m(k)), θk=-(r2(k)-D2(k)-d2(k)c¯0-1/2σ22(k)) and ϑk=-(r3(k)-1/2σ32(k)-e2(k)/m(k)). As the generator matrix Γ is irreducible, for each βk, θk, and ϑk, there exists η=(η1,η2,…,ηm)T, δ=(δ1,δ2,…,δm)T, and ξ=(ξ1,ξ2,…,ξm)T, respectively, which is a solution of the Poisson system (59)Γηk-βk=-∑j=1mπjβj,Γδk-θk=-∑j=1mπjθj,Γϑk-ϑk=-∑j=1mπjϑj.Therefore we have(60)-r1k-D1k-d1kc¯0-12σ12k-e1kmk-∑i≠kqkiηi-ηk=-∑k∈Sπkr1k-D1k-d1kc¯0-12σ12k-e1kmk,(61)-r2k-D2k-d2kc¯0-12σ22k-e2kmk-∑i≠kqkiδi-δk=-∑k∈Sπkr2k-D2k-d2kc¯0-12σ22k-e2kmk,(62)-r3k-12σ32k-∑i≠kqkiξi-ξk=-∑k∈Sπkr3k-12σ32k.Combining (60), (61) and (62) with (58), we obtain(63)LV~1x1,x2,y,k≤ν1-νηk-∑k∈Sπkr1k-D1k-d1kc¯0-12σ12k-e1kmk+12νσ12k-νηk1-νηk∑i≠kqkiηi-ηkx1-ν+ν1-νδk-∑k∈Sπkr2k-D2k-d2kc¯0-12σ22k-e2kmk+12νσ22k-νδk1-νδk∑i≠kqkiδi-δkx2-ν+ν1-νγk-∑k∈Sπkr3k-12σ32k+12νσ32k-νξk1-νξk∑i≠kqkiξi-ξky-ν+ν1-νηkb1kx11-ν+ν1-νδkb2kx21-ν+ν1-νξkb3ky1-ν.Applying Itô’s formula to V~2(x1,x2,y,k) yields(64)LV~2x1,x2,y,k=ε1ke1kx1r1k-d1kc¯0-b1kx1-e1kyx1+mky+ε1ke1kD1kx2-x1+ε2ke2kx2r2k-d2kc¯0-b2kx2-e1kyx2+mky+ε2ke2kD2kx1-x2+yr3k-b3ky+ε1kx1x1+mky+ε2kx2x2+mky=-b1kx12+ε1ke1kr1k-d1kc¯0-ε1ke1kD1k+ε2ke2kD2kx1-b2kx22+ε2ke2kr2k-d2kc¯0+ε1ke1kD1k-ε2ke2kD2kx2-b3ky2-r3ky≤-b^1x12+εˇ1e^1rˇ1+εˇ2e^2Dˇ2x1-b^2x22+εˇ2e^2rˇ2+εˇ2e^2Dˇ2x2-b^3y2-r^3y.Thus LV~=LV~1+LV~2 can be estimated as follows:(65)LV~≤ν1-νηk-∑k∈Sπkr1k-D1k-d1kc¯0-12σ12k-e1kmk+12νσ12k-νηk1-νηk∑i≠kqkiηi-ηkx1-ν+ν1-νηkb1kx11-ν-b^1x12+εˇ1e^1rˇ1+εˇ2e^2Dˇ2x1+ν1-νδk-∑k∈Sπkr2k-D2k-d2kc¯0-12σ22k-e2kmk+12νσ22k-νδk1-νδk∑i≠kqkiδi-δkx2-ν+ν1-νδkb2kx21-ν-b^2x22+εˇ2e^2rˇ2+εˇ2e^2Dˇ2x2+ν1-νξk-∑k∈Sπkr3k-12σ32k+12νσ32k-νξk1-νξk∑i≠kqkiξi-ξky-ν+ν1-νξkb3ky1-ν-b^3y2-r^3y. It is easy to see that (66)LV~≤f1x1+f2x2+f3x3≤f1x1+f2ux2+f3ux3,ifx1→0orx1→+∞,f1ux1+f2x2+f3ux3,ifx2→0orx2→+∞,f1ux1+f2ux2+f3x3,ifx3→0orx3→+∞.Consequently, we know that, for a sufficiently small ϵ, (67)LV~x1,x2,y,k≤-1,∀x1,x2,y,k∈E¯×S.
By Lemma2, we get that (x1(t),x2(t),y(t),γ(t)) is ergodic and positive recurrent; that is, system (4) has a unique stationary distribution.
This completes the proof of Theorem8.


## 5. Examples and Simulations

In this section, we introduce some examples and numerical simulations to illustrate our main results.Example 1.
We choose the parameter values in system (3) as follows: r1=0.4, r2=0.3, r3=0.2, d1=0.02, d2=0.01, b1=0.3, b2=0.3, b3=0.4, e1=0.01, e2=0.005, m=0.2, ε1=0.004, ε2=0.003, D1=0.1, D2=0.2, T=4, σ1=σ2=σ3=0.03. In addition, we let u=0.5, k=0.5, g=0.05, d=0.05, h=0.1 for the subsystem. By computation, we obtain that τ1=1/T∫0T(r1-d1c0∗(s)-1/2σ12)ds=0.27455>0 and τ2=1/T∫0Tr2-d2c0∗s-1/2σ22ds=0.23705>0, satisfying the conditions in Theorem 4.Figure1 verifies Theorem 4 from two-dimension and three-dimension. Therefore, there exists a T-periodic solution (x1∗(t),x2∗(t),y∗(t)) for system (3).Figure 1
(a) The phase portrait ofx1 and x2 in the stochastic dispersal model (3). (b) The phase portrait of x1, x2, and y.
(a) (b)Example 2.
Firstly, we choose the same parameter values in Example1, satisfying the conditions in Theorem 6. From the left figure of Figure 2, we find that the prey species and the predator are persistent in mean for system (3).Figure 2
The left figure is the permanence in mean ofx1, x2, and y for the stochastic dispersal model (3). The right figure stands for the extinction of the prey species for model (3).For the extinction of system (3), D1 decreases from 0.1 to 0.003 and D2 decreases from 0.2 to 0.002. Additionally, we choose d1=0.65, d2=0.7, e1=0.1, e2=0.05, ε1=0.04, and ε2=0.03. The other parameter values are the same as the values in Example 1. The results are shown in the right figure of Figure 2, in which the prey species both in patch 1 and 2 are extinct.We suppose that the limit system of system (4) has k regimes, i.e., S={1,2,…,k}. Then system (48) has k subsystems (68)dx1t=x1tr1i-d1ic¯0-b1ix1t-e1iytx1t+miyt+D1ix2t-x1tdt+σ1ix1tdB1t,dx2t=x2tr2i-d2ic¯0-b2ix2t-e2iytx2t+miyt+D2ix1t-x2tdt+σ2ix2tdB2t,dyt=ytr3i-b3iyt+ε1ix1tx1t+miyt+ε2ix2tx2t+miytdt+σ3iytdB3t,where i=1,2,…,k.Example 3.
Letγ(t) be a right-continuous Markov chain taking values in S={1,2} (i.e., k=2) and the generator Γ of the Markov chain is (69)Γ=-111-1.Hence,π=(π1,π2)=(0.5,0.5). The Markov chain is illustrated in Figure 3.Figure 3
Markov chain withπ=(0.5,0.5).In regime 1, we chooser1(1)=0.3, r2(1)=0.2, r3(1)=0.2, d1(1)=0.6, d2(1)=0.7, b1(1)=0.4, b2(1)=0.3, b3(1)=0.4, e1(1)=0.1, e2(1)=0.3, m(1)=0.2, ε1(1)=0.04, ε2(1)=0.03, D1(1)=0.05, D2(1)=0.07, σ1(1)=0.2, σ2(1)=0.15, σ3(1)=0.2.In regime 2, we chooser1(2)=0.2, r2(2)=0.3, r3(2)=0.1, d1(2)=0.7, d2(2)=0.6, b1(2)=0.3, b2(2)=0.2, b3(2)=0.2, e1(2)=0.2, e2(2)=0.1, m(2)=0.1, ε1(2)=0.03, ε2(2)=0.02, D1(2)=0.07, D2(2)=0.05, σ1(2)=0.1, σ2(2)=0.1, σ3(2)=0.1.Whenc¯0=0.6, we can observe that π1(Lmax(1)-Θmin2(1)/4)+π2(Lmax(2)-Θmin2(2)/4)=-0.0656<0. The conditions in Theorem 7 hold.From Figure4, the prey species in two patches are extinct, which is consistent with Theorem 7.Figure 4
The left figure is the prey in patch 1. The right figure stands for the prey in patch 2.Example 4.
The Markov chain is identical with that in Example3.In regime 1, we chooser1(1)=0.3, r2(1)=0.2, r3(1)=0.2, d1(1)=0.06, d2(1)=0.07, b1(1)=0.4, b2(1)=0.3, b3(1)=0.4, e1(1)=0.01, e2(1)=0.03, m(1)=0.2, ε1(1)=0.004, ε2(1)=0.003, D1(1)=0.01, D2(1)=0.02, σ1(1)=0.2, σ2(1)=0.15, σ3(1)=0.2.In regime 2, we chooser1(2)=0.2, r2(2)=0.3, r3(2)=0.1, d1(2)=0.07, d2(2)=0.06, b1(2)=0.3, b2(2)=0.2, b3(2)=0.2, e1(2)=0.02, e2(2)=0.01, m(2)=0.1, ε1(2)=0.003, ε2(2)=0.002, D1(2)=0.02, D2(2)=0.01, σ1(2)=0.1, σ2(2)=0.1, σ3(2)=0.1.Whenc¯0=0.3, we obtain that π1(r1(1)-D1(1)-d1(1)c¯0-1/2σ12(1)-e1(1)/m(1))+π2(r1(2)-D1(2)-d1(2)c¯0-1/2σ12(2)-e1(2)/m(2))=0.078>0 and π1(r2(1)-D2(1)-d2(1)c¯0-1/2σ22(1)-e2(1)/m(1))+π2(r2(2)-D2(2)-d2(2)c¯0-1/2σ22(2)-e2(2)/m(2))=0.0824>0. The conditions in Theorem 8 hold.The analytical results are consistent with the simulation results as illustrated in Figures5 and 6.Figure 5
The left figure is the persistence in mean of the prey species in patch 1. The right figure implies the stationary distribution of the prey species in patch 1.Figure 6
The left figure is the persistence in mean of the prey species in patch 2. The right figure implies the stationary distribution of the prey species in patch 2.From Figures5 and 6, we know that there exist two density functions which are displayed in Figure 7. Therefore, the prey species has a stationary distribution.Figure 7
The density functions of the prey species in patch 1 and patch 2, respectively.


## 6. Discussion

It is theoretically and practically valuable to investigate the dispersal population system. In this paper, we explore a predator-prey Lotka-Volterra system in two patches. Taking white noise, regime switching and impulsive toxicant input into account, two stochastic impulsive patchy systems are established. We discuss some properties of the two systems and the results are shown as follows.(i) Suppose thatτi>0  (i=1,2); then, there exists a positive T-periodic solution for system (3).(ii) For any initial value(x1(0),x2(0),y(0))∈R3 of system (3), if lmax-σmin2/4<0, then the prey will be extinct.(iii) For any initial value(x1(0),x2(0),y(0))∈R3 of system (3), if lmin-emax/m-σmax2/2>0, then the prey species will be persistent in mean.(iv) For system (4), given initial value (x1(0),x2(0),y(0),γ(0))∈R+3×S, if Υ=Σk∈SLmax(k)-Θmin2(k)/4<0, then the prey species will be extinct.(v) Assume∑k∈Sπk(r1(k)-D1(k)-d1(k)c¯0-1/2σ12(k)-e1(k)/m(k))>0, ∑k∈Sπk(r2(k)-D2(k)-d2(k)c¯0-1/2σ22(k)-e2(k)/m(k))>0 and ∑k∈Sπkr3k-1/2σ32k>0; then, the stochastic process (x1(t),x2(t),y(t),γ(t)) of model (4) is ergodic and has a unique stationary distribution in R+3×S.From (ii) and (iii), we can confirm that the extinction and persistence of the dispersal prey species are partly determined byΔ12 and Δ21, where Δ12=D1-D2 and Δ21=D2-D1. When max⁡{Δ12,Δ21} is sufficiently small such that lmax-σmin2/4<0, the prey species is destined to be extinct. In contrary, if min⁡{Δ12,Δ21} is large enough such that lmin-emax/m-σmax2/2>0, then the prey species will be persistent almost surely. Example 2 also further illustrates the conclusions.When the dispersal predator-prey system is perturbed by telegraph noise, from (iv), the condition for extinction of the prey species is similar to (ii). However, the parameters are relevant to the regimek. That is to say, the variations of dispersal rates Δ12(k)=D1(k)-D2(k) and Δ21(k)=D2(k)-D1(k) are also subject to the regimes. Different from (iii), the conditions for persistence of the prey in patch i(i=1,2) are mainly affected by the dispersal rate in the i-th patch in (v). If the dispersal rate in patch i is not large enough such that ∑k∈Sπkrik-Dik-dikc¯0-1/2σi2k-eik/mk>0, then the prey in patch i will be persistent. Examples 3 and 4 simulate the results (iv) and (v).With the construction of manufacturing industries and highways, refuge can be introduced into the Lotka-Volterra system to enrich the population models, which can simulate the natural reserves for animals. On the other hand, this also contributes to the evolution of the species, which can be found in [48, 49]. Also, Song et al. have explored the spatiotemporal dynamics of diffusive predator-prey models, such as steady-state bifurcations and Turing patterns (see [50–52]). Moreover, in predator-prey systems, time delay [53–56] also should be considered since the digestion process of the predator consumes time. These interesting issues deserve further investigation.



---

*Source: 8057153-2019-04-01.xml*