![MathJax Logo](/templates/jsp/_style2/_tandf/pb2/images/math-jax.gif)
Abstract
We propose a discrete-time host–parasitoid model with stage structure in both species. For this model, we establish conditions for the existence and global stability of the extinction and parasitoid-free equilibria as well as conditions for the existence and local stability of an interior equilibrium and system persistence. We study the model numerically to examine how pesticide spraying may interact with natural enemies (parasitoids) to control the pest (host) species. We then extend the model to an impulsive difference system that incorporates both periodic pesticide spraying and augmentation of the natural enemies to suppress the pest population. For this system, we determine when the pest-eradication periodic solution is globally attracting. We also examine how varying the control measures (pesticide concentration, natural enemy augmentation and the frequency of applications) may lead to different pest outbreak or persistence outcomes when eradication does not occur.
Dedicated to the memory of Abdul-Aziz Yakubu.
1. Introduction
Biological control refers to the use of a natural enemy, such as a predator, parasitoid or pathogen, to suppress a pest population. These natural enemies are released either periodically in small numbers, referred to as inoculative release, or all at once in large quantities, referred to as inundative release [Citation17]. For example, the periodic introduction of the parasitoid Encarsia formosa has been found to be effective at controlling the greenhouse whitefly, Trialeurodes vaporariorum, a serious pest of tomato and cucumber crops [Citation2]. Insect parasitoids are the most commonly utilized biological control agents [Citation13]. The larvae of parasitoids grow inside or on a host which they eventually kill [Citation25].
Integrative pest management (IPM) schemes employ a number of strategies, such as biological and chemical controls (e.g. pesticides or insecticides), to suppress a pest species [Citation20]. Since single chemical control strategies frequently fail as a result of evolution of resistance in the pest, these schemes often use chemical control in conjunction with biological control agents to reduce the abundance of a pest. However, the compatibility of these controls is critical to the success of an integrated pest management strategy [Citation12, Citation24, Citation36, Citation37]. Not only may different natural enemies vary widely in the response to toxicants, but variation in population structure may also impact results. For example, Diaeretiella rapae (M'Intosh), a common parasitoid of the cabbage aphid, was shown to experience differential stage susceptibility to the insecticide imidacloprid, with mortality in the adult stage being more than double that of the pupae (mummy) stage for various concentrations of the insecticide [Citation37].
There are many existing mathematical models of host–parasitoid systems (see [Citation27] and the references therein). One of the earliest discrete-time host–parasitoid models was developed by Thompson in 1924, while the more familiar Nicholson–Bailey model was developed later in 1935. A modification to the Nicholson–Bailey model with type-II functional responses was developed by Holling in 1959. This model is the first to assume that the parasitoid is search limited, in contrast to Thompson's model which assumed the parasitoid is egg-limited. The differential equation system developed independently by Lotka (1925) and Volterra (1926) for predator–prey interactions has also been used as a basis for modelling host–parasitoid interactions in continuous time [Citation27]. Ample extensions of these foundational host–parasitoid models have been developed to describe different host–parasitoid scenarios. These includes host–parasitoid models with different nonlinearities [Citation4, Citation16, Citation17, Citation30, Citation31], with different functional responses [Citation32, Citation33], with Allee effects in the host [Citation15, Citation21, Citation22], with stage structure [Citation3, Citation6, Citation15, Citation35, Citation42], age structure [Citation5], and co-evolution between host and parasitoid [Citation11, Citation29]. Numerous integrated pest management host–parasitoid systems have also been developed, including [Citation10, Citation20, Citation38–41, Citation43, Citation44]. These impulsive difference systems assume that control measures consisting of both chemical control and the release of natural enemies such as parasitoids are applied either periodically or when the pest (host) density reaches a designated economic threshold. However, we note that these systems often assume unstructured populations.
Motivated by the distinctive developmental stages often observed in both hosts and parasitoids as well as the observation that both species may exhibit stage-specific differential susceptibility to toxicants [Citation37], in this paper we develop a discrete-time host–parasitoid model with stage structure in both species. Specifically, we assume individuals in each species are classified according to two developmental stages, either juvenile or adult, resulting in a four-dimensional system. We first provide a rigorous analysis of the equilibria dynamics and persistence of this system. Then, viewing the host as a pest species, we use the system to consider integrated pest management interventions that consist of both biological control (i.e. the parasitoid) and chemical control via the release of a pesticide.
We first consider the case where control measures are continuous leading to pesticide-induced mortality in every time unit. We find that the effectiveness of the combined management strategies depends heavily on the timing of pesticide spray (e.g. at the beginning versus end of the time unit), the susceptibility of the parasitoid to the pesticide and the type of attractor to which the system converges. In particular, when time series solutions are attracted to equilibria, then both spraying at the end of the time unit and pesticide-induced parasitoid mortality can result in scenarios where, taken together, the two control measures may perform worse at controlling the pest species than when only a single control strategy is employed. However, this scenario may change for certain parameters when the system exhibits oscillatory dynamics. To examine the case where control measures are applied periodically rather than continuously, we further extend the system to an impulsive difference equation that describes the periodic application of pesticides as well as the periodic release of additional parasitoids. For this impulsive model, we study the existence and attractivity of the host-eradication cycle. We then numerically study how the control measures impact host dynamics when host eradication is not achieved.
This paper organized as follows. In Section 2, we introduce the discrete-time host–parasitoid model with stage structure in both the parasitoid and the host. In Section 3, we study the existence and stability of the model equilibria including the global stability of the two boundary equilibria, namely the extinction and parasitoid-free equilibria, and the local stability of the interior equilibrium near the bifurcation the occurs when the parasitoid-free equilibrium destabilizes. This latter result is obtained using a general bifurcation theorem, presented in Appendix A, that applies to structured species with host–parasitoid-type interactions. We also show in Section 3.2 that the system is persistent. Next in Section 3.3, we consider how the interaction of natural enemies (i.e. the parasitoids) and chemical control can impact control of a pest species. Finally, in Section 4, we extend the model to study an impulsive difference equation system which allows us to consider both periodic spraying and periodic parasitoid supplementations.
2. The host–parasitoid model
We consider a stage-structured, host–parasitoid model in which both hosts and parasitoids have two stages, namely juveniles and adults. Let and
denote the juvenile and adult host population and
and
denote the juvenile and adult parasitoid population. It is assumed that adult parasitoids
may parasitize either the juvenile
or adult
host while the juvenile parasitoid stage
occurs inside the host and, thus, does not contribute to parasitism. Specifically, the model is given by
(1)
(1) where
denotes a survival probability,
denotes a transition probability, and
denotes a stage-specific parasitoid searching efficiency. The nonlinearity f describes host fecundity, while
gives the average number of viable parasitoid eggs laid in a single host. The model assumes that parasitism occurs at the beginning of the time interval and that parasitized adult hosts cannot reproduce but do contribute to density dependence in fecundity. The nonlinearity f is assumed to satisfy the following conditions:
,
for
and
.
System (Equation1(1)
(1) ) can be written in the matrix form
where
and the projection matrix A is given by
3. Dynamics of the host–parasitoid model
An equilibrium of the model (Equation1(1)
(1) ) must satisfy the following system of equations:
(2)
(2) Solving the above system, we find that model (Equation1
(1)
(1) ) has two boundary equilibria, namely the extinction equilibrium
where both populations die out and the parasitoid-free equilibrium
where the host survives but the parasitoid goes extinct. Later, we also establish the existence of an interior equilibrium where both population densities are positive. In this section, we discuss the existence and stability of these equilibria.
In Lemma 3.1, we verify that solutions of system (Equation1(1)
(1) ) remain non-negative and bounded, and therefore system (Equation1
(1)
(1) ) is biologically sound.
Lemma 3.1
Solutions of system (Equation1(1)
(1) ) remain non-negative and are bounded for t>0.
Proof.
It is clear that solutions of the Equation (Equation1(1)
(1) ) with non-negative initial conditions remain non-negative for all forward time. Next, by the assumptions on f we note that
for some D>0. Therefore we have
It follows that
From the equation for
we have
which implies
Thus for any
, there exists a
, which depends on
and
, such that for
,
Define the compact set
, where
and
are defined above. Then every forward solution of (1) enters K in finite time and remains in K forever after.
Next, we consider the parasitoid equations. Notice that for any we have
Iterating the right-hand side of this inequality and taking
, we find
where
. Thus there exists a
, which depends on
such that for
It follows from an analogous argument that
. Hence, all solutions of system (Equation1
(1)
(1) ) are bounded.
3.1. Boundary dynamics
We first consider the two boundary equilibria of model (Equation1(1)
(1) ). In Theorem 3.2, we show that the extinction equilibrium is globally asymptotically stable when the inherent net reproductive number for the host, given by
(3)
(3) is less than 1. This quantity, which gives the expected number of offspring produced by a host individual over its lifetime in the absence of density effects and parasitism, is on the same side of one as the dominant eigenvalue of the inherent projection matrix (see, e.g., Theorem 1.1.3 of [Citation7]).
Theorem 3.2
The extinction equilibrium of model (Equation1
(1)
(1) ) is globally asymptotically stable if
and unstable if
.
Proof.
The Jacobian matrix of model (Equation1(1)
(1) ) evaluated at the extinction equilibrium
is given by the block diagonal matrix
From the lower
block, it is clear that two eigenvalues are
and
. Note that
for i = 1, 2.
Next, consider the upper block. Note that this matrix is the inherent projection matrix
of the host subsystem
(4)
(4) when the parasitoid is absent, that is
. Therefore, to determine when the eigenvalues of this matrix have magnitude less than 1, we calculate the inherent net reproductive number
for the host. To do this, we first decompose
as
, where the transition matrix
and the fertility matrix
are given by
Then
is the spectral radius of
, where
is the
identity matrix. A calculation shows that
is given by Equation (Equation3
(3)
(3) ). Thus the extinction equilibrium is locally asymptotically stable if
and unstable if
.
It remains to show is globally asymptotically stable when
. Consider the host subsystem given in (Equation4
(4)
(4) ). Note that
(5)
(5) Since the inherent projection matrix,
, of subsystem (Equation4
(4)
(4) ) is primitive, it has a positive, strictly dominant eigenvalue r. Furthermore, since
, it follows from Theorem 1.1.3 in [Citation7] that r<1 and, thus,
. Now since the projection matrix of the host subsystem is monotone, for any
we have that
by (Equation5
(5)
(5) ). By induction it follows that
where
converges to 0 as
Thus we have whenever
. As a result,
, and hence for
all solutions of (Equation1
(1)
(1) ) with non-negative initial conditions converge to the extinction equilibrium
.
Next in Theorem 3.3, we show that the parasitoid-free equilibrium exists when
and is locally asymptotically stable if
, where
(6)
(6) is the invasion net reproductive number of the parasitoid. This quantity gives the inherent net reproductive number of the parasitoid when the host is at its parasitoid-free equilibrium.
Theorem 3.3
The parasitoid-free equilibrium of model (1) exists if . It is locally asymptotically stable if
and unstable if
.
Proof.
Setting in the equilibrium Equation (Equation2
(2)
(2) ) and solving, we find that the host components of the parasitoid-free equilibrium are given by
(7)
(7) Note that from the properties of f we have
Hence the parasitoid-free equilibrium exists if and only if the
Next, we consider the Jacobian matrix of system (Equation1(1)
(1) ) evaluated at
given by
Since this matrix is block triangular, we first show that the eigenvalues of the upper
block have magnitude less than 1. Decomposing this block as
, where
is the transition matrix for
and
we observe that the assumptions on f mean
. Therefore, we can again apply Theorem 1.1.3 in [Citation7] to get that the dominant eigenvalue of
is on the same side of one as
where
denotes the spectral radius.
It follows that the local stability of the parasitoid-free equilibrium is determined by the dominant eigenvalue of the lower diagonal block. The dominant eigenvalue of this matrix is on the same side of one as the invasion net reproductive number
which is calculated in the same manner as
. Namely, we first decompose the matrix as
where
Then
is defined as
, which results in (Equation6
(6)
(6) ). Thus the parasitoid-free equilibrium is locally asymptotically stable when
and unstable when
.
Next, to establish the global stability of the parasitoid-free equilibrium, in Lemma 3.4 we first establish the global stability of the positive equilibrium of the host subsystem (Equation4(4)
(4) ) when parasitoids are absent.
Lemma 3.4
Consider the host subsystem (Equation4(4)
(4) ) when parasitoids are absent, that is
. If
, then there exists a unique positive equilibrium
defined by Equation (Equation7
(7)
(7) ) that is globally asymptotically stable.
Proof.
Note that the existence and local asymptotic stability of the unique positive equilibrium follow from the proof of Theorem 3.3. Thus it remains to establish the global attractivity of this equilibrium. Note that the map
is monotone and has a unique fixed point in
. Therefore, we make use of Lemma 3 in [Citation1] which states that if there is an ordered interval
containing the fixed point that satisfies
and
, then every solution sequence starting in
converges to the unique fixed point.
To apply Lemma 3, we first observe that every solution starting on the boundary of but not at
enters the positively invariant set int(
). Therefore, it is enough to consider solutions in int(
). Moreover, since Lemma 3.1 shows that every forward solution of (Equation4
(4)
(4) ) enters the compact set K in finite time and remains in K, it suffices to consider
int(
)
. The equilibrium
clearly belongs to K. Choose
(the maximal element in K). Then by Lemma 3.1, we have
. Since
is primitive, its spectral radius r (which we know is larger than one since
) has a corresponding positive eigenvector v,
In addition, the monotonicity of F together with r>1 mean that
for all
sufficiently small. Now, for a given
we can pick a sufficiently small
such that
Thus it follows from an application of Lemma 3 in [Citation1] that every solution of (Equation4
(4)
(4) ) starting in
converges to the equilibrium
.
Applying Lemma 3.4 together with the theory on asymptotically autonomous systems (see, e.g., Theorem 3.2 in [Citation28]), we now establish the global stability of the parasitoid-free equilibrium.
Theorem 3.5
Let . The parasitoid free equilibrium
of system (Equation1
(1)
(1) ) is globally asymptotically stable in
if
.
Proof.
We first improve upon the upper bound for the host that was obtained in Lemma 3.1. From the second equation in system (Equation1(1)
(1) ), we obtain
Moreover, from the first equation of the system (Equation1
(1)
(1) ), we have
Algebraic manipulations of this inequality in combination with the previous inequality result in
From this we have that for any
, there exists a
such that for
,
Iterating the right-hand side of this inequality and taking
we have
Substituting these improved bounds into the parasitoid equations in system (Equation1
(1)
(1) ) and using
we define the linear system
where
is sufficiently small so that
Note the existence of such an η is guaranteed since
. Then there exists a
such that when
and
we have for
Following the same steps as in the proof of Theorem 3.3, it is straightforward to show that
and
whenever
. Hence the non-autonomous host subsystem
is asymptotic to the limiting system,
Since
converges uniformly to F and
int
, it follows from Lemma 3.4 and Theorem 3.2 in [Citation28] that any solution of system (Equation1
(1)
(1) ) with
converges to
. Thus the parasitoid-free equilibrium is globally asymptotically stable.
3.2. Interior dynamics
In this section, we consider the interior dynamics of model (Equation1(1)
(1) ). We first show in Corollary 3.6 that, when the parasitoid-free equilibrium destabilizes at
, a branch of interior equilibria bifurcates from this equilibrium and the equilibria are locally asymptotically stable for
. This result follows from a more general bifurcation result presented in Appendix A which is an extension of the work developed in [Citation26].
Corollary 3.6
Assume . Then at
there exists a branch of interior equilibria bifurcating from parasitoid-free equilibrium
of model (Equation1
(1)
(1) ). The bifurcation is forward (supercritical) and the equilibria are locally asymptotically stable for
.
Proof.
Here we apply Theorem A.1 to model (Equation1(1)
(1) ) using
as the bifurcation parameter. Let τ denote the parasitoid-free equilibrium
and let
denote the value of
at which this equilibrium destabilizes, that is,
is defined such that
. The Jacobian matrix evaluated at
is given by
This matrix has a simple eigenvalue 1 at
, with right eigenvector
where
and left eigenvector
Next we calculate the direction of bifurcation κ and stability of the bifurcating equilibria, as determined by
, as given in Theorem A.1. For system (Equation1
(1)
(1) ), we have
,
Thus we obtain
Since
we have
and
. Thus it follows that, at
, there exists a forward bifurcating branch of stable equilibria.
Our next goal is to establish the uniform persistence of the system. We do this by first proving that the host population is uniformly persistent and then proving that the parasitoid population is uniformly persistent. To prove these results we utilize the theory established in [Citation34] and use similar notation as in that paper. First, observe that from Lemma 3.1 it follows that there exists a compact and positively invariant set that attracts all solutions in
(also, refer to Theorem 2.6 in [Citation23]).
Theorem 3.7
If , then there exists an
such that
for all initial conditions satisfying
.
Proof.
Using the notation of [Citation34], we let ,
and we take the continuous matrix
and write model (Equation1
(1)
(1) ) in the form presented in [Citation34] as follows:
(8)
(8) Let
and
. Clearly, Z and
are positively invariant sets. Define the set
. Clearly, M is nonempty compact and positively invariant. Furthermore, it is easy to see that
. Since
, from Theorem 3.2 we have that the spectral radius
Since,
is primitive, it follows from Proposition 3 and Corollary 1 in [Citation34] that M is a uniformly weak repeller. Thus, from Theorem 3.3 in [Citation34] we obtain that the total host population is uniformly persistent, i.e. there exist an
such that
for
, i.e. for
satisfying
. Now, utilizing the primitivity of
, we see that the assumptions of Proposition 4 in [Citation34] are satisfied using the compact set
and, hence, conclude that there exists an
such that
for all
.
Theorem 3.8
If , then there exists an
such that
for all initial conditions satisfying
and
Proof.
To establish persistence of , using the notation in [Citation34] we let
and
and define
We then take the continuous matrix
to be
and write model (Equation1
(1)
(1) ) in the form
(9)
(9) Now, let
,
and
. Then M is non-empty, compact and bounded away from zero. Thus from Lemma 3.4 it follows that
. Since
, from Theorem 3.3 we have that the spectral radius
Since,
is primitive, it follows from Proposition 3 and Corollary 1 in [Citation34] that M is a uniformly weak repeller. Thus, from Theorem 3.3 in [Citation34] we obtain that the total parasitoid population is uniformly persistent, i.e. there exist an
such that
for
, i.e. for
satisfying
and
. Now, utilizing the primitivity of
, we see that the assumptions of Proposition 4 in [Citation34] are satisfied using the compact set
and, hence, conclude that there exists an
such that
for all
.
Finally, the following corollary, which easily follows from combining Theorem 3.7 and Theorem 3.8 and letting , shows that system (Equation1
(1)
(1) ) is persistent.
Corollary 3.9
If , then there exists an
such that
for all initial conditions satisfying
and
Remark 3.1
Note that, by Corollary 3.9 and Lemma 3.1, the system is permanent. In particular, there exists m, M>0, such that for any initial condition satisfying and
, there exist a T (which depends on the initial condition) such that for t>T
Furthermore, by an application of Brouwer's fixed-point theorem it follows that there is at least one interior equilibrium when
(see [Citation14]).
3.3. Impact of continuous pesticide spraying
In this section, we consider the combined impact of pesticide spraying and natural enemies (i.e. the parasitoid) to control the pest (i.e. the host). For simplicity, we first assume that pesticide spraying is continuous, meaning that it occurs every time unit. Later, in Section 4 we extend this to the more realistic case where spraying occurs periodically.
We consider two cases, pesticide spraying at the beginning of the time interval or at the end of the time interval, and assume that pesticide spraying results in mortality in each stage with
. When pesticide spraying occurs at the end of the time interval, then the right hand side of each equation in system (Equation1
(1)
(1) ) is proportionally reduced by
resulting in the system:
(10)
(10) Meanwhile, if spraying occurs at the beginning of the time interval, then each density at time t is proportionally reduced by
resulting in the system:
(11)
(11) We note that the analysis of model (Equation1
(1)
(1) ) established in Section 3 applies to models (Equation10
(10)
(10) ) and (Equation11
(11)
(11) ) when
and
are given by
For convenience, we summarize the order-of-event assumptions underlying these two models in Tables and .
Table 1. Order of events within the time interval when spraying occurs at the end of the time interval.
Table 2. Order of events within the time interval when spraying occurs at the beginning of the time interval.
In Example 3.10, we consider the impact of pesticide spraying when spraying occurs at the beginning or end of the time interval. We consider two cases, when the only direct effects of the pesticide are on the host and when both species are directly impacted. We compare these two cases to a single control strategy, namely only parasitoids and only chemical control. We find that the timing of pesticide spraying can have a significant impact on the effectiveness of these combined control strategies. For all numerical simulations, we define host fecundity using the Beverton–Holt nonlinearity , where
.
Example 3.10
A comparison of spraying strategies
We consider the impact of varying the pesticide-induced mortality in the host assuming equal impacts on both stages, i.e. . We suppose that either pesticide spraying does not directly impact the parasitoid, that is
, or the two parasitoid stages are impacted by spraying but not as significantly as the host, namely,
. All other parameters are fixed at the following values:
For this choice of parameters, solutions are attracted to stable equilibria for all values of ϵ. Note that, in both cases, the parasitoid alone can never eradicate the host without additional parasitoid supplementation.
We observe that when spraying occurs at the end of the time interval (Figure A–B) and , then the combination of the two controls works are desired, with the total host abundance being lower than for chemical or biological control alone and the host going extinct for sufficiently large ϵ. Meanwhile, when the pesticide does result in direct mortality in the parasitoid,
, then the host density experiences a moderate increase as ϵ increases but does not exceed the host density with only chemical control. Conversely, when spraying occurs at the beginning of the time unit (Figure C–D), then the combination of chemical and biological controls may produce a worse scenario than biological control alone. This occurs even when the chemical control has no direct effects on the parasitoid due to increased indirect mortality (refer to the equation for
in model (Equation11
(11)
(11) )).
Figure 1. The bifurcation diagrams for pesticide spraying at the end of the time unit (A)–(B) and pesticide spraying at the beginning of the time unit (C)–(D). These diagrams were obtained by finding the time-series solution out to time 20,000 for initial conditions using parameter values provided in Example 3.10. We consider three cases: only chemical controls are applied (
) (magenta circle), both chemical and biological controls are applied with no direct effects on the parasitoid (
) (red square), and both chemical and biological controls are applied with direct effects on the parasitoid (
) (blue triangle).
![Figure 1. The bifurcation diagrams for pesticide spraying at the end of the time unit (A)–(B) and pesticide spraying at the beginning of the time unit (C)–(D). These diagrams were obtained by finding the time-series solution out to time 20,000 for initial conditions H1(0)=H2(0)=P1(0)=P2(0)=1 using parameter values provided in Example 3.10. We consider three cases: only chemical controls are applied (P1(t)≡P2(t)≡0) (magenta circle), both chemical and biological controls are applied with no direct effects on the parasitoid (ϵ3=ϵ4=0) (red square), and both chemical and biological controls are applied with direct effects on the parasitoid (ϵ3=ϵ4=12ϵ) (blue triangle).](/cms/asset/1fcd6af7-4acd-486a-a587-f7c1146ad75a/tjbd_a_2241483_f0001_oc.jpg)
However, when parameter values are chosen so that the model exhibits oscillatory behaviour, this pattern can change. This is illustrated in Figure where all parameters are the same as in Figure except for host fecundity which has been increase to . This increase in
results in stable oscillatory dynamics for smaller values of ϵ. Specifically, for
these oscillations occur in both cases when
and for
these oscillations occur in both cases when
. We observe that increasing ϵ within the oscillatory region may result in an increase in parasitoid abundance and a subsequent decrease in host abundance when the pesticide has direct mortality effects on the parasitoid.
Figure 2. The bifurcation diagrams for pesticide spraying at the end of the time unit (A)–(B) and pesticide spraying at the beginning of the time unit (C)–(D). These diagrams were obtained by finding the time-series solution out to time 20,000 for initial conditions using parameter values provided in Example 3.10 except with host fecundity increased to
. We then plot the average host and parasitoid densities using the last 500 time units. We consider three cases: only chemical controls are applied (
) (magenta circle), both chemical and biological controls are applied with no direct effects on the parasitoid (
) (red square), and both chemical and biological controls are applied with direct effects on the parasitoid (
) (blue triangle).
![Figure 2. The bifurcation diagrams for pesticide spraying at the end of the time unit (A)–(B) and pesticide spraying at the beginning of the time unit (C)–(D). These diagrams were obtained by finding the time-series solution out to time 20,000 for initial conditions H1(0)=H2(0)=P1(0)=P2(0)=1 using parameter values provided in Example 3.10 except with host fecundity increased to f0=4000. We then plot the average host and parasitoid densities using the last 500 time units. We consider three cases: only chemical controls are applied (P1(t)≡P2(t)≡0) (magenta circle), both chemical and biological controls are applied with no direct effects on the parasitoid (ϵ3=ϵ4=0) (red square), and both chemical and biological controls are applied with direct effects on the parasitoid (ϵ3=ϵ4=12ϵ) (blue triangle).](/cms/asset/4bd3fbd5-b23a-414c-a069-ddffb516e101/tjbd_a_2241483_f0002_oc.jpg)
4. Host–parasitoid model with periodic impulsive effects
In Section 3.3, we observed that with continuous control it is possible for the combination of chemical control and natural enemies to produce worse outcomes in terms of host control than natural enemies alone if the chemical control is not sufficiently effective. However, more realistically we may expect chemical controls to be applied periodically rather than continuously. In this section, we extend model (Equation1(1)
(1) ) to an impulsive difference equation that includes periodic pesticide spraying along with periodic releases of parasitoids. This extension is obtained by following a model formulation similar to that studied in [Citation39].
We assume that at every qth time step there is a perturbation which incorporates a proportional decrease of the host and parasitoid population due to the application of pesticides and an increase of the parasitoid population due to the release of a constant number of parasitoids which is independent of the current parasitoid density. This perturbation is assumed to occur at the end of the time unit. The model is given by following system of difference equations: (12)
(12) where q is a fixed integer and denotes the period of the impulsive effect, which means that an integrated control strategy is applied during time unit t = qk for
. Here
,
, for i = 1, 2, 3, 4 denotes the proportion of individuals in each stage that die as a result of pesticide spraying at time qk, and
for i = 1, 2 gives the amount of juvenile and adult parasitoids that are released at this time.
and
for i = 1, 2 are the densities of the host and parasitoid, respectively, at time qk prior to the impulsive perturbation, and
and
are the densities of host and parasitoid at time qk after an impulsive perturbation. The initial conditions are
. Here for convenience we denote the initial density after an impulsive perturbation at time t = 0.
Note that if q = 1 and , then system (Equation12
(12)
(12) ) is equivalent to continuous spraying at the end of the time interval, as given by system (Equation10
(10)
(10) ).
4.1. Host-eradication periodic solutions
Our goal is to determine when system (Equation12(12)
(12) ) results in host eradication. To this end, we first establish the existence of host-eradication periodic solutions. Consider the parasitoid-only subsystem of model (Equation12
(12)
(12) ):
(13)
(13) In Lemma 4.1, we show that all solutions of system (Equation13
(13)
(13) ) converge to a positive periodic solution.
Lemma 4.1
System (Equation13(13)
(13) ) has a periodic solution
defined as
where
and
Moreover, for every solution
of (Equation13
(13)
(13) ) we have
as
.
Proof.
It follows from the periodicity of the system (Equation13(13)
(13) ) that the solution can be defined at the impulsive sub-interval
with
, where
means that we take the density of the parasitoid after an impulsive perturbation as the initial value in this interval. Initiating the population at time
and applying system (Equation13
(13)
(13) ) we find that the solution at time
is given by
Applying the impulsive condition we find
(14)
(14) Therefore, the q-periodic solution of system (Equation13
(13)
(13) ) is determined by the solution to the system
Solving this system we find
Note that the restrictions on the model parameters ensure that these two quantities are always positive.
Next, we note that the convergence of solutions to the periodic solution follows from the linearity of the projection matrix in Equation (Equation14
(14)
(14) ) and the fact that the eigenvalues of this matrix have absolute values less than 1. Specifically, define
Then iterating Equation (Equation14
(14)
(14) ), we obtain
Since
, taking the limit as
we find that the solution converges to
which is the equilibrium solution to Equation (Equation14
(14)
(14) ).
Next, in Theorem 4.2 we establish sufficient conditions under which a given control strategy results in host eradication.
Theorem 4.2
The host-eradication cycle is a globally attracting solution of system (Equation12
(12)
(12) ) if
(15)
(15) In particular, a sufficient condition is
(16)
(16) where
Proof.
Note that
Consider the following impulsive difference equation:
(17)
(17) where
and
It follows from Lemma 4.1 that if
and
, then we have
and
and
and
as
Hence, for any
there exists
such that
and
for
. Therefore, for
we have
(18)
(18) It follows that
Thus, if inequality (Equation15
(15)
(15) ) is satisfied, then
,
as
Consequently,
and
as
The sufficient condition for (Equation15(15)
(15) ) can be obtained by noting that
4.2. Numerical study of the impulsive model
In this section, we further study the dynamics of the impulsive model (Equation12(12)
(12) ). In Example 4.3, we first consider the impact of combined control strategies, i.e. chemical and biological control, when controls are applied periodically. Next, in Example 4.4 consider how host density varies for differing control strategies in the case when host eradication is not achieved. We find that there are parameter ranges for which multiple attractors occur and, thus, the effectiveness of control is initial condition dependent.
Example 4.3
Interaction of control strategies in the impulsive model
We consider the impact of varying the pesticide-induced mortality in the host assuming equal impacts on both stages, i.e. . We consider three scenarios: pesticide spraying does not directly impact the parasitoid (
), has moderate impacts on the parasitoid (
) or has severe impacts on the parasitoid (
). All other parameters are fixed at the following values:
We observe that when there are no direct mortality effects on the parasitoid, then the average host density is a decreasing function of the toxicant effect ϵ (Figure A). Since the impulsive model corresponds to spraying at the end of the time unit, this is in agreement with the dynamics observed in Example 3.10 (Figure A). Meanwhile, when direct effects on the parasitoid are present we again observe that the average host density can increase with increasing pesticide concentration resulting in host densities almost doubling in certain situations. Note here that, due to parasitoid supplementation, the parasitoid population survives even when the hosts die out (Figure B). Also note that, similar to the case where models (Equation10
(10)
(10) ) and (Equation11
(11)
(11) ) exhibit oscillatory dynamics (Figure ), it is possible for parasitoid densities to be higher when they experience greater pesticide mortality.
Figure 3. The bifurcation diagrams for the impulsive system (Equation12(12)
(12) ). These diagrams were obtained by finding the time-series solution out to time 5000 for initial conditions
using parameter values provided in Example 4.3. We then plot the average host and parasitoid densities using the last 500 time units. We consider three cases: no direct effects on the parasitoid (
) (red square), moderate direct effects on the parasitoid (
) (blue triangle) and severe direct effects on the parasitoid (
) (cyan diamond).
![Figure 3. The bifurcation diagrams for the impulsive system (Equation12(12) H1(t+1)=f(H2(t))H2(t)e−a2P2(t)+s1(1−γ1)H1(t)e−a1P2(t),H2(t+1)=s1γ1H1(t)e−a1P2(t)+s2H2(t)e−a2P2(t),P1(t+1)=β1H1(t)(1−e−a1P2(t))+β2H2(t)(1−e−a2P2(t))+s3(1−γ3)P1(t),P2(t+1)=s3γ3P1(t)+s4P2(t),H1(qk+)=(1−ϵ1)H1(qk),H2(qk+)=(1−ϵ2)H2(qk),P1(qk+)=(1−ϵ3)P1(qk)+τ1,P2(qk+)=(1−ϵ4)P2(qk)+τ2,(12) ). These diagrams were obtained by finding the time-series solution out to time 5000 for initial conditions H1(0)=H2(0)=P1(0)=P2(0)=1 using parameter values provided in Example 4.3. We then plot the average host and parasitoid densities using the last 500 time units. We consider three cases: no direct effects on the parasitoid (ϵ3=ϵ4=0) (red square), moderate direct effects on the parasitoid (ϵ3=ϵ4=12ϵ) (blue triangle) and severe direct effects on the parasitoid (ϵ3=ϵ4=ϵ) (cyan diamond).](/cms/asset/31be3430-49e7-4914-8a13-a938d2fe7a9e/tjbd_a_2241483_f0003_oc.jpg)
Example 4.4
Impact of control strategies on host densities
In this example, we first consider the impact of varying the control parameters q and . We fix all other parameters at the following values:
In Figure (A), we give the average host density against the control frequency q and the number of adult parasitoids released
. We observe that the host is eradicated for sufficiently strong control measures (i.e. high frequency pesticide spraying and/or a large parasitoid release number), while its average density follows a generally increasing trend as these control measures are weakened. However, there are parameter regions for which this increasing trend does not hold. This is a result of the existence of multiple stable attractors with significant differences in average host densities in certain regions (Figure B). To examine this phenomenon further, we fix
and find the average host density against q and
(Figure A). We observe a region of intermediate q values for which these multi-attractors produce noticeable changes in host densities. Investigation of the larger q parameter space shows that though changes in host density are not as significant, the system dynamics are quite complicated and exhibit a large numbers of attractors (Figure B–C).
Figure 4. (A) The average host density against q and obtained from the last 100,000 values of a 500,000 time unit simulation. All parameters are as given in Example 4.4 and initial conditions are
. The red line corresponds to the host eradication condition (Equation16
(16)
(16) ). Note that q must be integer-valued. (B) The final time solutions of system (Equation12
(12)
(12) ) when q = 30 and
and all other parameters are provided in Example 4.4. For
, the system has four stable attractors, two with large amplitudes and two with small amplitudes.
![Figure 4. (A) The average host density against q and τ2 obtained from the last 100,000 values of a 500,000 time unit simulation. All parameters are as given in Example 4.4 and initial conditions are H1(0)=H2(0)=P1(0)=P2(0)=1. The red line corresponds to the host eradication condition (Equation16(16) λq(1−min{ϵ1,ϵ2})∏t=qsq(s+1)−1e−min{a1,a2}P2∗(t)<1,(16) ). Note that q must be integer-valued. (B) The final time solutions of system (Equation12(12) H1(t+1)=f(H2(t))H2(t)e−a2P2(t)+s1(1−γ1)H1(t)e−a1P2(t),H2(t+1)=s1γ1H1(t)e−a1P2(t)+s2H2(t)e−a2P2(t),P1(t+1)=β1H1(t)(1−e−a1P2(t))+β2H2(t)(1−e−a2P2(t))+s3(1−γ3)P1(t),P2(t+1)=s3γ3P1(t)+s4P2(t),H1(qk+)=(1−ϵ1)H1(qk),H2(qk+)=(1−ϵ2)H2(qk),P1(qk+)=(1−ϵ3)P1(qk)+τ1,P2(qk+)=(1−ϵ4)P2(qk)+τ2,(12) ) when q = 30 and τ2=4 and all other parameters are provided in Example 4.4. For 0.01≤P2(0)≤50, the system has four stable attractors, two with large amplitudes and two with small amplitudes.](/cms/asset/b4b17cd8-b4d1-4a65-8e90-ac462d7957dc/tjbd_a_2241483_f0004_oc.jpg)
5. Conclusion
We developed a discrete-time host–parasitoid model with stage structure in both species. This stage structure allows us to account for distinctive developmental stages in each species, which may impact species interactions, as well as differential susceptibility to pesticides which may be used in conjunction with a parasitoid species to control a pest. After establishing the equilibria dynamics and persistence of this system, we studied two forms of integrated pest management. For the first case, we assumed continuous pesticide spraying which means that pesticides are applied each time unit. Though this simple formulation makes analysis of the model mathematically tractable, it is less realistic particularly if the unit of time of the model is short, such as daily pesticide spraying. For the second case, we assumed that management occurs periodically, resulting in an impulsive difference equation system. For this latter case, we also included periodic parasitoid supplementation and established sufficient conditions for which the pest species is eradicated.
Figure 5. (A) The average host density against and q obtained from the last 100,000 values of a 500,000 time unit simulation. These graphs were obtained using
and
, with all other parameters given in Example 4.4. (B) The time series solutions when
for initial adult parasitoid density
. (C) A sample phase plane trajectory when q = 40 and
consisting of the last 3000 values of a 50,000 time unit simulation.
![Figure 5. (A) The average host density against P2(0) and q obtained from the last 100,000 values of a 500,000 time unit simulation. These graphs were obtained using τ2=4 and H1(0)=H2(0)=P1(0)=1, with all other parameters given in Example 4.4. (B) The time series solutions when q=40 for initial adult parasitoid density P2(0)=1:1:50. (C) A sample phase plane trajectory when q = 40 and P2(0)=1 consisting of the last 3000 values of a 50,000 time unit simulation.](/cms/asset/5ef20733-ea8d-4c26-b8e0-e70927b20b65/tjbd_a_2241483_f0005_oc.jpg)
We observed that the effectiveness of integrated pest management strategies is dependent on a number of factors. When spraying is continuous, pesticide spraying at the beginning of the time unit produces higher pest densities than spraying at the end of the time unit. Meanwhile, for both continuous and periodic spraying, if spraying directly impacts the parasitoid, then increases in pesticide-induced mortality can produce larger host densities if the pesticide is not sufficiently effective (that is, the
are not large enough). Interestingly, we also observe that when attractors are oscillatory, parasitoid density may be higher for larger direct pesticide induced mortality (larger
and
). However, this higher density does not necessarily translate into lower host density (see Figures and ).
We further found that, in certain parameter regimes, the effectiveness of control strategies may be initial condition dependent and, thus, both pest outbreak and low-density persistence may be obtained for the same control strategies (Figure A). This occurs as a result of multiple oscillatory attractors. Such scenarios are not unexpected when models contain Ricker-type nonlinearities and are a common observation for unstructured host-parasitoid systems [Citation19, Citation21], unstructured impulsive host-parasitoid systems [Citation39, Citation40, Citation43], as well as competing species models with Ricker-type nonlinearities [Citation8, Citation9]. Of note, however, is that a greater number of attractors does not necessarily mean greater variability in the mean host density as initial conditions are varied (compare q = 30 to q = 40 in Figure (B) and Figures (A)–(B).)
All together, our results demonstrate the need for integrated pest management strategies to account for interactions between different control strategies, particularly when natural enemies are also impacted by chemical controls. Here we observed that, if chemical control is not sufficiently effective, then the combination of chemical and biological control may produce worse outcomes than biological control alone. Of course, additional factors may further complicate outcomes, such as the development of resistance to the pesticide [Citation20] (which, in particular, may occur for long term and frequent application) or coevolution between host and parasitoid [Citation29].
Disclosure statement
No potential conflict of interest was reported by the author(s).
Additional information
Funding
References
- A.S. Ackleh and P. De Leenheer, Discrete three-stage population model: persistence and global stability results, J. Biol. Dyn. 2(4) (2008), pp. 415–427.
- T.S. Bellows and T.W. Fisher, Handbook of Biological Control: principles and applications of biological control, Academic Press, San Diego, CA, 1999.
- T.S. Bellows and M.P. Hassell, The dynamics of age-structured host–parasitoid interactions, J. Anim. Ecol. (57(1) (1988), pp. 259–268.
- C. Bernstein, Density dependence and the stability of host–parasitoid systems, Oikos (47(2) (1986), pp. 176–180.
- C.J. Briggs, Competition among parasitoid species on a stage-structured host and its effect on host suppression, Am. Nat. 141(3) (1993), pp. 372–397.
- K.M. Crowe and J.M. Cushing, Optimal instar parasitization in a stage structured host–parasitoid model, Nat. Resour. Model. 8(2) (1994), pp. 119–138.
- J.M. Cushing, An Introduction to Structured Population Dynamics, SIAM, 1998.
- J.M. Cushing, S.M. Henson, and C.C. Blackburn, Multiple mixed-type attractors in a competition model, J. Biol. Dyn. 1(4) (2007), pp. 347–362.
- J.L. Edmunds, Multiple attractors in a discrete competition model, Theor. Popul. Biol. 72(3) (2007), pp. 379–388.
- M. He, S. Tang, and R.A. Cheke, A holling type ii discrete switching host-parasitoid system with a nonlinear threshold policy for integrated pest management. Discrete Dynamics in Nature and Society, 2020, 2020.
- H.J. Henter and S. Via, The potential for coevolution in a host–parasitoid system. i. genetic variation within an aphid population in susceptibility to a parasitic wasp, Evolution 49(3) (1995), pp. 427–438.
- M.P. Hill, S. Macfadyen, and M.A. Nash, Broad spectrum pesticide application alters natural enemy communities and may facilitate secondary pest outbreaks, PeerJ 5 (2017), pp. e4179.
- M.E. Hochberg and R.D. Holt, The uniformity and density of pest exploitation as guides to success in biological control, in Theoretical approaches to biological control, Cambridge University Press, New York, NY, 1999. pp. 71–88.
- V. Hutson, The existence of an equilibrium for permanent systems, Rocky Mt. J. Math. 20(4) (1990), pp. 1033–1040.
- S.R.-J. Jang, Allee effects in a discrete-time host–parasitoid model with stage structure in the host, Discrete Contin. Dyn. Syst. B 8(1) (2007), pp. 145.
- S.R.-J. Jang and J.-L. Yu, Dynamics of a discrete-time host–parasitoid system with cannibalism, J. Biol. Dyn. 5(5) (2011), pp. 419–435.
- S.R.-J. Jang and J.-L. Yu, Discrete-time host–parasitoid models with pest control, J. Biol. Dyn. 6(2) (2012), pp. 718–739.
- T. Kato, A Short Introduction to Perturbation Theory for Linear Operators, Springer Science & Business Media, 2012.
- R. Kon, Multiple attractors in host–parasitoid interactions: coexistence and extinction, Math. Biosci. 201(1-2) (2006), pp. 172–183.
- J. Liang, S. Tang, and R.A. Cheke, A discrete host–parasitoid model with development of pesticide resistance and ipm strategies, J. Biol. Dyn. 12(1) (2018), pp. 1059–1078.
- H. Liu, Z. Li, M. Gao, H. Dai, and Z. Liu, Dynamics of a host–parasitoid model with allee effect for the host and parasitoid aggregation, Ecol. Complex. 6(3) (2009), pp. 337–345.
- G. Livadiotis, L. Assas, B. Dennis, S. Elaydi, and E. Kwessi, A discrete-time host–parasitoid model with an allee effect, J. Biol. Dyn. 9(1) (2015), pp. 34–51.
- Pierre Magal and Xiao-Qiang Zhao, Global attractors and steady states for uniformly persistent dynamical systems, SIAM J. Math. Anal. 37(1) (2005), pp. 251–275.
- A.R. Main, E.B. Webb, K.W. Goyne, and D. Mengel, Neonicotinoid insecticides negatively affect performance measures of non-target terrestrial arthropods: A meta-analysis, Ecol. Appl. 28(5) (2018), pp. 1232–1244.
- C.A. McKenzie, On the ecology of Diaeretiella rapae (M'Intosh). PhD thesis, Adelaide, 1977.
- E.P. Meissen, Invading a Structured Population: A Bifurcation Approach, PhD thesis, The University of Arizona, 2017.
- N.J. Mills and W.M. Getz, Modelling the biological control of insect pests: a review of host–parasitoid models, Ecol. Modell. 92(2–3) (1996), pp. 121–143.
- K. Mokni, S. Elaydi, M. Ch-Chaoui, and A. Eladdadi, Discrete evolutionary population models: a new approach, J. Biol. Dyn. 14(1) (2020), pp. 454–478.
- A. Sasaki and H.C.J. Godfray, A model for the coevolution of resistance and virulence in coupled host–parasitoid interactions, Proc. R. Soc. Lond. B. Biol. Sci. 266(1418) (1999), pp. 455–463.
- S.J. Schreiber, Host–parasitoid dynamics of a generalized Thompson model, J. Math. Biol. 52(6) (2006), pp. 719–732.
- S.J. Schreiber, Periodicity, persistence, and collapse in host–parasitoid systems with egg limitation, J. Biol. Dyn. 1(3) (2007), pp. 273–287.
- A. Singh, Analytical discrete-time models of insect population dynamics. 2020.
- A. Singh, Fluctuations in population densities inform stability mechanisms in host–parasitoid interactions. bioRxiv, 2021, pp. 2020–12.
- H.L. Smith and P.L. Salceanu, Lyapunov exponents and persistence in discrete dynamical systems, Discrete Contin. Dyn. Syst. B 12(1) (2009), pp. 187–203.
- T. Spataro and C. Bernstein, Combined effects of intraspecific competition and parasitoid attacks on the dynamics of a host population: a stage-structured model, Oikos 105(1) (2004), pp. 148–158.
- J.D. Stark, J.E. Banks, and S. Acheampong, Estimating susceptibility of biological control agents to pesticides: influence of life history strategies and population structure, Biol. Control. 29(3) (2004), pp. 392–398.
- J.D. Stark, J.K. McIntyre, and J.E. Banks, Population viability in a host–parasitoid system is mediated by interactions between population stage structure and life stage differential susceptibility to toxicants, Sci. Rep. 10(1) (2020), pp. 1–7.
- S. Tang and R.A. Cheke, Models for integrated pest control and their biological implications, Math. Biosci. 215(1) (2008), pp. 115–125.
- S. Tang, Y. Xiao, and R.A. Cheke, Multiple attractors of host–parasitoid models with integrated pest management strategies: eradication, persistence and outbreak, Theor. Popul. Biol. 73(2) (2008), pp. 181–197.
- S. Tang, G. Tang, and R.A. Cheke, Optimum timing for integrated pest management: modelling rates of pesticide application and natural enemy releases, J. Theor. Biol. 264(2) (2010), pp. 623–638.
- P. Wang, W. Qin, and G. Tang, Modelling and analysis of a host–parasitoid impulsive ecosystem under resource limitation. Complexity, 2019, 2019.
- S.M. White, S.M. Sait, and P. Rohani, Population dynamic consequences of parasitised-larval competition in stage-structured host–parasitoid systems, Oikos 116(7) (2007), pp. 1171–1185.
- C. Xiang, Z. Xiang, S. Tang, and J. Wu, Discrete switching host–parasitoid models with integrated pest control, Int. J. Bifurc. Chaos 24(09) (2014), pp. 1450114.
- C. Xiang, Z. Xiang, and Y. Yang, Dynamic complexity of a switched host–parasitoid model with Beverton–Holt growth concerning integrated pest management. Journal of Applied Mathematics, 2014, 2014.
Bifurcation analysis of a resident–invader system with host–parasitoid-type interactions
In this section, we extend the results obtained in [Citation26] to derive a general bifurcation result for a resident–invader system with host-parasitoid-type interactions. Consider a system of two species, a resident population with m stages and an invader population with n stages, given by the matrix equation
(A1)
(A1) where the projection matrix is defined as
(A2)
(A2) Here
is the resident projection matrix, matrices
and
determine the invader dynamics, and β is a scalar parameter appearing in these latter matrices. We require the entries of the projection matrix to be non-negative and that any general nonlinearities contained in these entries are sufficiently smooth. Specifically, we assume the following:
(H1): The entries of P are in , where
is an open interval and
is an open set, and are nonnegative for
and
.
The Jacobian matrix of system (EquationA2(A2)
(A2) ) is given by
where we use superscript to indicate whether the matrix corresponds to the resident (R) or the invader (I) equations and subscript to represent differentiation with respect to R or I. We assume that system (EquationA2
(A2)
(A2) ) possesses an invader-free non-trivial resident cycle of period
which we denote by
. That is, τ denotes a resident-only fixed point of the θ-composite system
We further assume that the stability of this resident-only cycle is not determined by the resident equations, namely:
(H2): the absolute value of the Jacobian of the dominant eigenvalue of is less than 1.
We consider what happens when the introduction of an invader species destabilizes this resident cycle. Specifically, we assume:
(H2): The matrix is non-negative and primitive and there exists a
such that the dominant eigenvalue of
equals 1. The corresponding right
and left
eigenvectors of
satisfy
, where superscript
denotes evaluation at the bifurcation point
. Moreover,
.
Remark A.1
In [Citation26], Meissen considered a general form of resident-invader interactions that encompasses cases such as two competing species or predator–prey interactions. However, that analysis assumes and, thus, does not cover host–parasitoid-type interactions. Here we extend the analysis presented in [Citation26] to allow for this matrix to be non-zero. Note that the assumption
means new invader individuals cannot be created if no invaders are present.
Theorem A.1 establishes the direction of bifurcation and stability of the branch of θ-cycles that bifurcate from the resident-only cycle at . This analysis closely follows that presented in [Citation26] (specifically, the proof of Theorem 5 in [Citation26]).
Theorem A.1
There exists a branch of solutions of the form (A3)
(A3) bifurcating from the resident-only cycle
of system (EquationA1
(A1)
(A1) ) at
. The bifurcation is forward if
and backward if
where
(A4)
(A4) Near the bifurcation, the cycle is locally asymptotically stable if
and unstable if
where
(A5)
(A5) Here
denote the right and left eigenvectors, respectively, of the Jacobian matrix
corresponding to the eigenvalue 1, and
and
denote the bottom left
block and bottom right
block, respectively, of the matrix
,
.
Proof.
The Jacobian of the θ-composite of the system (EquationA1(A1)
(A1) ) evaluated at τ can be written as
(A6)
(A6) From the block form of
, the right eigenvector of
corresponding to the eigenvalue 1 has the form
where the
vector
may contain negative components but the
vector
is the positive eigenvector of
corresponding to eigenvalue 1 which is guaranteed by the Perron–Frobenius theorem, i.e.
. The left eigenvector is given by
where the
vector
is the positive left eigenvector of
corresponding to eigenvalue 1, i.e.,
.
By Assumptions (H2) and (H3) the strictly dominant eigenvalue of is also the dominant eigenvalue of
. At
,
and hence
has a simple dominant eigenvalue of 1. If
, then following Lyapunov–Schmidt reduction which utilizes the Implicit Function Theorem and the Fredholm Alternative, we obtain the existence of a branch of solutions to
of the form
which bifurcates from
. (See the proof of this result in Theorem 1 of [Citation26] which is a generalization of a result provided in Appendix B.2 of [Citation7].) Note that since
, this branch exists for
. Therefore
(
) indicates a forward (backward) bifurcation.
To determine the direction of bifurcation we parameterize the variables in ϵ along the branch of bifurcating fixed points. Namely, we expand the state variable and the bifurcation parameter
, as given above, as well as the Jacobian
and the dominant eigenvalue
,
The eigenvalue expansion is guaranteed by Theorem 2 in [Citation26] (or see Theorem 5.13a in [Citation18]) and the expansion of
is guaranteed by the smoothness assumption in (H1).
We first compute the Taylor expansion of the right hand side of in terms of x and β near
. Recall that superscript 0 denotes evaluation at
We obtain
where
contains all terms of higher order than two in combination of x and
Note that terms 2 and 4 on the right-hand side evaluate to zero due to the structure of the zeros in
and
Next, we substitute the expansion of and
from (EquationA3
(A3)
(A3) ) into the Taylor expansion and equate orders of ϵ. At order
, we see that
. At order
, we obtain
which holds from the eigenvector equation. Meanwhile, at order
we have
We denote this as
where
and
.
By the Fredholm Alternative,
is solvable for u if
The structure of
along with that of
results in
. Solving for κ we find,
We use
to simplify the denominator. To simplify the numerator further, note that the form
means that only the bottom left
and bottom right
blocks of
, denoted
and
, respectively, affect κ. Notice that, since
, if the bifurcation parameter is strictly beneficial to the invader so that
, then the direction of bifurcation is determined by the numerator of κ.
Finally, we calculate to determine the stability of the branch near the bifurcation at
. Let
denote the eigenvector of
corresponding to
. The eigenvalue equation and its derivative with respect to ϵ are given by
Evaluated at
, these reduce to
The first equation holds since v is a right eigenvector of
corresponding to eigenvalue 1. The second equation is equivalent to
By the Fredholm alternative this equation has a solution if
Thus
We then compute
using
Now letting
, we write this as
Note that
and recall that
. Therefore, the numerator of
becomes
which results in
being given by Equation (EquationA5
(A5)
(A5) ).