# Dynamics of a Discrete Host-Parasitoid System with Stocking

Authors: Sophia R.-J. Jang
Journal: Discrete Dynamics in Nature and Society
 (2015)

---

## Abstract

Motivated by the biological control of
pests, we present discrete-time models of host-parasitoid
interactions to study the effects of external stocking upon the
systems. It is assumed that density dependence of the hosts occurs
first followed by parasitism. We prove that the constant stocking
can eliminate the pest population if the stocking is sufficiently
large. Furthermore, stocking can  simplify the dynamics of the
interaction by stabilizing the coexisting steady state.


---

## Body

## 1. Introduction

There are many natural host populations that are also pests. Biological control is the reduction of pest populations by natural enemies, also known as the biological control agents. Many species of wasps and some flies are parasitoids and most of the parasitoids have a narrow host range which can be used as biological control agents. Biological controls often involve supplemental release of natural enemies. Relatively few natural enemies may be released at a critical time of the season (inoculative release) or literally millions may be released in a single time (inundative release) [1].Several discrete-time mathematical models have been proposed to study the dynamical effects of external stocking or the inoculative release of the control agents. See, for example, AlSharawi and Rhouma [2], Chow and Jang [3], Elaydi and Yakubu [4], Jang and Yu [5], Kulenović and Nurkanović [6], and references cited therein. The classical Leslie-Gower competition model with stocking in one of the two populations is analyzed in [6] and the Ricker type competition system with stocking occurring in one of the two competing population is studied in [4]. The work in [2] investigates a multispecies population model with constant harvesting/stocking and a model of three interacting populations with stocking in one of the two competing populations is analyzed in [3]. The theory of optimal control is applied to study a host-parasitoid model in [5] where the constant stocking of parasitoids is used as a control strategy for the hosts.Insect populations frequently suffer from some density dependent effect in addition to mortality from insect parasitoids [7]. The ordering of density dependence and parasitism in the host life cycle can have significant impacts on the dynamics of the interactions [7]. In this investigation, we propose a discrete-time host-parasitoid model to study the effect of external stocking upon the host-parasitoid interaction. Unlike the work in [5], where parasitism occurs first followed by density dependence, it is assumed in the present study that density dependence of the hosts acts first followed by parasitism. Furthermore, the hosts are also pests and the parasitoids are used as biological control agents to control the pests. There is a constant level of the external parasitoids released into the interaction at each generation. It is shown that external stocking of the parasitoids can eliminate the hosts population and can also simplify the dynamics of the interaction.In the following section, the host-parasitoid model with no stocking is proposed and analyzed first. The corresponding model with stocking is then presented and studied. We derive sufficient conditions for the existence of interior steady state and discuss its magnitude relative to the model parameters for both systems. Numerical examples using MatLab are provided to illustrate our findings. The final section provides a brief summary.


## 2. The Models and Stability Analysis

In this section, we present two general host-parasitoid models, one with stocking and the other with no stocking.


### 2.1. The Model of No Stocking

Letx(t) and y(t) denote the host and parasitoid populations in generation t=0,1,…, respectively. The parasitoid has a very narrow range of hosts and is specialized to this particular host population. It is assumed that density dependence of the hosts occurs first followed by parasitism. Moreover, the number of encounter between hosts and parasitoids is distributed randomly and the probability of an individual host escaping from being parasitized is modeled by the zero term of a Poisson distribution (cf. [7]). The parameter a>0 denotes the average number of encounters per unit time per parasitoid and is also referred to as the searching efficiency of the parasitoid. The host-parasitoid interaction without external stocking of the parasitoids is described by the following system:(1)xt+1=λxtgxte-ayt,yt+1=βxtgxt1-e-ayt,where λ>0 is the intrinsic growth rate of the hosts and β>0 is the average number of parasitoids produced by each parasitized host. The per capita growth rate g of the host population satisfies the following assumptions; namely,(H1)
g∈C2[0,∞), g(0)=1, g(x)>0, g′(x)<0 for x≥0, limx→∞g⁡(x)=0, and sup⁡{xg(x):x≥0}=l<∞.It is clear that the classical Beverton-Holt and Ricker type growth rates [8] satisfy (H1). Furthermore, the growth rate of the Beverton-Holt model has the following monotone properties:(2)xgx′>0∀x≥0,(3)xgx′′<0∀x≥0.By defining new state variable y^=ay and new parameter β^=βa and by ignoring the hats, system (1) is converted into(4)xt+1=λxtgxte-yt,yt+1=βxtgxt1-e-ytwith nonnegative initial conditions.We first study system (4). Notice that system (4) always has a trivial steady state E00=(0,0), where both populations are extinct. Moreover, (4) has another boundary steady state E01=(x¯,0), if λ>1, where x¯ solves g(x¯)=1/λ. The Jacobian matrices evaluated at these two steady states E00 and E01 are given by(5)JE00=λ000,JE01=1+x¯g′x¯-ax¯0βx¯gx¯,respectively, where 1+x¯g′(x¯)<1 by (H1). In addition, if (2) is satisfied, then 1+x¯g′(x¯)=xgx′x=x¯>0 and thus E10 is asymptotically stable if βx¯g(x¯)<1.It is straightforward to show thatE00 is globally asymptotically stable for (4) if λ<1. Furthermore, E01 is globally asymptotically stable in {(x,y)∈R+2:x>0} if λ>1, βx¯g(x¯)<1, and (2) holds. Indeed, (2) implies that the map induced by z(t+1)=λz(t)g(z(t)) is strictly increasing and so that x¯ is globally asymptotically stable on (0,∞) for this scalar equation. It follows that limsupt→∞⁡xt≤x¯ for all solutions of (4). Using the assumption of βx¯g(x¯)<1, one can then prove that limt→∞⁡y(t)=0. Consequently, E10 is globally asymptotically stable in {(x,y)∈R+2:x>0}.Theorem 1.
Solutions of (4) remain nonnegative and are bounded for t>0. The following statements hold for (4).(a)
Ifλ<1, then E00=(0,0) is globally asymptotically stable in R+2.(b)
Ifλ>1, then E00 is unstable and E10=(x¯,0) exists. In addition, if (2) is satisfied, then E10 is globally asymptotically stable in {(x,y)∈R+2:x>0} if βx¯g(x¯)<1.If the monotone property (2) is not assumed, then dynamics of the host population in the absence of the parasitoids may be complicated. For example, the Ricker model undergoes a cascade of period-doubling bifurcations to chaos as the intrinsic growth rate of the population increases [8]. In such a case, dynamics of the host-parasitoid model (4) may not equilibrate even when βx¯g(x¯)<1.Letλ>1. Note that the x component of an interior steady state satisfies(6)H0x≔ln⁡λgx-βxgx+βxλ=0.Let L(x)=ln⁡(λg(x)) and R0(x)=βxg(x)-βx/λ. Then(7)L0=ln⁡λ>0,Lx¯=0,L′x=g′xgx<0,R00=0,R0x¯=0,R0′x=βxgx′-βλ.Therefore, H0(0)=ln⁡λ>0, H0(x¯)=0, H0′(0)=g′(0)+β(1/λ-1)<0, and H0′(x¯)=g′(x¯)/g(x¯)(1-βx¯g(x¯)). Consequently, if H0′(x¯)>0, or equivalently, if βx¯g(x¯)>1, then H0(x)=0 has at least one solution x∗ in (0,x¯). Therefore, system (4) has at least one interior steady state if βx¯g(x¯)>1. Next, R0′′(x)=βxgx′′. If (3) holds, then R0(x) is concave down and we can conclude that the positive solution x∗ of L(x)=R0(x) is unique if βx¯g(x¯)>1. Consequently, when (2) and (3) are satisfied, then, by Theorem 1, (4) has a unique interior steady state E∗0=(x∗,y∗) if and only if βx¯g(x¯)>1, where 0<x∗<x¯ and 0<y∗<ln⁡λ. We summarize the discussion as the following theorem.Theorem 2.
Ifλ>1 and βx¯g(x¯)>1, then (4) has at least one interior steady state E∗0=(x∗,y∗), where 0<x∗<x¯ and 0<y∗<ln⁡λ. If, in addition, (2) and (3) hold, then (4) has a unique interior steady state E∗0=(x∗,y∗) if and only if βx¯g(x¯)>1.Letλ>1, βx¯g(x¯)>1, and let E∗0=(x∗,y∗) denote an interior steady state of the model (4). Its magnitude relative to the parameter β can be obtained via(8)dx∗dβ=x∗g2x∗1-e-y∗g′x∗1-βx∗gx∗e-y∗-βgx∗1-e-y∗gx∗+x∗g′x∗,dy∗dβ=g′x∗gx∗·dx∗dβ.Since(9)βx∗gx∗=y∗e-y∗1-e-y∗=y∗ey∗-1<1,it follows from (8) that dx∗/dβ<0 and dy∗/dβ>0 if (2) holds.The stability ofE∗0 depends on the Jacobian matrix of (4) evaluated at E∗0 given by(10)JE∗0=λgx+xg′xe-y-λxβgx+xgx1-e-yβxgxe-y,where x=x∗ and y=y∗. For simplicity, (10) can be rewritten as(11)abcd=J∗0.According to the Jury conditions [8], E∗0 is asymptotically stable if tr⁡J∗0<1+det⁡J∗0<2. If (2) is assumed, then a>0 and hence tr⁡J∗0>0. Moreover, tr⁡J∗0=tr⁡J∗0<1+det⁡J∗0 is equivalent to (1-a)(1-d)-bc>0, which holds trivially since a<1, d<1, and bc<0. Therefore, E0∗ is asymptotically stable if(12)det⁡J∗0=xxgx′x=x∗<1.Differentiatingdet⁡J∗0 with respect to β, we have(13)ddet⁡J∗0dβ=x+βx′xgx′+βxx′xgx′′x=x∗,where x′=dx∗/dβ.Observe thatR0′(x)<0 and det⁡J∗0<1 are equivalent to xgx′<1/λ and xgx′<1/βx, respectively. Moreover, 1/λ=g(x)e-y and 1/βx=g(x)(1-e-y)/y hold at any interior steady state (x,y). A direct computation shows that g(x)e-y<g(x)(1-e-y)/y is equivalent to y<ey-1, which holds trivially. Therefore, if the interior steady state is at the intersection point of L(x)=R0(x) for which R0(x) is decreasing, then such a steady state is asymptotically stable. However, the converse is not true. Indeed, since 1/λ<1/βx at any interior steady state (x,y), it is possible that det⁡J∗0<1, while R0′(x)>0.Proposition 3.
Letλ>1 and let βx¯g(x¯)>1 and assume (2). Then an interior steady state E∗0=(x∗,y∗) of model (4) is asymptotically stable if (12) holds. In particular, E∗0 is asymptotically stable if R0′(x∗)<0. Moreover, x∗ and y∗ are decreasing and increasing functions of β, respectively.If the per capita growth rateg of the hosts follows the Beverton-Holt model, g(x)=1/1+kx, k>0, then x¯=λ-1/k, and (2) and (3) hold. Moreover,(14)det⁡J∗0=βx1+kx2,ddet⁡J∗0dx=β1-kx1+kx3,ddet⁡J∗0dβ=x1+kx+β1-kxx′1+kx3.The maximum of det⁡J∗0 occurs at x=xm≔1/k with max⁡(det⁡J∗0)=β/4k. Let β0=kλ/λ-1 and notice βx¯g(x¯)>1 is equivalent to β>β0. Therefore, if β0<β<4k, then E∗0 is asymptotically stable, where we implicitly assumed that λ>4/3. On the other hand, if either λ≤4/3 or λ>4/3 and β>4k hold, then max⁡(det⁡J∗0)>1 and E∗0 is unstable for some parameter β region. It is expected that a Neimark-Sacker bifurcation occurs when det⁡J∗0=1 [9, 10].We use a numerical example withλ=2 and k=1 to illustrate this observation. Then β0=2 and the unique interior steady state E∗0 exists if β>β0=2. Moreover, E∗0 is asymptotically stable if β<4k=4. Figures 1(a) and 1(b) plot bifurcation diagrams of model (4) using β as the bifurcation parameter. We see that E∗0 is asymptotically stable if β<6. Furthermore, the x component of E∗0 is decreasing, while the y component is increasing as the parameter β increases.Figure 1
Bifurcation diagrams for system (4) with a Beverton-Holt growth rate. Parameter values are λ=2 and k=1.
(a)
(b)


### 2.2. The Model with Stocking

Since the host population is regarded as a pest and the biological control of inoculative release of parasitoids is implemented, the model of host-parasitoid interaction is now given by(15)xt+1=λxtgxte-yt,yt+1=βxtgxt1-e-yt+u,where g satisfies (H1) and u>0 denotes the constant external stocking of the parasitoids. Notice that, similar to model (4), the parameter u has been rescaled with u^=au and the hat is dropped.System (15) always has a unique boundary steady state E0=(0,u), where the host population is extinct and the parasitoid population is stabilized at the constant stocking level u. The stability of E0 depends on the Jacobian matrix at E0,(16)JE0=λe-u0β1-e-u0.Since y(t)≥u for t≥1 implies x(t+1)≤λx(t)g(x(t))e-u≤λx(t)e-u for t≥1, we see that E0 is globally asymptotically stable in R+2 if λe-u<1.Letλe-u>1. Then E0 is a saddle point with its stable manifold lying on the nonnegative y-axis. The x component of an interior steady state of (15) can be shown to satisfy(17)Hx≔ln⁡λgx-βxgx+βxλ-u=0.Let L(x)=ln⁡(λg(x)) be defined as in model (4) and let R(x)=βxg(x)-βx/λ+u. Then(18)R0=u,Rx¯=u,R′x=βxgx′-βλ,R′′x=βxgx′′.It follows that H(0)=ln⁡λ-u>0 and H(x¯)=-u<0. Therefore, H(x)=0 has at least one positive solution x∗ in (0,x¯) and thus (15) has at least one interior steady state. Similar to model (4), if (3) is assumed, then the interior steady state is unique. Notice that the existence of an interior steady state only depends on λ and u but not on β.Theorem 4.
The following statements hold for model (15).(a)
Ifλe-u<1, then E0=(0,u) is globally asymptotically stable in R+2.(b)
Ifλe-u>1, then E0 is unstable and (15) has at least one interior steady state E∗=(x∗,y∗), where 0<x∗<x¯ and u<y∗<ln⁡λ. If in addition (3) is satisfied, then the interior steady state E∗=(x∗,y∗) is unique.The Jacobin matrix of (15) evaluated at an interior steady state E∗, J∗, is the same as that of (10) with x=x∗ and y=y∗. Observe that(19)βx∗gx∗e-y∗=e-y∗y∗-u1-e-y∗<1.If (2) is assumed, then tr⁡J∗=tr⁡J∗<1+det⁡J∗ holds and E∗ is locally asymptotically stable if det⁡J∗=βxxgx′x=x∗<1. The magnitude of the interior steady state with respect to the constant stocking u can be shown to satisfy(20)dx∗du=gx∗g′x∗1-βx∗gx∗e-y∗-βgx∗1-e-y∗gx∗+x∗g′x∗,dy∗du=g′x∗gx∗·dx∗du.Consequently, dx∗/du<0 and dy∗/du>0 if (2) is satisfied.Observe thatR′(x)<0 and det⁡J∗<1 are equivalent to xgx′<1/λ and xgx′<1/βx, respectively. At any interior steady state (x,y), we have 1/λ=g(x)e-y and 1/βx=g(x)(1-e-y)/y-u. It can be verified that 1/λ<1/βx for any interior steady state (x,y). We therefore conclude that R′(x∗)<0 implies det⁡J∗<1.Proposition 5.
Letλe-u>1. If (2) is satisfied, then an interior steady state E∗=(x∗,y∗) of model (15) is asymptotically stable if (12) holds. In particular, E∗ is asymptotically stable if R′(x∗)<0. Moreover, x∗ and y∗ are decreasing and increasing functions of u, respectively.Letλe-u>1. For the Beverton-Holt growth rate g(x)=1/1+kx, k>0, system (15) has a unique interior steady state E∗. We also have(21)det⁡J∗=βx1+kx2,ddet⁡J∗dx=β1-kx1+kx3,ddet⁡J∗du=β1-kx1+kx3·dxdu.The maximum of det⁡J∗ occurs at x=xm≔1/k with max⁡(det⁡J∗)=β/4k. In particular, if β<4k, then E∗ is asymptotically stable. Therefore, stocking has no effect on the stability of the interior steady state if β<4k. If β>4k and λ≤2, then 1/k≥x¯ and thus ddet⁡J∗/dx>0 on (0,x¯). It follows that ddet⁡J∗/du<0 and stocking can stabilize the interior steady state E∗. If β>4k and λ>2, then max⁡(det⁡J∗)>1. Since x∗ is a decreasing function of u, it is expected that stocking can stabilize the interior steady state by (21).Consider a numerical example withλ=2, k=1, and β=7>4k. Then the unique interior steady state E∗ exists if u<ln⁡λ≈0.6931. From Figure 1, it is known that the interior steady state of model (4) is unstable when there is no stocking. We then use u as the bifurcation parameter. Figure 2 illustrates that stocking can stabilize the host-parasitoid interaction. Moreover, x and y components of E∗ are clearly decreasing and increasing functions of u, respectively.Figure 2
Bifurcation diagrams for system (15) with a Beverton-Holt growth rate. Parameter values are λ=2, k=1, and β=7.
(a)
(b)


## 2.1. The Model of No Stocking

Letx(t) and y(t) denote the host and parasitoid populations in generation t=0,1,…, respectively. The parasitoid has a very narrow range of hosts and is specialized to this particular host population. It is assumed that density dependence of the hosts occurs first followed by parasitism. Moreover, the number of encounter between hosts and parasitoids is distributed randomly and the probability of an individual host escaping from being parasitized is modeled by the zero term of a Poisson distribution (cf. [7]). The parameter a>0 denotes the average number of encounters per unit time per parasitoid and is also referred to as the searching efficiency of the parasitoid. The host-parasitoid interaction without external stocking of the parasitoids is described by the following system:(1)xt+1=λxtgxte-ayt,yt+1=βxtgxt1-e-ayt,where λ>0 is the intrinsic growth rate of the hosts and β>0 is the average number of parasitoids produced by each parasitized host. The per capita growth rate g of the host population satisfies the following assumptions; namely,(H1)
g∈C2[0,∞), g(0)=1, g(x)>0, g′(x)<0 for x≥0, limx→∞g⁡(x)=0, and sup⁡{xg(x):x≥0}=l<∞.It is clear that the classical Beverton-Holt and Ricker type growth rates [8] satisfy (H1). Furthermore, the growth rate of the Beverton-Holt model has the following monotone properties:(2)xgx′>0∀x≥0,(3)xgx′′<0∀x≥0.By defining new state variable y^=ay and new parameter β^=βa and by ignoring the hats, system (1) is converted into(4)xt+1=λxtgxte-yt,yt+1=βxtgxt1-e-ytwith nonnegative initial conditions.We first study system (4). Notice that system (4) always has a trivial steady state E00=(0,0), where both populations are extinct. Moreover, (4) has another boundary steady state E01=(x¯,0), if λ>1, where x¯ solves g(x¯)=1/λ. The Jacobian matrices evaluated at these two steady states E00 and E01 are given by(5)JE00=λ000,JE01=1+x¯g′x¯-ax¯0βx¯gx¯,respectively, where 1+x¯g′(x¯)<1 by (H1). In addition, if (2) is satisfied, then 1+x¯g′(x¯)=xgx′x=x¯>0 and thus E10 is asymptotically stable if βx¯g(x¯)<1.It is straightforward to show thatE00 is globally asymptotically stable for (4) if λ<1. Furthermore, E01 is globally asymptotically stable in {(x,y)∈R+2:x>0} if λ>1, βx¯g(x¯)<1, and (2) holds. Indeed, (2) implies that the map induced by z(t+1)=λz(t)g(z(t)) is strictly increasing and so that x¯ is globally asymptotically stable on (0,∞) for this scalar equation. It follows that limsupt→∞⁡xt≤x¯ for all solutions of (4). Using the assumption of βx¯g(x¯)<1, one can then prove that limt→∞⁡y(t)=0. Consequently, E10 is globally asymptotically stable in {(x,y)∈R+2:x>0}.Theorem 1.
Solutions of (4) remain nonnegative and are bounded for t>0. The following statements hold for (4).(a)
Ifλ<1, then E00=(0,0) is globally asymptotically stable in R+2.(b)
Ifλ>1, then E00 is unstable and E10=(x¯,0) exists. In addition, if (2) is satisfied, then E10 is globally asymptotically stable in {(x,y)∈R+2:x>0} if βx¯g(x¯)<1.If the monotone property (2) is not assumed, then dynamics of the host population in the absence of the parasitoids may be complicated. For example, the Ricker model undergoes a cascade of period-doubling bifurcations to chaos as the intrinsic growth rate of the population increases [8]. In such a case, dynamics of the host-parasitoid model (4) may not equilibrate even when βx¯g(x¯)<1.Letλ>1. Note that the x component of an interior steady state satisfies(6)H0x≔ln⁡λgx-βxgx+βxλ=0.Let L(x)=ln⁡(λg(x)) and R0(x)=βxg(x)-βx/λ. Then(7)L0=ln⁡λ>0,Lx¯=0,L′x=g′xgx<0,R00=0,R0x¯=0,R0′x=βxgx′-βλ.Therefore, H0(0)=ln⁡λ>0, H0(x¯)=0, H0′(0)=g′(0)+β(1/λ-1)<0, and H0′(x¯)=g′(x¯)/g(x¯)(1-βx¯g(x¯)). Consequently, if H0′(x¯)>0, or equivalently, if βx¯g(x¯)>1, then H0(x)=0 has at least one solution x∗ in (0,x¯). Therefore, system (4) has at least one interior steady state if βx¯g(x¯)>1. Next, R0′′(x)=βxgx′′. If (3) holds, then R0(x) is concave down and we can conclude that the positive solution x∗ of L(x)=R0(x) is unique if βx¯g(x¯)>1. Consequently, when (2) and (3) are satisfied, then, by Theorem 1, (4) has a unique interior steady state E∗0=(x∗,y∗) if and only if βx¯g(x¯)>1, where 0<x∗<x¯ and 0<y∗<ln⁡λ. We summarize the discussion as the following theorem.Theorem 2.
Ifλ>1 and βx¯g(x¯)>1, then (4) has at least one interior steady state E∗0=(x∗,y∗), where 0<x∗<x¯ and 0<y∗<ln⁡λ. If, in addition, (2) and (3) hold, then (4) has a unique interior steady state E∗0=(x∗,y∗) if and only if βx¯g(x¯)>1.Letλ>1, βx¯g(x¯)>1, and let E∗0=(x∗,y∗) denote an interior steady state of the model (4). Its magnitude relative to the parameter β can be obtained via(8)dx∗dβ=x∗g2x∗1-e-y∗g′x∗1-βx∗gx∗e-y∗-βgx∗1-e-y∗gx∗+x∗g′x∗,dy∗dβ=g′x∗gx∗·dx∗dβ.Since(9)βx∗gx∗=y∗e-y∗1-e-y∗=y∗ey∗-1<1,it follows from (8) that dx∗/dβ<0 and dy∗/dβ>0 if (2) holds.The stability ofE∗0 depends on the Jacobian matrix of (4) evaluated at E∗0 given by(10)JE∗0=λgx+xg′xe-y-λxβgx+xgx1-e-yβxgxe-y,where x=x∗ and y=y∗. For simplicity, (10) can be rewritten as(11)abcd=J∗0.According to the Jury conditions [8], E∗0 is asymptotically stable if tr⁡J∗0<1+det⁡J∗0<2. If (2) is assumed, then a>0 and hence tr⁡J∗0>0. Moreover, tr⁡J∗0=tr⁡J∗0<1+det⁡J∗0 is equivalent to (1-a)(1-d)-bc>0, which holds trivially since a<1, d<1, and bc<0. Therefore, E0∗ is asymptotically stable if(12)det⁡J∗0=xxgx′x=x∗<1.Differentiatingdet⁡J∗0 with respect to β, we have(13)ddet⁡J∗0dβ=x+βx′xgx′+βxx′xgx′′x=x∗,where x′=dx∗/dβ.Observe thatR0′(x)<0 and det⁡J∗0<1 are equivalent to xgx′<1/λ and xgx′<1/βx, respectively. Moreover, 1/λ=g(x)e-y and 1/βx=g(x)(1-e-y)/y hold at any interior steady state (x,y). A direct computation shows that g(x)e-y<g(x)(1-e-y)/y is equivalent to y<ey-1, which holds trivially. Therefore, if the interior steady state is at the intersection point of L(x)=R0(x) for which R0(x) is decreasing, then such a steady state is asymptotically stable. However, the converse is not true. Indeed, since 1/λ<1/βx at any interior steady state (x,y), it is possible that det⁡J∗0<1, while R0′(x)>0.Proposition 3.
Letλ>1 and let βx¯g(x¯)>1 and assume (2). Then an interior steady state E∗0=(x∗,y∗) of model (4) is asymptotically stable if (12) holds. In particular, E∗0 is asymptotically stable if R0′(x∗)<0. Moreover, x∗ and y∗ are decreasing and increasing functions of β, respectively.If the per capita growth rateg of the hosts follows the Beverton-Holt model, g(x)=1/1+kx, k>0, then x¯=λ-1/k, and (2) and (3) hold. Moreover,(14)det⁡J∗0=βx1+kx2,ddet⁡J∗0dx=β1-kx1+kx3,ddet⁡J∗0dβ=x1+kx+β1-kxx′1+kx3.The maximum of det⁡J∗0 occurs at x=xm≔1/k with max⁡(det⁡J∗0)=β/4k. Let β0=kλ/λ-1 and notice βx¯g(x¯)>1 is equivalent to β>β0. Therefore, if β0<β<4k, then E∗0 is asymptotically stable, where we implicitly assumed that λ>4/3. On the other hand, if either λ≤4/3 or λ>4/3 and β>4k hold, then max⁡(det⁡J∗0)>1 and E∗0 is unstable for some parameter β region. It is expected that a Neimark-Sacker bifurcation occurs when det⁡J∗0=1 [9, 10].We use a numerical example withλ=2 and k=1 to illustrate this observation. Then β0=2 and the unique interior steady state E∗0 exists if β>β0=2. Moreover, E∗0 is asymptotically stable if β<4k=4. Figures 1(a) and 1(b) plot bifurcation diagrams of model (4) using β as the bifurcation parameter. We see that E∗0 is asymptotically stable if β<6. Furthermore, the x component of E∗0 is decreasing, while the y component is increasing as the parameter β increases.Figure 1
Bifurcation diagrams for system (4) with a Beverton-Holt growth rate. Parameter values are λ=2 and k=1.
(a)
(b)


## 2.2. The Model with Stocking

Since the host population is regarded as a pest and the biological control of inoculative release of parasitoids is implemented, the model of host-parasitoid interaction is now given by(15)xt+1=λxtgxte-yt,yt+1=βxtgxt1-e-yt+u,where g satisfies (H1) and u>0 denotes the constant external stocking of the parasitoids. Notice that, similar to model (4), the parameter u has been rescaled with u^=au and the hat is dropped.System (15) always has a unique boundary steady state E0=(0,u), where the host population is extinct and the parasitoid population is stabilized at the constant stocking level u. The stability of E0 depends on the Jacobian matrix at E0,(16)JE0=λe-u0β1-e-u0.Since y(t)≥u for t≥1 implies x(t+1)≤λx(t)g(x(t))e-u≤λx(t)e-u for t≥1, we see that E0 is globally asymptotically stable in R+2 if λe-u<1.Letλe-u>1. Then E0 is a saddle point with its stable manifold lying on the nonnegative y-axis. The x component of an interior steady state of (15) can be shown to satisfy(17)Hx≔ln⁡λgx-βxgx+βxλ-u=0.Let L(x)=ln⁡(λg(x)) be defined as in model (4) and let R(x)=βxg(x)-βx/λ+u. Then(18)R0=u,Rx¯=u,R′x=βxgx′-βλ,R′′x=βxgx′′.It follows that H(0)=ln⁡λ-u>0 and H(x¯)=-u<0. Therefore, H(x)=0 has at least one positive solution x∗ in (0,x¯) and thus (15) has at least one interior steady state. Similar to model (4), if (3) is assumed, then the interior steady state is unique. Notice that the existence of an interior steady state only depends on λ and u but not on β.Theorem 4.
The following statements hold for model (15).(a)
Ifλe-u<1, then E0=(0,u) is globally asymptotically stable in R+2.(b)
Ifλe-u>1, then E0 is unstable and (15) has at least one interior steady state E∗=(x∗,y∗), where 0<x∗<x¯ and u<y∗<ln⁡λ. If in addition (3) is satisfied, then the interior steady state E∗=(x∗,y∗) is unique.The Jacobin matrix of (15) evaluated at an interior steady state E∗, J∗, is the same as that of (10) with x=x∗ and y=y∗. Observe that(19)βx∗gx∗e-y∗=e-y∗y∗-u1-e-y∗<1.If (2) is assumed, then tr⁡J∗=tr⁡J∗<1+det⁡J∗ holds and E∗ is locally asymptotically stable if det⁡J∗=βxxgx′x=x∗<1. The magnitude of the interior steady state with respect to the constant stocking u can be shown to satisfy(20)dx∗du=gx∗g′x∗1-βx∗gx∗e-y∗-βgx∗1-e-y∗gx∗+x∗g′x∗,dy∗du=g′x∗gx∗·dx∗du.Consequently, dx∗/du<0 and dy∗/du>0 if (2) is satisfied.Observe thatR′(x)<0 and det⁡J∗<1 are equivalent to xgx′<1/λ and xgx′<1/βx, respectively. At any interior steady state (x,y), we have 1/λ=g(x)e-y and 1/βx=g(x)(1-e-y)/y-u. It can be verified that 1/λ<1/βx for any interior steady state (x,y). We therefore conclude that R′(x∗)<0 implies det⁡J∗<1.Proposition 5.
Letλe-u>1. If (2) is satisfied, then an interior steady state E∗=(x∗,y∗) of model (15) is asymptotically stable if (12) holds. In particular, E∗ is asymptotically stable if R′(x∗)<0. Moreover, x∗ and y∗ are decreasing and increasing functions of u, respectively.Letλe-u>1. For the Beverton-Holt growth rate g(x)=1/1+kx, k>0, system (15) has a unique interior steady state E∗. We also have(21)det⁡J∗=βx1+kx2,ddet⁡J∗dx=β1-kx1+kx3,ddet⁡J∗du=β1-kx1+kx3·dxdu.The maximum of det⁡J∗ occurs at x=xm≔1/k with max⁡(det⁡J∗)=β/4k. In particular, if β<4k, then E∗ is asymptotically stable. Therefore, stocking has no effect on the stability of the interior steady state if β<4k. If β>4k and λ≤2, then 1/k≥x¯ and thus ddet⁡J∗/dx>0 on (0,x¯). It follows that ddet⁡J∗/du<0 and stocking can stabilize the interior steady state E∗. If β>4k and λ>2, then max⁡(det⁡J∗)>1. Since x∗ is a decreasing function of u, it is expected that stocking can stabilize the interior steady state by (21).Consider a numerical example withλ=2, k=1, and β=7>4k. Then the unique interior steady state E∗ exists if u<ln⁡λ≈0.6931. From Figure 1, it is known that the interior steady state of model (4) is unstable when there is no stocking. We then use u as the bifurcation parameter. Figure 2 illustrates that stocking can stabilize the host-parasitoid interaction. Moreover, x and y components of E∗ are clearly decreasing and increasing functions of u, respectively.Figure 2
Bifurcation diagrams for system (15) with a Beverton-Holt growth rate. Parameter values are λ=2, k=1, and β=7.
(a)
(b)


## 3. Discussion

Pests are animals that are detrimental to humans or to human concerns. They cause damage to agriculture by feeding on crops or parasitising livestock. Biological control is the reduction of pest populations by natural enemies that typically involves an active human role. Very often biological control requires supplemental release of natural enemies. Relatively few natural enemies may be released at a critical time of the season. Such an augmentation is termed as aninoculative release or a stocking [1].Natural enemies of insect pests are often parasitoids. Since insects have distinctive life stages, we propose discrete-time host-parasitoid models to study the effects of constant stocking of parasitoids upon the host-parasitoid interaction. The novelty of this work is that density dependence of the hosts occurs first followed by parasitism. This assumption is revealed in the parasitoid equation of models (5) and (15). It is shown that the stocking clearly eliminates the steady state for which the host survives and the parasitoid is extinct. Therefore, with the biological control strategy of stocking, it is impossible for the hosts to survive alone without the parasitoids. On the other hand, it is possible to eliminate the host population if the stocking is large. Moreover, stocking can simplify the host-parasitoid interaction by stabilizing the coexisting steady state. In such case, the host and the parasitoid densities at the steady state are decreasing and increasing functions of the stocking, respectively.We conclude from this study that stocking of the parasitoids can drive the pests to extinction even if parasitism occurs after density dependence. This is consistent with the model studied in [5], where parasitism occurs prior to density dependence. However, it may take a longer or a shorter period of time to eliminate the pests depending on whether density dependence occurs first or after the parasitism.



---

*Source: 905957-2015-06-10.xml*