Nonlinear Pulse Vaccination in an SIR Epidemic Model with Resource Limitation

Abstract

Mathematical models can assist in the design and understanding of vaccination strategies when resources are limited. Here we propose and analyse an SIR epidemic model with a nonlinear pulse vaccination to examine how a limited vaccine resource affects the transmission and control of infectious diseases, in particular emerging infectious diseases. The threshold condition for the stability of the disease free steady state is given. Latin Hypercube Sampling/Partial Rank Correlation Coefficient uncertainty and sensitivity analysis techniques were employed to determine the key factors which are most significantly related to the threshold value. Comparing this threshold value with that without resource limitation, our results indicate that if resources become limited pulse vaccination should be carried out more frequently than when sufficient resources are available to eradicate an infectious disease. Once the threshold value exceeds a critical level, both susceptible and infected populations can oscillate periodically. Furthermore, when the pulse vaccination period is chosen as a bifurcation parameter, the SIR model with nonlinear pulse vaccination reveals complex dynamics including period doubling, chaotic solutions, and coexistence of multiple attractors. The implications of our findings with respect to disease control are discussed.

1. Introduction

Epidemiology is the study of the spread of diseases with the objective of tracing factors that are responsible for or contribute to their occurrence and serves as the foundation and logic of interventions made in the interest of public health and preventive medicine. Mathematical models describing the population dynamics of infectious diseases have played an important role in better understanding epidemiological patterns and disease control for a long time. Various epidemic models have been proposed and explored extensively and considerable progress has been achieved in the studies of disease control and prevention (see [1–3] and the references therein).

Outbreaks of infectious diseases have not only caused the loss of billions of lives but have often also rapidly damaged social economic systems, bringing about much human misery. Consequently, the focus of our research has been on how to prevent and cure infectious diseases effectively. It is well known that one of the most important concerns in the analysis of epidemic logical models is the efficacy of vaccination programmes. This subject gained prominence as a result of highly successful application of vaccinations for the worldwide eradication of smallpox and the restriction of diseases such as poliomyelitis, measles, tetanus, diphtheria, pertussis, and tuberculosis. Vaccination is one of the most cost effective of health investments as it can prevent or ameliorate morbidity due to infections. The effectiveness of vaccination has been widely studied and verified for the influenza vaccine [4], the HPV vaccine [5], the chicken pox vaccine [6], and others.

In particular, Agur et al. [7] first proposed a pulse vaccination strategy (PVS), which consists of periodical repetitions of impulsive vaccinations of all the age cohorts in a population, which has been confirmed as an important and effective strategy for the elimination of infectious diseases. At each vaccination time a constant fraction 𝑝 of susceptible people is vaccinated. This kind of vaccination is called impulsive since all the vaccine doses are applied within a very short time, with respect to the dynamics of the target disease. PVS allows the eradication of a disease with some practical advantages, as discussed in [7–9]. Recently, epidemic models with pulse vaccination have been the subject of intense theoretical analysis. For example, the dynamical behaviors of various epidemic models with PVS are studied in [8, 10–13]. The theoretical results show that PVS can be distinguished from the conventional strategies in leading to disease eradication at relatively low values of vaccination (see [7]).

Traditional epidemic models with PVS of population dynamics have assumed that the pulse vaccination proportion 𝑝 is constant, which implies that medical resources such as drugs, vaccines, hospital beds, and isolation units are sufficient for the infectious disease in question. However, in reality, every community or country has an appropriate or limited capacity for treatment, especially for emerging infectious diseases, and so understanding resource limitation is critical to effective management.

During the last decade, the outbreaks of SARS in 2003 and avian influenza among humans (H5N1 in 2003, H1N1 in 2009, and H7N9 in 2013) emphasized the need to enhance the capacity to fight emerging infectious diseases, which remains a challenge for public health. Several different vaccines have been developed specifically targeting avian influenza among humans, but limited production has exposed weaknesses [14–16]. A lack of sufficient vaccine at the outset of a new pandemic is to be expected, given the crucial lag between the timing of the emergence of a new strain and when a vaccine has been developed and is ready for distribution. Although most vaccination programmes aim to vaccinate as many susceptible individuals as possible, this is not always the case for emerging infectious diseases due to limited vaccine availability, particularly in developing countries.

However, resource limitation is difficult to quantify partly because it is a dynamic process (see [17]). Recently, continuous SIR models concerning limited medical resources have been developed and investigated (see [18, 19]). Chow et al. [18] explored the dynamics of an SIR epidemic model to understand how limited medical resources and their supply efficiency affect the transmission of infectious diseases. Zhou and Fan [19] studied a multigroup SIR epidemiological model to explore the effects of limited medical resources and group-targeted vaccination strategies on disease control and prevention. To the best of our knowledge, no work has been done for the effects of resource limitation on an SIR model with PVS. In order to investigate the effect of limited vaccine availability on the spread of infectious disease, a saturation phenomenon of limited medical resources is considered. That is, we will study the dynamic behavior of an SIR epidemic model with nonlinear pulse vaccination. Such mathematical models are suitable for simulating processes with short duration perturbations during their development (see [20]).

In this paper, the dynamical behavior of an SIR epidemic model with nonlinear pulse vaccination is proposed and analyzed. The main purpose is to address and understand how a limited vaccine resource affects the transmission and control of infectious diseases. The rest of this paper is organized as follows. In the next section, an SIR epidemic model with nonlinear PVS is introduced. In Section 3, we investigate the existence and stability of the disease-free periodic solution by using the method of differential inequality, qualitative analysis, a discrete dynamical system determined by the stroboscopic map, and comparison theorem. Meanwhile, Latin Hypercube Sampling (LHS)/Partial Rank Correlation Coefficients (PRCCs) uncertainty and sensitivity analysis techniques are employed to investigate the key control parameters which are most significantly related to threshold values (see [21–23]). Section 4 focuses on finding the sufficient condition under which the endemic periodic solution exists by using the bifurcation theorem. The paper ends with some interesting biological conclusions and numerical bifurcation analyses, which complement the theoretical findings.

2. Model Formulation

In the classical epidemiological model (see [24–26]), the total population 𝑁 is composed of three groups of individuals: susceptible (𝑆), infective (𝐼), and recovered (𝑅). Let 𝑆(𝑡) denote the number of members of a population susceptible to the disease at time 𝑡; 𝐼(𝑡) represents the number of infective members; and 𝑅(𝑡) denotes the number of members who have been removed from the population. This leads to the following SIR model: d𝑆(𝑡)(𝑡)𝐼(𝑡)(𝑡),(𝑡)d𝑡=𝑏−𝛽𝑆−𝑏𝑆d𝐼(𝑡)𝐼(𝑡)(𝑡)(𝑡),(𝑡)d𝑡=𝛽𝑆−𝛾𝐼−𝑏𝐼d𝑅(𝑡)(𝑡).d𝑡=𝛾𝐼−𝑏𝑅 (1) In model (1), we assume that the three classes of subpopulations have the same constant birth and death rates 𝑏 (i.e., the total population has a constant size), which can be normalized to unity; that is, 𝑁(𝑡)=𝑆(𝑡)+𝐼(𝑡)+𝑅(𝑡)=1, and the infectious individuals recover at a rate 𝛾>0, so that 1/𝛾 is the mean infectious period. Susceptibles become infected at a rate 𝛽𝐼, where 𝛽>0 is the transmission rate. In practice, the equation for d𝑅(𝑡)/d𝑡 is redundant because 𝑅(𝑡) can be obtained from 𝑁(𝑡)=1. A detailed description of model (1) and its dynamics may be found in [9, 27, 28].

Let 𝑇>0 be the time between two consecutive pulse vaccinations and let 𝑝(𝑡) (0≤𝑝(𝑡)<1) be the fraction of susceptible subjects who are inoculated with vaccine, which depends on the number of susceptibles in the population at time 𝑡. This yields the following model with PVS [8, 9]: d𝑆(𝑡)(𝑡)𝐼(𝑡)(𝑡),(𝑡)d𝑡=𝑏−𝛽𝑆−𝑏𝑆d𝐼(𝑡)𝐼(𝑡)(𝑡)(𝑡),𝑆𝑡d𝑡=𝛽𝑆−𝛾𝐼−𝑏𝐼𝑡≠𝑛𝑇,𝑛∈𝒩,₊=(𝑡)𝑆(𝑡),𝐼𝑡1−𝑝₊=𝐼(𝑡),𝑡=𝑛𝑇,𝑛∈𝒩, (2) here 𝒩={0,1,2,…}.

The formula for 𝑝(𝑡) plays a key role in investigating the effect of the limited capacity for treatment on the spread of infectious disease. To do this we assume, without loss of generality, that the PVS is implemented in a developing country but one in which vaccines are relatively abundant. Note that there are often major shortages of medical personnel in rural areas and in certain specialities such as surgery, paediatrics, and obstetrics [14]. Let vaccination pulses occur every 𝑇 time units. When a pulse day occurs, the medical personnel have to go out into the rural population and vaccinate as many susceptibles as they can find on that particular day. If there are very few susceptibles, then it will be hard to find them, so the vaccination coverage will be low. The more susceptibles there are, the easier it is to find them, so the coverage per pulse will be greater. But there is a limit to how fast the medical team can find and handle each susceptible, and this causes the vaccination rate 𝑝(𝑡) to saturate, that is, a saturation phenomenon of limited medical resources. By employing the Holling type II functional response for a predator-prey model, we can define 𝑝(𝑡) as follows: 𝑝(𝑡)=𝑝_max𝑆(𝑡)𝑆(𝑡)+𝜃,0≤𝑝_max<1, (3) where 𝑝_max denotes the maximal vaccination proportion and 𝜃 is the half saturation constant. Note that the vaccination coverage is less than 100% and in practice a vaccine efficacy is also less than 100%, which lead to the parameter 𝑝_max<1. Furthermore generalized formulae such as the Hill function [29] can be employed to characterize the saturation phenomenon of the limited medical resources.

Taking into account the resource limitation and saturation effects, for convenience and simplification, then model (2) can be rewritten as d𝑆(𝑡)(𝑡)𝐼(𝑡)(𝑡),(𝑡)d𝑡=𝑏−𝛽𝑆−𝑏𝑆d𝐼(𝑡)𝐼(𝑡)(𝑡)(𝑡),𝑆𝑡d𝑡=𝛽𝑆−𝛾𝐼−𝑏𝐼𝑡≠𝑛𝑇,𝑛∈𝒩,₊=1−𝑝_max𝑆²(𝑡)+𝜃𝑆(𝑡)𝑆(𝑡),𝐼𝑡+𝜃₊=𝐼(𝑡),𝑡=𝑛𝑇,𝑛∈𝒩, (4) where the parameters 𝑏,𝛽,𝛾, and 𝑝_max,𝜃 are defined in models (1) and (3), respectively. Note that model (4) is reduced to the classical SIR model with constant pulse vaccination when 𝜃=0, which was investigated by Agur et al. [7], Stone et al. [8], and Shulgin et al. [9]. The parameter definitions of model (4) are summarized in Table 1.

Table 1 

Parameter definitions, values and source for model (4).

3. Existence of the Disease-Free Periodic Solution

We first demonstrate the existence of a disease-free periodic solution of model (4), in which infectious individuals are entirely absent from the population permanently; that is, 𝐼(𝑡)=0 for a sufficiently long time. Under some threshold conditions, we show below that the susceptible population oscillates with period 𝑇 in synchronization with the periodic pulse vaccination and the infectious class will die out eventually. Thus, in this section, we focus on determining the global attraction of this disease-free periodic solution. To address this, we first consider the following submodel over the time interval 𝑛𝑇<𝑡≤(𝑛+1)𝑇: d𝑆(𝑡)(𝑡)𝑆𝑡d𝑡=𝑏−𝑏𝑆,𝑡≠𝑛𝑇,₊=1−𝑝_max𝑆²(𝑡)+𝜃𝑆(𝑡)𝑆(𝑡)𝑆0+𝜃,𝑡=𝑛𝑇,₊=𝑆₀. (5) Integrating the first equation of model (5) between pulses yields 𝑆(𝑡)=𝑆𝑛𝑇⁺−1e^{−𝑏(𝑡−𝑛𝑇)}+1 (6) for 𝑛𝑇<𝑡≤(𝑛+1)𝑇. It follows from the second equation of model (5) that 𝑆(𝑛+1)𝑇⁺=1−𝑝_max𝑆²((𝑛+1)𝑇)+𝜃𝑆((𝑛+1)𝑇)𝑆)𝑇).((𝑛+1+𝜃 (7) Substitution of (6) into (7) gives 𝑆(𝑛+1)𝑇⁺=1−𝑝_max𝑆𝑛𝑇⁺−1e^−𝑏𝑇+1²+𝜃𝑆𝑛𝑇⁺−1e^−𝑏𝑇×𝑆+1𝑛𝑇⁺−1e^−𝑏𝑇+1+𝜃⁻¹. (8) Denote 𝑆(𝑛𝑇⁺)=𝑆_𝑛 then the previous equation can be rewritten as the following difference equation: 𝑆_𝑛+1=1−𝑝_max𝑆_𝑛−1e^−𝑏𝑇+1²+𝜃𝑆_𝑛−1e^−𝑏𝑇+1𝑆_𝑛−1e^−𝑏𝑇𝑆+1+𝜃≐𝑓_𝑛, (9) which is the so-called stroboscopic map of model (5), and it describes the relations of the number of susceptibles in the population between any two successive pulse vaccinations. Consequently, the existence of the positive steady state of model (9) implies the existence of a positive periodic solution of model (5). Therefore, we first discuss the conditions which guarantee the existence of a positive steady state of (9). Taking the derivative of 𝑓(𝑆_𝑛) with respect to 𝑆_𝑛 yields 𝑓^′𝑆_𝑛=1−𝑝_max+𝑝_max𝜃²𝑆_𝑛−1e^−𝑏𝑇+1+𝜃²e^−𝑏𝑇, (10) and it is easy to see that 0<𝑓^′(𝑆_𝑛)<1 holds true.

For convenience, we denote the positive fixed point of the stroboscopic map (9) by 𝑆, which satisfies the following equation: 𝐴₁𝑆²+𝐵₁𝑆+𝐶₁=0, (11) where 𝐴₁=e^−𝑏𝑇1−1−𝑝_maxe^−𝑏𝑇𝐵>0,₁=e^−𝑏𝑇−121−𝑝_maxe^−𝑏𝑇−(,𝐶𝜃+1)₁=−1−𝑝_max1−e^−𝑏𝑇2+𝜃1−e^−𝑏𝑇<0. (12) Obviously, (11) has a unique positive root; that is, 𝑆=−𝐵₁+𝐵²₁−4𝐴₁𝐶₁2𝐴₁, (13) which is stable due to 0<𝑓^′(𝑆)<1.

According to the above discussion on the relations between a fixed point of the stroboscopic map (9) and the periodic solution of model (5), we conclude that the submodel (5) has a unique nontrivial positive periodic solution denoted by 𝑆^∗(𝑡), and 𝑆^∗(𝑡)=𝑆−1e^{−𝑏(𝑡−𝑛𝑇)}+1,for𝑛𝑇<𝑡≤(𝑛+1)𝑇 (14) and it follows from 0<𝑓^′(𝑆_𝑛)<1 and [30] that we have the following lemma.

Lemma 1. Model (5) has a positive periodic solution 𝑆^∗(𝑡), and for any solution 𝑆(𝑡) of (5), we have |𝑆(𝑡)−𝑆^∗(𝑡)|→0 as 𝑡→+∞.

By using Lemma 1, we have thus shown that the unique fixed point 𝑆 of (9) is globally stable.

In the following, we will present the sufficient condition for the global attraction of the disease-free periodic solution (𝑆^∗(𝑡),0) of model (4).

Theorem 2. The disease-free periodic solution (𝑆^∗(𝑡),0) of model (4) is globally asymptotically stable in the first quadrant provided that 𝑅¹₀<1, where 𝑅¹₀=𝛽𝑆−11−e^−𝑏𝑇+𝑏𝑇,𝑏𝑇𝛾+𝑏 (15) here 𝑆 is defined in (13).

The proof of Theorem 2 is given in Appendix A.

Note that although the threshold value 𝑅¹₀ depends on all parameters of model (4), the most interesting parameters here are pulse vaccination period 𝑇, maximum proportion of pulse vaccination 𝑝_max, and the parameter related to resource limitation, that is, 𝜃. In order to address how those factors (𝑇 and 𝜃) affect the threshold value 𝑅¹₀, we first note that lim_𝑇→0𝑅¹₀=lim_𝑇→0𝛽𝑆−11−e^−𝑏𝑇+𝑏𝑇𝑏𝑇𝛾+𝑏=0, (16) which indicates the importance of the frequency of vaccination in eradicating an infectious disease. Then letting 𝜃 and 𝑇 vary and fixing all other parameters as those in Figure 1, we see that the threshold value 𝑅¹₀ is a monotonically increasing function of 𝑇 for fixed 𝜃, and it is a monotonically increasing function of 𝜃 for fixed period 𝑇. Moreover, there are some critical values of 𝑇 and 𝜃 such that 𝑅¹₀=1, which indicate that there exists one maximum allowable period of pulse vaccination 𝑇^𝜃_max such that 𝑅¹₀<1 for all 𝑇<𝑇^𝜃_max. Moreover, the larger 𝜃 the smaller the value that 𝑇^𝜃_max has. Unfortunately, we cannot get its analytical expression due to the complexity of 𝑅¹₀. So we turn to find some approximations for 𝑇^𝜃_max in the following.

Figure 1 

The effects of resource limitation parameter 𝜃 and pulse vaccination period 𝑇 on the threshold value 𝑅10 with 𝑝max=0.75. All other parameter values are given in Table 1.

If the PVS is applied frequently enough such that d𝐼(𝑡)/d𝑡≤0 for all 𝑡>0 in model (4), then the number of infectious individuals is a decreasing function of time. It is possible to satisfy this condition if pulsing ensures that 𝑆(𝑡)≤𝑆_𝑐≐(𝑏+𝛾)/𝛽 for all 𝑡; that is, pulse vaccination is applied every time once the 𝑆(𝑡) approaches the threshold value 𝑆_𝑐 (see [7, 27] and the references therein).

It follows from the first equation of model (4) that d𝑆(𝑡)(𝑡).d𝑡≤𝑏−𝑏𝑆 (17) Consider the following comparison equation: d𝑥(𝑡)(𝑡)(0)=d𝑡=𝑏−𝑏𝑥,𝑥1−𝑝_max𝑆²_𝑐+𝜃𝑆_𝑐𝑆_𝑐,+𝜃 (18) which gives 𝑥(𝑡)=1+1−𝑝_max𝑆²_𝑐+𝜃𝑆_𝑐𝑆_𝑐e+𝜃−1^−𝑏𝑡. (19) Hence 𝑆(𝑡)≤𝑥(𝑡)=1+1−𝑝_max𝑆²_𝑐+𝜃𝑆_𝑐𝑆_𝑐e+𝜃−1^−𝑏𝑡≤𝑆_𝑐, (20) and solving the above inequality with respect to 𝑡 yields 10≤𝑡≤𝑏𝑝ln1+^max𝑆²_𝑐𝑆_𝑐+𝜃1−𝑆_𝑐, (21) which indicates that 𝑇^𝜃_max=1𝑏𝑝ln1+^max𝑆²_𝑐𝑆_𝑐+𝜃1−𝑆_𝑐. (22) Thus, if the nonlinear pulse vaccination strategy is applied periodically with a period less than 𝑇^𝜃_max, then the number of infectious in the population will decrease forever and eventually die out.

We also note that if 𝜃=0, the maximum allowable period 𝑇^𝜃_max becomes 𝑇⁰_max=(1/𝑏)ln(1+((𝑝_max𝑆_𝑐)/(1−𝑆_𝑐))), and this is the case studied in [7, 9] for a constant proportion pulse vaccination strategy. Obviously, 𝑇^𝜃_max is a monotonically decreasing function with respect to 𝜃, which gives that 𝑇^𝜃_max<𝑇⁰_max. This theoretical result indicates that if we aim to eradicate an infectious disease, then it is necessary to carry out pulse vaccination programs more frequently under resource limitation than when necessary resources are available. When various nonlinear factors including 𝑝_max,𝜃 are seriously considered, our results show that it is getting more and more difficult to eliminate the endemic diseases. It also confirms that getting vaccinated is by far the most effective action that susceptibles can take to protect themselves and their family from infectious disease. Other options available to health practitioners in response to limited medical resources include the release of prescriptive health education, training vaccinators, and promoting the use of vaccine-delivery patch with dissolving microneedles (see [14–16] for details).

A typical solution of the SIR model with a nonlinear pulse vaccination strategy is shown in Figure 2, where we observe how the variable 𝑆(𝑡) oscillates in a stable cycle (as shown in Figure 2(a)). In contrast, the infection 𝐼(𝑡) rapidly decreases to zero (Figure 2(b)) if we set 𝜃=0. That is, without resource limitation the infectious disease is eventually eradicated (Figure 2(b)). However, if there is resource limitation, that is, the parameter 𝜃 is larger than zero, then both susceptible and infected populations oscillate periodically, as shown in Figures 2(c) and 2(d) with 𝜃=0.2.

(a)

(b)

(c)

(d)

Figure 2 

The effects of resource limitation on dynamics of model (4) with 𝑝max=0.5, 𝑇=2, and all other parameter values are given in Table 1. (a and b) 𝜃=0; (c and d) 𝜃=0.2.

Therefore, if we have enough vaccine and vaccinators to combat the infectious disease, then it is reasonable that the susceptible population vaccinated each time is assumed to be proportional to the number of susceptible individuals (here 𝑝_max). However, this approximation cannot reflect the real case if we aim to mitigate the emerging infectious disease. Usually, the number of susceptible individuals needing to be vaccinated may exceed the capacity of local medical conditions due to shortages of vaccine and doctors, especially in rural areas in many developing countries, where reaching all of the target population may be difficult. For instance in 2010 in Lesotho, diphtheria, tetanus toxoid and, pertussis (DTP3) immunization coverage among 1-year-olds was only 83% (but an improvement on 78% in 2001, see [31]) and both the measles and tuberculosis (BCG) vaccination programs were hampered by supplies of vaccine running out in some parts of the country and then allowing vaccinations only for children born in hospitals and not for those born elsewhere in the communities. Thus, it is important to design pulse vaccination campaigns carefully. For example, how frequently should the vaccination strategy be applied? And how can the size of the susceptible population be evaluated? Those are key issues for evaluating vaccination coverages and the distribution of vaccines. Therefore, with the help of mathematical models the vaccination strategy can be properly designed.

It is well known that in order to eradicate an infectious disease, the vaccination strategy must be designed such that the threshold value 𝑅¹₀ is less than 1, which is equivalent to 𝑅²₀<0, where 𝑅²₀=𝛽𝑆−11−e^−𝑏𝑇𝑏+𝛽−𝛾−𝑏𝑇. (23) In the following we employ the formula 𝑅²₀ to investigate the important factors which affect the threshold values most significantly by using uncertainty and sensitivity analysis (see [21, 22, 32]). Sensitivity analysis of the most significant parameters including the birth or death rate (𝑏), the contact rate (𝛽), the maximum proportion vaccinated (𝑝_max), and the pulse vaccination period (𝑇) was performed by evaluating the PRCCs for various input parameters against the threshold condition 𝑅²₀ with the LHS method [21, 22, 32]. By using LHS with 3,000 samples, uncertainty and sensitivity analyses for all parameters in model (4) were determined. A uniform distribution function was used and tested for significant PRCCs for all parameters with wide ranges, such as 𝑏∼𝑈(0.001,0.1), 𝛽∼𝑈(1000,2500), 𝑇∼𝑈(0.1,20), and 𝜃∼𝑈(0,1), and the baseline values of all parameters are given in the figure legend of Figure 3. Note that for comparative purposes, we chose the same baseline parameter values as those used in literatures [8, 9, 33]. Unfortunately, we do not have the realistic ranges of parameter values and we chose those parameter ranges just for illustrative purposes only and then provided a scientific basis for disease control.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

Figure 3 

PRCC results and PRCC scatter plots of the baseline parameters as those shown in Table 1. Here we fix 𝑝max=0.75,𝜃=0.2,𝑇=5. (a) PRCC results. (b–g) PRCC scatter plots. The title of each plot represents the PRCC value with the corresponding P value.

Figure 3(a) shows the PRCC results which illustrate the dependence of 𝑅²₀ on each parameter, and PRCC scatter plots for each parameter are given in Figures 3(b)–3(g), respectively. We considered absolute values of PRCC values which belong to (0.4,1), (0.2,0.4), and (0,0.2) as indicating very important, moderate, or correlations of little significance between input parameters and output variables, respectively. The positive signs of their PRCCs indicate that if the parameters are increased, the value of 𝑅²₀ increases (vice versa). The negative sign suggests that if increased, the value of 𝑅²₀ decreases (vice versa).

Therefore, the parameters 𝑏, 𝛽, 𝜃, and 𝑇 are responsible for increasing the values of 𝑅²₀, so decreasing all those four parameters can reduce the threshold value and consequently are beneficial for infectious disease control. As the parameter maximum proportion of pulse vaccination 𝑝_max changes, it is negatively and significantly correlated with large PRCCs. Thus, increasing this parameter will result in significantly reducing the threshold value. The most significant control parameters are 𝑏, 𝛽, 𝑝_max, and 𝑇 (as shown in Figures 3(b), 3(c), 3(e), and 3(g)). It is interesting to note that the period, 𝑇, and the maximum proportion vaccinated, 𝑝_max, are strongly correlated with the threshold value 𝑅²₀, which indicates that carefully designing a pulse vaccination strategy is crucial for infectious disease control.

4. Existence of Endemic Periodic Solution and Complex Dynamics

An important issue is how the infectious disease develops if the threshold condition 𝑅¹₀ is larger than one. In particular, what types of dynamic behavior may model (4) have once the period of pulse vaccination exceeds certain threshold levels. Therefore, in this section we first focus on the existence of an endemic periodic solution, that is, to investigate under what conditions will the infected population oscillate periodically with small amplitude, and then we choose the period of pulse vaccination as a bifurcation parameter to investigate the complex patterns of the infectious disease dynamics.

Next, we analyze the existence of an endemic periodic solution of model (4) near the disease-free periodic solution (𝑆^∗(𝑡),0) using bifurcation theory [34] and set the impulsive period 𝑇 as bifurcation parameter. For ease of narrative, we rewrite the model (4) as d𝑆(𝑡)(𝑡)𝐼(𝑡)(𝑡)d𝑡=𝑏−𝛽𝑆−𝑏𝑆≐𝐹₁(𝑆(𝑡)(𝑡,(𝑡),𝐼))d𝐼(𝑡)𝐼(𝑡)(𝑡)(𝑡)d𝑡=𝛽𝑆−𝛾𝐼−𝑏𝐼≐𝐹₂(𝑆(𝑡)(𝑡,𝑆𝑡,𝐼))𝑡≠𝑛𝑇,𝑛∈𝒩,₊=1−𝑝_max𝑆²(𝑡)+𝜃𝑆(𝑡)𝑆(𝑡)+𝜃≐Θ₁(𝑆(𝑡)(𝑡,𝐼𝑡,𝐼))₊=𝐼(𝑡)≐Θ₂(𝑆(𝑡),𝐼(𝑡)),𝑡=𝑛𝑇,𝑛∈𝒩. (24)

According to the bifurcation theorem in [34], we can obtain the following result.

Theorem 3. If there exists a 𝑇₀>0 satisfying 𝑅¹₀(𝑇₀)=1 and all the parameters in model (24) satisfy 𝛽𝑆𝑇₀−1e^−𝑏𝑇0𝑏𝑇₀𝐴₃1−𝐴₃e^−𝑏𝑇0−1<𝛽−𝑟−𝑏, (25) where 𝐴₃=1−𝑝_max+((𝜃²𝑝_max)/(𝜃+𝑆^∗(𝑇₀))²), then there exists a supercritical branch at the point 𝑇₀.

The existence of 𝑇₀ has been indicated in Figure 1 and the proof of Theorem 3 is given in Appendix B. Theorem 3 shows that there exists an endemic periodic solution under some conditions, providing that the disease-free periodic solution becomes unstable. It follows from Theorem 3 that if 𝑇>𝑇₀ and is close to 𝑇₀, the endemic periodic solution of model (4) is stable.

To investigate the complex dynamics that model (4) can have, we chose the pulse period 𝑇 as a bifurcation parameter and fixed all other parameters as those in Figure 4 for two different 𝜃 values. Figure 4(a) is a bifurcation diagram without resource limitation; that is, 𝜃=0, which was obtained by Shulgin et al. [9]. Figure 4(b) is a bifurcation diagram with resource limitation; that is, 𝜃=0.2. Comparing those two bifurcation diagrams we conclude that the nonlinear pulse vaccination can produce more complex dynamics than those for the model with linear pulse vaccination. Once the pulse period 𝑇 exceeds 𝑇^𝜃_max defined in (22), both susceptible and infected populations can oscillate periodically with a large amplitude that corresponds to periodic outbreaks of epidemics. As the pulse period is further increased, a series of complex and interesting bifurcation phenomena are observed (see Figures 4(a) and 4(b)). Figures 4(b)–4(d) indicate that the dynamical behavior of model (4) is very complex, including period doubling bifurcation, chaotic attractors, multistability, periodic-adding, chaos crisis, and periodic windows. Meanwhile, bifurcation analyses (i.e., Figures 4(b)–4(d)) also indicate that the model (4) has several different attractors which can coexist for a wide range of parameters. For example, Figure 5 provides an example of two attractor coexistence when 𝑇=7.8. Multiple attractor coexistence indicates that the infectious outbreak patterns depend on the initial values, which may cause difficulties for infectious disease control.

(a)

(b)

(c)

(d)

Figure 4 

Bifurcation diagrams of model (4) with respect to bifurcation parameter 𝑇. Here 𝑝max=0.75 and all other parameter values are given in Table 1. (a) 𝜃=0 which are reproduced from [9]; (b) 𝜃=0.2; (c and d) a magnified part of Figure 4(b).

(a)

(b)

(c)

(d)

Figure 5 

Coexistence of two stable attractors of model (4). Here 𝑝max=0.75,  𝜃=0.2,  𝑇=7.8, and all other parameter values are given in Table 1. (a and b) The initial value (𝑆0,𝐼0)=(0.12,0.05); (c and d) the initial value (𝑆0,𝐼0)=(0.1,0.02).

5. Discussion and Biological Conclusions

In order to understand the effect of resource limitation, in particular lack of vaccine, on the transmission of infectious disease, we deliberately investigated the dynamical behavior of an SIR epidemic model which incorporates a nonlinear pulse vaccination strategy. To this end, we introduced a nonlinear form 𝑝(𝑡)=(𝑝_max𝑆(𝑡))/(𝑆(𝑡)+𝜃) as vaccination proportion. This resulted in interesting and dramatic changes in the dynamical behavior of solutions and they became more and more complicated, which means that it is a very difficult task to control infectious diseases, in particular emerging infectious diseases, under resource limitation. We discussed the control strategy based on the threshold value in more detail through theoretical analysis, numerical studies, and sensitivity analysis. As pointed out by Agur et al. [7], the costs, risks, and the efficacy of the pulse vaccinations are key factors in eliminating transmission of the measles virus. In order to depict those factors, the nonlinear resource limitations should be taken into account [9], and the main purpose of the present work is to formulate the disease model to address the effects of resource limitations on disease control.

The results indicate that vaccination as an approach to controlling epidemics must be committed to a long-term strategy, especially when resources are limited. To understand the interactions amongst 𝑆(𝑡) and 𝐼(𝑡) which are important for successful control of the spread of infectious disease, we have proposed detailed modeling methods involving vaccination on 𝑆(𝑡) and determining the most significant parameters for the basic reproductive number 𝑅²₀ by using LHS/PRCC uncertainty and sensitivity analysis techniques (see Figure 3). Our results indicate that we should pay more attention to the vaccination period and vaccination proportion which reduce the threshold value 𝑅²₀ or 𝑅¹₀ to prevent the outbreak of disease. Our results clarify that it is more significant to improve immunization programs (i.e., vaccinations as mentioned in [15, 16]) under conditions of resource limitation.

The results also indicate that the dynamic behavior of model (4) may be dramatically affected by small changes in the value of initial densities of susceptible and infected with resource limitation. Bifurcation diagrams shown in Figure 4 indicate that there are many forms of complexities in model (4), which are related to chaotic bands with periodic windows and attractor crises, and Figure 5 reveals that the occurrence of multiple attractors is a common property of the SIR model with resource limitation, which can help us to further understand the application of a nonlinear pulse vaccination strategy in an SIR model (for more information about bifurcation diagrams of impulsive control strategies, see [35, 36]). Some complexities are related to the long-term behavior of population dynamics, characterized either by well behaving relatively simple or very complicated strange attractors. We also found that the routes to chaos are very complicated. For instance, with resource limitation there are several hidden factors that can adversely affect the control strategy. The increasing number of potential complexities predicted by the theory does not seem to make this task any easier. Nevertheless, identifying complicated, possibly chaotic, dynamics in population data has remained a major challenge in epidemiology studies. In particular, the suggestion that the dynamics of measles are chaotic has been much debated (e.g., see [37]) and, when seeking to establish whether measles dynamics in New York City were chaotic or not, Sugihara and May [38] omitted data from after 1963 when vaccinations were introduced, as the immunisations had altered the intrinsic dynamics. However, our results imply that even with immunisations, disease dynamics could be chaotic under some circumstances and so a reanalysis of the New York data including the vaccination period might be worthwhile.

In the present paper, we have focused on studying the effects of dynamical behaviors on the SIR model with nonlinear pulse vaccination. For our simple model, we have assumed that the maximal vaccination proportion 𝑝_max is a constant and the pulse period 𝑇 can be varied. However, due to limited vaccine availability, the maximal vaccination proportion 𝑝_max should depend on the vaccine obtained in each pulse period. Thus, to describe the impact of a limited vaccine stockpile, a separate differential equation for the amount of vaccine, denoted by 𝑉(𝑡), could be included in the model, and then 𝑝(𝑡) should be a function of 𝑉(𝑡) and 𝑆(𝑡). Meanwhile, it is essential to link the costs of developing and implementing controls to population dynamic modeling of disease epidemics in order to consider other resource limitations such as quantities of drugs, availability of isolation units, numbers of hospital beds, and equipment. These topics will be considered in further work in the future.

Appendices

A. Proof of Theorem 2

Firstly, we prove the local stability of a periodic solution (𝑆^∗(𝑡),0), which may be determined by considering the behavior of small amplitude perturbations of the solution.

Define 𝑆(𝑡)=𝑈(𝑡)+𝑆^∗(𝑡),𝐼(𝑡)=𝑉(𝑡), (A.1) here 𝑈(𝑡),𝑉(𝑡) are small perturbations, which may be written as 𝑈(𝑡)𝑉(𝑡)(𝑡)𝑈(0)𝑉(0)=Φ,0≤𝑡<𝑇, (A.2) where the fundamental matrix Φ(𝑡) satisfies dΦ(𝑡)=d𝑡−𝑏−𝛽𝑆^∗(𝑡)0𝛽𝑆^∗(𝑡)Φ(𝑡)−𝛾−𝑏, (A.3) and Φ(0)=𝐼 is the identity matrix. So Φ(𝑇)=e^−𝑏𝑇★0e^{∫𝑇0(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}, (A.4) where ★ is not necessarily computed in detail as it is not required in the following analysis.

The resetting impulsive conditions of model (4) become 𝑈𝑛𝑇⁺𝑉𝑛𝑇⁺=1−𝑝_max+𝜃²𝑝_max(𝜃+𝑆^∗(𝑛𝑇))²0𝑈()𝑉()()𝑈()𝑉().01𝑛𝑇𝑛𝑇≐𝐵𝑛𝑇𝑛𝑇𝑛𝑇 (A.5) Hence, according to Floquet theory, if the moduli of both eigenvalues of the matrix (𝑇)Φ(𝑇)=𝑀=𝐵1−𝑝_max+𝜃²𝑝_max(𝜃+𝑆^∗(𝑛𝑇))²0×e01_−𝑏𝑇★0e^{∫𝑇0(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡} (A.6) are less than one, then the periodic solution (𝑆^∗(𝑡),0) is locally stable. In fact, two Floquet expressions multiplied are thus 𝜆₁=1−𝑝_max+𝜃²𝑝_max(𝜃+𝑆^∗(𝑇))²e^−𝑏𝑇,𝜆₂=e^{∫𝑇0(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}. (A.7) It is obvious that 0<𝜆₁<1. Thus, the stability of (𝑆^∗(𝑡),0) is decided by whether 𝜆₂<1 or not. So we conclude that the disease-free periodic solution (𝑆^∗(𝑡),0) is locally stable if ∫^𝑇₀(𝛽𝑆^∗(𝑡)−𝛾−𝑏)d𝑡<0; that is, 1𝑇^𝑇₀𝑆^∗(𝑡)d𝑡=𝑆−11−e^−𝑏𝑇+𝑏𝑇<𝑏𝑇𝛾+𝑏𝛽, (A.8) and by simple calculation we have 𝛽𝑆−11−e^−𝑏𝑇+𝑏𝑇𝑏𝑇𝛾+𝑏<1, (A.9) which holds true due to 𝑅¹₀<1. Therefore, we can draw a conclusion from the above analysis that if 𝑅¹₀<1, the disease-free periodic solution (𝑆^∗(𝑡),0) of model (4) is locally asymptotically stable.

Next, we will prove the global attraction of the disease-free periodic solution (𝑆^∗(𝑡),0) of model (4). Since 𝑅¹₀<1, we can choose 𝜀₁>0 sufficiently small such that 1𝑇^𝑇₀𝑆^∗(𝑡)+𝜀₁d𝑡<𝛾+𝑏𝛽. (A.10) It follows from model (4) that d𝑆(𝑡)(𝑡)𝑆𝑡d𝑡≤𝑏−𝑏𝑆,𝑡≠𝑛𝑇,₊=1−𝑝_max𝑆²(𝑡)+𝜃𝑆(𝑡)𝑆(𝑡)+𝜃,𝑡=𝑛𝑇. (A.11) Consider the comparison equation d𝑍(𝑡)(𝑡)𝑍𝑡d𝑡=𝑏−𝑏𝑍,𝑡≠𝑛𝑇,₊=1−𝑝_max𝑍²(𝑡)+𝜃𝑍(𝑡)𝑍(𝑡)+𝜃,𝑡=𝑛𝑇. (A.12)

From Lemma 1 and the comparison theorem on impulsive differential equations, we have 𝑆(𝑡)≤𝑍(𝑡) and 𝑍(𝑡)→𝑆^∗(𝑡) as 𝑡→+∞. Hence, there exists a 𝑡₁>0 such that 𝑆(𝑡)≤𝑍(𝑡)<𝑆^∗(𝑡)+𝜀₁, (A.13) for all 𝑡≥𝑡₁; here 𝜀₁>0 is small enough.

From the second equation of model (4) and (A.13) we get d𝐼(𝑡)≤𝛽𝑆d𝑡^∗(𝑡)+𝜀₁−𝛾−𝑏𝐼(𝑡), (A.14) for 𝑡>𝑡₁; then we consider the following comparison equation with pulses: d𝑦₁(𝑡)=𝛽𝑆d𝑡^∗(𝑡)+𝜀₁−𝛾−𝑏𝑦₁(𝑡)𝑦,𝑡≠𝑛𝑇,₁𝑡⁺=𝑦₁(𝑡),𝑡=𝑛𝑇. (A.15) Integrating model (A.15) between pulses (𝑛𝑇,(𝑛+1)𝑇], we have 𝑦₁((𝑛+1)𝑇)=𝑦₁(𝑛𝑇)e^{∫(𝑛+1)𝑇𝑛𝑇[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡}=𝑦₁(𝑛𝑇)e^{∫𝑇0[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡}. (A.16) Then by using step by step iterations 𝑦₁(𝑛𝑇)=𝑦₁((𝑛−1)𝑇)e^{∫𝑇0[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡}=𝑦₁((𝑛−2)𝑇)e^{2∫𝑇0[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡}=⋯=𝑦₁(0)e^{𝑛∫𝑇0[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡}, (A.17) where 𝑦₁(0)=𝑦₁(0⁺)>0, we obtain lim_𝑛→∞𝑦₁(𝑛𝑇)=0 according to (A.10). On the other hand, integrating and solving the first equation of model (A.15) between pulses, it gives 𝑦₁(𝑡)=𝑦₁(𝑛𝑇)e^{∫𝑡𝑛𝑇[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡},𝑡∈(𝑛𝑇,(𝑛+1)𝑇], (A.18) incorporating into the boundedness of e^{∫𝑡𝑛𝑇[𝛽(𝑆∗(𝑡)+𝜀1)−𝛾−𝑏]d𝑡} due to 𝑅¹₀<1, we have lim_𝑡→∞𝑦₁(𝑡)=0.

Let (𝑆(𝑡),𝐼(𝑡)) be any solution of model (4) with initial value (𝑆₀,𝐼₀), and 𝑆₀=𝑆(0⁺)>0, 𝐼₀=𝐼(0⁺)=𝑦₁(0)>0, according to the comparison theorem we have lim_𝑡→∞sup𝐼(𝑡)≤lim_𝑡→∞sup𝑦₁(𝑡)=0. Incorporating into the positivity of 𝐼(𝑡), we know that lim_𝑡→∞𝐼(𝑡)=0.

Consequently, there exists a 𝑡₂>𝑡₁ such that 0<𝐼(𝑡)≤𝜀₂ for 𝑡≥𝑡₂, and 𝜀₂ is small enough. Further, we have 𝑏−𝛽𝜀₂𝑆(𝑡)≤(𝑡)+𝑏d𝑆(𝑡),d𝑡≤𝑏−𝑏𝑆 (A.19) for all 𝑡>𝑡₂, from which we can obtain the following equation: d𝑦₂(𝑡)d𝑡=𝑏−𝛽𝜀₂+𝑏𝑦₂(𝑡)𝑦,𝑡≠𝑛𝑇,₂𝑡⁺=1−𝑝_max𝑦²₂(𝑡)+𝜃𝑦₂(𝑡)𝑦₂(𝑡)+𝜃,𝑡=𝑛𝑇. (A.20)

By employing the same methods as used for the proof of Lemma 1 we get that the model (A.20) has a positive periodic solution 𝑦^∗₂(𝑡), which is globally attractive, where 𝑦^∗₂(𝑡)=𝑦₂−𝑏𝛽𝜀₂e+𝑏^{−(𝛽𝜀2+𝑏)(𝑡−𝑛𝑇)}+𝑏𝛽𝜀₂,(()𝑇],𝑦+𝑏𝑡∈𝑛𝑇,𝑛+1₂=𝑦₂𝑛𝑇⁺=𝑦^∗₂0⁺=−𝐵₂+𝐵²₂−4𝐴₂𝐶₂2𝐴₂,𝐴₂=e^{−(𝛽𝜀2+𝑏)𝑇}𝑏𝛽𝜀₂−+𝑏1−𝑝_maxe^{−(𝛽𝜀2+𝑏)𝑇},𝐵₂=e^{−(𝛽𝜀2+𝑏)𝑇}−𝑏𝛽𝜀₂×2+𝑏1−𝑝_maxe^{−(𝛽𝜀2+𝑏)𝑇}−𝑏𝜃+𝛽𝜀₂,𝐶+𝑏₂=−1−𝑝_max𝑏𝛽𝜀₂+𝑏−e^{−(𝛽𝜀2+𝑏)𝑇2}𝑏+𝜃𝛽𝜀₂+𝑏−e^{−(𝛽𝜀2+𝑏)𝑇}. (A.21) It follows from the comparison theorem on impulsive differential equations that 𝑦₂(𝑡)≤𝑆(𝑡)≤𝑍(𝑡). Moreover, 𝑦₂(𝑡)→𝑦^∗₂(𝑡) and 𝑍(𝑡)→𝑆^∗(𝑡) as 𝑡→+∞. Consequently, there exists a 𝑡₃ for 𝜀₃ small enough such that 𝑡₃≥𝑡₂ and 𝑦^∗₂(𝑡)−𝜀₃<𝑆(𝑡)<𝑆^∗(𝑡)+𝜀₃ (A.22) for 𝑡≥𝑡₃. Let 𝜀₂→0. Then 𝑆^∗(𝑡)−𝜀₃<𝑆(𝑡)<𝑆^∗(𝑡)+𝜀₃. (A.23)

Therefore, 𝑆(𝑡)→𝑆^∗(𝑡) as 𝑡→+∞. Thus, we have proved the global stability of the disease-free periodic solution (𝑆^∗(𝑡),0) of model (4). This completes the proof.

B. Proof of Theorem 3

In order to apply the bifurcation theory of [34] in the main text, we make the following calculations: 𝑑^′₀=1−e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}. (B.1) If 𝑑^′₀=0, then 𝑇₀ satisfies the following equation: 𝛽𝑆−11−e^−𝑏𝑇0+𝑏𝑇₀𝑏𝛾+𝑏𝑇₀=1, (B.2) which indicates that there exists a 𝑇₀ such that 𝑅¹₀=1 and the disease-free periodic solution 𝛿=(𝑆^∗(𝑡),0) loses its stability.

Further, 𝜕Φ₁𝑇₀,𝑋₀𝜕𝑆=e^−𝑏𝑇0>0,𝜕Φ₂𝑇₀,𝑋₀𝜕𝐼=e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}>0,𝜕Φ₁𝑇₀,𝑋₀𝜕𝐼=−^𝑇0₀𝛽𝑆^∗(𝑣)e^{−𝑏(𝑇0−𝑣)}e^{∫𝑣0(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}𝑎d𝑣≐𝐴<0,_′0=1−𝐴₃e^−𝑏𝑇0,𝑏^′₀=1−𝐴₃𝜕𝐴,₂Φ₂𝑇₀,𝑋₀𝜕𝐼𝜕𝑆=𝛽𝑇₀e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}𝜕>0,₂Φ₂𝑇₀,𝑋₀𝜕𝐼²=−𝛽^2𝑇0₀e^{∫𝑇0𝑣(𝛽𝑆∗(𝜉)−𝛾−𝑏)d𝜉}×^𝑣₀𝑆^∗(𝜃)e^{∫𝜃0(𝛽𝑆∗(𝜉)−𝛾−𝑏)d𝜉−𝑏(𝑣−𝜃)}𝜕d𝜃d𝑣<0,₂Φ₂𝑇₀,𝑋₀𝜕𝐼𝜕𝑇=𝛽𝑆^∗𝑇₀−𝛾−𝑏e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡},𝜕Φ₁𝑇₀,𝑋₀𝜕𝑇=−𝑏𝑆−1e^−𝑏𝑇0,𝐵=𝛽𝑏𝑇₀𝐴₃e^−𝑏𝑇01−𝐴₃e^−𝑏𝑇0𝑆𝑇₀−1−𝛽𝑆^∗𝑇₀−𝛾−𝑏×e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡},𝐶=𝛽^2𝑇0₀e^{∫𝑇0𝑣(𝛽𝑆∗(𝜉)−𝛾−𝑏)d𝜉}×^𝑣₀e^{−𝑏(𝑡−𝜃)}𝑆^∗(𝜃)e^{∫𝜃0(𝛽𝑆∗(𝜉)−𝛾−𝑏)d𝜉}d𝜃d𝑣+2𝛽𝑇₀1−𝐴₃𝐴1−𝐴₃e^−𝑏𝑇0e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}>0. (B.3) It follows that 𝐵𝐶<0 and we get the condition, that is, if the parameters satisfy condition 𝐵<0, then the model (4) has a supercritical branch at 𝑇₀.

In fact, e^{∫𝑇00(𝛽𝑆∗(𝑡)−𝛾−𝑏)d𝑡}>0. (B.4) Therefore, 𝐵<0 is equivalent to 𝛽𝑏𝑇₀𝐴₃e^−𝑏𝑇01−𝐴₃e^−𝑏𝑇0𝑆𝑇₀−1−𝛽𝑆^∗𝑇₀−𝛾−𝑏<0. (B.5) Obviously, 𝐵<0 due to (25). The proof is completed.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgments

This work is supported by the National Natural Science Foundation of China (NSFCs, 11171199, 11371030, 11301320) and the Fundamental Research Funds for the Central Universities (GK201305010). The authors would like to thank the editor and the referees for their careful reading of the original paper and many valuable comments and suggestions that greatly improved the presentation of this paper.