Optimal Control Policies of Pests for Hybrid Dynamical Systems

Abstract

We improve the traditional integrated pest management (IPM) control strategies and formulate three specific management strategies, which can be described by hybrid dynamical systems. These strategies can not only effectively control pests but also reduce the abuse of pesticides and protect the natural enemies. The aim of this work is to study how the factors, such as natural enemies optimum choice in the two kinds of different pests, timings of natural enemy releases, dosages and timings of insecticide applications, and instantaneous killing rates of pesticides on both pests and natural enemies, can affect the success of IPM control programmes. The results indicate that the pests outbreak period or frequency largely depends on the optimal selective feeding of the natural enemy between one of the pests and the control tactics. Ultimately, we obtain the only pest 𝑥₂ needs to be controlled below a certain threshold while not supervising pest 𝑥₁.

1. Introduction

A pest is an insect which is detrimental to humans or human concerns (as agriculture or livestock production). In its broadest sense, a pest is a competitor of humanity. Often insects are regarded as pests as they cause damage to agriculture by feeding on crops or parasitizing livestock, such as codling moth on apples or boll weevil on cotton. An animal could also be a pest when it causes damage to a wild ecosystem or carries germs within human habitats. Examples of these include those organisms which are vector-borne human diseases, such as rats and fleas which carry the plague disease, mosquitoes which are vector-borne malaria, and ticks which carry Lyme disease. The most serious pests (in the order of economic importance) are insects. Pesticides are chemicals and other agents (e.g., beneficial microorganisms) that are used to control or protect other organisms from insect pests. To control these insect pests, farmers rely strongly on intervention with chemical pesticides, which remain a significant component of the cost of production and ecological problems from pesticide resistance in key pests. In order to address the issue, researchers are increasingly embracing more components of the integrated pest management (IPM) [1–6] systems approach that is always applied in ecology. With the rise of interdisciplinary research, the mathematical ecology has also emerged and developed rapidly. A variety of mathematical methods can be used in ecological science. There have been numerous publications [7–15] over the last ten years using ecological mathematical model to research IPM strategy (spraying pesticides and introducing additional natural enemy into a pest-natural enemy system). When we study the dynamic property between the pest and its natural enemy (predator-prey), one of the most important components of the predator-prey relationship is the so-called functional responses. In [7–15], the Holling functional responses and Beddington-DeAngelis functional response are introduced. The Holling type extends the range of values of 𝑥 and 𝑦 over which the feeding term is realistic. However, in some situations, the increase of the feeding rate is not proportional to the increase of the predator density, as a result of mutual interference between predators, which decreases the efficiency of predation [16]. In addition, it is shown that most of the mathematical models on IPM include only one pest and one natural enemy. In [1], Finch and Collier’s study concerning IPM strategies in field vegetable crops focuses on two key pests, the cabbage root fly (Delia radicum) and the carrot fly (Psila rosae), the two major root feeding pests. Thus, we consider a continuous three-level food web model with Beddington-DeAngelis functional response, which consists of two competing pests (𝑥₁ and 𝑥₂) and a natural enemy (𝑦), and it can be represented as follows [17]: 𝑑𝑥₁𝑑𝑡=𝑎₁𝑥₁𝑥1−¹𝐾₁−𝑐₁𝑥₂𝐾₁−𝐹₁𝑥₁,𝑥₂,𝑦𝑦,𝑑𝑥₂𝑑𝑡=𝑎₂𝑥₂𝑥1−²𝐾₂−𝑐₂𝑥₁𝐾₂−𝐹₂𝑥₁,𝑥₂,𝑦𝑦,𝑑𝑦𝑑𝑡=𝜔₁𝐹₁𝑥₁,𝑥₂,𝑦𝑦+𝜔₂𝐹₂𝑥₁,𝑥₂𝑥,𝑦𝑦−𝑑𝑦,_𝑖(0)≥0,𝑦(0)≥0,𝑖=1,2, (1) where 𝐹_𝑖(𝑥₁,𝑥₂,𝑦)=𝛼_𝑖𝑥_𝑖/(𝑟+𝑥₁+𝑏𝑥₂+𝑐𝑦),  𝑖=1,2. 𝑎_𝑖, 𝐾_𝑖, 𝜔_𝑖, 𝑐_𝑖, and 𝑑 are positive (𝑖=1,2). The prey 𝑥_𝑖 grows with intrinsic growth rate 𝑎_𝑖 and carries capacity 𝐾_𝑖 in the absence of the predator. The constant 𝑐_𝑖 is the interspecific competition rate between the two prey species. The predator 𝑦 consumes the prey 𝑥_𝑖 with a functional response of Beddington-DeAngelis type 𝐹_𝑖(𝑥₁,𝑥₂,𝑦) and contributes to its growth with rate 𝜔_𝑖𝐹_𝑖(𝑥₁,𝑥₂,𝑦). The constant 𝑑 is the death rate of the predator, and the term, 𝑐𝑦, measures the mutual interference between predators. The constants 𝛼_𝑖(𝑖=1,2) show the predation capacity of predators, and 𝑟 is the saturating functional response parameters. The parameter 𝑏 is the predator’s relative preference on 𝑥₂ with respect to 𝑥₁.

There are numbers of diagnostic tools to detect the qualitative behavior of the system (1). We give a simple investigation of the dynamics of the system (1) through constructing a bifurcation diagram using the software XPPAUT [18]. To construct the bifurcation diagram, we use a numerical integration with varying a key parameter and keeping other parameters fixed; our choice of parameters is guided by two assumptions: first, the system has to be biologically feasible, and second, natural enemies have food preference phenomenon [19]; that is, 𝛼₁>𝛼₂ or 𝛼₁<𝛼₂ in this paper. The bifurcation diagram as a function of 𝛼₁ in the range 0.3<𝛼₁<0.7 is drawn in Figure 1. As shown in Figure 1, the thick black lines represent the stable limit cycle, the black solid curves represent the stable equilibrium state, and the black dashed curves represent the unstable equilibrium state. When 𝛼₁=0.4601, system (1) undergoes a Hopf bifurcation, where the stable focus becomes unstable and a stable limit cycle emanates from it. When the parameter 𝛼₁ increases further, system (1) undergoes a Hopf bifurcation again at 𝛼₁=0.5644, where an unstable focus becomes stable and the other unstable focus is still unstable. When 𝛼₁=0.5698, a transcritical bifurcation appears, where the stability of two focuses occurs in exchange. In Figure 2, we fixed the parameter 𝛼₁=0.6000, and the other parameters are the same as those in Figure 1. As we can see in Figure 2, the pest 𝑥₁ tends to be extinct eventually due to the impact of the food preference phenomenon of natural enemies before recurrence of the next generation of pests, while the pest 𝑥₂ and natural enemies can coexist through stabilizing the boundary equilibrium.

Figure 1 

Bifurcation diagram of system (1) for 0.3<𝛼1<0.7, where 𝑋=𝑥1, 𝑌=𝑥2, 𝑍=𝑦, and 𝐴1=𝛼1 and the initial value 𝑥1=1.3;𝑥2=1.5; and 𝑦=5 with parameters as follows:𝑎1=0.9,  𝑎2=0.8,  𝛼2=0.5,  𝑐1=0.55,  𝐾1=20,  𝑟=1.2,  𝑏=0.25,  𝑐=0.3,  𝛼2=0.5,  𝐾2=30,  𝑐2=0.55,  𝑊1=0.85,  𝑊2=0.75, and  𝑑=0.25.

Figure 2 

Time series and phase portrait of system (1) with 𝛼1=0.6000 and other parameters the same as Figure 1.

In order to manage pests through spraying pesticides and introducing additional natural enemies (IPM strategy), we apply the impulsive differential equations (IDE) [20–22] to integrate system (1). In [13], the authors constructed a predator-prey impulsive system to show the process of releasing natural enemies periodically and spraying pesticides twice at different fixed times in a period. The authors, obviously, avoid the side effects of pesticides on natural enemies existing in [7–9]. But only two pulses in a period cannot effectively control pests. And most of the research is single species of pests while few papers have discussed multispecies of pests. However, most real pests communities are more complex than the community previously analysed by them. In the present paper, we make the following improvements: (a) in the modeling, we introduce a one-natural enemy and two-pest model, and the natural enemy shows optimal foraging between pest 𝑥₁ and pest 𝑥₂, which is a well-known behavior of many predators [23]; (b) in the control strategy, we can control two pests more selectively by controlling pulse frequencies appropriately.

2. Model Formulation and Auxiliary Lemmas

Considering the previous factors, firstly, two models of different control strategies are discussed as follows.

Case 1. We suppose pesticides are sprayed several times in a release period, and the kill rates of pesticides to pests (𝛿^𝑖_𝑘, 𝑖=1,2) and natural enemies (𝛿³_𝑘) are different in different impulsive moments in a release period. That is, we consider the following model: 𝑑𝑥₁𝑑𝑡=𝑎₁𝑥₁𝑥1−¹𝐾₁−𝑐₁𝑥₂𝐾₁−𝛼₁𝑥₁𝑦𝑟+𝑥₁+𝑏𝑥₂,+𝑐𝑦𝑑𝑥₂𝑑𝑡=𝑎₂𝑥₂𝑥1−²𝐾₂−𝑐₂𝑥₁𝐾₂−𝛼₂𝑥₂𝑦𝑟+𝑥₁+𝑏𝑥₂,+𝑐𝑦𝑑𝑦=𝜔𝑑𝑡¹𝛼₁𝑥₁𝑦𝑟+𝑥₁+𝑏𝑥₂+𝜔+𝑐𝑦₂𝛼₂𝑥₂𝑦𝑟+𝑥₁+𝑏𝑥₂+𝑐𝑦−𝑑𝑦,𝑡≠𝜏_𝑛𝑘𝑥,𝑡≠𝑛𝑇,₁𝜏⁺_𝑛𝑘=1−𝛿¹_𝑘𝑥₁𝜏_𝑛𝑘,𝑥₂𝜏⁺_𝑛𝑘=1−𝛿²_𝑘𝑥₂𝜏_𝑛𝑘,𝑦𝜏⁺_𝑛𝑘=1−𝛿³_𝑘𝑦𝜏_𝑛𝑘,𝑡=𝜏_𝑛𝑘,𝑥₁𝑛𝑇⁺𝑥=0,₂𝑛𝑇⁺𝑦=0,𝑛𝑇⁺=𝑦(𝑛𝑇)+𝜃,𝑡=𝑛𝑇, (2) where 𝑛𝑇=𝜏_𝑛0<𝜏_𝑛1<𝜏_𝑛2<⋯<𝜏_𝑛𝑞<𝜏_𝑛(𝑞+1)=(𝑛+1)𝑇,𝑛,𝑞∈𝑍⁺; 𝑘=0,1,2,…,𝑞; 𝜏_𝑛(𝑘+1)−𝜏_𝑛𝑘=𝑇/(𝑞+1); 0≤𝛿^𝑖_𝑘<1(𝑖=1,2,3) represents an effective kill rate to pest 𝑥₁, pest 𝑥₂, and natural enemy 𝑦 at time 𝑡=𝜏_𝑛𝑘, respectively. 𝜃 represents natural enemy number of additional release at 𝑡=𝑛𝑇; 𝑇 is a release period; 𝜏_𝑛𝑘 is the 𝑘th spraying moment in the 𝑛th release period. Other parameters are the same as those in system (1).

Case 2. We suppose that the natural enemies are released several times in a spraying period. The release number of natural enemies (𝜃_𝑘) and an effective kill rate to natural enemies (𝜇_𝑘) are different in different impulsive moments in a spraying period. That is, we consider the following model: 𝑑𝑥₁𝑑𝑡=𝑎₁𝑥₁𝑥1−¹𝐾₁−𝑐₁𝑥₂𝐾₁−𝛼₁𝑥₁𝑦𝑟+𝑏₁𝑥₁+𝑏₂𝑥₂,+𝑐𝑦𝑑𝑥₂𝑑𝑡=𝑎₂𝑥₂𝑥1−²𝐾₂−𝑐₂𝑥₁𝐾₂−𝛼₂𝑥₂𝑦𝑟+𝑏₁𝑥₁+𝑏₂𝑥₂,+𝑐𝑦𝑑𝑦=𝜔𝑑𝑡¹𝛼₁𝑥₁𝑦𝑟+𝑏₁𝑥₁+𝑏₂𝑥₂+𝜔+𝑒𝑦₂𝛼₂𝑥₂𝑦𝑟+𝑏₁𝑥₁+𝑏₂𝑥₂+𝑒𝑦−𝑑𝑦,𝑡≠𝜆_𝑛𝑘,𝑡≠𝑛𝑇₁,𝑥₁𝜆⁺_𝑛𝑘𝑥=0,₂𝜆⁺_𝑛𝑘𝑦𝜆=0,_+𝑛𝑘=1−𝜇_𝑘𝑦𝜆_𝑛𝑘+𝜃_𝑘,𝑡=𝜆_𝑛𝑘,𝑥₁𝑛𝑇⁺₁=1−𝛿¹₁𝑥₁𝑛𝑇₁,𝑥₂𝑛𝑇⁺₁=1−𝛿²₁𝑥₂𝑛𝑇₁,𝑦𝑛𝑇⁺₁=1−𝛿³₁𝑦𝑛𝑇₁,𝑡=𝑛𝑇₁, (3) where 𝑛𝑇₁=𝜆_𝑛0<𝜆_𝑛1<𝜆_𝑛2<⋯<𝜆_𝑛𝑝<𝜆_𝑛(𝑝+1)=(𝑛+1)𝑇₁,𝑛,𝑝∈𝑍⁺; 𝑘=0,1,2,…,𝑝;𝜆_𝑛(𝑘+1)−𝜆_𝑛𝑘=𝑇₁/(𝑝+1); 0≤𝛿^𝑖₁<1(𝑖=1,2,3) represents an effective kill rate to pest 𝑥₁, pest 𝑥₂ and natural enemy 𝑦 at time 𝑡=𝜆_𝑛𝑘, respectively. 𝜇_𝑘 represents reduced proportion of natural enemies owing to the delay effect of pesticides and eating the deleterious pests; 𝜃_𝑘 represents natural enemy number of additional release at 𝑡=𝜆_𝑛𝑘; 𝑇₁ is a spraying period; 𝜆_𝑛𝑘 is the 𝑘th releasing moment in the spraying 𝑛th period. Other parameters are the same as those in system (1).
Furthermore, some essential notations, definitions and lemmas are given as follows.
Let 𝑅₊=[0,∝) and 𝑅³₊={𝐱=(𝑥₁(𝑡),𝑥₂(𝑡),𝑦(𝑡)∈𝑅³∶𝑥₁(𝑡),𝑥₂(𝑡),𝑦(𝑡)≥0)}. Denote 𝐟=(𝑓₁,𝑓₂,𝑓₃) the map defined by the right-hand side of the first, second, and third equations of system (2). The solution of system (2), denoted by 𝐱∶𝑅₊→𝑅³₊, is piecewise continuous, and it is continuous on ((𝑛−1)𝑇,𝜏_𝑛1],(𝜏_𝑛1,𝜏_𝑛2],…,(𝜏_𝑛𝑞,𝑛𝑇]. Note that 𝐱((𝜏⁺_𝑛𝑘))=lim_𝑡→𝜏⁺_𝑛𝑘𝑥(𝑡) and 𝑥(𝑛𝑇⁺)=lim_{𝑡→𝑛𝑇}⁺𝑥(𝑡) exist. The existence and uniqueness of solutions of system (2) are guaranteed by the smoothness of 𝐟 (see [22]).

Lemma 1. Assume that 𝐱(𝑡) is a solution of system (2) with 𝐱(0⁺)≥0, then 𝐱(𝑡)≥0 for 𝑡≥0.

Lemma 2 (see Lakshmikantham et al. [22]). Consider the following impulsive differential inequalities: 𝑚^′(𝑡)≤(≥)𝑝(𝑡)𝑚(𝑡)+𝑞(𝑡),𝑡≠𝑡_𝑘,𝑚𝑡⁺_𝑘≤(≥)𝑑_𝑘𝑚𝑡_𝑘+𝑏_𝑘,𝑡=𝑡_𝑘,𝑘∈𝑁, (4) where 𝑝(𝑡),𝑞(𝑡)∈𝐶(𝑅₊,𝑅), 𝑑_𝑘≥0, and 𝑏_𝑘 are constants.
Assume that:(A₀) the sequence {𝑡_𝑘} satisfies 0≤𝑡₀<𝑡₁<𝑡₂<⋯, with lim_𝑡→∝𝑡_𝑘=∝;(A₁)𝑚(𝑡)∈𝑃𝐶^′(𝑅₊,𝑅) and 𝑚(𝑡) is left continuous at 𝑡_𝑘, 𝑘=1,2,…;(A₂) for 𝑘=1,2,…,𝑡≥𝑡₀.
Then, 𝑚(𝑡)≤(≥)𝑚𝑡_{0𝑡0<𝑡𝑘<𝑡}𝑑_𝑘exp^𝑡_𝑡0𝑝(𝑠)𝑑𝑠+_{𝑡0<𝑡𝑘<𝑡𝑡𝑘<𝑡𝑗<𝑡}𝑑_𝑗exp^𝑡_𝑡𝑘𝑝(𝑠)𝑑𝑠𝑏_𝑘+^𝑡_{𝑡0𝑠<𝑡𝑘<𝑡}𝑑_𝑘exp^𝑡_𝑠𝑝(𝜎)𝑑𝜎𝑞(𝑠)𝑑𝑠,𝑡≥𝑡₀. (5)

Remark 3. The previous definitions and lemmas may similarly be applied in system (3).

3. Dynamical Analysis of Case 1 and Its Biological Implications

For Case 1, the basic properties of the following subsystem: 𝑑𝑦(𝑡)𝑑𝑡=−𝑑𝑦,𝑡≠𝑛𝑇,𝑡≠𝜏_𝑛𝑘,𝑦𝑡⁺=1−𝛿³_𝑘𝑦(𝑡),𝑡=𝜏_𝑛𝑘𝑦𝑡,𝑘=1,2,…,𝑞,₊(𝑡)𝑦0=𝑦+𝜃,𝑡=𝑛𝑇,₊=𝑦₀ (6) play a key role in analyzing the pest control.

It is shown in the Appendix that there exists a globally stable periodic solution 𝑦_𝑇(𝑡) for the subsystem (6). Therefore, the complete expression for the pest-eradication periodic solution of system (6) over the (𝑛−1)th time interval (𝑛−1)𝑇≤𝑡≤𝑛𝑇 is given by (0,0,𝑦_𝑇(𝑡)). Furthermore, if the following threshold condition: 𝜆_𝑖=^𝑞_𝑘=11−𝛿^𝑖_𝑘exp𝑎_𝑖𝛼𝑇+^𝑖𝐵₁()𝑐𝑑<1𝑖=1,2 (7) is satisfied, then the pest-eradication periodic solution (0,0,𝑦_𝑇(𝑡)) is globally asymptotically stable, where 𝑌^∗=𝜃/[1−(1−𝛿³₁)⋯(1−𝛿³_𝑞)𝑒^−𝑑𝑇] and 𝐵₁=ln𝑟+𝑐𝑌^∗𝑒^{−𝑑𝑇/(𝑞+1)}𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−2𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^−𝑑𝑇×𝑟+𝑐𝑌^∗𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)−1}. (8)

Denote 1−𝛿¹_𝑘=𝛿₁, 1−𝛿²_𝑘=𝛿₂, and 1−𝛿³_𝑘=𝛿₃, where 𝑘=1,2,…,𝑞. What we want to address in the following is how control tactics including residual rates 𝛿₁, 𝛿₂, and 𝛿₃ release constant 𝜃, timing of pesticide application 𝑞, and timing of release period 𝑇 affect the threshold condition 𝜆_𝑖. Firstly, in the following section, we take Figure 2 as an example to control pests by the control strategy of Case 1. It implies that only the pest 𝑥₂ needs be eradicated. We firstly use the traditional control method, spraying insecticides, in order to know well the impact of insecticides on pests and natural enemy species. In general, pesticides tend to be harmful to most natural enemies [24], which may be associated with the acute poisoning. It is significant to understand the acute poisoning of insecticides to natural enemies for the research of IPM strategy.

We only apply insecticides but we do not release natural enemies as shown in Figure 3. The simulation results indicate that pesticide applications (number of pesticide applications 𝑞=1 in a period shown in Figure 3(a)) do not lead to the extinction of pest 𝑥₂, and on the contrary, they can result in the recurrence of pest 𝑥₁, and with the increase of the number 𝑞, the quantityof both pest 𝑥₁ and pest 𝑥₂ increases (see Figure 3(b)). Only when the dosages of pesticides are increased enough can the pests become extinct, but natural enemies also become extinct at the same time (see Figure 3(c)). This shows that the extinction of pests needs plenty of pesticides. Nevertheless, pesticide abuse can bring about environmental contamination, which can also result in human exposure through consumption of residues of pesticides in food and, possibly, drinking water. While developed countries have systems already in place to register pesticides and control their trade and usage, this is not always the case elsewhere, especially in China. Moreover, the pesticides have a serious impact on the natural enemies (Figures 3(a)–3(c)), and the repeated use of the same pesticides can result in one or more population pest outburst (Figures 3(a)-3(b)). The previous results show that pest control of multispecies is much more complicated than single pests [7–9].

(a)

(b)

(c)

Figure 3 

The impact of the number and dosage of application insecticides on pests and natural enemy with parameters of (a) 𝛿1=0.9,  𝛿2=0.85,  𝛿3=0.85,  𝜃=0,  𝑇=6, and  𝑞=1, (b)  𝛿1=0.95,  𝛿2=0.85,  𝛿3=0.85,  𝜃=0,  𝑇=6, and  𝑞=6, (c) 𝛿1=0.58,  𝛿2=0.62,  𝛿3=0.85,  𝜃=0,  𝑇=6, and  𝑞=6, and other parameters the same as Figure 1.

By the previous analysis, the additional release of natural enemies is an indispensable part for pest control. Without loss of generality, we assume that natural enemies have food preference phenomenon with pest 𝑥₁; that is, 𝛼₁>𝛼₂, the intrinsic growth rate 𝑎₁≤𝑎₂, and the residual ratio of pests 𝛿₁≤𝛿₂ (here the reason is explained in the following section) after spraying; thus, according to Theorem A.2, 𝜆₁<𝜆₂ and the condition of pest extinction only needs 𝜆₂<1. In Figure 4, we fix the other parameters of 𝜆₂ and let the residual ratio 𝛿₃ vary. The simulation results indicate that if the pesticide poisons the natural enemies with a relatively low residual ratio 𝛿₃ (e.g., 𝛿₃=0.9), the threshold value 𝜆₂ is a monotonically increasing function with respect to the number of pesticide applications 𝑞 (Figure 4(a)). This further explains that if the pesticide has a severe impact on the natural enemies, repeated use of the same pesticides can result in pest resurgence. If the residual ratio 𝛿₃ is slightly increased from 0.9 to 0.92, Figure 4(b) shows that the threshold value 𝜆₂ is not monotonic with respect to the number of pesticide applications 𝑞. So in this case we must carefully select the number of pesticide applications (one to three events in this case). If the pesticides do not kill the natural enemies so much, Figures 4(c) and 4(d) clarify that the threshold value 𝜆₂ is a monotonically decreasing function with respect to the number of pesticide applications 𝑞. In Figures 5(a)–5(d), similarly, we fix the other parameters of 𝜆₂ and let the release period 𝑇 vary. The simulation results indicate that the small change of release period 𝑇 can lead to the change of the number of pesticide applications 𝑞. All these simulations show that the releasing period, the number of times of spraying pesticides within this period, and the residual ratio of natural enemies are crucial to eradicate pest 𝑥₁ and pest 𝑥₂. Figure 6 shows the relationship between the controllable parameters and the threshold condition 𝜆₂. All simulation results demonstrate that 𝜆₂ seems to be quite sensitive to small changes in residual ratio 𝛿₁ and 𝛿₃, release constant 𝜃, and release period 𝑇. By the simulation results, we obtain that the optimum time and frequency of pesticide application; the protection of natural enemies and pesticides choice are the key factors to pest outburst or eradication. The results may provide a theoretical basis for agricultural practitioners to guide them to spray pesticides and release natural enemies more efficiently.

(a)

(b)

(c)

(d)

Figure 4 

The effects of number of pesticide applications and the residual ratio of natural enemy on the threshold level 𝜆2. The parameter values are as follows:  𝛿1=0.85,  𝛿2=0.85,  𝑇=10.7,  𝜃=6.3,  𝑑=0.25,  𝑟=1.2,  𝑐=0.3,  𝑎1=0.75,  𝑎2=0.8, and  𝛼2=0.7.

(a)

(b)

(c)

(d)

Figure 5 

The effects of number of pesticide applications and release period 𝑇 on the threshold level 𝜆2 with the parameter value 𝛿3=0.93 and other parameters the same as Figure 4.

(a)

(b)

(c)

(d)

Figure 6 

The effects of times of spraying pesticides and parameter sets on the threshold level 𝜆2. The parameter values are as follows:  𝑎1=0.8,  𝑑=0.25,  𝑟=1.2,  𝑐=0.3, and  𝛼2=0.7. ((a), (b))  𝛿2=0.75,  𝛿3=0.94,  𝜎=6.3, and  𝑇=10; ((c), (d))  𝛿2=0.7,  𝛿3=0.91,  𝜎=6.3, and  𝑇=10.

4. Dynamical Analysis of Case 2 and Its Biological Implications

In this section, similar to the study method of system (2), we will give the dynamic property of system (3).

According to Floquet theory [16], if the following threshold condition: 𝜆_𝑖=𝛿_𝑖exp𝑎_𝑖𝑇₁+𝛼_𝑖𝐵₂()𝑐𝑑<1𝑖=1,2 (9) is satisfied, then the pest-eradication periodic solution (0,0,𝑦_𝑇1(𝑡)) of system (3) is globally asymptotically stable, where 𝐵₂=ln𝑟+𝑐𝑌^∗₁𝑒^{−𝑑𝑇1/(𝑝+1)}×𝑟+𝑐𝑌^∗₁𝜂₁𝑒^{−2𝑑𝑇1/(𝑝+1)}+𝜃₁𝑒^{−𝑑𝑇1/(𝑝+1)}⋯𝑟+𝑐𝑌^∗₁^𝑝_𝑘=1𝜂_𝑘𝑒^−𝑑𝑇1+𝜃₁^𝑝_𝑘=2𝜂_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)}+⋯+𝜃_𝑝𝑒^{−𝑑𝑇1/(𝑝+1)}×𝑟+𝑐𝑌^∗₁𝑟+𝑐𝑌^∗₁𝜂₁𝑒^{−𝑑𝑇1/(𝑝+1)}+𝜃₁⋯𝑟+𝑐𝑌^∗₁^𝑝_𝑘=1𝜂_𝑘𝑒^{−𝑝𝑑𝑇1/(𝑝+1)}+𝜃₁^𝑝_𝑘=2𝜂_𝑘𝑒^{(−𝑑(𝑝−1)𝑇)/(𝑝+1)}+⋯+𝜃_𝑝⁻¹,𝑌^∗₁=𝑆1−1−𝜇₁1−𝜇₂⋯1−𝜇_𝑝𝛿₃𝑒^−𝑑𝑇1,𝜂_𝑘=1−𝜇_𝑘,𝛿₁=1−𝛿¹₁,𝛿₂=1−𝛿²₁,𝛿₃=1−𝛿³₁. (10)

For Case 2, there are 𝑝 times releasing natural enemy during spraying period 𝑇₁. Denote 𝜂_𝑘=𝜇, 𝜃_𝑘=𝜃, where 𝑘=1,2,…,𝑝. Since the release of natural enemies in this case is more frequent than spraying pesticides, the side effects of pesticides on the natural enemy population are largely reduced. Moreover, the threshold condition 𝜆₂ can be strongly affected by the additional release of natural enemies. In Figure 7, we fix all the other parameters and choose a different releasing constant 𝜃 and different release times 𝑝. The simulation results indicate that slight perturbation of the release constant 𝜃 can rapidly reduce the threshold value 𝜆₂ (Figure 7), while increasing the number of natural enemy releases as well. This shows that repeated releases of a small number of natural enemies in key time of the season can effectively control the pest outburst. In practice, an example of Liriomyza sativae and Trialeurodes vaporariorum occurs in heliogreenhouse. The parasitic rate of parasitic wasps which is their natural enemy can amount to over fifty percent without drugs. Thence, periodic releases of the parasitic wasps have been used to control the Liriomyza sativae and Trialeurodes vaporariorum in Anshan city in Liaoning Province, China, where greenhouse agriculture is developing rapidly.

Figure 7 

The effects of times of releasing natural enemies and releasing constant 𝜃 on the threshold level 𝜆2. The parameter values are as follows:  𝛿1=0.92,  𝛿2=0.93,  𝛿3=0.92,  𝜇=0.99,  𝑇1=5,  𝑑=0.25,  𝑟=1.2,  𝑐=0.3,  𝑎1=0.8, and  𝛼2=0.7.

Remark 4. In this paper, we suppose that the natural enemies can first select which are their favorite prey between pest 𝑥₁ and pest 𝑥₂. It means that the favorite natural enemies may be a profitable pest to them. Thence, as seen in [19], the profitable pest is classified as palatable and the other as unpalatable. In the following section, by numerical simulation, we will explain under what condition the natural enemies can prey on pest 𝑥₁ or pest 𝑥₂.

5. Hybrid Impulsive Model with Economic Threshold

As the previous simulation indicates, pesticides may seriously influence the survival of natural enemies. They may impact natural enemies indirectly by killing or contaminating their hosts or prey. It is essential to avoid pesticide abuse when biological control is feasible, as shown in systems (2) and (3). Probably the best method for reducing the side effects of pesticides on natural enemies is to apply pesticides only when the sum of density of two pest populations reaches the economic threshold (ET), since a small number of insect pests may have compensation effect on crops [25]. Thence, we formulate the model as follows: 𝑑𝑥₁𝑑𝑡=𝑎₁𝑥₁𝑥1−¹𝐾₁−𝑐₁𝑥₂𝐾₁−𝛼₁𝑥₁𝑦𝑟+𝑥₁+𝑏𝑥₂,+𝑐𝑦𝑑𝑥₂𝑑𝑡=𝑎₂𝑥₂𝑥1−²𝐾₂−𝑐₂𝑥₁𝐾₂−𝛼₂𝑥₂𝑦𝑟+𝑥₁+𝑏𝑥₂,+𝑐𝑦𝑑𝑦=𝜔𝑑𝑡¹𝛼₁𝑥₁𝑦𝑟+𝑥₁+𝑏𝑥₂+𝜔+𝑐𝑦₂𝛼₂𝑥₂𝑦𝑟+𝑥₁+𝑏𝑥₂𝑥+𝑐𝑦−𝑑𝑦,₁+𝑥₂𝑥<ET,𝑡≠𝑛𝑇,₁𝜏⁺_𝑛𝑘=𝛿₁𝑥₁𝜏_𝑛𝑘,𝑥₂𝜏⁺_𝑛𝑘=𝛿₂𝑥₂𝜏_𝑛𝑘,𝑦𝜏⁺_𝑛𝑘=𝛿₃𝑦𝜏_𝑛𝑘,𝑥₁+𝑥₂𝑥=ET,₁𝑛𝑇⁺𝑥=0,₂𝑛𝑇⁺𝑦=0,𝑛𝑇⁺=𝑦(𝑛𝑇)+𝜃,𝑡=𝑛𝑇, (11) where 𝛿₁, 𝛿₂, and 𝛿₃ are the same as in Section 5 and ET is the economic threshold. 𝑇 is the releasing period of natural enemies.

The previously mentioned facts show that the effect among the intrinsic growth rate of pests, 𝑎₁ and 𝑎₂, the predation capacity of natural enemies (or functional response parameter), 𝛼₁ and 𝛼₂, residual ratio 𝛿₁, 𝛿₂, and 𝛿₃, releasing quantity 𝜃, and other factors (such as ET) may determine dynamic behavior of pests and natural enemy species. How do these key parameters affect the control strategies? In particular, what we want to achieve is to study how the ET, 𝛼₁ or 𝛼₂, and controllable parameters (such as 𝑇) affect the control strategies.

For a fixed ET, by simulation, we obtain the result that the successful control strategies are affected by the predation capacity of natural enemy to different types of pests. To show this, we vary the key parameter 𝛼₁ while the other parameters are fixed as those in Figure 8. In Figure 8(a), for the predation capacity of natural enemy 𝛼₁=0.74, the simulation result indicates that the sum of density of the two pests population never reaches the given ET, which implies that 𝛼₁≥0.74 is free from spraying. If we set 𝛼₁=0.737, Figure 8(b) indicates that the system is free from chemical control after spraying pesticides. If we set 𝛼₁=0.736 or 𝛼₁=0.732 or 𝛼₁=0.73, Figures 8(c)–8(e) indicate that the system is free from chemical control after three, four, or five pesticide applications. If we further reduce the predation capacity of natural enemy and set 𝛼₁=0.71, the pest outbreak frequency is sharply increased, as shown in Figure 8(f). As mentioned in Remark 4, in Figure 8, we will expound the interdependent relationship among the natural enemy, pest 𝑥₁, and pest 𝑥₂. In the beginning, the density of pest 𝑥₂ is larger, and it is regarded as palatable for the natural enemy,which implies that 𝛼₂≥𝛼₁, as shown in Figure 8(f), the pests can break out. However, with the release of natural enemies, an expanded population of natural enemies causes the reduction of pest 2, but when pest 𝑥₂ falls below a certain critical value 𝑥^∗₂ [19], the natural enemy begins to eat not only pest 𝑥₂ but also pest 𝑥₁, which largens 𝛼₁ and makes pests gradually no longer break out as shown in Figures 8(e)–8(a). Meanwhile, it causes an immediate recovery of pest 𝑥₂, and when pest 𝑥₂ raises above the certain critical value 𝑥^∗₂, the natural enemy begins to eat palatable pest 𝑥₂ again, which forms a cycle (it is also called switching between pest 𝑥₁ and pest 𝑥₂ by natural enemy [26]). The previous analysis illustrates that if the natural enemy and the two pests meet the above relationship, we only need to control pest 𝑥₂ falling below the certain critical value without supervising pest 𝑥₁.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 8 

The effects of ET and 𝛼1 of model (11) on the control strategies. The initial values are 𝑥1=1.3;  𝑥2=1.5; and  𝑦=5, and parameters are fixed as follows: 𝛿1=0.65,  𝛿2=0.7,  𝛿3=0.95,  𝜎=1,  ET=8,  𝑇=4,  𝑎1=0.75,  𝑎2=0.85,  𝛼2=0.85,  𝑐1=0.2,  𝐾1=20,  𝑏=0.25,  𝑐=0.3,  𝐴2=0.85,  𝐾2=30,  𝑐2=0.65,  𝑊1=0.45,  𝑊2=0.75,  𝑑=0.25, and  𝑟=1.2.

Finally, we will introduce the definition of pest outbreak duration (period) and analyse the relationship between pest outbreak period and the controllable parameters. We denote the time points at which the solution reaches ET as 𝑡_𝑛(𝑛=1,2,…). If mod(𝑡_𝑛,𝑇)=0, a chemical control is applied at 𝑡_𝑛, and after that a biological control is also applied at the same time. If mod(𝑡_𝑛,𝑇)≠0, only a chemical control programme is applied. Further, denote 𝑇_𝑛=𝑡_𝑛−𝑡_𝑛−1 (12) with 𝑡₀=0 as pests-outbreak duration (or period), where 𝑛 may be finite or infinite which depends on the solutions of the models. The relationship among 𝑇, ET, 𝛿₁, 𝛿₂, 𝛿₃, or 𝜃 and mean outbreak period of pests can be calculated from model (11) and formula (12) numerically (Figure 9). Mean pest outbreaks period is an average over several pest outbreak (here outbreaks indicate that the sum of densities of the pest 𝑥₁ and pest 𝑥₂ reaches the given ET). Model (11) predicts that the pests do not break out if the natural enemies are released more transitorily (𝑇≤2, Figure 9(d)) and the mean outbreak period is decreasing as the release period 𝑇 or residual ratios of the pests 𝛿₁ and 𝛿₂ increase (Figures 9(a), 9(b), and 9(d)). Conversely, model (11) predicts that with the increase of ET or residual ratios of the natural enemy 𝛿₃, the mean outbreak period becomes longer (Figures 9(c) and 9(e)). In Figure 9(f), let the release period 𝑇 and the other parameters fixed, and let the release constant 𝜃 vary. This indicates that when the release period 𝑇 is smaller (here 𝑇=4), with the increase of release quantity, the pests will not break out (here 𝜃≥1.5). And, more remarkably, the mean outbreak period can suddenly jump from a small value to a larger value at some critical points of 𝛿₃, 𝑇, and 𝜃, which implies that the protection of natural enemies, the selection of releasing time, and the quantity may be crucial in prolonging the pest outbreak period. Moreover, the different 𝛿₂ or 𝛿₃ or different values of the release constant 𝜃 may have the same mean outbreak period (Figures 9(b), 9(c) and 9(f)). For system (11), the relationship between pest outbreak period and other parameters such as 𝑎₁, 𝑎₂, 𝛼₁, and 𝛼₂ can be researched similarly.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 9 

The relationship between mean outbreak period and the parameters 𝛿1, 𝛿2, 𝛿3, 𝑇, ET, and 𝜃 of model (11). The other parameters are the same as Figure 8.

6. Discussion

The agricultural pests management plays a decisive role on the survival of people all over the world especially that the impacts of extreme climate change are severer for pest control. For example, the armyworm which is the typical pest threatening corn growth in fall has been widely seen in North China Plain and Northeast China producing regions in August 2012. The leaves of corn stalks in portions of the above regions have been eaten up, cutting corn harvest prospect. The pests, common pests, became the most serious threat to the production of corn this summer in the country’s major grain-producing regions. This is mainly due to frequent cyclonic activities since mid-July provided favorable conditions to the migration of the armyworms, and then heavy rainfall forced them to stay in the north and northeastern parts of the country. The pest outbreak occurred at the same time as severe droughts in the United States, where the driest conditions in more than half a century have battered corn and soybean crops, causing an upsurge in global grain prices [27, 28].

Thence, it is particularly significant to explore an effective control strategy. In this paper, based on the IPM strategies, we give three different control strategies, which improve traditional IPM control strategies. For system (2), in a releasing natural enemies period, we spray pesticides several times. By the theoretical derivation to system (2), the critical value of pests eradication is figured out. When only using the traditional control method, insecticides, by the numerical simulation, we obtain that the pest 𝑥₁ and pest 𝑥₂ may break out with the increase of spray times, which is different from the control of a single pest. To better understand how the controllable parameters (here 𝑇, 𝑞, 𝛿₁, 𝛿₂, and 𝛿₃) impact the pest control, by the numerical simulation, we give the relationship between the critical value of pests eradication and the controllable parameters. All these results express that the selection of spray times and the protection of the natural enemies are of vital importance for pests eradication. From system (3), we know the real embodiment of the significance of protecting natural enemies. By the analysis of systems (2) and (3), in summary, when insecticides are used excessively, the pests are killed and the natural enemies of the pests are wiped out. In the absence of natural enemies, the surviving population of insect pests multiply rapidly and reach epidemic proportions. Indiscriminate use of pesticides also leads to the development of resistance in pests. This occurs as a result of killing the susceptible genotypes and selecting the more resistant genotypes at every pesticide application event. After several years of using the same pesticide, there would come a time when that particular pesticide will have no effect on the pests because they have developed resistance to the pesticide. Considering the factors, we formulate model (11), which is to apply pesticides only when necessary. By model (11), we also obtain some important conclusions.

Most real communities are more complex than the community analysed here. Therefore, in the future, the factors on the pests of more species and natural enemies, resistance to the pesticide of pests and so on should be considered in the model to depict the dynamic behavior between populations much more accurately.

Appendix

In any time interval ((𝑛−1)𝑇,𝑛𝑇], we investigate the dynamical behavior of model (6). In fact, integrating the first equation of model (6) from (𝑛−1)𝑇 to 𝜏_𝑛1 yields 𝑦(𝑡)=𝑦(𝑛−1)𝑇⁺𝑒^{−𝑑(𝑡−(𝑛−1)𝑇)},𝑡∈(𝑛−1)𝑇,𝜏_𝑛1. (A.1) At time 𝜏_𝑛1, one pesticide application occurs and 𝑦𝜏⁺_𝑛1=1−𝛿³₁𝑦(𝑛−1)𝑇⁺𝑒^{−𝑑𝑇/(𝑞+1)}. (A.2) Again, integrating the first equation of model (6) from 𝜏_𝑛1 to 𝜏_𝑛2 yields 𝑦(𝑡)=𝑦𝜏⁺_𝑛1𝑒^{−𝑑(𝑡−𝜏𝑛1)},𝑡∈𝜏_𝑛1,𝜏_𝑛2. (A.3) At time 𝜏_𝑛2, a single pesticide application occurs and 𝑦𝜏⁺_𝑛2=1−𝛿³₂𝑦𝜏⁺_𝑛1𝑒^{−𝑑𝑇/(𝑞+1)}=1−𝛿³₁1−𝛿³₂𝑦(𝑛−1)𝑇⁺𝑒^{−2𝑑𝑇/(𝑞+1)}. (A.4) By induction, we can see that 𝑦(𝑡)=1−𝛿³₁1−𝛿³₂⋯1−𝛿³_𝑞−1𝑒^{(−𝑑(𝑞−1)𝑇)/(𝑞+1)}×𝑦(𝑛−1)𝑇⁺𝑒^{−𝑑(𝑡−𝜏𝑛(𝑞−1))},𝑡∈𝜏_{𝑛(𝑞−1)},𝜏_𝑛𝑞. (A.5) At time 𝜏_𝑛𝑞, the last time pesticide is applied in the 𝑛th period and 𝑦𝜏⁺_𝑛𝑞=1−𝛿³₁⋯1−𝛿³_𝑞𝑦(𝑛−1)𝑇⁺𝑒^{−𝑑𝑞𝑇/(𝑞+1)}. (A.6) Finally, integrating the first equation of model (6) from 𝜏_𝑛𝑞 to 𝑛𝑇 yields 𝑦(𝑡)=𝑦𝜏⁺_𝑛𝑞𝑒^{−𝑑(𝑡−𝜏𝑛𝑞)}=1−𝛿³₁⋯1−𝛿³_𝑞×𝑦(𝑛−1)𝑇⁺𝑒^{−𝑑𝑞𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏𝑛𝑞)},𝑡∈𝜏_𝑛𝑞,𝑛𝑇. (A.7) At time 𝑛𝑇, release of natural enemies occurs once and 𝑦𝑛𝑇⁺=1−𝛿³₁⋯1−𝛿³_𝑞𝑦(𝑛−1)𝑇⁺𝑒^−𝑑𝑇+𝜃. (A.8) Denote 𝑌_𝑛=𝑦(𝑛𝑇⁺), then we have the following difference equation: 𝑌_𝑛+1=1−𝛿³₁⋯1−𝛿³_𝑞𝑌_𝑛𝑒^−𝑑𝑇+𝜃, (A.9) which has a unique steady state 𝑌^∗=𝜃1−1−𝛿³₁⋯1−𝛿³_𝑞𝑒^−𝑑𝑇. (A.10) Let 𝐹=(1−𝛿³₁)⋯(1−𝛿³_𝑞)𝑌𝑒^−𝑑𝑇+𝜃, since |𝑑𝐹/𝑑𝑌|=1−(1−𝛿³₁)⋯(1−𝛿³_𝑞)𝑒^−𝑑𝑇<1. Therefore, 𝑌^∗ is a globally asymptotically stable equilibrium of model (A.8), then system (6) has a globally stable 𝑇 periodic solution 𝑦_𝑇(𝑡), which can be calculated as follows:𝑦_𝑇(𝑡)=𝑌^∗𝑒^{−𝑑(𝑡−(𝑛−1)𝑇)},𝑡∈(𝑛−1)𝑇,𝜏_𝑛1,𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏𝑛1)},𝑡∈𝜏_𝑛1,𝜏_𝑛2,⋮𝑌^∗1−𝛿³_𝑞−1⋯1−𝛿³₁𝑒^{(−𝑑(𝑞−1)𝑇)/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏𝑛𝑞−1)},𝑡∈𝜏_𝑛𝑞−1,𝜏_𝑛𝑞,𝑌^∗1−𝛿³_𝑞⋯1−𝛿³₁𝑒^{−𝑑𝑞𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏𝑛𝑞)},𝑡∈𝜏_𝑛𝑞,𝑛𝑇. (A.11) We have the following.

Lemma A.1. System (6) has a positive periodic solution 𝑦_𝑇(𝑡) and for every solution 𝑦(𝑡) of system (6) one has 𝑦(𝑡)→𝑦_𝑇(𝑡) as 𝑡→∝.

Furthermore, we can obtain the complete expression for the prey-free periodic solution of system (2), (0,0,𝑦_𝑇(𝑡)), for 𝑡∈((𝑛−1)𝑇,𝑛𝑇]. Now, we give the conditions which assure the globally asymptotical stability of the pest-eradication (0,0,𝑦_𝑇(𝑡)).

Theorem A.2. Let (𝑥₁(𝑡),𝑥₂(𝑡),𝑦(𝑡)) be any solution of (2); then, (0,0,𝑦_𝑇(𝑡)) is globally asymptotically stable provided 1ln∏^𝑞_𝑘=11−𝛿^𝑖_𝑘>𝑎_𝑖𝛼𝑇+^𝑖𝐵₁𝑐𝑑,𝑖=1,2, (A.12) where 𝐵₁=ln𝑟+𝑐𝑌^∗𝑒^{−𝑑𝑇/(𝑞+1)}𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−2𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^−𝑑𝑇×𝑟+𝑐𝑌^∗𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)−1}. (A.13)

Proof. The local stability of periodic solution (0,0,𝑦_𝑇(𝑡)) may be determined by considering the behavior of small amplitude perturbations of the solution. Defining 𝑥₁(𝑡)=𝑈(𝑡),  𝑥₂(𝑡)=𝑉(𝑡), and  𝑦(𝑡)=𝑦_𝑇(𝑡)+𝑊(𝑡), it may be written 𝑈(𝑡)𝑉(𝑡)𝑊(𝑡)(𝑡)𝑈(0)𝑉(0)𝑊(0)=𝜙, (A.14) where 𝜙 satisfies 𝑑𝜙=𝑎𝑑𝑡¹−𝛼₁𝑦_𝑇(𝑡)𝑟+𝑐𝑦_𝑇(𝑡)000𝑎₂−𝛼₂𝑦_𝑇(𝑡)𝑟+𝑐𝑦_𝑇(𝑡)0𝜙(𝑡)00−𝑑, (A.15) and 𝜙(0)=𝐼, the identity matrix. Hence, the fundamental solution matrix is 𝜙(𝑡)=𝑒^{∫𝑡0(𝑎1−(𝛼1𝑦𝑇(𝑠)/(𝑟+𝑐𝑦𝑇(𝑠))))𝑑𝑠}000𝑒^{∫𝑡0(𝑎2−(𝛼2𝑦𝑇(𝑠)/(𝑟+𝑐𝑦𝑇(𝑠))))𝑑𝑠}000𝑒^−𝑑𝑡. (A.16) The impulsive conditions of system (2) become 𝑈𝜏⁺_𝑛𝑘𝑉𝜏⁺_𝑛𝑘𝑊𝜏⁺_𝑛𝑘=1−𝛿¹_𝑘0001−𝛿²_𝑘0001−𝛿³_𝑘𝑈𝜏_𝑛𝑘𝑉𝜏_𝑛𝑘𝑊𝜏_𝑛𝑘, (A.17) where 𝑘=1,2,…,𝑞, and 𝑈𝑛𝑇⁺𝑉𝑛𝑇⁺𝑊𝑛𝑇⁺=𝑈()𝑉()𝑊()100010001𝑛𝑇𝑛𝑇𝑛𝑇. (A.18) The stability of the periodic solution (0,0,𝑦_𝑇(𝑡)) is determined by the eigenvalues of the matrix 𝑀=^𝑞_𝑘=11−𝛿¹_𝑘0001−𝛿²_𝑘0001−𝛿³_𝑘𝜙(𝑇)100010001, (A.19) which are 𝜆₁=^𝑞_𝑘=11−𝛿¹_𝑘𝑒^{∫𝑇0(𝑎1−(𝛼1𝑦𝑇(𝑡)/(𝑟+𝑐𝑦𝑇(𝑡))))𝑑𝑡},𝜆₂=^𝑞_𝑘=11−𝛿²_𝑘𝑒^{∫𝑇0(𝑎2−(𝛼2𝑦𝑇(𝑡)/(𝑟+𝑐𝑦𝑇(𝑡))))𝑑𝑡},𝜆₃=^𝑞_𝑘=11−𝛿³_𝑘𝑒^−𝑑𝑇<1. (A.20) According to Floquet theory [16], (0,0,𝑦_𝑇(𝑡)) is locally stable if |𝜆₁|,|𝜆₂|<1.
Next, integrating the 𝑦_𝑇(𝑡)/(𝑟+𝑐𝑦_𝑇(𝑡)) from 0 to 𝑇 yields ^𝑇₀𝑦_𝑇(𝑡)𝑟+𝑐𝑦_𝑇(𝑡)=𝑑𝑡_𝜏110𝑌^∗𝑒^−𝑑𝑡𝑟+𝑐𝑌^∗𝑒^−𝑑𝑡+𝑑𝑡_{𝜏12𝜏+11}𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏+11)}𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏+11)}+𝑑𝑡+⋯_{𝜏1𝑞𝜏+1𝑞−1}𝑌^∗∏^𝑞−1_𝑘=11−𝛿³_𝑘𝑒^{(−𝑑(𝑞−1)𝑇)/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏+1𝑞−1)}𝑟+𝑐𝑌^∗∏^𝑞−1_𝑘=11−𝛿³_𝑘𝑒^{(−𝑑(𝑞−1)𝑇)/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏+1𝑞−1)}×𝑑𝑡+^𝑇_𝜏1𝑞𝑌^∗∏^𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏+1𝑞)}𝑟+𝑐𝑌^∗∏^𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)}𝑒^{−𝑑(𝑡−𝜏+1𝑞)}=𝑑𝑡−1𝑐𝑑ln𝑟+𝑐𝑌^∗𝑒^{−𝑑𝑇/(𝑞+1)}𝑟+𝑐𝑌^∗+−1𝑐𝑑ln𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−2𝑑𝑇/(𝑞+1)}𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}+−1𝑐𝑑ln𝑟+𝑐𝑌^∗∏^𝑞−1_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)}𝑟+𝑐𝑌^∗∏^𝑞−1_𝑘=11−𝛿³_𝑘𝑒^{(−𝑑(𝑞−1)𝑇)/(𝑞+1)}+−1𝑐𝑑ln𝑟+𝑐𝑌^∗∏^𝑞_𝑘=11−𝛿³_𝑘𝑒^−𝑑𝑇𝑟+𝑐𝑌^∗∏^𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)}. (A.21) Thus, we can obtain 𝜆₁=^𝑞_𝑘=11−𝛿¹_𝑘×exp𝑎₁𝛼𝑇+¹𝑐𝑑×ln𝑟+𝑐𝑌^∗𝑒^{−𝑑𝑇/(𝑞+1)}×𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−2𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^−𝑑𝑇×𝑟+𝑐𝑌^∗×𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)−1},𝜆₂=^𝑞_𝑘=11−𝛿²_𝑘×exp𝑎₂𝛼𝑇+²𝑐𝑑×ln𝑟+𝑐𝑌^∗𝑒^{−𝑑𝑇/(𝑞+1)}×𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−2𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^−𝑑𝑇×𝑟+𝑐𝑌^∗×𝑟+𝑐𝑌^∗1−𝛿³₁𝑒^{−𝑑𝑇/(𝑞+1)}⋯𝑟+𝑐𝑌^∗𝑞_𝑘=11−𝛿³_𝑘𝑒^{−𝑑𝑞𝑇/(𝑞+1)−1}. (A.22)(A.23) Clearly, the condition (A.12) of Theorem A.2 may be obtained, if we set |𝜆₁|,|𝜆₂|<1. This completes the proof of local stability of periodic solution (0,0,𝑦_𝑇(𝑡)).
We now need to prove the global attractiveness. Choose an 𝜀>0 such that 𝜉_𝑖≜^𝑞_𝑘=11−𝛿^𝑖_𝑘𝑒^{∫𝑇0(𝑎𝑖−(𝛼𝑖(𝑦𝑇(𝑡)−𝜀)/(𝑟+𝑐(𝑦𝑇(𝑡)−𝜀))))𝑑𝑡}<1(𝑖=1,2). (A.24) Note that 𝑑𝑦(𝑡)/𝑑𝑡≥−𝑑𝑦(𝑡); from Lemma A.1 and comparison theorem of impulsive equation, we have 𝑦(𝑡)>𝑦_𝑇(𝑡)−𝜀 (A.25) for all sufficiently large 𝑡. For simplification, we assume that the inequality (A.25) holds for all 𝑡≥0. Consider the following impulse differential inequalities: 𝑑𝑥_𝑖≤𝑎𝑑𝑡^𝑖−𝛼_𝑖𝑦_𝑇(𝑡)−𝜀𝑟+𝑐𝑦_𝑇(𝑡)𝑥−𝜀^𝑖(𝑡),𝑡≠𝜏_𝑛𝑘,𝑥_𝑖𝜏⁺_𝑛𝑘=1−𝛿^𝑖_𝑘𝑥₁𝜏_𝑛𝑘,𝑡=𝜏_𝑛𝑘, (A.26) where 𝑘=1,2,…,𝑞 and 𝑖=1,2. By using Lemma 2, we have 𝑥_𝑖(𝑡)≤𝑥_𝑖((𝑛−1)𝑇)×_{()𝑛−1𝑇<𝜏𝑛𝑘<𝑡}1−𝛿^𝑖_𝑘𝑒^{∫𝑡(𝑛−1)𝑇𝑒𝑎𝑖𝑖𝑇𝑇−(𝛼(𝑦(𝑠)−𝜀)/(𝑟+𝑐(𝑦(𝑠)−𝜀)))𝑑𝑠}; (A.27) then, 𝑥_𝑖(𝑛𝑇)≤𝑥_𝑖((𝑛−1)𝑇)×_{()𝑛−1𝑇<𝜏𝑛𝑘<𝑛𝑇}1−𝛿^𝑖_𝑘𝑒^{∫𝑛𝑇(𝑛−1)𝑇𝑒𝑎𝑖𝑖𝑦𝑇𝑇−(𝛼(𝑠)/(𝑟+𝑐𝑦(𝑠)))𝑑𝑠}=𝑥_𝑖((𝑛−1)𝑇)𝜉_𝑖. (A.28) Hence, 𝑥_𝑖(𝑛𝑇)≤𝑥_𝑖(0⁺)𝜉^𝑛_𝑖 and 𝑥_𝑖(𝑛𝑇)→0 as 𝑛→∞. Therefore, 𝑥_𝑖(𝑡)→0 as 𝑡→∞, since 0<𝑥_𝑖(𝑡)<𝑥_𝑖((𝑛−1)𝑇)∏^𝑞_𝑘=1(1−𝛿^𝑖_𝑘)𝑒^𝑎𝑖𝑇 for (𝑛−1)𝑇<𝑡<𝑛𝑇.
Next, we prove that 𝑦(𝑡)→𝑦_𝑇(𝑡) as 𝑡→∞. For 𝜀>0 sufficiently small, there exists a 𝑇^′ such that 0<𝑥₁(𝑡)<𝜀 and 0<𝑥₂(𝑡)<𝜀,𝑡≥𝑇^′. Without any loss of generality, we assume that 0<𝑥₁(𝑡)<𝜀 and 0<𝑥₂(𝑡)<𝜀 for all 𝑡≥0. Then, from the system (2) we obtain (𝑡)≤(𝑡)−𝑑𝑦𝑑𝑦≤𝑦(𝑡)𝜔𝑑𝑡−𝑑+𝛾𝜀,𝛾=¹𝛼₁𝑟+𝜔₂𝛼₂𝑟. (A.29) From Lemmas 2 and A.1, we have 𝑢(𝑡)≤𝑦(𝑡)≤𝑣(𝑡) and 𝑢(𝑡)→𝑦_𝑇(𝑡),𝑣(𝑡)→𝑦_𝑇(𝑡) as 𝑡→∝, where 𝑢(𝑡) and 𝑣(𝑡) are solutions of the equations 𝑑𝑢(𝑡)(𝑡)𝑑𝑡=−𝑑𝑢,𝑡≠𝑛𝑇,𝑡≠𝜏_𝑛𝑘,𝑢𝑡⁺=1−𝛿³_𝑘𝑢(𝑡),𝑡=𝜏_𝑛𝑘𝑢𝑡,𝑘=1,2,…,𝑞,₊(𝑡)𝑢0=𝑢+𝜃,𝑡=𝑛𝑇,₊=𝑦₀,𝑑𝑣(𝑡)=𝑣(𝑡)𝑑𝑡−𝑑+𝛾𝜀,𝑡≠𝑛𝑇,𝑡≠𝜏_𝑛𝑘,𝑣𝑡⁺=1−𝛿³_𝑘𝑣(𝑡),𝑡=𝜏_𝑛𝑘𝑣𝑡,𝑘=1,2,…,𝑞,₊(𝑡)𝑣0=𝑣+𝜃,𝑡=𝑛𝑇,₊=𝑦₀, (A.30) respectively.
We also have 𝑣_𝑇(𝑡)=𝑌^∗𝑒^{(−𝑑+𝛾𝜀)(𝑡−(𝑛−1)𝑇)},𝑡∈(𝑛−1)𝑇,𝜏_𝑛1,𝑌^∗1−𝛿³₁𝑒^{(−𝑑+𝛾𝜀)𝑇/(𝑞+1)}𝑒^{(−𝑑+𝛾𝜀)(𝑡−𝜏𝑛1)},𝑡∈𝜏_𝑛1,𝜏_𝑛2,⋮𝑌^∗1−𝛿³_𝑞−1⋯1−𝛿³₁𝑒^{((−𝑑+𝛾𝜀)(𝑞−1)𝑇)/(𝑞+1)}×𝑒^{(−𝑑+𝛾𝜀)(𝑡−𝜏𝑛𝑞−1)},𝑡∈𝜏_{𝑛(𝑞−1)},𝜏_𝑛𝑞,𝑌^∗1−𝛿³_𝑞⋯1−𝛿³₁𝑒^{((−𝑑+𝛾𝜀)𝑞𝑇)/(𝑞+1)}×𝑒^{(−𝑑+𝛾𝜀)(𝑡−𝜏𝑛𝑞)},𝑡∈𝜏_𝑛𝑞,𝑛𝑇. (A.31) Then, for any 𝜀₁ there exists a 𝑇₁ such that 𝑢_𝑇(𝑡)−𝜀₁<𝑦(𝑡)<𝑣_𝑇(𝑡)+𝜀₁,𝑡>𝑇₁. (A.32) Let 𝜀→0, we have 𝑦_𝑇(𝑡)−𝜀<𝑦(𝑡)<𝑦_𝑇(𝑡)+𝜀 (A.33) for 𝑡 large enough, which implies that 𝑦(𝑡)→𝑦_𝑇(𝑡) as 𝑡→∝. This completes the proof.

Acknowledgments

This work is supported by the National Natural Science Foundation of China (10971001) and Excellent Talents Support Project of Universities and Colleges in Liaoning.