Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Fundamental limits of parasitoid-driven host population suppression: Implications for biological control

  • Abhyudai Singh

    Roles Conceptualization, Formal analysis, Methodology, Visualization, Writing – original draft, Writing – review & editing

    absingh@udel.edu

    Affiliation Departments of Electrical and Computer Engineering, Biomedical Engineering and Mathematical Sciences, University of Delaware, Newark, DE, United States of America

Abstract

Parasitoid wasps are increasingly being used to control insect pest populations, where the pest is the host species parasitized by the wasp. Here we use the discrete-time formalism of the Nicholson-Bailey model to investigate a fundamental question—are there limits to parasitoid-driven suppression of the host population density while still ensuring a stable coexistence of both species? Our model formulation imposes an intrinsic self-limitation in the host’s growth resulting in a carrying capacity in the absence of the parasitoid. Different versions of the model are considered with parasitism occurring at a developmental stage that is before, during, or after the growth-limiting stage. For example, the host’s growth limitation may occur at its larval stage due to intraspecific competition, while the wasps attack either the host egg, larval or pupal stage. For slow-growing hosts, models with parasitism occurring at different life stages are identical in terms of their host suppression dynamics but have contrasting differences for fast-growing hosts. In the latter case, our analysis reveals that wasp parasitism occurring after host growth limitation yields the lowest pest population density conditioned on stable host-parasitoid coexistence. For ecologically relevant parameter regimes we estimate this host suppression to be roughly 10-20% of the parasitoid-free carrying capacity. We further expand the models to consider a fraction of hosts protected from parasitism (i.e., a host refuge). Our results show that for a given host reproduction rate there exists a critical value of protected host fraction beyond which, the system dynamics are stable even for high levels of parasitism that drive the host to arbitrary low population densities. In summary, our systematic analysis sheds key insights into the combined effects of density-dependence in host growth and parasitism refuge in stabilizing the host-parasitoid population dynamics with important implications for biological control.

1 Introduction

Parasitoids are often used as agents of biological control to manage insect pets in forest and agricultural ecosystems. The primary goal of such biological control is to suppress the pest population density through the introduction of natural enemies, such as parasitoids, and hence minimize or eliminate the usage of insecticides. This goal sets up an interesting tradeoff on the level of parasitism—while a certain minimal level is needed for parasitoid establishment, high levels of parasitism can destabilize the system eventually leading to parasitoid extinction [1]. This contribution uses population dynamic models to rigorously quantify these tradeoffs and determine the optimal parasitism levels that yield the lowest host population density while still ensuring the stable coexistence of both species.

There is a long-standing tradition of modeling host-parasitoid population dynamics using discrete-time models [16]. This is primarily motivated by populations living in the temperate regions of the world where annual insect life stages are synchronized by season and reproduction occurs at specific times in the year. A typical life cycle is illustrated in Fig 1 where adult hosts in a given year emerge and lay eggs that hatch into larvae. Host larvae feed and grow on resources and then pupate to metamorphosize as adults the following year. At one of its life stages, the host is parasitized by a parasitoid wasp (for example, in Fig 1 the wasps attack host eggs). Adult female parasitoids locate and oviposit an egg into the host that hatches into a juvenile parasitoid. The juvenile develops using the host’s body as a food source and goes on to transform into an adult parasitoid the next year, killing the parasitized host in the process [7, 8]. It is important to point out that the host is only needed for juvenile development, and adult parasitoids are free-living insects.

thumbnail
Fig 1. Schematic of the insect life cycles, with parasitoids attacking the host egg stage.

Self-limitation in host growth due to resource competition is assumed to occur at the larval stage and is captured by a density-dependent mortality rate. Bottom-left: Dynamical behaviors of the model (18) for k = K = 1, where for a given host reproduction rate R, increasing levels of parasitism results in—parasitoid extinction (low attack rate), stable coexistence of both the host and its parasitoid (intermediate attack rate), and bounded oscillations in population densities (high attack rate). Right: Representative population trajectories for the different dynamical scenarios with parameters taken as K = k = 1, R = 2, and cp = 0.49, 0.8, 1.7 in (18) from bottom to top.

https://doi.org/10.1371/journal.pone.0295980.g001

A simple model capturing the ecological dynamics of this interaction is given by (1a) (1b) where Ht and Pt are the adult host and parasitoid population densities in year t, respectively. If R > 1 is the number of eggs produced by each adult host, then RHt is the host density vulnerable to parasitism. In the text, R is referred to as the host reproduction rate. The escape response f(Pt) is the fraction of hosts escaping parasitism to become next year’s adult hosts with f(0) = 1. Following this model, RHt[1 − f(Pt)] is the density of parasitized hosts with k ≥ 1 adult female parasitoids developing per infected host. A classical form of (1) in the Nicholson-Bailey model [9] (2) where T is the duration of the host vulnerable stage, and cp > 0 is the parasitoid attack rate. One can interpret the attack rate as 1/cp being the average time taken by a parasitoid to forage, locate, and parasitize a given host. Key assumptions underlying the Nicholson-Bailey model is that the parasitoids are search-limited (but not egg-limited) and have fast handling times.

The Nicholson-Bailey model has a single non-trivial fixed point (3) that is unstable resulting in diverging population density oscillations [9]. Here and throughout the manuscript log represents the natural logarithm. A Type II functional response can be incorporated in (3) by replacing cp with (4) where Th is the handling time, and this has been shown to further destabilize the population dynamics [1, 10, 11]. A multitude of mechanisms are known to stabilize population dynamics and they can be classified into two types [10, 12]:

  1. The first class of mechanisms includes those stabilizing factors, where the escape response f(Pt) only depends on the parasitoid density. One such example is a proportional host refuge, where a fraction μ of hosts are protected from parasitism due to elevated host defenses or being inaccessible to parasitoids. In this case, the Nicholson-Bailey escape response is modified to (5) [13]. Other stabilizing mechanisms that fall in this category of a parasitoid-dependent f(Pt) include interference between parasitoids [14], host-to-host differences in susceptibility to parasitism [1520], and aggregation of parasitoid attacks on high-risk hosts [2123]. In all these cases, a necessary and sufficient condition for stable host-parasitoid coexistence is (6) i.e., the adult host equilibrium density (in the parasitoid’s presence) is an increasing function of R [24]. This is also reflected in the instability of the Nicholson-Bailey fixed point (3) where H* is a decreasing function of R.
  2. The second class of mechanisms includes a Type III functional response where the parasitoid attack rate increases (or accelerates) in response to higher host density [25, 26]. Such responses have been reported for several parasitoids [27, 28], and involve a change in behavior where the consumer is able to exploit the resource more efficiently at higher densities [29]. For example, parasitoids can have much faster handling times [30], or spend more time searching for hosts at higher host densities [31]. A key difference with the earlier case is that here the escape response depends on both the host and parasitoid densities and stability arises with the adult host density being a decreasing function of R [32].

In this contribution, we primarily focus on the case of a constant attack rate cp (as in the Nicholson-Bailey model), and later in the manuscript we investigate the impact of proportional host refuge as outlined in the first class of mechanisms.

The Nicholson-Bailey model has no intrinsic self-limitation in host growth with geometric expansion in host numbers (7) in the parasitoid’s absence. Previous works have shown the stabilizing effects of including density-dependent self-limitation in host growth [3338]. Motivated by this, we model the parasitoid-free host’s population dynamics as per the Beverton-Holt model (8) that has been previously reported in the context of intraspecific competitions [3941]. The parameter ch > 0 quantifies the strength of this competition and (8) has a stable fixed point (9) for all values of R, ch and K denotes the carrying capacity. Without loss of any generality, we assume that this growth limitation acts at the host larval stage and consider different parasitism scenarios:

  • Parasitoids attack host eggs—parasitism occurs at a stage before host growth limitation.
  • Parasitoids attack host larvae. This leads to two different models depending on whether intraspecific competition acts only on the unparasitized larvae or on all (unparasitized and parasitized) larvae.
  • Parasitoids attack host pupae—parasitism occurs at a stage after host growth limitation.

A key question driving this investigation is how much suppression in host density below its parasitoid-free carrying capacity is possible while still ensuring the stable coexistence of both species. Are there parameter regimes where the population dynamics of host suppression is invariant to the relative timing of parasitism with respect to host growth limitation? For regimes where this relative timing is critically important, which scenario provides the most efficient suppression of host density? Finally, we expand the study by coupling host growth limitation with an additional stabilizing factor of host refuge to investigate their combined effects on stable coexistence.

2 General model formulation and analysis

We start by reviewing a mechanistic derivation of the Beverton-Holt model using the semi-discrete approach, where update rules in discrete-time models are obtained by solving a system of continuous-time differential equations [4249]. Let L(t, τ) denote the host larval density in year t at time τ ∈ [0, 1] within the larval stage, where τ = 0 and τ = 1 correspond to the start and end of the stage, respectively. We consider a per capita density-dependent larval mortality rate chL(t, τ) that scales linearly with the population density and acts continuously throughout the stage. This mortality could be a result of intraspecific competition for resources or predation by natural enemies other than the parasitoid in consideration. Then, the larval density decays continuously as per the ordinary differential equation (10)

Solving (10) yields (11) and using the initial condition L(t, 0) = RHt together with Ht+1 = L(t, 1) (i.e., surviving larvae at the end of the stage become next year’s adults) results in the Beverton-Holt model (8).

In the parasitoid’s presence, (8) transforms to (12a) (12b) where f and g are continuously differentiable functions that depend on the host and parasitoid population densities with (13)

Apart from the trivial fixed point that excludes both species, the model’s fixed points are given by simultaneously solving (14) where H* and P* represent the host and parasitoid equilibrium densities, respectively. One of these fixed points corresponds to parasitoid extinction (P* = 0) and the host at its carrying capacity (H* = K). We are primarily interested in the existence of alternative fixed points that allow for the stable coexistence of both the host and the parasitoid—and directly connected to it—what is the lowest possible value of H*/K that quantifies the limit of parasitoid-mediated suppression of pest population density. We next present stability analysis tools for nonlinear systems of the form (12).

The stability of the fixed point can be assessed by linearizing the nonlinearities in (12) for small perturbations around the fixed point. This process yields the linear discrete-time dynamical system (15a) (15b) (15c) (15d) (15e) (15f)

The fixed point is stable, if and only if, all the following three conditions hold (16) [50], where (17) are the trace and determinant of the 2 × 2 Jacobian matrix A, respectively.

3 Parasitoids attack the host egg stage

We start with the scenario where the parasitoids parasitize the host eggs. As per the Nicholson-Bailey model, the unparasitized and parasitized host densities at the end of the egg stage are RHt exp(−cpPt) and RHt(1 − exp(−cpPt)), respectively. While the parasitized eggs become next year’s adult parasitoids, the unparasitized eggs hatch into larvae to face intraspecific competition and later develop into next year’s adult hosts. This results in the following discrete-time model describing the host-parasitoid population dynamics (18a) (18b) where (18a) corresponds to (11) with initial condition L(t, 0) = RHt exp(−cpPt). Here and in other models, for the sake of convenience, we assume the duration of the host vulnerable stage T = 1.

Standard stability analysis shows that the no-parasitoid fixed point (H* = K, P* = 0) is stable for (19)

One can also see from (18b) that a sufficient large attack rate is needed for parasitoid establishment ensuring population number increase from low densities (i.e., Pt+1/Pt > 1 when Ht = K and Pt → 0). When , there exists a unique fixed point corresponding to the coexistence of both species that is given as the solution to (20a) (20b)

Our analysis reveals that this fixed point is stable for (21) with attack rates above a critical value destabilizing the population dynamics. This critical value is obtained by solving (22) where P* in (22) is given by (20a) with . Beyond the system exhibits bounded oscillations in population densities whose amplitude and time period increases with increasing attack rate. The range (21) of attack rates allowing for stable coexistence is illustrated in Fig 1 along with the representative population trajectories corresponding to different dynamical outcomes:

  • Parasitoid extinction and the host at its carrying capacity for .
  • Stable host-parasitoid coexistence for .
  • Bounded oscillations in population densities for .

The plot in Fig 1 showing the stability region was generated in Wolfram Mathematica, where for a given parameter set, was obtained by numerically solving (22) (see S1 File that has the Wolfram Mathematica code used for generating the stability region). The population density trajectories were plotted in Microsoft Excel by iteratively solving the discrete-time model (18). The lowest possible stable suppression of host density occurs when and is given by (23) where z* is the unique solution to (24)

It is interesting to note that this limit of biological control only depends on R, and in this case, it increases with R varying from 40%(R = 2) to 55%(R = 10). We next contrast this result with parasitism occurring at the host larval or pupal stage.

4 Parasitoids attack the host larval stage

When parasitoids parasitize the host larvae, both parasitism and density-dependent mortality from intraspecific competition occurs concurrently and continuously throughout the stage. Our previous work analyzed this case assuming that intraspecific competition only acts on the unparasitized larvae [42]. We review these results and also consider the scenario where both unparasitized/parasitized larvae compete for resources.

4.1 Density-dependent mortality on unparasitized larvae

Using the semi-discrete formalism for the mechanistic derivation of discrete-time models, the continuous changes in population densities are described by the ordinary differential equations (25a) (25b) where L(t, τ) and I(t, τ) are the unparasitized and parasitized larval densities in year t at time τ ∈ [0, 1] within the stage. Solving (25) with L(t, 0) = RHt, I(t, 0) = RHt yields the model (26a) (26b)

Here, a minimum attack rate (27) is needed for parasitoid establishment. Given this sufficiently large attack rate, there exists a unique fixed point (28) where both species are present. This coexistence equilibrium is stable for attack rates in the range (29) where is the solution to [42]. This range (29) is shown in Fig 2, with lower and higher values of cp resulting in parasitoid extinction and bounded oscillation, respectively.

thumbnail
Fig 2. Difference in population dynamics for larval vs. pupal parasitoids.

The top plots show the life cycles of a larval (left) and pupal (right) parasitoid with host growth limitation occurring at the larval stage. The bottom plots show the corresponding range of stabilizing attack rates given by (29) and (35) as a function of the host reproduction rate R. While the region of stability contracts with increasing R for a larval parasitoid (left), it expands for a pupal parasitoid (right).

https://doi.org/10.1371/journal.pone.0295980.g002

The largest attack rate allowing a stable coexistence yields the lowest host density (30) corresponding to a larval parasitoid. In contrast to egg parasitoids, this limit (30) slightly increases close to R ≈ 1 (Fig 3), and then decreases with R varying from 35%(R = 2) to 30%(R = 10).

thumbnail
Fig 3. Minimal possible host population density conditioned on stable host-parasitoid coexistence.

Plots of H*/K (host density with parasitoid normalized to its parasitoid-free carrying capacity) as a function of the host reproduction rate R as given by (23) for an egg parasitoid (host growth limitation after parasitism), (30) for a larval parasitoid (host growth limitation occurring concurrently with parasitism), and (37) for a pupal parasitoid (host growth limitation before parasitism). All cases consider Beverton-Holt type host population dynamics (8) in the absence of the parasitoid with density-dependent self-limitation in growth assumed to occur at the host larval stage.

https://doi.org/10.1371/journal.pone.0295980.g003

4.2 Density-dependent mortality on all larvae

When both unparasitized and parasitized larvae have density-dependent mortality then (25) is altered to (31a) (31b) yielding the model (32a) (32b)

Incidentally, this is the exact same model as that for a pupal parasitoid. To see this, note that after intraspecific competition RHt/(1 + chRHt) is the host density at the start of the pupal stage. While fraction exp(−cpPt) of host pupae escape parasitism to metamorphose into adult hosts, the other fraction 1 − exp(−cpPt) are parasitized.

5 Parasitoids attack the host pupal stage

Population dynamics of a pupal parasitoid as described by the model (32) shares qualitatively similar dynamical outcomes to its counterparts (18) and (26), but with contrasting quantitative parameter ranges. The establishment of a pupal parasitoid requires a minimal attack rate (33) that ensures a unique equilibrium (34a) (34b) where both species are present. This fixed point is stable for attack rates in the range (35) where is given by (36)

Plotting (35) as a function of R in Fig 3 reveals key differences with previous cases:

  • The range of stabilizing attack rates is much broader for a pupal parasitoid (Fig 2; right) as compared to an egg (Fig 1) or larval parasitoid (Fig 2; left).
  • The range (35) expands with R—a fixed lower bound and an upper bound that increases with increasing R (Fig 2; right).
  • This expanding stability region is in sharp contrast to when parasitism occurs in earlier life stages, where with increasing R both the lower/upper bound of stabilizing attack rate decreases with a contracting region of stability (Figs 1 & 2).

The lowest host population density corresponding to is given by (37) where z* is the unique solution to (24), and this limit is found to sharply decrease from 20%(R = 2) to 10%(R = 10) (Fig 3).

6 Inclusion of host refuge

As mentioned in the Introduction, diverse ecological mechanisms are known to stabilize the Nicholson-Bailey model. Here we consider one such mechanism, where a fixed fraction 0 ≤ μ < 1 of hosts are protected from parasitism [13]. For example, these hosts could be in specific locations that are less accessible to parasitoids or have an elevated immune response against parasitism [5154]. For example, the host Bactrocera dorsalis (oriental fruit fly) is parasitized by the pupal parasitoid Dirhinus giffardii. B. dorsalis larvae pupate below ground, and data shows that pupation depth determines the risk of parasitism—pupae further underground experience much lower rates of parasitism compared to pupae closer to the surface [55]. We focus on how such refuges work concertedly with the host’s growth limitation to stabilize the host-parasitoid population dynamics.

6.1 Without host self-limitation

We first consider a host refuge fraction μ in the Nicholson-Bailey model resulting in (38a) (38b) where (39) is the fraction of hosts escaping parasitism [13]. Note that here μ > 1/R would result in unbounded population growth, and this is prevented in the next section by including a carrying capacity. When μ < 1/R, (38) has a non-trivial fixed point given by (40)

Since in (39) the fraction of host escaping parasitism only depends on Pt, (40) is stable iff (6) holds. Substituting H* from (40) in (6) it can be seen that there exists a critical refuge μ* that only depends on R and is given as the unique solution to (41) such that the fixed point (40) is stable when (42) for all values of cp and k. While μ* cannot be analytically solved from (41), when R ≈ 1, then using (43) we can approximate (44) in (41) to solve for μ* yielding (45)

The actual value of μ* as obtained by numerically solving (41) is shown in Fig 4 (bottom orange line in the top-left plot). While as predicted by (45), μ* → 0.5 as R → 1, μ* decreases slower with R than as predicted by (45) with μ* = 0.3 (R = 2) and μ* = 0.075 (R = 10). Outside the stability region a weak refuge (μ < μ*) results in bounded oscillations, while a strong refuge (μ > 1/R) leads to unbounded growth in population densities.

thumbnail
Fig 4. The impact of host refuge on host-parasitoid population dynamics.

In the absence of any host growth limitation, stable host-parasitoid coexistence occurs for a range of refuge fractions as given by (42). Outside this range, a weak refuge results in bounded oscillations in the population densities, and a strong refuge leads to unbounded population growth. Top-right to bottom-left to bottom-right: In the presence of host growth limitation, the range of refuge fractions (49) allowing stable coexistence is shown with decreasing carrying capacity K. Here a strong refuge causes parasitoid extinction and this region expands as K is lowered. The parameter space corresponding to bounded oscillation shrinks and vanishes for low-carrying capacities. For this plot, other parameters are fixed as cp = k = 1 and K = 10, 4.5, 2.

https://doi.org/10.1371/journal.pone.0295980.g004

6.2 With host self-limitation

To prevent unbounded population growth we now include a carrying capacity. For this, we consider the previously analyzed case of the pupal parasitoid and modify model (32) to include a fixed host refuge fraction μ (46a) (46b)

Not surprisingly, with some hosts protected from wasps a higher attack rate (47) is now needed for parasitoid establishment as compared with (33). With parasitoid establishment, there exists a unique equilibrium (48a) (48b) that is stable for a range of host refuge (49)

The lower limit of host refuge need for stability satisfies , and for sufficiently large K. These results imply that μ > μ* and a large enough attack rate (47) ensuring parasitoid establishment is sufficient for stable host-parasitoid population dynamics.

The stability region (49) is illustrated in Fig 4 with varying carrying capacity. Since strong limitation in the host’s growth can by itself stabilize the host-parasitoid interaction in the absence of refuge, the lower limit with decreasing K (Fig 4; bottom right). These results have important consequences for biological control, where in the absence of host refuge a high attack rate destabilizes the coexistence equilibrium (Fig 5; left), but in the presence of host refuge μ > μ*, the system remains stable for all values of cp with high parasitism levels driving the host population to arbitrarily low levels (Fig 5; right).

thumbnail
Fig 5. Increasing attack rate destabilizes host-parasitoid population dynamics.

Trajectories of host population densities with increasing parasitoid attack rates cp (represented by darker shades of blue) as obtained from model (32) (without host refuge; left) and (46) (with host refuge; right). While in the absence of any refuge, high attack rates destabilize the population dynamics (left), with host refuge the population dynamics remain stable providing stably low suppression of the host population. Other parameters are taken as K = k = 1, R = 2, μ = 0 (left) and μ = 0.2 (right).

https://doi.org/10.1371/journal.pone.0295980.g005

7 Discussion

We have analyzed a class of models where host larval density decreases throughout the stage as a result of density-dependent mortality due to intraspecific competition or parasitism. The per capita mortality rate is proportional to the host density resulting in the Beverton-Holt model describing the host population dynamics in the absence of the parasitoid [42, 56]. This is a key difference from previous work, where host population dynamics follows the Ricker model [33]. Coupled with the Beverton-holt model is the parasitoid population dynamics leading to three different models—(18) for an egg parasitoid, (26) for a larval parasitoid, (32) for a pupal parasitoid. The larval parasitoid here refers to the scenario where density-dependent mortality acts only on the unparasitized larvae, as when mortality acts on both parasitized and unparasitized larvae, then the ecological population dynamics are identical for larval and pupal parasitoids. This latter scenario is likely when larval mortality is a result of predation by other natural enemies that do not differentiate between the two types of larvae.

We systematically analyzed these host-parasitoid models in the context of biological control of insect pests, where the primary goal is to suppress their population density via introductions of natural enemies, such as parasitoids [5763]. Such pests can be attacked by multiple parasitoid species that could parasitize different host life-cycle stages. For example, the European corn borer Ostrinia nubilalis (a major pest of grains) has parasitoids that attack the egg stage (Trichogramma ostriniae) [64], and the larval stage (Macrocentrus grandii and Lydella thompsoni) [65]. Motivated by this we specifically contrasted population dynamic models where hosts are vulnerable to parasitoids in different developmental stages. All models (18), (26) and (32) (corresponding to parasitism of egg, larval, and pupal host stage, respectively), share a common feature of destabilized population dynamics for sufficiently large parasitoid attack rates. This can be intuitively understood from the fact that high parasitism levels drive the host density significantly below the carrying capacity, reducing (18), (26) and (32) to the unstable Nicholson-Bailey model. We quantify the limit of host density suppression across these models using H*/K—the ratio of host equilibrium density (just before stability is lost for high parasitoid attack rates) and the host’s parasitoid-free carrying capacity. Before summarizing this limit we discuss the minimal levels of parasitism needed for parasitoid establishment.

To assess the potential for parasitoid establishment, one should consider the product of the four dimensionless parameters cpkKT. This product combines the following terms: the parasitoid attack rate (per unit time per host per parasitoid), T (duration of host vulnerable stage that we assumed to be one time units earlier), k (number of parasitoids per parasitized host), K (host carrying capacity). The later the parasitism occurs in the host’s life cycle, the higher this parameter needs to be for parasitoids to grow from small densities and establish. More specifically, parasitoid establishment requires (50a) (50b) (50c)

Since k = K = T = 1 in Figs 1 & 2, the y-axis on the stability region plots can be interpreted as this dimensionless parameter cpkKT. An interesting finding from our analysis is that when R ≈ 1, the range of cpkKT allowing stable host-parasitoid coexistence is identical in all three models and given by (51)

Here the lower limit corresponds to parasitoid establishment and is obtained by taking R → 1 in (50). Crossing the upper limit cpkKT > 3 destabilizes the coexistence resulting in limit cycles (Figs 1 & 2). In terms of host suppression, for R ≈ 1 the lowest possible value of H*/K = 1/3 occurs when cpkKT = 3 (Fig 3). Thus, for slow-growing host populations, the timing of parasitism may not have an appreciable impact on the population dynamics with a limit of host suppression that is 33% of the carrying capacity.

With increasing R, we see contrasting differences in the range of stabilizing values of cpkKT, with this range contracting for an egg and larval parasitoids, but expanding for a pupal parasitoid (Figs 1 and 2). These differences directly impact the host density just before stability is again lost for high attack rates with H*/K ≈ 0.55 (egg parasitoid), H*/K ≈ 0.3 (larval parasitoid), and H*/K ≈ 0.1 (pupal parasitoid) for R = 10 (Fig 3). It is important to point out that for a fixed attack rate, the host density is the lowest for an egg parasitoid. However, for a pupal parasitoid, the coexistence equilibrium remains stable for a much broader range and larger values of attack rates leading to a lower stable host density as compared to an egg parasitoid.

How does the lower limit of H*/K depend on the form of density-dependent self-limitation in host growth? To see this, we consider a different model for host population dynamics (52) where q = 1 is the Beverton-Holt model and we consider q = 2. For values of q ∈ {1, 2}, (52) has a stable equilibrium H* = K for all R > 1. When q = 1, the lowest stable suppression of host density varies from ≈20% (R = 2) to ≈10% (R = 10) in Fig 3. Our simulation results show that these limits increase to ≈40% (R = 2) to ≈20% (R = 10) for q = 2 suggesting quantitative differences in host suppression capabilities depending on the form of host growth limitation, even though the parasitoid attack rates needed for establishment are the same for q = 1 & 2.

We next expanded these results to consider a fraction of host refuge 0 ≤ μ < 1. Consistent with previous analysis [13], we find that in the absence of host growth limitation, stable coexistence arises in a small range of refuge fractions (Fig 4; top-left). Our contribution here is to show that this range is approximated by (53) implying close to 50% protection for slow-growing hosts is needed for stability. The range of stabilizing μ shrinks with increasing R and is given by 0.075 < μ < 0.1 for R = 10 (Fig 4). Our results show that with the inclusion of host growth limitation, this range as given by (49) dramatically increases for medium/high carrying capacities (Fig 4) and expands with increasing R. The dimensionless parameters cpkKT needed for parasitoid establishment is also higher by a factor of 1/(1 − μ) with respect to (50). Since host refuge can stabilize the population dynamics even in the absence of any host growth limitation, above the critical refuge μ* < μ the system remains stable even for high attack rates that drive the host to arbitrary low densities (Fig 5). Hence, the lower limit of H*/K decreases with increasing μ and reaches zero for μ > μ*.

In summary, our investigation reveals key insights into how the relative sequence of parasitoid attack and density-dependent host growth limitation impact the overall population dynamics. For the ecologically relevant scenario of slow-growing hosts, the system dynamics become invariant to the specific timing of parasitism. Our results coupling host refuge with the Beverton-Holt model reveal large regions of parameter space allowing stable coexistence (Fig 4) consistent with field studies implicating proportional refuges in stabilizing host-parasitoid interactions [6670]. Future work will extend this study in several directions, such as incorporating a Type II functional response implemented using the semi-discrete approach, exploring how spatial effects alter the fundamental limits of biological control [7175], and considering multiple parasitoid species attacking different host stages. Related to the last point, recent work has provided novel conditions for the coexistence of multiple parasitoids on a single host [76], and in many cases the requirement (6) of adult equilibrium host density increasing with R emerges as a universal criterion for population stability in these complex ecological communities.

Another direction of future work would be to consider Allee effects in both the host and the parasitoid as done in recent population dynamics models [77, 78]. Inclusion of an Allee effect in the host could be interesting as it could lead to two different scenarios at high parasitism levels- one where the host is driven to extinction and the other where the population dynamics is destabilized before the Allee effect comes into play leading to parasitoid extinction. Inclusion of an Allee effect in the parasitoid should result in higher attack rates needed for parasitoid establishment than as predicted by (50), but should not impact the fundamental limits of host density suppression where the parasitoid density is relatively high. Finally, it is important to point out that while here we have used the traditional discrete-time formalism of the Nicholson-Bailey model to investigate the limits of host density suppression, it will be interesting to extend these studies with the continuous-time framework of Lotka-Volterra which is more appropriate for modeling population dynamics of insects in the tropics [1, 11].

Supporting information

S1 File. Mathematica file.

Wolfram Mathematica code used for generating the stability regions shown in Figs 1 and 2, and the limit of host suppression in Fig 3.

https://doi.org/10.1371/journal.pone.0295980.s001

(PDF)

References

  1. 1. Murdoch WW, Briggs CJ, Nisbet RM. Consumer-Resouse Dynamics. Princeton,NJ: Princeton University Press; 2003.
  2. 2. Singh A. Stochastic dynamics of predator-prey interactions. PLoS One. 2021;16(8):e0255880. pmid:34383811
  3. 3. Hassell MP. The Dynamics of Arthopod Predator-Prey Systems.(MPB-13), Volume 13. vol. 111. Princeton University Press; 2020.
  4. 4. Jervis MA, Kidd NAC, Mills NJ, van Nouhuys S, Singh A, Yazdani M. In: Hardy ICW, Wajnberg E, editors. Population Dynamics. Cham: Springer International Publishing; 2023. p. 591–667. Available from: https://doi.org/10.1007/978-3-031-23880-2_7.
  5. 5. Hassell MP. The Spatial and Temporal Dynamics of Host Parasitoid Interactions. New York, NY: Oxford University Press; 2000.
  6. 6. Schreiber SJ. Host-parasitoid dynamics of a generalized Thompson model. Journal of Mathematical Biology. 2006;52(6):719–732. pmid:16699833
  7. 7. Godfray HCJ. Parasitoids; Behavioral and Evolutionary Ecology. 41 William St, Princeton, NJ 08540: Princeton University Press; 1994.
  8. 8. Waage JK, Greathead D. Insect Parasitoids. Academic Press; 1986.
  9. 9. Nicholson AJ, Bailey VA. The balance of animal populations. Part 1. Proc of Zoological Society of London. 1935;3:551–598.
  10. 10. Singh A, Emerick B. Generalized stability conditions for host–parasitoid population dynamics: Implications for biological control. Ecological Modelling. 2021;456:109656.
  11. 11. Singh A. A comparative approach to stabilizing mechanisms between discrete-and continuous-time consumer-resource models. Plos one. 2022;17(4):e0265825. pmid:35413067
  12. 12. Singh A. Stochasticity in host-parasitoid models informs mechanisms regulating population dynamics. Scientific Reports. 2021;11(1):16749. pmid:34408234
  13. 13. Hassell MP, May RM. Stability in insect host-parasite models. Journal of Animal Ecology. 1973;42(3):693–726.
  14. 14. Hassell MP, Varley GC. New Inductive Population Model for Insect Parasites and its Bearing on Biological Control. Nature. 1969;223:1133–1137. pmid:5810684
  15. 15. Taylor AD. Heterogeneity in host-parasitoid interactions: ‘aggregation of risk’ and the ‘CV2 > 1 rule’. Trends in Ecology and Evolution. 1993;8:400–405. pmid:21236211
  16. 16. Hassell MP, May RM, Pacala SW, Chesson PL. The persistence of host–parasitoid associations in patchy environments. I. A general criterion. American Naturalist. 1991;138:568–583.
  17. 17. Pacala SW, Hassell MP. The persistence of host– parasitoid associations in patchy environments. II. Evaluation of field data. American Naturalist. 1991;138:584–605.
  18. 18. Cobbold CA, Roland J, Lewis MA. The impact of parasitoid emergence time on host-parastioid population dynamics. Theoretical Population Biology. 2009;75(2):201–215. pmid:19275911
  19. 19. Liere H, Jackson D, Vandermeer J. Ecological complexity in a coffee agroecosystem: spatial heterogeneity, popoulation persistence and biological control. PLoS One. 2012;7(9).
  20. 20. Zoroa N, Lesigne E, Fernandez-Saez MJ, Zoroa P, Casas J. The coupon collector urn model with unequal probabilities in ecology and evolution. Journal of The Royal Society Interface. 2017;14(127). pmid:28179550
  21. 21. May RM. Host–parasitoid systems in patchy environments: a phenomenological model. Journal of Animal Ecology. 1978;47:833–844.
  22. 22. Hassell MP, May RM. Aggregation of predators and insect parasites and its effect on stability. Journal of Animal Ecology. 1974;43(2):567–594.
  23. 23. Rohani P, Godfray H, Hassell M. Aggregation and the dynamics of host-parasitoid systems: a discrete-generation model with within-generation redistribution. The American Naturalist. 1994;144(3):491–509.
  24. 24. Singh A, Murdoch WW, Nisbet RM. Skewed attacks, stability, and host suppression. Ecology. 2009;90(6):1679–1686. pmid:19569382
  25. 25. Singh A. Attack by a common parasitoid stabilizes population dynamics of multi-host communities. Journal of Theoretical Biology. 2021;531:110897. pmid:34506808
  26. 26. Hassell MP, Comins HN. Sigmoid functional responses and population stability. Theoretical Population Biology. 1978;14:62–66. pmid:741397
  27. 27. Hassell M, Lawton J, Beddington J. Sigmoid functional responses by invertebrate predators and parasitoids. The Journal of Animal Ecology. 1977; p. 249–262.
  28. 28. Fernández-arhex V, Corley JC. The functional response of parasitoids and its implications for biological control. Biocontrol Science and Technology. 2003;13(4):403–413.
  29. 29. Kalinkat G, Rall BC, Uiterwaal SF, Uszko W. Empirical evidence of type III functional responses and why it remains rare. Frontiers in Ecology and Evolution. 2023;11:1033818.
  30. 30. Collins M, Ward S, Dixon A. Handling time and the functional response of Aphelinus thomsoni, a predator and parasite of the aphid Drepanosiphum platanoidis. The Journal of Animal Ecology. 1981; p. 479–487.
  31. 31. Takahashi F. Functional response to host density in a parasitic wasp, with reference to population regulation. Population Ecology. 1968;10(1):54–68.
  32. 32. Singh A, Nisbet RM. Semi-discrete host-parasitoid models. Journal of Theoretical Biology. 2007;247(4):733–742. pmid:17524428
  33. 33. May RM, Hassell MP, Anderson RM, Tonkyn DW. Density dependence in host-parasitoid models. Journal of Animal Ecology. 1981;50:855–865.
  34. 34. Smith JM, Slatkin M. The stability of predator-prey systems. Ecology. 1973;54(2):384–391.
  35. 35. Beddington J, Free C, Lawton J. Dynamic complexity in predator-prey models framed in difference equations. Nature. 1975;255(5503):58–60.
  36. 36. Beddington J, Free C, Lawton J. Concepts of stability and resilience in predator-prey models. The Journal of Animal Ecology. 1976; p. 791–816.
  37. 37. Marcinko K, Kot M. A comparative analysis of host-parasitoid models with density dependence preceding parasitism. Journal of Biological Dynamics. 2020;14(1):479–514. pmid:32603259
  38. 38. Jang SRJ, Yu JL. Discrete-time host–parasitoid models with pest control. Journal of Biological Dynamics. 2012;6(2):718–739. pmid:22873614
  39. 39. Geritz SA, Kisdi E. On the mechanistic underpinning of discrete-time population models with complex dynamics. Journal of Theoretical Biology. 2004;228(2):261–269. pmid:15094020
  40. 40. Brännström Å, Sumpter DJ. The role of competition and clustering in population dynamics. Proceedings of the Royal Society B: Biological Sciences. 2005;272(1576):2065–2072. pmid:16191618
  41. 41. Beverton RJ, Holt SJ. On the dynamics of exploited fish populations. vol. 11. Springer Science & Business Media; 2012.
  42. 42. Singh A, Nisbet RM. Variation in risk in single-species discrete-time models. Mathematical Biosciences and Engineering. 2008;5:859–875. pmid:19278287
  43. 43. Pachepsky E, Nisbet RM, Murdoch WW. Between discrete and continuous: Consumer-resource dynamics with synchronized reproduction. Ecology. 2007;89(1):280–288.
  44. 44. Eskola TM, Geritz SA. On the mechanistic derivation of various discrete-time population models. Bulletin of Mathematical Biology. 2007;69:329–346. pmid:16838083
  45. 45. Dugaw CJ, Hastings A, Preisser EL, Strong DR. Seasonally limited host supply generates microparasite population cycles. Bulletin of Mathematical Biology. 2004;66(3):583–594. pmid:15006450
  46. 46. Bonsall MB, Hassell MP. Parasitoid-mediated effects: apparent competition and the persistence of host-parasitoid assemblages. Researches on Population Ecology. 1999;41(1):59–68.
  47. 47. Briggs CJ, Godfray HCJ. The dynamics of insect-pathogen interactions in seasonal environments. Theoretical Population Biology. 1996;50:149–177. pmid:8955031
  48. 48. Emerick BK, Singh A. The effects of host-feeding on stability of discrete-time host-parasitoid population dynamic models. Mathematical Biosciences. 2016;272:54–63. pmid:26686008
  49. 49. Singh A, Emerick B. Hybrid systems framework for modeling host-parasitoid population dynamics. In: 2020 59th IEEE Conference on Decision and Control (CDC). IEEE; 2020. p. 4628–4633.
  50. 50. Elaydi S. An Introduction to Difference Equations. New York: Springer; 1996.
  51. 51. Beckage NE. Modulation of immune responses to parasitoids by polydnaviruses. Parasitology. 1998;116(S1):S57–S64. pmid:9695110
  52. 52. Edson KM, Vinson SB, Stoltz DB, Summers MD. Virus in a parasitoid wasp: suppression of the cellular immune response in the parasitoid’s host. Science. 1981;211(4482):582–583. pmid:7455695
  53. 53. Smilanich AM, Dyer LA, Gentry GL. The insect immune response and other putative defenses as effective predictors of parasitism. Ecology. 2009;90(6):1434–1440. pmid:19569356
  54. 54. Strand MR, Pech LL. Immunological basis for compatibility in parasitoid-host relationships. Annual review of entomology. 1995;40(1):31–56. pmid:7810989
  55. 55. Okuyama T. Density-dependent distribution of parasitism risk among underground hosts. Bulletin of Entomological Research. 2019;109(4):528–533. pmid:30457061
  56. 56. Gurney WSC, Nisbet RM. Ecological Dynamics. Oxford University Press; 1998.
  57. 57. Abram PK, Brodeur J, Burte V, Boivin G. Parasitoid-induced host egg abortion; An underappreciated component of biological control services provided by egg parasitoids. Biological Control. 2016;(98):52–60.
  58. 58. Jervis MA, Hawkin BA, Kidd NAC. The usefulness of destructive host-feeding parasitoids in classical biological control: theory and observation conflict. Ecological Entomology. 1996;21(1):41–46.
  59. 59. Ueno T. Selective host-feeding on parasitized hosts by the parasitoid Itoplectis naranyae (Hymenoptera: Ichneumonidae) and its implication for biological control. Bulletin of Entomological Research. 1998;88(4):461–466.
  60. 60. Reeve JD, Murdoch WW. Aggregation by parasitoids in the successful control of the california red scale: a test of theory. Journal of Animal Ecology. 1985;54(3):797–816.
  61. 61. Jang SR, Yu JL. Discrete-time host-parasitoid models with pest control. Journal of Biological Dynamics. 2012;6(2):718–739. pmid:22873614
  62. 62. Hassell MP, Varley GC. New inductive population model for insect and its bearing on biological control. Nature. 1969;223(1):1133–1137. pmid:5810684
  63. 63. Kaser JM, Nielsen AL, Abram PK. Biological control effects of non-reproductive host mortality caused by insect parasitoids. Ecological Applications. 2018;28(4):1081–1092. pmid:29485221
  64. 64. Hoffmann M. Early-Season Establishment of Trichogramma ostriniae for Season-Long Suppression of European Corn Borer in Sweet Corn. Progress Report, New York State Integrated Pest Management Program, 1998.
  65. 65. Romig R, Mason C, Burbutis P. Parasitism of European corn borer by Lydella thompsoni and Macrocentrus grandii in southeast Pennsylvania and Delaware. Entomol News. 1985;96:121–128.
  66. 66. Vinson SB, Iwantsch G. Host suitability for insect parasitoids. Annual Review of Entomology. 1980;25(1):397–419.
  67. 67. Hawkins B, Browning H, Smith J. Field evaluation of Allorhogas pyralophagus [Hym.: Braconidae], imported into Texas for biological control of the stalkborer Eoreuma loftini [Lep.: Pyralidae] in sugar cane. Entomophaga. 1987;32:483–491.
  68. 68. Murdoch WW, Luck RF, Walde SJ, Reeve JD, Yu DS. A refuge for red scale under control by Aphytis: structural aspects. Ecology. 1989;70(6):1707–1714.
  69. 69. Hochberg ME, Hawkins BA. Refuges as a predictor of parasitoid diversity. Science. 1992;255(5047):973–976. pmid:17793160
  70. 70. Reeve JD, Cronin JT, Strong DR. Parasitoid aggregation and the stabilization of a salt marsh host– parasitoid system. Ecology. 1994;75:288–295.
  71. 71. Rohani P, Miramontes O. Host-parasitoid metapopulations: the consequences of parasitoid aggregation on spatial dynamics and searching efficiency. Proceedings of the Royal Society B: Biological Sciences. 1995;260:335–342.
  72. 72. Cronin JT, Reeve JD. Host-parasitoid spatial ecology: A plea for a landscape-level synthesis. Proceedings: Biological Sciences. 2005;272(1578):2225–2235. pmid:16191634
  73. 73. Adler FR. Migration alone can produce persistence of host-parasitoid models. The American Naturalist. 1993;141(4):642–650.
  74. 74. Comins HN, Hassell MP, May RM. The spatial dynamics of host-parasitoid systems. Journal of Animal Ecology. 1992;61(3):735–748.
  75. 75. Emerick B, Singh A, Chhetri SR. Global redistribution and local migration in semi-discrete host–parasitoid population dynamic models. Mathematical Biosciences. 2020;327:108409. pmid:32615211
  76. 76. Singh A, Emerick B. Coexistence conditions in generalized discrete-time models of insect population dynamics. Ecological Modelling. 2022;474:110148.
  77. 77. Bompard A, Amat I, Fauvergue X, Spataro T. Host-parasitoid dynamics and the success of biological control when parasitoids are prone to Allee effects. PLoS One. 2013;8(10):233–253. pmid:24116153
  78. 78. Livadiotis G, Assas L, Dennis B, Elaydi S, Kwessi E. A discrete-time host-parasitoid model with Allee effect. Journal of Biological Dynamics. 2015;9:34–51. pmid:25431970