Global Stability of Two-Group Epidemic Models with Distributed Delays and Random Perturbation

Abstract

We discuss a two-group SEIR epidemic model with distributed delays, incorporating random fluctuation around the endemic equilibrium. Our research shows that the endemic equilibrium of the model with distributed delays and random perturbation is stochastically asymptotically stable in the large. In addition, a sufficient stability condition is obtained by constructing suitable Lyapunov function.

1. Introduction

For the research of control of disease in populations, significant progress has been made in the theory and application of epidemiology modeling by mathematical research [1–8]. One of the main problems for the theory of differential equations and their applications is connected with stability. Most traditional compartmental models in mathematical epidemiology descend from the classical SIR model of Kermack and McKendrick [9], where the population is divided into the classes of susceptible, infected, and recovered individuals. For some diseases, such as influenza and tuberculosis, on adequate contact with an infectious individual, a susceptible becomes exposed for a while, that is, infected but not yet infectious. Thus it is realistic to introduce a latent compartment; the total population can be partitioned into four compartments: susceptible, latent or exposed, infectious, and recovered, with sizes denoted by S, E, I, and R, respectively. The resulting models are of SEI, SEIR, or SEIRS type, respectively. SEIR model has been widely discussed in the literature. Local and global stability analysis of the disease-free and endemic equilibria has been carried out using different assumptions and contact rates in [6, 7]. Greenhalgh [5] considered SEIR models that incorporate density dependence in the death rate. Korobeinikov [6] considers the global properties for SEIR and SEIS by means of the Lyapunov functions. In fact, there are real benefits to be gained in using stochastic models because real life is full of randomness and stochasticity. Recently, some stochastic epidemic models have been studied by many authors, see [3, 10]. Dalal et al. [4] showed that stochastic models had nonnegative solutions and carried out analysis on the asymptotic stability of models. Tornatore et al. [10] studied the stability of disease-free equilibrium of a stochastic SIR model with or without distributed time delay. On the other hand, taking into account environmental variability, white noise stochastic perturbations around the positive endemic equilibrium of epidemic models was considered in [3, 11]. Beretta et al. proved the stability of epidemic model with stochastic time delays influenced by probability under certain conditions [3]. Such type of stochastic perturbations firstly was proposed in [3, 12] and later was successfully used in many other papers for many other different systems (see, for instance, [13–20]). A more general multigroup epidemic model is proposed to describe the disease spread in a heterogeneous host population with general age structure and varying infectivity by Li et al. [1]. They investigated a class of multigroup epidemic models with distributed delays and established the global dynamics determined by the basic reproduction number 𝑅₀. More specifically, they proved that, if 𝑅₀≤1, then the disease-free equilibrium is globally asymptotically stable; if 𝑅₀>1, then there exists a unique endemic equilibrium, and it is globally asymptotically stable. However, to the best of the authors’ knowledge, no literature exists regarding SEIR model with random perturbation. Thus, the current study hopes to serve such a need and is inspired by the report of [1]. In this paper, based on the SEIR model of [1], we consider the white noise stochastic perturbations around its endemic equilibrium and use the methods, which is similar to [3]. We construct a class of the Lyapunov functions, as it is useful to study the global properties of stochastic models. By means of it, we prove the SEIR model is stochastically asymptotically stable in the large under certain condition.

The paper is organized as follows. In Section 2 we recall the deterministic SEIR model and its main results by Li et al. [1]. We introduce the model with stochastic perturbations around the endemic equilibrium in Section 3. In Section 4 the global stability of the endemic equilibrium is proved by the method of the Lyapunov functions.

2. Preliminaries

We briefly review the following results obtained by Li et al. [1]. Let 𝑆_𝑘,𝐸_𝑘,𝐼_𝑘, and 𝑅_𝑘 denote the susceptible, infected but noninfectious, infectious, and recovered populations in the 𝑘th group, respectively. Let 𝑖_𝑘(𝑡,𝑟) denote the population of infectious individuals in the 𝑘th group with respect to the age of infection 𝑟 at time 𝑡, and 𝐼_𝑘(𝑡)=∫^∞_𝑟=0𝑖_𝑘(𝑡,𝑟)𝑑𝑟. Let ℎ_𝑘(𝑟)⩾0 be a continuous kernel function that represents the infectivity at the age of infection 𝑟. The disease incidence in the 𝑘-th group, assuming a bilinear incidence form, can be calculated as ∑^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)∫^∞_𝑟=0ℎ_𝑗(𝑟)𝑖_𝑗(𝑡,𝑟)𝑑𝑟, where the sum takes into account cross-infections from all groups and 𝛽_𝑘𝑗 represents the transmission coefficient between compartments 𝑆_𝑘 and 𝐼_𝑗. In the special case ℎ_𝑘(𝑟)≡1, the incidence becomes ∑^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)𝐼_𝑗(𝑡) as in [2]. Therefore, the model in [2] can be generalized to the following system of differential equations𝑆^′_𝑘=Λ_𝑘−^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)^∞_𝑟=0ℎ_𝑗(𝑟)𝑖_𝑗(𝑡,𝑟)𝑑𝑟−𝑑^𝑆_𝑘𝑆_𝑘,𝐸^′_𝑘=^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)^∞_𝑟=0ℎ_𝑗(𝑟)𝑖_𝑗(𝑡,𝑟)𝑑𝑟−𝑑^𝐸_𝑘+𝜖_𝑘𝐸_𝑘,𝐼^′_𝑘=𝜖_𝑘𝐸_𝑘−𝑑^𝐼_𝑘+𝛾_𝑘𝐼_𝑘,𝑅^′_𝑘=𝛾_𝑘𝐼_𝑘−𝑑^𝑅_𝑘𝑅_𝑘,𝑘=1,2,…,𝑛. (2) Here Λ_𝑘 represents influx of individuals into the 𝑘th group, 𝑑^𝑆_𝑘,𝑑^𝐸_𝑘,𝑑^𝐼_𝑘, and 𝑑^𝑅_𝑘 represent death rates of 𝑆,𝐸,𝐼, and 𝑅 populations in the 𝑘th group, respectively, 𝜖_𝑘 represents the rate of becoming infectious after a latent period in the 𝑘-th group, and 𝛾_𝑘 represents the recovery rate of infectious individuals in the 𝑘-th group. All parameter values are assumed to be nonnegative and Λ_𝑘,𝑑^𝑆_𝑘,𝑑^𝐸_𝑘>0 for all 𝑘. Note that𝜕+𝜕𝜕𝑡𝑖𝜕𝑟^𝑘(𝑡,𝑟)=−𝑑^𝐼_𝑘+𝛾_𝑘𝑖_𝑘(),𝑖𝑡,𝑟_𝑘(𝑡,0)=𝜖_𝑘𝐸_𝑘(𝑡), (2.2) whose solution is𝑖_𝑘(𝑡,𝑟)=𝑖_𝑘(𝑡−𝑟,0)𝑒^{−(𝑑𝐼𝑘+𝛾𝑘)𝑟}=𝜖_𝑘𝐸_𝑘(𝑡−𝑟)𝑒^{−(𝑑𝐼𝑘+𝛾𝑘)𝑟}. (2.3) Substituting (2.3) into (2), we obtain𝑆^′_𝑘=Λ_𝑘−^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)^∞_𝑟=0ℎ_𝑗(𝑟)𝜖_𝑗𝐸_𝑗(𝑡−𝑟)𝑒^{−(𝑑𝐼𝑗+𝛾𝑗)𝑟}𝑑𝑟−𝑑^𝑆_𝑘𝑆_𝑘,𝐸^′_𝑘=^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)^∞_𝑟=0ℎ_𝑗(𝑟)𝜖_𝑗𝐸_𝑗(𝑡−𝑟)𝑒^{−(𝑑𝐼𝑗+𝛾𝑗)𝑟}𝑑𝑟−𝑑^𝐸_𝑘+𝜖_𝑘𝐸_𝑘,𝐼^′_𝑘=𝜖_𝑘𝐸_𝑘−𝑑^𝐼_𝑘+𝛾_𝑘𝐼_𝑘,𝑅^′_𝑘=𝛾_𝑘𝐼_𝑘−𝑑^𝑅_𝑘𝑅_𝑘,𝑘=1,2,…,𝑛. (2.4)

Since the variables 𝐼_𝑘 and 𝑅_𝑘 do not appear in the first two equations of (2.4), Li et al. consider the following reduced system with distributed time delays and general kernel functions [1]:𝑆^′_𝑘=Λ_𝑘−^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)^∞_𝑟=0𝑓_𝑗(𝑟)𝐸_𝑗(𝑡−𝑟)𝑑𝑟−𝑑^𝑆_𝑘𝑆_𝑘,𝐸^′_𝑘=^𝑛_𝑗=1𝛽_𝑘𝑗𝑆_𝑘(𝑡)^∞_𝑟=0𝑓_𝑗(𝑟)𝐸_𝑗(𝑡−𝑟)𝑑𝑟−𝑑^𝐸_𝑘+𝜖_𝑘𝐸_𝑘. (2.5)

Here the kernel function 𝑓_𝑘(𝑟)⩾0 is continuous and 1⩾∫^∞_𝑟=0𝑓_k(𝑟)𝑑𝑟=ℎ_𝑘>0. System (2.5) can be interpreted as a multigroup model for an infectious disease whose latent period 𝑟 in hosts has a general probability density function (1/ℎ_𝑘)𝑓_𝑘(𝑟)𝑑𝑟, for the 𝑘-th group. Let 𝑆⁰_𝑘=Λ_𝑘/𝑑^𝑆_𝑘,ℎ_𝑘=∫^∞_𝑟=0𝑓_𝑘(𝑟)𝑑𝑟. The next-generation matrix for system (2.5) is𝑀₀=𝛽_𝑘𝑗𝑆⁰_𝑘ℎ_𝑘𝑑^𝐸_𝑘+𝜖_{𝑘𝑛×𝑛}. (2.6) Define the basic reproduction number as the spectral radius of 𝑀₀,𝑅₀=𝜌𝑀₀. (2.7) In the special case when 𝑓_𝑘(𝑟) is an exponential function, 𝑅₀ reduces to that for the resulting ODE models. Make the following assumption on the kernel function 𝑓_𝑘(𝑟) in (2.5):^∞_𝑟=0𝑓_𝑘(𝑟)𝑒^{−𝜆𝑘𝑟}𝑑𝑟<∞, (2.8) where 𝜆_𝑘 is a positive number, 𝑘=1,2,…,𝑛. Define the following Banach space of fading memory type:𝐶_𝑘=𝜙∈((−∞,0],ℝ)𝜙(𝑠)𝑒^𝜆𝑘𝑠isuniformlycontinuouson(−∞,0]andsup_𝑠⩽0𝜙(𝑠)𝑒^𝜆𝑘𝑠<∞ (2.9) with norm ‖𝜙‖_𝑘=sup_𝑠⩽0|𝜙(𝑠)|𝑒^𝜆𝑘𝑠. For 𝜙∈𝐶_𝑘, let 𝜙_𝑡∈𝐶_𝑘 be such that 𝜙_𝑡(𝑠)=𝜙(𝑡+𝑠),𝑠∈(−∞,0]. Let 𝑆_𝑘,0∈ℝ₊ and 𝜙_𝑘∈𝐶_𝑘 such that 𝜙_𝑘(𝑠)⩾0,𝑠∈(−∞,0]. We consider solutions (𝑆₁(𝑡),𝐸_1𝑡,𝑆₂(𝑡),𝐸_2𝑡,…,𝑆_𝑛(𝑡),𝐸_𝑛𝑡) of system (2.5) with initial conditions𝑆_𝑘(0)=𝑆_𝑘,0,𝐸_𝑘,0=𝜙_𝑘,𝑘=1,2,…,𝑛. (2.10)

Standard theory of functional differential equations implies 𝐸_𝑘𝑡∈𝐶_𝑘 for 𝑡>0. We consider system (2.5) in the phase space𝑋=^𝑛_𝑘=1ℝ×𝐶_𝑘. (2.11) It can be verified that solutions of (2.5) in 𝑋 initial conditions (2.10) remain nonnegative. In particular, 𝑆_𝑘(𝑡)>0 for t>0. The following set is positively invariant for system (2.5):Θ=𝑆₁,𝐸₂(⋅),𝑆₂,𝐸₂(⋅),…,𝑆_𝑛,𝐸_𝑛(⋅)∈𝑋∣0⩽𝑆_𝑘⩽Λ_𝑘𝑑^𝑆_𝑘,0⩽𝑆_𝑘+𝐸_𝑘(0)⩽Λ_𝑘𝑑^∗_𝑘,𝐸_𝑘(𝑠)(].⩾0,𝑠∈−∞,0,𝑘=1,2,…,𝑛 (2.12) All positive semiorbits in Θ are precompact in 𝑋 and thus have nonempty 𝜔-limit sets. We have the following results [1].

Lemma 2.1. All positive semi-orbits in Θ have non-empty 𝜔-limit sets. Let Θ^∘=𝑆₁,𝐸₂(⋅),𝑆₂,𝐸₂(⋅),…,𝑆_𝑛,𝐸_𝑛(⋅)∈𝑋∣0<𝑆_𝑘<Λ_𝑘𝑑^𝑆_𝑘,0<𝑆_𝑘+𝐸_𝑘(0)<Λ_𝑘𝑑^∗_𝑘,𝐸_𝑘(𝑠)(].>0,𝑠∈−∞,0,𝑘=1,2,…,𝑛 (2.13) It can be shown that Θ^∘ is the interior of Θ.

Lemma 2.2. Assume that 𝐵=(𝛽_𝑘𝑗) is irreducible.(1)If 𝑅₀⩽1, then 𝑃₀ is the only equilibrium for system (2.5) in Θ.(2)If 𝑅₀>1, then there exist two equilibria for system (2.5) in Θ: the disease-free equilibrium 𝑃₀ and a unique endemic equilibrium 𝑃^∗.

Lemma 2.3. Assume that 𝐵=(𝛽_𝑘𝑗) is irreducible.(1)If 𝑅₀⩽1, then the disease-free equilibrium 𝑃₀ of system (2.5) is globally asymptotically stable in Θ. If 𝑅₀>1, then 𝑃₀ is unstable.(2)If 𝑅₀>1, then the endemic equilibrium 𝑃^∗ of system (2.5) is globally asymptotically stable in Θ^∘.

Biologically, Lemma 2.3 implies that if the basic reproduction number 𝑅₀⩽1, then the disease always dies out from all groups; if 𝑅₀>1, then the disease always persists in all groups at the unique endemic equilibrium level, irrespective of the initial conditions.

3. Stochastic Model Derivation

In this paper, based on system (2.5), we consider the case of 𝑘=1,2 in the following system(3.1):𝑆_𝑘̇=Λ_𝑘−²_𝑗=1𝛽_𝑘𝑗𝑆_k^∞_𝑟=0𝑓_𝑗(𝑟)𝐸_𝑗(𝑡−𝑟)𝑑𝑟−𝑑^𝑆_𝑘𝑆_𝑘,𝐸̇_𝑘=²_𝑗=1𝛽_𝑘𝑗𝑆_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝐸_𝑗(𝑡−𝑟)𝑑𝑟−𝑑^𝐸_𝑘+𝜖_𝑘𝐸_𝑘,𝐼_𝑘̇=𝜖_𝑘𝐸_𝑘−𝑑^𝐼_𝑘+𝛾_𝑘𝐼_𝑘,𝑅̇_𝑘=𝛾_𝑘𝐼_𝑘−𝑑^𝑅_𝑘𝑅_𝑘. (3.1) It is easy to see that equilibrium for system (3.1) is given by 𝑃^∗=(𝑆^∗;𝐸^∗;𝐼^∗;𝑅^∗),𝑆^∗_𝑘=𝑑^𝐸_𝑘+𝜖_𝑘∑²_𝑗=1𝛽_𝑘𝑗ℎ_𝑗,𝐸^∗_𝑘=Λ_𝑘−𝑑^𝑆_𝑘∑²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗,𝐼^∗_𝑘=𝜖_𝑘𝐸^∗_𝑘𝑑^𝐼_𝑘+𝛾_𝑘,𝑅^∗_𝑘=𝛾_𝑘𝐼^∗_𝑘𝑑^𝑅_𝑘. (3.2)

We assume stochastic perturbations are of white noise type, which are directly proportional to distances 𝑆_𝑘(𝑡),𝐸_𝑘(𝑡),𝐼_𝑘(𝑡),𝑅_𝑘(𝑡) from values of 𝑆^∗_𝑘,𝐸^∗_𝑘,𝐼^∗_𝑘,𝑅^∗_𝑘, influence the 𝑆_𝑘̇(𝑡),𝐸̇_𝑘(𝑡),𝐼_𝑘̇(𝑡),𝑅̇_𝑘(𝑡) respectively. So system (2.4) results in𝑆_𝑘̇=Λ_𝑘−²_𝑗=1𝛽_𝑘𝑗𝑆_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝐸_𝑗(𝑡−𝑟)𝑑𝑟−𝑑^𝑆_𝑘𝑆_𝑘+𝜎_1𝑘𝑆_𝑘−𝑆^∗_𝑘𝐵̇_1𝑘,𝐸̇_𝑘=²_𝑗=1𝛽_𝑘𝑗𝑆_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑗𝐸_𝑗(𝑡−𝑟)𝑑𝑟−𝑑^𝐸_𝑘+𝜖_𝑘𝐸_𝑘+𝜎_2𝑘𝐸_𝑘−𝐸^∗_𝑘𝐵̇_2𝑘,𝐼_𝑘̇=𝜖_𝑘𝐸_𝑘−𝑑^𝐼_𝑘+𝛾_𝑘𝐼_𝑘+𝜎_3𝑘𝐼_𝑘−𝐼^∗_𝑘𝐵̇_3𝑘,𝑅̇_𝑘=𝛾_𝑘𝐼_𝑘−𝑑^𝑅_𝑘𝑅_𝑘+𝜎_4𝑘𝑅_𝑘−𝑅^∗_𝑘𝐵̇_4𝑘,𝑘=1,2, (3.3) where 𝐵_1𝑘(𝑡),𝐵_2𝑘(𝑡),𝐵_3𝑘(𝑡),𝐵_4𝑘(𝑡) are independent standard Brownian motions and 𝜎²_𝑖𝑘>0 represent the intensities of 𝐵_𝑖𝑘(𝑡)(𝑖=1,2,3,4), respectively. Obviously, stochastic system (3.3) has the same equilibrium points as system (3.1). In the next section, we will investigate the stability of the equilibrium 𝑃^∗ of system (3.3). Below we will construct a class of different Lyapunov functions to achieve our proof under certain conditions.

4. Stochastic Stability of the Endemic Equilibrium

In this paper, unless otherwise specified, let (Ω,ℱ,{ℱ_𝑡}_𝑡⩾𝑡0,𝑃) be a complete probability space with a filtration {ℱ_𝑡}_𝑡⩾𝑡0 satisfying the usual conditions (i.e., it is increasing and right continuous while ℱ₀ contains all 𝑃-null sets). Let 𝛽_𝑖(𝑡) be the Brownian motions defined on this probability space. If 𝑅₀>1, then the stochastic system (3.3) can be centered at its endemic equilibrium 𝑃^∗=(𝑆^∗_𝑘,𝐸^∗_𝑘,𝐼^∗_𝑘,𝑅^∗_𝑘), by the change of variables𝑢_𝑘=𝑆_𝑘−𝑆^∗_𝑘,𝑣_𝑘=𝐸_𝑘−𝐸^∗_𝑘,𝑤_𝑘=𝐼_𝑘−𝐼^∗_𝑘,𝑧_𝑘=𝑅_𝑘−𝑅^∗_𝑘. (4.1) we obtain 𝑢̇_𝑘=−²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟−²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟−²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝐸^∗_𝑗ℎ_𝑗−𝑑^𝑆_𝑘𝑢_𝑘+𝜎_1𝑘𝑢_𝑘𝐵̇_1𝑘,𝑣_𝑘̇=²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝐸^∗_𝑗ℎ_𝑗−𝑑^𝐸_𝑘+𝜖_𝑘𝑣_𝑘+𝜎_2𝑘𝑣_𝑘𝐵̇_2𝑘,𝑤̇_𝑘=𝜖_𝑘𝑣_𝑘−𝑑^𝐼_𝑘+𝛾_𝑘𝑤_𝑘+𝜎_3𝑘𝑤_𝑘𝐵̇_3𝑘,𝑧_𝑘̇=𝛾_𝑘𝑤_𝑘−𝑑^𝑅_𝑘𝑧_𝑘+𝜎_4𝑘𝑧_𝑘𝐵̇_4𝑘,𝑘=1,2. (4.2)

It is easy to see that the stability of the equilibrium of the system (3.3) is equivalent to the stability of zero solution of system (4.2). Before proving the main theorem we put forward a lemma in [21]. Consider the 𝑑-dimensional stochastic differential equation𝑑𝐱(𝑡)=𝑓(𝐱(𝑡),𝑡)𝑑𝑡+𝑔(𝐱(𝑡),𝑡)𝑑𝐵(𝑡),𝑡⩾𝑡₀. (4.3) Assume that the assumptions of the existence-and-uniqueness theorem are fulfilled. Hence, for any given initial value 𝐱(𝑡₀)=𝐱₀∈𝑅^𝑑, (4.3) has a unique global solution that is denoted by 𝐱(𝑡;𝑡₀,𝑥₀). Assume furthermore that 𝑓(𝟎,𝑡)=0 and 𝑔(𝟎,𝑡)=0 for all 𝑡⩾𝑡₀. So (4.3) has the solution 𝐱(𝑡)≡0 corresponding to the initial value 𝐱(𝑡₀)=0. This solution is called the trivial solution or equilibrium position. Denote by 𝐶^2,1(𝑅^𝑑×[𝑡₀,∞);𝑅₊) the family of all nonnegative functions 𝑉(𝐱,𝑡) defined on 𝑅^𝑑×[𝑡₀,∞) such that they are continuously twice differentiable in 𝐱 and once in 𝑡. Define the differential operator 𝐿 associated with (4.3) by𝜕𝐿=+𝜕𝑡^𝑑_𝑖=1𝑓_𝑖()𝜕𝐱,𝑡𝜕𝑥_𝑖+12^𝑑_𝑖,𝑗=1𝑔^𝑇(𝐱,𝑡)𝑔(𝐱,𝑡)_𝑖𝑗𝜕²𝜕𝑥_𝑖𝑥_𝑗. (4.4) If 𝐿 acts on a function 𝑉∈𝐶^2,1(𝑅^𝑑×[𝑡₀,∞);𝑅₊), then𝐿𝑉(𝐱,𝑡)=𝑉_𝑡(𝐱,𝑡)+𝑉_𝑡()𝑓()+1𝐱,𝑡𝐱,𝑡2Tr𝑔^𝑇(𝐱,𝑡)𝑉_𝑥𝑥(𝐱,𝑡)𝑔(𝐱,𝑡). (4.5)

Definition 4.1. (1) The trivial solution of (4.3) is said to be stochastically stable or stable in probability if for every pair of 𝜀∈(0,1) and 𝑟>0, there exists a 𝛿=𝛿(𝜀,𝑟,𝑡₀)>0 such that 𝑃𝐱𝑡;𝑡₀,𝑥₀<𝑟∀𝑡⩾𝑡₀⩾1−𝜀 (4.6) whenever |𝐱₀|<𝛿. Otherwise, it is said to be stochastically unstable.
(2) The trivial solution is said to be stochastically asymptotically stable if it is stochastically stable, and, moreover, for every 𝜀∈(0,1), there exists a 𝛿₀=𝛿₀(𝜀,𝑡₀)>0 such that 𝑃lim_𝑡→∞𝐱𝑡;𝑡₀,𝑥₀=0⩾1−𝜀 (4.7) whenever |𝐱₀|<𝛿₀.
(3)The trivial solution is said to be stochastically asymptotically stable in the large if it is stochastically asymptotically stable and, moreover, for all 𝐱₀∈𝑅^𝑑, 𝑃lim_𝑡→∞𝐱𝑡;𝑡₀,𝑥₀=0=1. (4.8)

Lemma 4.2 (see [21]). If there exists a positive-definite decrescent radially unbounded function 𝑉(𝐱,𝑡)∈𝐶^2,1(𝑅^𝑑×[𝑡₀,∞);𝑅₊) such that 𝐿𝑉(𝐱,𝑡) is negative definite, then the trivial solution of (4.3) is stochastically asymptotically stable in the large.

From the above lemma, we can obtain the stochastically asymptotically stability of equilibrium as follows.

Theorem 4.3. Assume that 𝐵=(𝛽_𝑘𝑗) is irreducible and 𝑅₀>1; then, if the following condition is satisfied 𝜎²_1𝑘<2𝑑^𝑆_𝑘,𝜎²_2𝑘<2𝑑^𝐸_𝑘+𝜖_𝑘,𝜎²_3𝑘<2𝑑^𝐼_𝑘+𝛾_𝑘,𝜎²_4𝑘<2𝑑^𝑅_𝑘, (4.9) the endemic equilibrium 𝑃^∗ of system (3.3) is stochastically asymptotically stable in the large.

Proof. It is easy to see that we only need to prove the zero solution of (4.2) is stochastically asymptotically stable in the large. Let 𝐱(𝑡)=(𝑢(𝑡),𝑣(𝑡),𝑤(𝑡),𝑧(𝑡))^𝑇. We define the Lyapunov function 𝑉(𝐱(𝑡)) as follows: 𝑉(𝐱)=²_𝑘=1𝑎_𝑘𝑢²_𝑘+𝑣²_𝑘+𝑏_𝑘𝑢_𝑘+𝑣_𝑘²+𝑎²_𝑘𝑤²_𝑘+𝑎⁴_𝑘𝑧²_𝑘+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)^𝑡_𝑡−𝑟𝑣²_𝑘(𝜏)𝑑𝜏𝑑𝑟, (4.10) where 𝑎_𝑘>0,𝑏_𝑘>0 are real positive constants to be chosen later. So it is obvious that 𝑉(𝐱) is positive definite and decrescent.
Using Itô’s formula, we compute =𝐿𝑉_2𝑘=12𝑎_𝑘𝑢_𝑘+𝑏_𝑘𝑢_𝑘+𝑣_𝑘×−²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟−²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟−²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝐸^∗_𝑗ℎ_𝑗−𝑑^𝑆_𝑘𝑢_𝑘+2𝑣_𝑘+𝑏_𝑘𝑢_𝑘+𝑣_𝑘²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗()+𝑡−𝑟𝑑𝑟_2𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝐸^∗_𝑗ℎ_𝑗−𝑑^𝐸_𝑘+𝜖_𝑘𝑣_𝑘+2𝑎²_𝑘𝑤_𝑘𝜖_𝑘𝑣_𝑘−𝑑^𝐼_𝑘+𝛾_𝑘𝑤_𝑘+2𝑎⁴_𝑘𝑧_𝑘𝛾_𝑘𝑤_𝑘−𝑑^𝑅_𝑘𝑧_𝑘+𝑎_𝑘+𝑏_𝑘𝜎²_1𝑘𝑢²_𝑘+1+𝑏_𝑘𝜎²_2𝑘𝑣²_𝑘+𝑎²_𝑘𝜎²_3𝑘𝑤²_𝑘+𝑎⁴_𝑘𝜎²_4𝑘𝑧²_𝑘+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘ℎ_𝑗𝑣²_𝑗−²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣²_𝑗(𝑡−𝑟)𝑑𝑟. (4.11) Using (3.2), we obtain 𝐿𝑉=²_𝑘=1−𝑢²_𝑘𝑎_𝑘+𝑏_𝑘2𝑑^𝑆_𝑘−𝜎²_1𝑘+2²_𝑗=1𝑎_𝑘𝛽_𝑘𝑗𝐸^∗_𝑗ℎ_𝑗−𝑣²_𝑘1+𝑏_𝑘2²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗−𝜎²_2𝑘−𝑎²_𝑘𝑤²_𝑘2𝑑^𝐼_𝑘+𝛾_𝑘−𝜎²_3𝑘−𝑎⁴_𝑘𝑧²_𝑘2𝑑^𝑅_𝑘−𝜎²_4𝑘+2²_𝑗=1𝛽_𝑘𝑗𝐸^∗_𝑗ℎ_𝑗−𝑏_𝑘𝑑^𝑆_𝑘+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗×𝑢_𝑘𝑣_𝑘+2𝑎²_𝑘𝜖_𝑘𝑣_𝑘𝑤_𝑘+𝛾_𝑘𝑤_𝑘𝑧_𝑘+2²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+2²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_kℎ_𝑗𝑣²_𝑗−²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣²_𝑗(𝑡−𝑟)𝑑𝑟. (4.12) In (4.12), we choose 𝑏_𝑘=∑²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗𝑑^𝑆_𝑘+∑²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗. (4.13) Then 𝐿𝑉=²_𝑘=1−𝑢²_𝑘𝑎_𝑘+𝑏_𝑘2𝑑^𝑆_𝑘−𝜎²_1𝑘+2²_𝑗=1𝑎_𝑘𝛽_𝑘𝑗𝐸^∗_𝑗ℎ_𝑗−𝑣²_𝑘1+𝑏_𝑘2²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗−𝜎²_2𝑘−𝑎²_𝑘𝑤²_𝑘2𝑑^𝐼_𝑘+𝛾_𝑘−𝜎²_3𝑘−𝑎⁴_𝑘𝑧²_𝑘2𝑑^𝑅_𝑘−𝜎²_4𝑘+2𝑎²_𝑘𝜖_𝑘𝑣_𝑘𝑤_𝑘+𝑎²_𝑘𝛾_𝑘𝑤_𝑘𝑧_𝑘+2²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+2²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗()+𝑡−𝑟𝑑𝑟_2𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘ℎ_𝑗𝑣²_𝑗−²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣²_𝑗(𝑡−𝑟)𝑑𝑟. (4.14) Moreover, using the Cauchy inequality to 2𝑎²_𝑘𝜖_𝑘𝑣_𝑘𝑤_𝑘 and 2𝑎⁴_𝑘𝛾_𝑘𝑤_𝑘𝑧_𝑘, we can obtain 2𝑎²_𝑘𝜖_𝑘𝑣_𝑘𝑤_𝑘⩽𝑎²_𝑘𝜖_𝑘𝑣²_𝑘𝑎_𝑘+𝑎_𝑘𝑤²_𝑘,2𝑎⁴_𝑘𝛾_𝑘𝑤_𝑘𝑧_𝑘⩽𝑎⁴_𝑘𝛾_𝑘𝑤²_𝑘𝑎_𝑘+𝑎_𝑘𝑧²_𝑘,2²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟⩽𝑎_𝑘²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘𝑢²_𝑘+^∞_𝑟=0𝑓_𝑗(𝑟)𝑣²_𝑗()+𝑡−𝑟𝑑𝑟_2𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘𝑣²_𝑘+^∞_𝑟=0𝑓_𝑗(𝑟)𝑣²_𝑗(𝑡−𝑟)𝑑𝑟. (4.15) Substituting (4.15) into (4.14) as well as using 1⩾∫^∞_𝑟=0𝑓_𝑘(𝑟)𝑑𝑟=ℎ_𝑘>0, yields 𝐿𝑉⩽²_𝑘=1−𝑢²_𝑘𝑏_𝑘2𝑑^𝑆_𝑘−𝜎²_1𝑘−²_𝑗=1𝑎_𝑘𝛽_𝑘𝑗𝑆^∗_𝑘−𝑣²_𝑘1+𝑏_𝑘2²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘−𝜎²_2𝑘−𝑎_𝑘𝜖_𝑘+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘−𝑎²_𝑘𝑤²_𝑘2𝑑^𝐼_𝑘+𝛾_𝑘−𝜎²_3𝑘−𝑎_𝑘𝜖_𝑘+𝛾_𝑘−𝑎⁴_𝑘𝑧²_𝑘2𝑑^𝑅_𝑘−𝜎²_4𝑘−𝑎_𝑘𝛾_𝑘+2²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘1+𝑎_𝑘𝑣²_𝑗=𝐿₀𝑉+2²²_𝑘=1𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟, (4.16) where 𝐿₀𝑉=²_𝑘=1−𝑢²_𝑘𝑏_𝑘2𝑑^𝑆_𝑘−𝜎²_1𝑘−²_𝑗=1𝑎_𝑘𝛽_𝑘𝑗𝑆^∗_𝑘−𝑣²_𝑘2𝑏_𝑘²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘−1+𝑏_𝑘𝜎²_2𝑘−𝑎_𝑘𝜖_𝑘+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘−𝑎²_𝑘𝑤²_𝑘2𝑑^𝐼_𝑘+𝛾_𝑘−𝜎²_3𝑘−𝑎_𝑘𝜖_𝑘+𝛾_𝑘−𝑎⁴_𝑘𝑧²_𝑘2𝑑^𝑅_𝑘−𝜎²_4𝑘−𝑎_𝑘𝛾_𝑘. (4.17) From (4.9) it follows that there exists 𝑎_𝑘>0 such that 𝑎_𝑘𝑏<min^𝑘2𝑑^𝑆_𝑘−𝜎²_1𝑘∑²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘,2𝑏_𝑘∑²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘−1+𝑏_𝑘𝜎²_2𝑘𝜖_𝑘+∑²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘,2𝑑^𝐼_𝑘+𝛾_𝑘−𝜎²_3𝑘𝜖_𝑘+𝛾_𝑘,2𝑑^𝑅_𝑘−𝜎²_4𝑘𝛾_𝑘. (4.18) Therefore, there exists 𝒞>0 such that 𝐿𝑉₀⩽−𝒞|(𝑢_𝑘,𝑣_𝑘,𝑤_𝑘,𝑧_𝑘)|².
Let us suppose that 𝑃{|𝑣_𝑘(𝑡)|<𝛿_𝑘}=1. Then 2²_𝑗=1𝛽_𝑘𝑗𝑢_𝑘𝑣_𝑘−𝑎_𝑘𝑢_𝑘^∞_𝑟=0𝑓_𝑗(𝑟)𝑣_𝑗(𝑡−𝑟)𝑑𝑟⩽²_𝑗=1𝛽_𝑘𝑗𝛿_𝑘𝑢²_𝑘1+2𝑎_𝑘+𝑣²_𝑘. (4.19) Therefore 𝐿𝑉⩽²_𝑘=1−𝑢²_𝑘𝑏_𝑘2𝑑^𝑆_𝑘−𝜎²_1𝑘−²_𝑗=1𝑎_𝑘𝛽_𝑘𝑗𝑆^∗_𝑘−²_𝑗=1𝛽_𝑘𝑗𝛿_𝑘1+2𝑎_𝑘−𝑣²_𝑘1+𝑏_𝑘2²_𝑗=1𝑎_𝑘𝛽_𝑘𝑗𝑆^∗_𝑘ℎ_𝑗−𝜎²_2𝑘−𝑎_𝑘𝜖_𝑘+²_𝑗=1𝛽_𝑘𝑗𝑆^∗_𝑘−²_𝑗=1𝛽_𝑘𝑗𝛿_𝑘−𝑎²_𝑘𝑤²_𝑘2𝑑^𝐼_𝑘+𝛾_𝑘−𝜎²_3𝑘−𝑎_𝑘𝜖_𝑘+𝛾_𝑘−𝑎⁴_𝑘𝑧²_𝑘2𝑑^𝑅_𝑘−𝜎²_4𝑘−𝑎_𝑘𝛾_𝑘. (4.20) Hence for sufficiently small 𝛿_𝑘>0, 𝐿𝑉(𝑥,𝑡) is negative definite in a sufficiently small neighborhood of 𝑥=0 for 𝑡⩾0. According to Lemma 4.2, we therefore conclude that the zero solution of (4.2) is stochastically asymptotically stable in the large. The proof is complete.

Acknowledgment

The first author was partially supported by the NSF of China no. 10926106.