A comparative approach to stabilizing mechanisms between discrete- and continuous-time consumer-resource models

There is rich literature on using continuous-time and discrete-time models for studying population dynamics of consumer-resource interactions. A key focus of this contribution is to systematically compare between the two modeling formalisms the stabilizing/destabilizing impacts of diverse ecological processes that result in a density-dependent consumer attack rate. Inspired by the Nicholson-Bailey/Lotka-Volterra models in discrete-time/continuous-time, respectively, we consider host-parasitoid interactions with an arbitrary parasitoid attack rate that is a function of both the host/parasitoid population densities. Our analysis shows that a Type II functional response is stabilizing in both modeling frameworks only when combined with other mechanisms, such as mutual interference between parasitoids. A Type III functional response is by itself stabilizing, but the extent of attack-rate acceleration needed is much higher in the discrete-time framework, and its stability regime expands with increasing host reproduction. Finally, our results show that while mutual parasitoid interference can stabilize population dynamics, cooperation between parasitoids to handle hosts is destabilizing in both frameworks. In summary, our comparative analysis systematically characterizes diverse ecological processes driving stable population dynamics in discrete-time and continuous-time consumer-resource models.


I. Introduction
Interaction between a resource (such as, a prey or host) and a consumer (such as, a predator or parasitoid) forms a core motif in ecological networks. Population dynamics of consumerresource interactions has been extensively studied using two different approaches: continuoustime and discrete-time models. The continuous-time framework is generally used to model populations with overlapping generations and all year-round reproduction. In contrast, discrete-time models are more suited for populations with non-overlapping generations that reproduce in a discrete pulse determined by season. Perhaps, the classical example of this in a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 continuous time is the Lotka-Volterra model dhðtÞ dt ¼ rhðtÞ À chðtÞpðtÞ ð1aÞ dpðtÞ dt ¼ chðtÞpðtÞ À gpðtÞ ð1bÞ [1][2][3][4][5][6][7][8][9][10]. Here h(t) and p(t) denote the population densities of the host and the parasitoid at time t, respectively. Parasitoids attack hosts with rate c with each parasitized host developing into a new parasitoid. Finally, 1/γ is the average lifespan of an individual parasitoid. It is wellknown that the steady-state equilibrium of the Lotka-Volterra model is neutrally stable resulting in cycling population densities with a period of � 2p= ffi ffi ffi ffi rg p [1].
The analogous counterpart of the Lotka-Volterra model in discrete-time is the Nicholson-Bailey model where we now use capital letters H t and P t to denote the adult host, and the adult parasitoid densities, respectively, in year t [4,[11][12][13][14][15][16]. The model reflects the synchronized annual life cycles of these insects living in the temperate regions of the world [17][18][19][20]. More specifically, if R > 1 denotes the number of viable eggs produced per adult host, then RH t is the host larval density that becomes exposed to parasitoid attacks. Parasitoids attack and parasitize hosts at a constant rate c, resulting in the fraction exp(−cP t ) that escapes parasitism to become the adult hosts for the next year. Similarly, the fraction 1 − exp(−cP t ) of parasitized larvae give rise to adult parasitoids in the next generation. The Nicholson-Bailey model is characterized by diverging oscillations in population densities resulting in an unstable population dynamics [11]. The neutrally stable Lokta-Volterra equilibrium and the unstable Nicholson-Bailey equilibrium catalyzed rich theoretical work in understanding how individual ecological mechanisms (such as a Type II functional response, interference between parasitoid in attacking hosts, etc.) promote stability of host-parasitoid interactions, especially in the context of biological control of pest species [1,[21][22][23][24]. A key contribution of this work is to study the combined effect of different forms of density-dependence in the parasitoid attack rate and quantify their overall stabilizing effect in a single stability criterion. This gives us a holistic picture of how a mixture of processes impact population dynamics. For example, while a Type II functional response is destabilizing in both the continuous-and discrete-time frameworks, combining it with some form of parasitoid interference can lead to stable population dynamics. Another novelty of this work is the comparison of stability regimes between the modeling formalisms, and how these stability regimes change with relevant ecological parameters, such as the extent of host reproduction. This comparison is specifically important to elucidate the general forms of densitydependence in the parasitoid attack rate that stabilizes population dynamics irrespective of the modeling framework, and hence are robust to model choice. Such general features are particularly relevant for actual populations that may not exactly follow the idealistic modeling assumption of either Lotka-Volterra or Nicholson-Bailey formalisms and may be best described by an intermediate hybrid modeling framework.
The manuscript is organized as follows: in Section II we formulate the generalized Lotka-Volterra model and derive criteria for stable population dynamics. The same process is repeated in Section III for an analogously formulated Nicholson-Bailey model. We compare stability regions in Section IV summarizing the findings in Fig 1, and put these results in the context of known literature and highlight new insights arising from this study.

II. Stability condition for a generalized Lotka-Volterra model
Consider a generalized Lotka-Volterra model with an attack rate f(h, p) that is an arbitrary continuously-differentiable function of both the host and parasitoid population densities. Such an attack rate captures diverse ecological processes: • Function f decreasing with host density corresponds to a Type II functional response. In contrast, f varies non-monotonically with host density for a Type III response with the attack rate increasing (accelerating) with host density at low densities and then decreasing at high
We assume that model (3) has a unique non-trivial equilibrium given as the solution to where h � and p � represent the equilibrium host and parasitoid densities, respectively. Note that the model also has a trivial equilibrium h � = p � = 0, but we are primarily interested in the properties of the non-trivial equilibrium. In Appendix A, we provide explicit conditions for a wide range of attack rates f that result in a unique non-trivial equilibrium.
Our analysis in Appendix B shows that this equilibrium is asymptotically stable, if and only if, are the dimensionless log sensitivities of the attack rate to the host and parasitoid densities, respectively. It is reasonable to have the net attack rate per host f(h, p)p be a non-decreasing function of the parasitoid density that constrains f p � −1. Similarly, we constrain f h � −1 that ensures the net attack rate per parasitoid f(h, p)h be a non-decreasing function of the host density. If any one of the inequalities in (5) does not hold, then the equilibrium is locally unstable and this instability often manifests in a stable limit cycle [42,43].

III. Stability conditions for a generalized Nicholson-Bailey model
We next consider a potential discrete-time counterpart of the generalized Lotka-Volterra model. One phenomenological approach to obtain this model is to simply substitute c in (2) with f(H t , P t ) yielding [44]. Note that during the time window when hosts are vulnerable to parasitoid attacks, the density of unparasitized host continuously decreases. To capture such effects of continuously changing populations, a semi-discrete or hybrid formalism has been proposed to mechanistically formulate the corresponding discrete-time model [45][46][47][48][49][50][51]. Briefly, an ordinary-differential equation model describing population interaction during the attack time window is solved to derive the update functions of the discrete-time model. Importantly, this mechanistic approach can yield contrasting results from the phenomenological approach-a Type III functional response is stabilizing in the semi-discrete approach [48,52], but always destabilizing if one considers (7) with an analogous attack rate [48].
Referring the reader to Appendix C for details, the semi-discrete approach yields the following model in the discrete-time formalism that is analogous to the continuous-time counterpart Recall from (6) that f h and f p are the log sensitivities of the parasitoid attack rate to the host and parasitoid densities, respectively. This model has a unique non-trivial equilibrium which is asymptotically stable, if and only if, Note that the second inequality is the same as in (5), and hence 1 + f h + f p > 0 is a necessary condition for stability in both modeling frameworks. One can consider small and large values of R for which the stability condition reduces to ð11aÞ respectively.

IV. Comparison of stability regimes
To understand the stabilizing/destabilizing effects of diverse ecological processes we plot the stability regions as predicted by inequalities (5) and (10) in Fig 1. For an analogous comparison, the hosts growth rate r in the continuous-time framework is related to R in the discretetime framework by r = log(R). In Fig 1, f h = f p = 0 represents the classical Lotka-Volterra/Nicholson-Bailey models and • Going right along the positive x-axis (f h > 0) corresponds to an attack rate that accelerates with increasing host density as in a Type III functional response.
• Going left along the negative x-axis (f h < 0) corresponds to an attack rate that decreases with increasing host density as in a Type II functional response.
• Going up along the positive y-axis (f p > 0) corresponds to an attack rate that increases with parasitoid density capturing cooperation between parasitoids to handle hosts.
A straightforward observation from Fig 1 is that the region of stability is larger in continuous time as compared with discrete time, but in some limits, they turn out to be identical.

A. Comparison of stability conditions with respect to host reproduction
An important result that emerges from the inequalities (5) and (10) is that for low levels of host reproduction (i.e., R ! 1 and r ! 0), both modeling frameworks have exactly the same stability criterion showing that an attack rate decreasing with increasing parasitoid density is a necessary condition for stability. For high levels of host reproduction (i.e., R ! 1 and r ! 1), f h > 0 becomes the necessary and sufficient condition for stability in the Lotka-Volterra framework. In the Nicholson-Bailey framework, the stability condition for R ! 1 is given by (11b) that corresponds to the two lines meeting at f h = 0 and f p = −1. This implies that f h > 0 is only a necessary condition for stability in the discrete-time framework for large values of R and one further requires f p < f h − 1 for stability. It is important to point that the limit R ! 1 is less ecologically relevant, as for most natural systems R is expected to be less than 10 [1].

B. Impact of functional responses
We first consider the effect of a Type II functional response where c 1 is the attack rate at low-host density and T h is the parasitoid handling time. In this case, the dimensionless log sensitivities of the attack rate to the host density is negative, and f h approaches −1 as handling times becomes longer and longer. In the context of Fig 1, this corresponds to the neutrally stable Lotka-Volterra equilibrium moving left on the negative x-axis, and this destabilizing effect of handling times is consistent with known literature [61]. Interestingly, a Type II functional response (f h < 0) can provide stability if combined with other mechanisms, such as, interference in parasitoid attack where f p < 0 (grey-shaded area in the third quadrant). Note that this grey-shaded area decreases in size with increasing reproduction, suggesting that the stability arising from the mixture of handling times and parasitoid interference is more likely to operate at low host proliferation. The meeting point of two lines in the third quadrant gives the necessary condition for stability that reduces to f h > −1 in the limit R ! 1 & r ! 0. These limits on f h have important implications suggesting that long handling times that cause f h ! −1 can drive instability even in the presence of parasitoid interference in both modeling frameworks.
A Type III functional response (f h > 0 and f p = 0) is stabilizing in both framework but the degree of attack-rate acceleration needed for stability is much higher in the discrete-time framework: f h > 1 in the discrete-time framework compared to f h > 0 in the continuous-time framework.
Recall that one of the inequalities needed for stability 1 + f h + f p > 0 is the same for both frameworks. The other inequality in the continuous-time framework is anchored at the origin and rotates anticlockwise with increasing R. In contrast, the other inequality in the discretetime framework is anchored at f h = 1, f h = 0 and also rotates anticlockwise with increasing R. This leads to the stability region related to a Type II functional response shrinking with increasing R. However, the stability region related to Type III responses expands with increasing R in the continuous-time framework. The non-origin anchoring in the discrete-time framework leads to the stability region related to Type III responses expanding with increasing R for f h > 1, but shrinking for 0 < f h < 1.

C. Impact of a parasitoid-dependent attack rate
For a parasitoid-dependent attack rate (f h = 0) one can see that f p > 0 (i.e., cooperation between parasitoids) is not stabilizing in both frameworks. In contrast, mutual interference between parasitoids is stabilizing with the stability criterion reducing to in continuous time ð16aÞ An interesting observation from the fourth-quadrant of the discrete-time stability region is that while values of and are by themselves not stabilizing, their combination can lead to stability. Thus, moderate levels of parasitoid interference together with attack-rate acceleration to host density can stabilize population dynamics in the Nicholson-Bailey formulation of host-parasitoid population dynamics. Finally, we point out that the stability resulting from a combination of cooperation between parasitoids and a Type III functional response (grey-shaded region in the first quadrant), but increasing parasitoid cooperation also requires a higher degree of attack-rate acceleration to keep the system dynamics stable.

V. Conclusion
In summary, this contribution combined the stabilizing effects of diverse ecological processes into a single stability criterion in both the Lotka-Volterra ad Nicholson-Bailey modeling formalisms. These results highlight general features that are stabilizing, irrespective of model choice. For example, while destabilizing effects of Type II functional response are well known in the literature, combining it with parasitoid interference can lead to stability (the third quadrant in Fig 1), and stability arising from this mixture is more likely to be prevalent for low values of R and r. Eq (16) quantifies the limits to such a strategy and suggests that long handling times cannot be stabilized irrespective of parasitoid interference. Recent fieldwork has documented cooperation between parasitoids, where groups of female Sclerodermus harmandi wasps are shown to be more successful in exploiting largesized host individuals as compared to individual wasps [40]. Our analysis shows that such cooperation between parasitoids can have a destabilizing effect, but this can be countered by the stabilizing effect of a Type III functional response. Indeed, given cooperation between parasitoids, a Type III response is necessary for stability. For high values of R and r, stability is much more likely to arise from a Type III response with either parasitoid interference or cooperation.
It is perhaps intuitive that the overall stability regime is larger in the continuous-time framework, but interestingly, our analysis reveals that the stability region is identical in the limit of low host reproduction. An important restriction of this work is that there is no limit to host growth in the absence of the parasitoid, and these results could be further generalized to consider a host carrying capacity, which can be stabilizing in both frameworks [1,46]. Another interesting direction of future work would be to expand these stability criteria in the context of two parasitoid species attacking the same host, or two different host species attacked by a common parasitoid providing a holistic impact of diverse mechanisms on the stability of complex consumer-resource models.

Appendix A: Uniqueness of non-trivial equilibrium in the generalized Lotka-Volterra model
A functional response in the parasitoid attack rate can be modeled using a function f that only depends on the host density with parameter c 1 > 0, T h is the parasitoid handling time, q = 0 corresponds to a Type II functional response, and q > 0 captures a Type III functional response. The non-trivial equilibrium is given as the solution to (4), which can be rewritten as Substituting (19) in (20), h � is the solution to Assuming the parasitoid handling time is smaller that the average parasitoid lifespan (i.e., T h < 1/γ), (21) will have a unique solution as the left-hand-side of (21) is a monotonically increasing function that starts at zero and saturates at 1/T h > γ, and hence, h � f(h � ) will intersect with γ only once. Having obtained the unique h � , the corresponding parasitoid density is given by p � ¼ rh � g . Let us consider a parasitoid-dependent attack rate where c 1 > 0 and −1 < α < 1, with negative (positive) values of α denoting parasitoid interference (cooperation). In this case the unique non-trivial equilibrium will be Appendix B: Stability criterion for the generalized Lotka-Volterra model Linearizing the right-hand-side of (3) around the equilibrium yields the following Jacobian matrix and stability requires a Hurwitz matrix whose eigenvalues have negative real parts [62,63]. For a two-dimensional system, the equilibrium is asymptotically stable, if and only if, the determinant of the A matrix is positive and its trace is negative [62,63]. This implies that the equilibrium obtained as the solution to (4) is asymptotically stable, if and only if, both inequalities in (5) hold.

Appendix C: Stability criterion for the generalized Nicholson-Bailey model
The semi-discrete approach models the host-parasitoid interaction during the host's vulnerable stage as an ordinary differential equation. Let τ denote the time within the host vulnerable stage that varies from 0 to T corresponding to the start and end of the vulnerable stage. The densities of parasitoids, un-parasitized and parasitized host larvae at time τ within the vulnerable stage of year t are represented by P(τ, t), L(τ, t), I(τ, t), respectively. These densities evolve as per where c represents the parasitoid's attack rate per host, and γ P , γ L , γ I are the death rates of the respective species. Assuming P t parasitoids, RH t host larvae, and no parasitized larvae at the start of the vulnerable period (τ = 0), solving the above differential equations with initial conditions Lð0; tÞ ¼ RH t ; Pð0; tÞ ¼ P t ; Ið0; tÞ ¼ 0 ð26Þ predicts the parasitized and unparasitized larval populations at the end of the season (τ = T). This leads to a more general discrete-time model where update functions are obtained by setting GðH t ; P t Þ ¼ IðT; tÞ: Solving (25) for a constant attack rate c with no mortalities (γ P = γ L = γ I = 0), and assuming T = 1 without loss of any generality, yields the Nicholson-Bailey model (2).
To capture a generalized parasitoid attack with log sensitivities f h (with respect to the host) and f p (with respect to the parasitoid), we replace c in (25) with the monomial Lðt; tÞ f h P f p t . In the absence of any mortalities (γ P = γ L = γ I = 0) and again assuming T = 1, the above semi-discrete approach results in the model (8). Linearizing the right-hand-side of (8) around the equilibrium (9) results in the Jacobian matrix and stability in the discrete-time formalism requires all eigenvalues of A to have an absolute value less than one [63,64]. For a 2 × 2 matrix, the stability criterion can be written in terms of the determinant and the trace of A-the equilibrium (9) is stable, if and only if, 1 þ TrðAÞ À DetðAÞ > 0 1 þ TrðAÞ þ DetðAÞ > 0: It turns out that the last inequality always holds, and the first two inequalities yield the stability conditions (10).