# Dynamics of a Multigroup SIR Epidemic Model with Nonlinear Incidence and Stochastic Perturbation

Authors: Yuguo Lin; Daqing Jiang
Journal: Abstract and Applied Analysis
 (2013)

---

## Abstract

We introduce stochasticity into a multigroup SIR model with nonlinear incidence. We prove that when the intensity of white noise is small, the solution of stochastic system converges weakly to a singular measure (i.e., a distribution) ifℛ0≤1 and there exists an invariant
distribution which is ergodic if ℛ0>1. This is the same situation as the corresponding deterministic
case. When the intensity of white noise is large, white noise controls this system. This means that
the disease will extinct exponentially regardless of the magnitude of ℛ0.


---

## Body

## 1. Introduction

Considering different contact patterns, distinct number of sexual partners, or different geography, and so forth, individual hosts are often divided into groups in modeling epidemic diseases. Multigroup models have been proposed in the literature to describe the transmission dynamics of infectious diseases in heterogeneous host populations. One of the earliest work on multigroup models is the seminal paper by Lajmanovich and Yorke [1] on a class of SIS multigroup models for the transmission dynamics of Gonorrhea. From then on, much research has been done on various forms of multi-group models, see, for example, [2–6]. It is well known that the global stability of the endemic equilibrium of multigroup models is a very challenging problem. Recently, Guo et al. [7, 8] and Li and Shuai [9] proposed a graph-theoretic approach to the method of global Lyapunov functions and completely solved this problem for some multigroup epidemic models. Subsequently, this approach or ideas of [7–9] were applied to the investigation into the dynamics of several classes of multigroup epidemic models (e.g., [10–13]).Now we consider a deterministic multigroup SIR (susceptible, infective, and recovered) epidemic model with nonlinear incidence:(1)S˙k=Λk-∑j=1n‍βkjSkϕj(Ij)-dkSk,I˙k=∑j=1n‍βkjSkϕj(Ij)-(ϵk+γk)Ik,R˙k=γkIk-δkRk,k=1,2,…,n.
Here Sk(t), Ik(t), and Rk(t) denote the susceptible, infective, and recovered population at time t in the kth group, respectively, k=1,2,…,n. Suppose the death rates of Sk, Ik, Rk in the kth group are different. The parameters in the model (1) are assumed to be positive and summarized in the following list: Λk
: influx of individuals into thekth group;βkj: transmission coefficient between compartments Sk and Ij;dk: death rate of S compartment in the kth group;ϵk: death rate of I compartment in the kth group;δk: death rate of R compartment in the kth group;γk: recovery rate of infectious individuals in the kth group.Considering that the death rates of infective and recovered population are usually no less than the susceptible’s, we assumedk≤min{ϵk,δk} for all k.Throughout this paper, we consider the following basic assumptions on functionsϕj(x),  x∈(0,+∞), j=1,2,…,n: (C1)
ϕj(x) is a C1-function on (0,+∞), ϕj(x)>0 for x>0 and ϕj(0)=0;(C2)
there exist positive constantscj such that ϕj(x)≤cjx, x∈(0,+∞);(C3)
limx→0+(ϕj(x)/x)=cj;(C4)
[ϕj(x)-ϕj(y)][ϕj(x)/x-ϕj(y)/y]≤0, for any x>0 and y>0.It is easy to check that classes ofϕj(x) satisfying (C1)–(C3) include common incidence functions such as ϕj(Ij)=Ij [7], ϕj(Ij)=Ij/(1+αjIj) [14], ϕj(Ij)=Ij/(1+αIj2) [15]. Furthermore, ϕj(Ij)=Ij, ϕj(Ij)=Ij/(1+αjIj) satisfy (C4).In Sun [13], the author considered a general multigroup SIR models with nonlinear incidence by using the same method as in [9]. Noting that model (1) is a special case of multigroup SIR models appeared in [13], according to the results of [13], the following results hold for system (1).Proposition 1.
Assume thatB=(βkj) is irreducible and ϕj(x) satisfies (C1)–(C3).(1)
Ifℛ0≤1, then for system (1),P0=(S10,0,0,S20,0,0,…,Sn0,0,0)is the unique equilibrium and it is globally asymptotically stable inΓ;(2)
ifℛ0>1and(C4)is satisfied, then there exists an endemic equilibriumP*=(S1*,I1*,R1*,…,Sn*,In*,Rn*)for system (1) andP*is globally asymptotically stable inΓ,
whereSk0=Λk/dk,  k=1,2,…,n, ℛ0=ρ(M0) (the spectral radius of M0), M0 = M(S10,S20,…,Sn0) = (βkjSk0cj/(ϵk+γk))n×n,  and Γ = {(S1,I1,R1,…,Sn,In,Rn)∈ℝ+3n:Sk≤Λk/dk,Sk+Ik≤Λk/dk,k= 1,2,…,n}.The aim of this paper is to evaluate the effect of stochastic parameter perturbation on system (1). Considering the existence of environmental noise, we introduce randomness into the model (1) by replacing the parameters dk, ϵk, and δk by
(2)dk⟶dk+αkB˙1k(t),ϵk⟶ϵk+βkB˙2k(t),δk⟶δk+σkB˙3k(t),
where B1k(t), B2k(t), B3k(t), k=1,2,…,n are mutual independent standard Brownian motions with B1k(0)=0,  B2k(0)=0,  B3k(0)=0, and the intensities of white noises αk2>0, βk2>0,  σk2>0, respectively. Then the stochastic version corresponding to the deterministic model (1) takes the following form:
(3)dSk(t)=[Λk-∑j=1n‍βkjSk(t)ϕj(Ij(t))-dkSk(t)]dt-αkSk(t)dB1k(t),dIk(t)=[∑j=1n‍βkjSk(t)ϕj(Ij(t))-(ϵk+γk)Ik(t)]dt-βkIk(t)dB2k(t),dRk(t)=[γkIk(t)-δkRk(t)]dt-σkRk(t)dB3k(t),k=1,2,…,n.
System (3) with bilinear incidence has been researched by Ji et al. [16]. They proved that if ℛ0≤1, the solution of the model is fluctuating around a steady state, which is the disease-free equilibrium of the corresponding deterministic model; while if ℛ0>1, there is a stationary distribution, which means the disease will prevail.In this paper, we will prove that when the intensity of white noise is small, the solution of system (3) converges weakly to a singular measure (i.e., a distribution) if ℛ0≤1 and there exists an invariant distribution which is ergodic if ℛ0>1. This is the same situation as in Proposition 1. When the intensity of white noise is large, the disease will extinct exponentially regardless of the magnitude of ℛ0. Compared to the results of Ji et al. [16] our results can provide a deep insight into the dynamics of corresponding multigroup model.It is worth mentioning that by now there are many excellent works about stochastic single-group SIR model. Beretta et al. [17] considered a stochastic SIR model with time delays and obtained asymptotic mean square stability conditions for positive equilibrium. In this paper, the authors assumed that stochastic perturbations are of white noise type, which are directly proportional to distances S(t),  I(t),  R(t) from values of S*,  I*,  R*, influence on the S˙(t),  I˙(t),  R˙(t), respectively. Tornatore et al. [18] discuss a single-group case of (3). They proved that 0<β<min{λ+μ-σ2/2,2μ} is a sufficient condition for the asymptotic stability of the disease-free equilibrium. And only by computer simulations they showed that if min{λ+μ-σ2/2,2μ}<β<λ+μ+σ2/2, the disease-free equilibrium is stable and the disease does not occur; if (λ+μ)(R0-1)>σ2/2, the disease-free equilibrium is unstable. Ji et al. [19] considered the same model as in [18]. They deduce the globally asymptotical stability and exponential mean square stability of the disease-free equilibrium under some conditions and investigate the asymptotic behavior of the solution around the endemic equilibrium of the deterministic model.The rest of this paper is organized as follows. In Section2, we show there is a unique nonnegative solution of system (3). In Section 3, if ℛ0≤1 combined with small or large enough intensity of white noise, we show that the solution converges weakly to a singular measure. Section 4 focuses on the persistence of the disease. By choosing appropriate Lyapunov function, we show that there is a stationary distribution for system (3) and it is ergodic, if ℛ0>1. Finally, for the self-contained, we present an Appendix which contains some results used in the previous sections.Throughout this paper, let(Ω,ℱ,{ℱt}t≥0,P) be a complete probability space with a filtration {ℱt}t≥0 satisfying the usual conditions (i.e., it is right continuous and ℱ0 contains all P-null sets). Denote
(4)ℝ+3n={x∈ℝ3n:xk>0∀1≤k≤3n},ℝ-+3n={x∈ℝ3n:xk≥0∀1≤k≤3n}.
In general, consider a d-dimensional stochastic differential equation
(5)dx(t)=f(x(t),t)dt+g(x(t),t)dB(t)ont≥t0
with initial value x(t0)=x0∈ℝd. B(t) denotes d-dimensional standard Brownian motion defined on the above probability space. Define the differential operator L associated with (5) by
(6)L=∂∂t+∑k=1d‍fk(x,t)∂∂xk+12∑k,j=1d‍[gT(x,t)g(x,t)]kj∂2∂xk∂xj.
If L acts on a function V∈C2,1(Rd×ℝ-+;ℝ-+), then
(7)LV(x,t)=Vt(x,t)+Vx(x,t)f(x,t)+12trace[gT(x,t)Vxx(x,t)g(x,t)],
where Vt=∂V/∂t,  Vx=(∂V/∂x1,…,∂V/∂xd), Vxx=(∂2V/∂xk∂xj)d×d. By Itô's formula, if x(t)∈Rd, then
(8)dV(x(t),t)=LV(x(t),t)dt+Vx(x(t),t)g(x(t),t)dB(t).Consider (5), assume f(0,t)=0 and g(0,t)=0 for all t≥t0. So x(t)≡0 is a solution of (5), called the trivial solution or equilibrium position.


## 2. Existence and Uniqueness of the Nonnegative Solution

For a population model, one is interested in whether the solution is nonnegative and global. Hence in this section we show that the solution of system (3) is global and nonnegative. By making the change of variables and Lyapunov analysis method (see [20]), we will show global existence and uniqueness of the positive solution. From now on, we denote the solution (S1(t),I1(t),R1(t),…,Sn(t),In(t),Rn(t)) of system (3) as Y(t).Theorem 2.
If assumption(C1) holds and B=(βkj)n×n is irreducible, then for any initial value Y(0)∈R+3n, there is a unique solution Y(t) of system (3) on t≥0, and the solution will remain in R+3n with probability 1.Proof.
First consider system(9)duk(t)=[Λke-uk(t)-∑j=1n‍βkjϕj(evj(t))-dk-αk22]dt-αkdB1k(t),dvk(t)=[∑j=1n‍βkjeuk(t)-vk(t)ϕj(evj(t))-(ϵk+γk)-βk22]dt-βkdB2k(t),dwk(t)=[γkevk(t)-δkwk(t)]dt-σkwk(t)dB3k(t),k=1,2,…,n.
With initial value uk(0)=lnSk(0), vk(0)=lnIk(0), wk(0)=Rk(0), 1≤k≤n. Noting that (C1) holds, the coefficients of system (9) satisfy the local Lipschitz condition; then there is a unique local solution (u1(t),v1(t),w1(t),…,un(t),vn(t),wn(t)) on t∈[0,τe), where τe is the explosion time. Therefore, by Itô’s formula, it is easy to see (eu1(t),ev1(t),w1(t),…,eun(t),evn(t),wn(t)) is the unique positive local solution to system (3) with initial value Y(0)∈R+3n.
Next, we will prove that this solution is global. The following proof is almost the same as in the proof of Theorem 3.1 in [16]. We do not alter any words except replacing Ij by ϕj(Ij) and therefore we can obtain that LV≤M~. Here we omit the details.


## 3. Exponential Stability of Infectious Disease

It is clear thatP0=(Λ1/d1,0,0,Λ2/d2,0,0,…,Λn/dn,0,0) is the disease-free equilibrium of system (1) but not (3). For system (1), which has been mentioned in Introduction, P0 is globally stable if ℛ0≤1, which means the disease will die out after some period of time. Hence, it is interesting to study the disease-free equilibrium for controlling infectious disease. Although there is none of disease-free equilibrium of system (3), in this section, we can still present sufficient conditions for the disease to extinct exponentially.Theorem 3.
Assume(C1),(C2) hold and B=(βkj)n×n is irreducible. If 2dk>αk2,  k=1,2,…,n, then the solution Y(t) of system (3) with initial value Y(0)∈R+3n has the property limt→+∞(log(∑k=1n‍akIk)/t)≤c,a.s. If c<0, the disease will extinct exponentially almost surely. Here
(10)M0=M(S0)=(βkjSk0cjϵk+γk)n×n,Sk0=Λkdk,ak=ωk(ϵk+γk),k=1,2,…,n,(ω1,ω2,…,ωn)isalefteigenvectorofM0correspondingtoρ(M0),c≔max0≤k,j≤n{βkjcjaj}∑k=1n‍akαkSk02dk-αk2+min1≤k≤n{ϵk+γk}(ℛ0-1)1ℛ0≤1+max1≤k≤n{ϵk+γk}(ℛ0-1)1ℛ0≥1-12∑k=1n(1/βk2).Proof.
LetXk, 1≤k≤n be the solution of the equation:
(11)dXk(t)=[Λk-dkXk(t)]dt-αkXk(t)dB1k(t),Xk(0)=Sk(0).
By comparison theorem for stochastic equations, we have
(12)Sk(t)≤Xk(t),a.s.
It is easy to know that Xk(t) has the following property:
(13)limt→∞∫0t‍Xkdtt=Sk0a.s.,(14)limt→∞∫0t‍|Xk-Sk0|dtt≤αkSk02dk-αk2,if2dk>αk2,k=1,2,…,n.
As for the proof of the above two properties, the reader may refer to the proof of Theorems 3.1 and 4.1 in [21] for details.
SinceB=(βkj)n×n is irreducible, βkj≥0,  k,j=1,2,…,n and Sk0>0,  ϵk+γk>0,  ck>0,   k=1,2,…,n, then M0 is also nonnegative and irreducible. Hence by Lemma A.1, there is a left eigenvector of M0 corresponding to ρ(M0), which is denoted as (ω1,ω2,…,ωn) and ωk>0,  k=1,2,…,n; that is,
(15)(ω1,ω2,…,ωn)ρ(M0)=(ω1,ω2,…,ωn)M0.
Let ak=ωk/(ϵk+γk), k=1,2,…,n. Define C2-function V:R+n→R+ as follows:
(16)V(I1,I2,…,In)=∑k=1n‍akIk.
By Itô’s formula, we compute
(17)d(logV)=L(logV)dt-1V∑k=1n‍akβkIkdB2k,
where
(18)L(logV)=1V∑k=1n‍ak[∑j=1n‍βkjSkϕj(Ij)-(ϵk+γk)Ik]-∑k=1nak2βk2Ik22V2≤1V∑k=1n‍ak[∑j=1n‍βkjXkcjIj-(ϵk+γk)Ik]-∑k=1nak2βk2Ik22V2=1V∑k=1n‍ak[∑j=1n‍βkjcj(Xk-Sk0)Ij+∑j=1n‍βkjcjSk0Ij-(ϵk+γk)Ik]-∑k=1nak2βk2Ik22V2
according to (C2) and (12). Let
(19)H1=1V∑k=1n‍ak∑j=1n‍βkjcj(Xk-Sk0)Ij,H2=1V∑k=1n‍ak[∑j=1n‍βkjcjSk0Ij-(ϵk+γk)Ik],H3=-∑k=1nak2βk2Ik22V2.
Then
(20)H1≤1V∑k=1n‍ak∑j=1n‍βkjcj|Xk-Sk0|Ij=1V∑k=1n‍(ak|Xk-Sk0|∑j=1n‍βkjcjIj)≤max0≤k,j≤n{βkjcjaj}1V∑k=1n‍(ak|Xk-Sk0|∑j=1n‍ajIj)=max0≤k,j≤n{βkjcjaj}∑k=1n‍(ak|Xk-Sk0|).
In view of the definition of ak and (15), we have
(21)H2=1V[∑k=1n‍∑j=1nakβkjcjSk0Ij‍-∑k=1n‍ak(ϵk+γk)Ik]=1V[∑k=1n‍∑j=1n‍ωkβkjcjSk0ϵk+γkIj-∑k=1n‍ωkIk]=1V(ω1,ω2,…,ωn)×(M0(I1,I2,…,In)T-(I1,I2,…,In)T)=1V(ℛ0-1)∑k=1n‍ωkIk=1V(ℛ0-1)∑k=1n‍(ϵk+γk)akIk≤min1≤k≤n{ϵk+γk}(ℛ0-1)1ℛ0≤1+max1≤k≤n{ϵk+γk}(ℛ0-1)1ℛ0≥1:=m1,H3=-∑k=1nak2βk2Ik22V2=-∑k=1nak2βk2Ik22(∑k=1n(ak/βk)βkIk)2≤-∑k=1nak2βk2Ik22(∑k=1nak2βk2Ik2)(∑k=1n(1/βk2))=-12∑k=1n(1/βk2):=m2.
Hence,
(22)L(logV)≤max0≤k,j≤n{βkjcjaj}×∑k=1n‍(ak|Xk-Sk0|)+m1+m2,
which combined with (17) yields
(23)logV(t)t≤logV(0)t+max0≤k,j≤n{βkjcjaj}∑k=1n‍(ak∫0t‍|Xk-Sk0|dtt)+m1+m2-∑k=1n‍akβk∫0t‍Ik/VdB2kt.
By [22, Lemma 2.6], we get
(24)limt→+∞∫0t‍Ik/VdB2kt=0,a.s.
This together with (14) implies
(25)limt→∞logVt≤max0≤k,j≤n{βkjcjaj}∑k=1n‍akαkSk02dk-αk2+m1+m2,a.s.
Thus, the proof of Theorem 3 is completed.Remark 4.
From Theorem3 we know that if ℛ0<1 and αk, 1≤k≤n are small, the disease will extinct exponentially; that is, Ik→0 as t→∞. By using the same arguments as in [20, 21], we get Sk(t)→wνk and Rk(t)→0 (t→+∞), where →w means the weak convergence and νk is a distribution in R+ such that ∫0∞‍xνk(dx)=Sk0 and its density is (Mkx2p(x))-1, where Mk is a normal constant and p(x)=x2dk/αk2exp(2dkSk0/αk2x), x>0. That is, the solution of system (3) converges weakly to a singular measure (i.e., a distribution) when ℛ0<1 and αk, 1≤k≤n are small.Theorem 5.
Assume(C1), (C2) hold and B=(βkj)n×n is irreducible. Then the solution Y(t) of system (3) with initial value Y(0)∈R+3n has the property limt→+∞(log(∑k=1n‍Ik)/t)≤c,a.s. If c<0, the disease will extinct exponentially almost surely. Here c:=max1≤k,j≤n{βkjcj}∑k=1n‍Sk0 − min1≤k≤n{ϵk+γk}-1/2∑k=1n‍(1/βk2).Proof.
DefineC2-function V:R+n→R+ by V(I1,I2,…,In)=∑k=1n‍Ik. By Itô's formula, we compute
(26)d(logV)=L(logV)dt-1V∑k=1n‍βkIkdB2k,
where
(27)L(logV)=1V∑k=1n‍[∑j=1n‍βkjSkϕj(Ij)-(ϵk+γk)Ik]-∑k=1nβk2Ik22V2≤1V∑k=1n‍[∑j=1n‍βkjXkcjIj-(ϵk+γk)Ik]-∑k=1nβk2Ik22(∑k=1n(1/βk)βkIk)2≤max1≤k,j≤n{βkjcj}∑k=1n‍Xk-min1≤k≤n{ϵk+γk}-12∑k=1n(1/βk2).
Then
(28)logV(t)t≤logV(0)t+max1≤k,j≤n{βkjcj}∑k=1n‍∫0t‍Xkdtt-min1≤k≤n{ϵk+γk}-12∑k=1n(1/βk2)-∑k=1n‍βk∫0t‍Ik/VdB2kt.
Taking t→∞, in view of (24) and (13), we get
(29)limt→+∞logV(t)t≤max1≤k,j≤n{βkjcj}∑k=1n‍Sk0-min1≤k≤n{ϵk+γk}-12∑k=1n(1/βk2),a.s.Remark 6.
Theorem5 tells us that the large perturbation forces the infective to expire regardless of the magnitude of ℛ0.


## 4. Ergodicity of System (3)

When studying epidemic dynamical systems, we are also interested in when the disease will prevail and persist in a population. For a deterministic model, this is usually solved by showing that the endemic equilibrium is a global attractor or is globally asymptotically stable. But, for the stochastic system (3), there is no endemic equilibrium. In this section, we explore a weak stability. We show that there is a stationary distribution based on the theory of [23] (see Appendix), which reveals the disease will prevail.From the proof of Theorem2, we obtain LV≤M~. Let V~=V+M~. Then LV~≤V~ and it is clear that infY∈R+3n∖DkV~(Y)→∞ as k→∞, where Dk=(1/k,k)×(1/k,k)×⋯×(1/k,k). Hence, by Remark 2 of Theorem 4.1 of Hasminskii ([23], p. 86), we obtain that the solution Y(t) is a homogeneous Markov process in R+3n.Theorem 7.
Assume(C1)–(C4) hold, B=(βkj)n×n is irreducible, and ℛ0=ρ(M0)>1. If 0<αk2<dk,  0<βk2<γk+ϵk,  0<σk2<δk,  k=1,2,…,n such that
(30)δ<min{c-k2(dk-αk2)Sk*,bk4(ϵk+γk-βk2)(Ik*)2,mk2(δk-σk2)(Rk*)2,k=1,2,…,n}.
Then, for any initial value Y(0)∈R+3n, there is a stationary distribution μ(·) for system (3) and it has ergodic property, where M0=M(S0)=(βkjSk0cj/(ϵk+γk))n×n, Sk0=Λk/dk,
(31)δ=∑k=1n‍[(a+22c-k+bkSk*)Sk*αk2+(a+12c-k+bkIk*)Ik*βk2+mk(Rk*)2σk2(a+22c-k+bkSk*)],a,  bk,  mk,  k=1,2,…,n are defined as in the proof, P*=(S1*,I1*,R1*,…,Sn*,In*,Rn*) is the endemic equilibrium of system (1), and c-k,  k=1,2,…,n denote the cofactor of the kth diagonal element of LB-, B-=(β-kj)n×n=(βkjSk*ϕj(Ij*))n×n.Proof.
Sinceℛ0>1 and (C1)–(C4) hold, from Proposition 1 there is an endemic equilibrium P* of system (1). Then
(32)∑j=1n‍βkjSk*ϕj(Ij*)+dkSk*=Λk,∑j=1n‍βkjSk*ϕj(Ij*)=(ϵk+γk)Ik*.
Define
(33)V(Y)=a∑k=1n‍c-k(Sk-Sk*-Sk*logSkSk*+Ik-Ik*-Ik*logIkIk*)+∑k=1n‍c-k(Ik-Ik*-logIkIk*)+12∑k=1n‍bk(Sk-Sk*+Ik-Ik*)2+12∑k=1n‍lk(Sk-Sk*)2+12∑k=1n‍mk(Rk-Rk*)2:=aV1+V2+V3+V4+V5,
where a,  bk,  mk,  lk,  k=1,2,…,n are positive constants to be determined later, and c-k>0,  k=1,2,…,n according to Lemma A.1. Then V is positive definite. In view of (32) and the inequality (x+y)2≤2(x2+y2), by direct calculation, we get LVi  (i=1,2,…,5) as follows:
(34)LV1=∑k=1n‍c-k(1-Sk*Sk)(Λk-∑j=1n‍βkjSkϕj(Ij)-dkSk)+∑k=1n‍c-kSk*αk22+∑k=1n‍c-k(1-Ik*Ik)(∑j=1n‍βkjSkϕj(Ij)-(ϵk+γk)Ik)+∑k=1n‍c-kIk*βk22=∑k=1n‍c-k(1-Sk*Sk)×(∑j=1n‍β-kj(1-Skϕj(Ij)Sk*ϕj(Ij*))-dk(Sk-Sk*))+∑k=1n‍c-kSk*αk22+∑k=1n‍c-k(1-Ik*Ik)×(∑j=1n‍β-kjSkϕj(Ij)Sk*ϕj(Ij*)-∑j=1n‍β-kjIkIk*)+∑k=1n‍c-kIk*βk22=∑k=1n‍c-k[∑j=1n‍β-kj(1-Sk*Sk+ϕj(Ij)ϕj(Ij*)-Skϕj(Ij)Sk*ϕj(Ij*))-dk(Sk-Sk*)2Sk(1-Sk*Sk+ϕj(Ij)ϕj(Ij*)-Skϕj(Ij)Sk*ϕj(Ij*))∑j=1n‍β-kj]+∑k=1n‍c-kSk*αk22+∑k=1n‍c-k[∑j=1n‍β-kj(1-IkIk*+Skϕj(Ij)Sk*ϕj(Ij*)-Skϕj(Ij)Ik*Sk*ϕj(Ij*)Ik)]+∑k=1n‍c-kIk*βk22=∑k=1n‍c-k[-dk(Sk-Sk*)2Sk+∑j=1n‍β-kj(2-Sk*Sk-IkIk*+ϕj(Ij)ϕj(Ij*)-Skϕj(Ij)Ik*Sk*ϕj(Ij*)Ik)+Sk*αk22+Ik*βk22(Sk-Sk*)2Sk].
From the calculation of LV1, we directly get
(35)LV2=∑k=1n‍c-k[∑j=1n‍β-kj(1-IkIk*+Skϕj(Ij)Sk*ϕj(Ij*)-Skϕj(Ij)Ik*Sk*ϕj(Ij*)Ik)]+∑k=1n‍c-kIk*βk22(36)LV3=∑k=1n‍bk{Ik2βk22(Sk-Sk*+Ik-Ik*)×[Λk-dkSk-(ϵk+γk)Ik]+Sk2αk22+Ik2βk22}=∑k=1n‍bk{Ik2βk22(Sk-Sk*+Ik-Ik*)×[-dk(Sk-Sk*)-(ϵk+γk)(Ik-Ik*)]+Sk2αk22+Ik2βk22}≤-∑k=1n‍bk[dk(Sk-Sk*)2+(ϵk+γk)(Ik-Ik*)2-(γk+dk+ϵk)(Sk-Sk*)(Ik-Ik*)]+∑k=1n‍bkαk2[(Sk-Sk*)2+(Sk*)2]+∑k=1n‍bkβk2[(Ik-Ik*)2+(Ik*)2]=-∑k=1n‍bk[(dk-αk2)(Sk-Sk*)2+(ϵk+γk-βk2)(Ik-Ik*)2]-∑k=1n‍bk(γk+dk+ϵk)(Sk-Sk*)(Ik-Ik*)+∑k=1n‍bk[αk2(Sk*)2+βk2(Ik*)2]≤-∑k=1n‍bk[(dk-αk2-(dk+γk+ϵk)22(γk+dk+ϵk-βk2))×(Sk-Sk*)2+12(ϵk+γk-βk2)(Ik-Ik*)2(dk-αk2-(dk+γk+ϵk)22(γk+dk+ϵk-βk2))(dk+γk+ϵk)22(γk+dk+ϵk-βk2)]+∑k=1n‍bk[αk2(Sk*)2+βk2(Ik*)2],LV4=∑k=1n‍lk(Sk-Sk*)×(Λk-∑j=1n‍βkjSkϕj(Ij)-dkSk)+∑k=1n‍lkSk2αk22=∑k=1n‍lk(Sk-Sk*)×(∑j=1n‍βkj(Sk*ϕj(Ij*)-Skϕj(Ij))-dk(Sk-Sk*)∑j=1n‍)+∑k=1n‍lkSk2αk22=-∑k=1n‍∑j=1n‍lkβkjSk*×(Sk-Sk*)(ϕj(Ij)-ϕj(Ij*))-∑k=1n‍∑j=1n‍lkβkj(Sk-Sk*)2ϕj(Ij)-∑k=1n‍lkdk(Sk-Sk*)2+∑k=1n‍lkSk2αk22≤-∑k=1n‍∑j=1n‍lkβkjSk*(Sk-Sk*)(ϕj(Ij)-ϕj(Ij*))-∑k=1n‍lk(dk-αk2)(Sk-Sk*)2+∑k=1n‍lk(Sk*)2αk2,LV5=∑k=1n‍mk(Rk-Rk*)(γkIk-δkRk)+∑k=1n‍mkRk2σk22=∑k=1n‍mk(Rk-Rk*)×[γk(Ik-Ik*)-δk(Rk-Rk*)]+∑k=1n‍mkRk2σk22≤∑k=1n‍mkγk(Ik-Ik*)(Rk-Rk*)-∑k=1n‍mk(δk-σk2)(Rk-Rk*)2+∑k=1n‍mk(Rk*)2σk2≤∑k=1n‍mkγk22(δk-σk2)(Ik-Ik*)2-12∑k=1n‍mk(δk-σk2)(Rk-Rk*)2+∑k=1n‍mk(Rk*)2σk2.
By Lemma A.3, we know
(37)∑k=1n‍c-k(∑j=1n‍β-kjϕj(Ij)ϕj(Ij*)-∑j=1n‍β-kjϕk(Ik)ϕk(Ik*))=0,(38)∑k=1n‍c-k(∑j=1n‍β-kjlogϕj(Ij)ϕj(Ij*)-∑j=1n‍β-kjlogϕk(Ik)ϕk(Ik*))=0.
Besides note that a-1-loga≥0 for a>0; then
(39)1-Sk*Sk≤-logSk*Sk,(40)1-Skϕj(Ij)Ik*Sk*ϕj(Ij*)Ik≤-logSkϕj(Ij)Ik*Sk*ϕj(Ij*)Ik=-logSkSk*-logϕj(Ij)ϕj(Ij*)-logIk*Ik,(41)-logϕk(Ik)ϕk(Ik*)-logIk*Ik=-logϕk(Ik)Ik*ϕk(Ik*)Ik≤ϕk(Ik*)Ikϕk(Ik)Ik*-1.
Successively substituting (39), (40), (37), (38), and (41) into (34) yields
(42)LV1≤-∑k=1n‍c-kdk(Sk-Sk*)2Sk+∑k=1n‍c-k∑j=1n‍β-kj(ϕj(Ij)ϕj(Ij*)-logSk*Sk-logSkSk*-logϕj(Ij)ϕj(Ij*)-logIk*Ik-IkIk*+ϕj(Ij)ϕj(Ij*))+∑k=1n‍c-k[Sk*αk22+Ik*βk22]=-∑k=1n‍c-kdk(Sk-Sk*)2Sk+∑k=1n‍c-k∑j=1n‍β-kj(-logϕk(Ik)ϕk(Ik*)-logIk*Ik-IkIk*+ϕk(Ik)ϕk(Ik*))+∑k=1n‍c-k[Sk*αk22+Ik*βk22]≤-∑k=1n‍c-kdk(Sk-Sk*)2Sk+∑k=1n‍c-k∑j=1n‍β-kj(ϕk(Ik*)Ikϕk(Ik)Ik*-1-IkIk*+ϕk(Ik)ϕk(Ik*))+∑k=1n‍c-k[Sk*αk22+Ik*βk22].
In view of (C4), we obtain
(43)ϕk(Ik*)Ikϕk(Ik)Ik*-1-IkIk*+ϕk(Ik)ϕk(Ik*)=1ϕk(Ik*)Ikϕk(Ik)(ϕk(Ik*)-ϕk(Ik))×(ϕk(Ik*)Ik*-ϕk(Ik)Ik)≤0.
Then
(44)LV1≤-∑k=1n‍c-kdk(Sk-Sk*)2Sk+∑k=1n‍c-k[Sk*αk22+Ik*βk22].
By using inequality a-1-loga≥0 for a>0 again, it follows that
(45)-logSkSk*≤Sk*Sk-1.
Successively substituting (40), (37), (45), (41), and (43) into (35), we obtain
(46)LV2≤∑k=1n‍c-k[∑j=1n‍β-kj(-IkIk*+Skϕj(Ij)Sk*ϕj(Ij*)-logSkSk*-logϕj(Ij)ϕj(Ij*)-logIk*Ik)∑j=1n‍β-kj]+∑k=1n‍c-kIk*βk22=∑k=1n‍c-k[∑j=1n‍β-kj(-IkIk*+Skϕj(Ij)Sk*ϕj(Ij*)-logSkSk*-logϕk(Ik)ϕk(Ik*)-logIk*IkSkϕj(Ij)Sk*ϕj(Ij*)Skϕj(Ij)Sk*ϕj(Ij*))∑j=1n‍β-kj]+∑k=1n‍c-kIk*βk22≤∑k=1nc-k‍[∑j=1n‍β-kj(-IkIk*+Skϕj(Ij)Sk*ϕj(Ij*)+Sk*Sk-1+Ikϕk(Ik*)Ik*ϕk(Ik)-1Skϕj(Ij)Sk*ϕj(Ij*))∑j=1n‍β-kj]+∑k=1n‍c-kIk*βk22=∑k=1n‍c-k(∑j=1n‍β-kj[(SkSk*-1)(ϕj(Ij)ϕj(Ij*)-1)+(SkSk*+Sk*Sk-2)+(ϕk(Ik*)Ikϕk(Ik)Ik*-1-IkIk*+ϕk(Ik)ϕk(Ik*))(ϕj(Ij)ϕj(Ij*)-1)])+∑k=1n‍c-kIk*βk22≤∑k=1n‍c-k(∑j=1n‍β-kj[(SkSk*-1)(ϕj(Ij)ϕj(Ij*)-1)+(SkSk*+Sk*Sk-2)ϕj(Ij)ϕj(Ij*)]∑j=1n‍β-kj)+∑k=1n‍c-kIk*βk22=∑k=1n‍c-k((Sk-Sk*)2Sk∑j=1n‍βkj(Sk-Sk*)(ϕj(Ij)-ϕj(Ij*))+∑j=1n‍βkjϕj(Ij*)(Sk-Sk*)2Sk)+∑k=1n‍c-kIk*βk22.
Hence
(47)aLV1+LV2+LV4≤-∑k=1n‍[(adk-∑j=1n‍βkjϕj(Ij*))∑j=1n‍βkj(lkSk*-c-k)(Sk-Sk*)(ϕj(Ij)-ϕj(Ij*))+c-k(adk-∑j=1n‍βkjϕj(Ij*))(Sk-Sk*)2Sk]-∑k=1n‍lk(dk-αk2)(Sk-Sk*)2+∑k=1n‍c-k[βk22(a+1)Ik*+αk22(a+2)Sk*].
Choose lk=c-k/Sk* and a=max{∑j=1n‍βkjIj*/dk,k=1,2,…,n} such that lkSk*-c-k=0 and adk − ∑j=1n‍βkjIj* ≥ 0,  k=1,2,…,n. Then
(48)aLV1+LV2+LV4≤-∑k=1n‍c-kSk*(dk-αk2)(Sk-Sk*)2+∑k=1n‍c-k[βk22(a+1)Ik*+αk22(a+2)Sk*].
Furthermore,
(49)aLV1+LV2+LV3+LV4+LV5≤-∑k=1n‍[bk(dk-αk2-(dk+γk+ϵk)22(γk+dk+ϵk-βk2))+c-kSk*(dk-αk2)(dk+γk+ϵk)22(γk+dk+ϵk-βk2)](Sk-Sk*)2-12∑k=1n‍(bk(ϵk+γk-βk2)-mkγk2δk-σk2)(Ik-Ik*)2-12∑k=1n‍mk(δk-σk2)(Rk-Rk*)2+∑k=1n‍[(a+22c-k+bkSk*)Sk*αk2+(a+12c-k+bkIk*)×Ik*βk2+mk(Rk*)2σk2(a+22c-k+bkSk*)].
Choose bk>0 such that bk(dk-αk2-(dk+γk+ϵk)2/2(γk+dk+ϵk-βk2)) + (c-k/Sk*)(dk-αk2) ≥ (c-k/2Sk*)(dk-αk2) and choose mk = (bk/2)(ϵk+γk-βk2)((δk-σk2)/γk2) such that bk(ϵk+γk-βk2) − mkγk2/(δk-σk2) = (bk/2)(ϵk+γk-βk2); then
(50)aLV1+LV2+LV3+LV4+LV5≤-∑k=1n‍c-k2Sk*(dk-αk2)(Sk-Sk*)2-14∑k=1n‍bk(ϵk+γk-βk2)(Ik-Ik*)2-12∑k=1n‍mk(δk-σk2)(Rk-Rk*)2+∑k=1n‍[(a+22c-k+bkSk*)Sk*αk2+(a+12c-k+bkIk*)Ik*βk2+mk(Rk*)2σk2(a+22c-k+bkSk*)]:=-∑k=1n‍c-k2Sk*(dk-αk2)(Sk-Sk*)2-∑k=1n‍bk4(ϵk+γk-βk2)(Ik-Ik*)2-∑k=1n‍mk2(δk-σk2)(Rk-Rk*)2+δ.
Note that
(51)δ<min{c-k2(dk-αk2)Sk*,bk4(ϵk+γk-βk2)(Ik*)2,mk2(δk-σk2)(Rk*)2,k=1,2,…,n}.
Then the ellipsoid
(52)-∑k=1n‍c-k2Sk*(dk-αk2)(Sk-Sk*)2-∑k=1n‍bk4(ϵk+γk-βk2)(Ik-Ik*)2-∑k=1n‍mk2(δk-σk2)(Rk-Rk*)2+δ=0
lies entirely in R+3n. We can take U to be a neighborhood of the ellipsoid with U-⊆E3n=R+3n, so for x∈E3n∖U,  LV≤-K(K is a positive constant), which implies condition (B.2) in Lemma A.5 is satisfied.
Note that the diffusion matrix of system (3) is
(53)A=(aij)n×n=diag(α12S12,β12I12,σ12R12,…,αn2Sn2,βn2In2,σn2Rn2).
For any bounded domain D⊂R+3n there is
(54)M=min{αi2Si2,βi2Ii2,σi2Ri2,i=1,2,…,n,Y∈D-}>0
such that
(55)∑i,j=13n‍aijξiξj=∑i=1n‍αi2Si2ξ3i-22+∑i=1n‍βi2Ii2ξ3i-12+∑i=1n‍σi2Ri2ξ3i2≥M∥ξ2∥
for all Y∈D-,  ξ∈R3n, which implies condition (B.1) is also satisfied. Therefore, by Lemma A.5, the stochastic system (3) has a stationary distribution μ(·) and it is ergodic.



---

*Source: 917389-2013-07-16.xml*