Abstract

We introduce stochasticity into an SIS epidemic model with vaccination. The stochasticity in the model is a standard technique in stochastic population modeling. In the deterministic models, the basic reproduction number is a threshold which determines the persistence or extinction of the disease. When the perturbation and the disease-related death rate are small, we carry out a detailed analysis on the dynamical behavior of the stochastic model, also regarding of the value of . If , the solution of the model is oscillating around a steady state, which is the disease-free equilibrium of the corresponding deterministic model, whereas, if , there is a stationary distribution, which means that the disease will prevail. The results are illustrated by computer simulations.

1. Introduction

Epidemiology is the study of the spread of diseases with the objective to trace factors that are responsible for or contribute to their occurrence. Significant progress has been made in the theory and application of epidemiology modeling by mathematical research. Controlling infectious diseases has been an increasingly complex issue in recent years. Vaccination is an important strategy for the elimination of infectious diseases [13]. The vaccination enables the vaccinated to acquire a permanent or temporary immunity. When the immunity is temporary, the immunity can be lost after a period of time. This temporary immunity is used in many references [47] that assume the process of losing immunity is in the exponential form.

In [7], Li and Ma discussed an SIS model with vaccination. The system has following form: Here denotes the number of members of a population who are susceptible to an infection at time . denotes the number of members who are infective at time . denotes the number of members who are immune to an infection at time as the result of vaccination. The derivatives , , and are denoted by , , and , respectively.

The parameters in the model are summarized in the following list:: a constant input of new members into the population per unit time;: a fraction of vaccinated for newborns;: transmission coefficient between compartments and ;: natural death rate of , , and compartments;: the proportional coefficient of vaccinated for the susceptible;: recovery rate of infectious individuals;: the rate of losing their immunity for vaccinated individuals;: disease-caused death rate of infectious individuals.

All parameter values are assumed to be nonnegative and .

In Li and Ma [7], system (1) always has the disease-free equilibrium . If , then is the unique equilibrium of (1), and it is globally stable in the invariant set , where . If , then is unstable and there is an endemic equilibrium of (1), which is globally asymptotically stable under the sufficient condition in the invariant set . The basic reproduction number of system (1) is as follows:

Considering real life is full of randomness and stochasticity, using stochastic models can gain more real benefits. Consequently, there are many authors who studied stochastic biological systems and stochastic epidemic models; see [824]. There may be different approaches which result in different effects on the population system to include random perturbations in the models. In [2023], the situation of the parameter perturbation was considered. Gray et al. [23] discuss the following stochastic SIS model: Here and denote the susceptible and infected numbers of individuals, respectively. With , they simplify system (3) into a single equation and show that if and , then the disease will die out with probability one. If , then the disease will be persistent, where and is the threshold of the corresponding deterministic system. can be considered as the threshold of system (3), which is less than the value of .

On the other hand, white noise stochastic perturbations around the positive endemic equilibrium of epidemic models were considered in [19, 25]. Beretta et al. [25] considered the following system differential equation: where , , and are the positive points of equilibrium for the corresponding deterministic system and are constants, are independent from other standard Wiener processes . In their paper, they proved the above stochastic system was stable in probability by using the general method of Lyapunov functionals construction and obtained the stability condition immediately in terms of the parameters of the system under consideration.

Another different approach to include stochastic perturbations in a biological model was considered by Imhof and Walcher in [11]. They introduced and analyzed a variant of the deterministic single-substrate chemostat model. Also they found all species were persistent in the multiple species by using a comparison principle. Subsequently, they turned to modeling the influence of random fluctuations on the deterministic chemostat model. They set up the following stochastic chemostat model: where , , and are independent Brownian motions. By analyzing this stochastic differential equation model, they proved that the stochastic model led to extinction even though the deterministic counterpart predicts persistence.

In this paper, our approach to include stochastic perturbation is analogous to that of Imhof and Walcher [11]. We are devoted to studying the following stochastic system: where are independent Brownian motions and are their intensities.

This paper is organized as follows. In Section 2, we show there is a unique positive solution of system (6) by the way mentioned in [20, 23]. In Section 3, when , we derive that the solution of the system (6) oscillates around the disease-free equilibrium of system (1). In Section 4, when , although the solution of system (6) does not converge to the endemic proportion equilibrium of system (1), we conclude there is a unique stationary distribution for system (6), and it has ergodic property under some conditions. In Section 5, we make simulations to confirm our analytical results. Finally, in order to be self-contained, we have the appendix.

2. Existence and Uniqueness of Positive Solution

To investigate the dynamical behavior, the first concern is whether the solution has a global existence. Moreover, for a population dynamics model, whether the value is nonnegative is also considered. Hence in this section, we first show that the solution of system (6) is global and nonnegative.

Theorem 1. There is a unique solution of system (6) on for any initial value , and the solution will remain in with probability 1, namely,   for all almost surely.

Proof. Our proof is motivated by the works of Mao et al. [17]. Since the coefficients of (6) are locally Lipschitz continuous for any given initial value , there is a unique local solution on , where is the explosion time (see [20]). To show that this solution is global, we need to show that a.s. Let be sufficiently large so that , , and all lie within the interval . For each integer , define the stopping time: where throughout this paper, we set (as usual denotes the empty set). According to the definition, is increasing as . Set , whence a.s. If we can show that a.s., then and a.s. for all . In other words, to complete the proof, all we need to show is that a.s. If this statement is false, then there exist a pair of constants and such that Hence there is an integer such that Define a -function (see the appendix) by The nonnegativity of this function can be seen from , for all . Let and be arbitrary. Applying the Itô formula, we obtain where is defined by Let , s.t. . Substituting this into (12), we get The remainder of the proof follows that in Mao et al. [17].

3. Asymptotic Behavior around

It is clear that is the disease-free equilibrium of system (1). If , then is globally stable, 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. But there is no disease-free equilibrium in system (6). In this section, we show the average oscillation around in time to exhibit whether the disease will die out.

Theorem 2. If , , and the following condition is satisfied then the solution of system (6) with initial value has the property where is the disease-free equilibrium of system (1), , and is positive constant that satisfied .

Proof. Let , , and . Then Define a -function by where . By Itô's formula, we compute If , then Note that Therefore, we can obtain where is positive constant to be specified later and Young’s inequality is used. If we choose such that . This implies Then we have Integrating this from 0 to and taking the expectation, we have Hence, Consequently, This finishes the proof of Theorem 2.

Remark 3. From Theorem 2, we show the solution of system (6) will oscillate around the disease-free equilibrium of system (1) under some conditions and the disturbance intensity is proportional to the intensity of the white noise. In a biological view, as the intensity of stochastic perturbations is small, we consider the disease will die out.
Besides, if and , then is also the disease-free equilibrium of system (6). From the proof of Theorem 2, we get
Thus, the solution of system (6) is stochastically asymptotically stable in the large if .

4. Ergodicity of System (6)

In studying epidemic dynamical system, we are also interested in when the disease will prevail and persist in a population. In the deterministic models, the problem is solved by showing that the endemic equilibrium of corresponding model is a global attractor or is globally asymptotically stable. But there is no endemic equilibrium in system (6). In this section, based on the theory of Has’minskiĭ [26] (see the appendix), we show that there is a stationary distribution subjected to some conditions on and the parameters of the model which reveals that the disease will prevail also.

System (6) can be written as a form of system (see the appendix).

Consider and the diffusion matrix is .

Remark 4. Theorem 1 shows that there exists a unique positive solution of system (6). Besides, from the proof of Theorem 1, we obtain . Define , then , and it is clear that as , where . Hence, by Remark 2 of Theorem 4.1 of Has’minskiĭ (1980, page 86) [26], we obtain that the solution is a homogeneous Markov process in .

Theorem 5. If , , and such that then, for any initial value , there is a stationary distribution for system (6) and it has an ergodic property, where , is the endemic proportion equilibrium of system (1), and and are the same constants as in Theorem 2.

Proof. If , there is an endemic proportion equilibrium of system (1); then
Define a -function by where . The nonnegativity of the function can be observed from on . Let be the generating operator of system (6). Then we get Then Note that Then where Young’s inequality is used and If , choose such that Note that Thus the ellipsoid lies entirely in . We can take to be a neighborhood of the ellipsoid with , so, for , ( is a positive constant), which implies that condition () in Lemma A.2 is satisfied. Hence, the solution is recurrent in the domain , which, together with Lemma A.4 and Remark 4, implies that is recurrent in any bounded domain . Besides, for all , there is such that for all , , which implies that condition () is also satisfied. Therefore, the stochastic system (6) has a stationary distribution , and it is ergodic.

5. Numerical Simulations

In order to confirm the results above, we numerically simulate the solution of system (6) with given initial value and the parameters. Using Milstein’s higher order method mentioned in [27], we get the discretization equation of system (6): where , , and , , are the independent Gaussian random variables . Choosing suitable parameters in the system, by Matlab we get the simulation figures with initial value and time step .

Example 6. Throughout the paper we will assume that the unit of time is one day, and the population sizes are measured in units of 1 million. Choose the parameters in system (6) as follows: , , , , , , , , , , and . Note that and , , and satisfy condition (14), then by Theorem 2, we show the solution of system (6) will oscillate around the disease-free equilibrium of system (1) in time which is globally asymptotically stable. Besides Theorem 2 tells us that the difference between the perturbed solution and is only related with white noises and . Numerical simulations are shown in Figure 1.

Example 7. In Figures 2 and 3, we choose the parameters in system (6) as follows: , , , , , , , , , , and . Note that and ; then the equilibrium of system (1) is globally asymptotically stable. But the white noise may make system (1) appear as different phenomena. In detail, the conditions (30) and (31) are also satisfied. Therefore, as Theorem 5 stated, there is a stationary distribution (see the histograms in the middle in Figure 2). In addition, the left pictures in Figure 2 show that the solution of system (6) has fluctuating in a small neighborhood. Moreover, from Figure 3, we find that 95% or more of the population distribution lies within a neighborhood, which can be imagined a circular or elliptic region centered at (see the red point in Figure 3). All of these imply system (6) has stochastic stability.
Testing these data for normality, all tests used were highly significant, conclusively rejecting normality in all cases. This is not surprising in view of the very large sample sizes (10,000), as even moderate deviations from the tested distribution will be significant; however the normal quantile-quantile plots (see the plots in the right in Figure 2) suggest that these data are not far from being normally distributed for smaller values of ().

6. Conclusions

As most real world problems are not deterministic, incorporating stochastic effects into the model give us a more realistic way of modeling epidemic models. In this paper, we have considered stochastic SIS epidemic models with vaccination. We first proved the positivity of the solutions. When the perturbation and the disease-related death rate are small, we illustrated the dynamical behavior of the SDE model according to or . When , we proved that the asymptotic behavior around the disease-free equilibrium of the deterministic model. When , we also proved that the SDE model has the ergodic property as the fluctuation is small, where the positive solution converges weakly to the unique stationary distribution. Simulations are also carried out to verify our analytical results.

Appendix

Here, we give some basic theory in stochastic differential equations (see [18]).

Throughout this paper, unless otherwise specified, let be a complete probability space with a filtration satisfying the usual conditions (i.e., it is right continuous, and contains all P-null sets). Denote In general, consider -dimensional stochastic differential equation: with initial value . denotes -dimensional standard Brownian motion defined on the above probability space. Define the differential operator associated with (A.2) by If acts on a function , then where , , and . By Itô's formula, if , then

Consider (A.2) assume and for all . So is a solution of (A.2) called the trivial solution or equilibrium position.

Lemma A.1 (strong law of large numbers). Let  be a real-value continuous local martingale vanishing at . Then and also

Next, there is some theory about stationary distributions (see Has’minskiĭ [26], 1980).

Let be a homogeneous Markov process in ( denotes Euclidean -space) described by The diffusion matrix is , .

Assumption B. There exists a bounded domain with regular boundary , having the following properties.() In the domain and some neighborhood thereof, the smallest eigenvalue of the diffusion matrix is bounded away from zero.() If , the mean time at which a path issuing from reaches the set is finite and for every compact subset .

Lemma A.2 (see [26]). If (B) holds, then the Markov process has a stationary distribution . Let be a function integrable with respect to the measure . Then

Remark A.3. The proof is given in Has’minskiĭ [26]. The existence of a stationary distribution with a density is given in Theorem 4.1, page 119, and Lemma 9.4, page 138. The weak convergence and the ergodicity are obtained in Theorem 5.1, page 121, and Theorem 7.1, page 130.
To validate (), it suffices to prove that is uniformly elliptical in any bounded domain , where ; that is, there is a positive number such that , , (see chapter 3, page 103 of Gard [28] and Rayleigh's principle in Strang [29], chapter 6, page 349). To verify (), it is sufficient to show that there exist some neighborhood and a nonnegative -function such that, for any , is negative (for details refer to Zhu and Yin [30], page 1163).

Lemma A.4. Let be a regular temporally homogeneous Markov process in . If is recurrent relative to some bounded domain , then it is recurrent relative to any nonempty domain in .

Acknowledgments

The work was supported by Program for Changjiang Scholars and Innovative Research Team in University, NSFC of China (no. 10971021), the Ph.D. Programs Foundation of Ministry of China (no. 200918), and the Scientific and Technological Research Project of Jilin Province's Education Department (no. 2012244).