The effects of allelochemical transfer on the dynamics of hosts, parasitoids, and competing hyperparasitoids

Allelochemicals produced by plants may be ingested by herbivorous insects and transferred to higher trophic levels with potentially deleterious effects. We develop a system of differential equations to investigate the effect of the transfer of allelochemicals, such as nicotine, on the population dynamics of a system of hosts, parasitoids, and two competing hyperparasitoids that attack different life stages of the parasitoids. We find both somewhat deleterious effects of nicotine on the larvae‐attacking hyperparasitoids and increased attack rates for the pupae‐attacking hyperparasitoids can promote coexistence. We also use an evolutionary game‐theoretic approach to determine the optimal distribution of hyperparasitoid attacks among nicotine‐producing and nicotine‐free plants. With strong deleterious effects of nicotine and increased attack rates for the pupae‐attacking hyperparasitoid, we find both species attack parasitoids on the nicotine‐free plant but only pupae‐attacking hyperparasitoids attack parasitoids on the nicotine‐producing plant.


| INTRODUCTION
The population dynamics of systems of interacting herbivorous hosts, parasitoids, and hyperparasitoids are affected by a lower trophic level-the host plants on which the herbivorous host feeds. Some plants produce allelochemicals as a defense mechanism against herbivory, and in turn, some herbivorous insects have evolved behavioral and physiological adaptations to overcome and specialize on plants that produce these specific allelochemicals. Plant allelochemicals are then transferred from specialist hosts to parasitoids and hyperparasitoids, with potentially deleterious effects.
Endophagous parasitoids deposit eggs and complete larval development within their herbivorous hosts, resulting in the eventual death of their hosts (Godfray et al., 1994). These primary parasitoids are themselves vulnerable to higher trophic-level parasitoids known as "hyperparasitoids," or "secondary parasitoids" (Sullivan, 1987). Both primary and secondary parasitoids vary in their ovipositional strategies; for example, some species attack the larval stage of the host or primary parasitoid, whereas others attack the prepupal or pupal stages (Godfray et al., 1994).
The tobacco hornworm, Manduca sexta (L.) (Lepidoptera: Sphingidae), is a specialist on the Solanaceae, which includes many alkaloid-producing species (Chowański et al., 2016). In North America, it is a common pest of cultivated tobacco (Nicotiana tobacum L.), which produces nicotine, and tomato (Lycopersicon esculentum Mill.), which does not produce nicotine. Manduca sexta is highly tolerant to alkaloids and particularly insensitive to nicotine, due to the rapid and efficient degradation of nicotine by inducible enzymes and excretion of nicotine metabolites (Self et al., 1964a(Self et al., , 1964bSnyder et al., 1994;Wink & Theile, 2002). Although most of the ingested nicotine is metabolized and excreted, enough remains in the host hemolymph to affect parasitoids and hyperparasitoids (Barbosa et al., 1986;Harvey et al., 2007).
Cotesia congregata (Say) (Hymenoptera: Braconidae) is the only hymenopterous endoparasitoid that attacks larvae of M. sexta. This parasitoid is gregarious, that is, females deposit multiple eggs into a host during a single ovipositional event; the preferred host stage is the third of five instars. Parasitoid larvae feed on nutrients dissolved in the host hemolymph and are thus exposed to low levels of undegraded nicotine that can deleteriously affect their development and survival (Barbosa et al., 1986). However, populations of C. congregata associated historically with M. sexta on tobacco have evolved behavioral and physiological adaptations, resulting in greater tolerance to nicotine (Kester & Barbosa, 1994).
While undergoing development in M. sexta, C. congregata is vulnerable to parasitism by several hymenopterous hyperparasitoids, including species that attack larvae, such as Mesochorus americanus Cresson (Ichenumonidae), and several species, including Hypopteromalus tabacum (Fitch) (Pteromalidae), that attack prepupae within cocoons (McNeil & Rabb, 1973). Allelochemicals ingested by the herbivorous host can be transferred not only to the third trophic level (parasitoids) but also, subsequently, to the fourth trophic level (hyperparasitoids) (Harvey et al., 2007). Since most of the alkaloids present in C. congregata larvae are excreted in the meconium or shed into cocoon silk (Barbosa et al., 1986), it is likely that hyperparasitoid species that attack larval parasitoids are exposed to higher levels of nicotine than hyperparasitoids that attack prepupal or pupal parasitoids. Interestingly, the diversity of hyperparasitoid species has been observed to be lower on nicotine-producing tobacco plants compared to nicotine-free tomato plants; no larvaeattacking species of hyperparasitoid have been observed emerging from C. congregata cocoons on tobacco plants ( Karen Kester, Personal Communication).
Many mathematical models have been created and analyzed for host-parasitoid systems (Hassell, 2000;Ives, 1992;Mills & Getz, 1996;Nicholson & Bailey, 1935), including systems with multiple parasitoids attacking a single host (Comins, 1996;Hassell & May, 1986;Hogarth & May, 1984;Kakehashi et al., 1984;May & Hassell, 1981). Bonsall et al. (2002) studied a stagestructured system with a single host and two competing parasitoids, where juveniles of each parasitoid species are vulnerable to attack by adults of the other parasitoid species. Briggs (1993) explored the result of competition between two species of parasitoids attacking different life stages of the host, egg and larval. Other models have incorporated a third trophic level (Beddington & Hammond, 1977;Comins, 1996;Holt & Hochberg, 1998;May & Hassell, 1981), showing the impact of the hyperparasitoid on the host and parasitoid dynamics. The effects of hyperparasitoids can have important implications for biological control, since parasitoids are frequently used as a means of pest suppression (Tougeron & Tena, 2019).
Motivated by interactions among M. sexta, C. congregata, M. americanus, and H. tabacum, we develop a system of differential equations to investigate the effect of allelochemical transfer on the population dynamics of a system of hosts, parasitoids, and both larval and prepupal or pupal hyperparasitoids. Our continuous-time model incorporates stage-structure for the hosts, parasitoids, and two species of hyperparasitoids within a habitat of a single plant, which may be either a nicotine-producing plant, such as tobacco, or a nicotine-free plant, such as tomato. Nicotine is assumed to affect the survival of the hyperparasitoid M. americanus, which attacks the larval stage of the primary parasitoid, C. congregata, but not the survival of H. tabacum, which attacks the prepupal stage of the primary parasitoid. We analyze the stability of the model's equilibria and determine invasion criteria for the competing hyperparasitoid species. We explore how increased levels of nicotine in the system, as well as differences in oviposition rates between the two hyperparasitoid species, affect the outcome of competition between the hyperparasitoids, and discuss the implications for hyperparasitoid diversity. We focus on mechanisms that promote coexistence of the two hyperparasitoid species.
We next extend our model to a more realistic two-plant system to investigate how plantspecific oviposition rates of competing larvae-attacking and pupae-attacking hyperparasitoid species might evolve when able to attack parasitoids on both nicotine-free plants and nicotineproducing plants. We use an evolutionary game theoretic approach to determine optimal oviposition strategies, as well as the resulting equilibrium spatial distributions, for each hyperparasitoid species. We explore how results change as the effects of nicotine on the nicotine-producing plant increase, and again discuss the implications for hyperparasitoid diversity.

| MODEL DEVELOPMENT
We construct an ordinary differential equation model to describe the dynamics of a system of hosts, parasitoids, and two competing hyperparasitoids on a single plant species. We use this model to explore the population dynamics under varying levels of plant nicotine production.

| Host-only model
We first consider a differential equation model for the host, M. sexta, in the absence of parasitoids and hyperparasitoids. The life-cycle of M. sexta has four stages: egg, larva, pupa, and adult. All stages are modeled except for eggs. Eggs are laid by adults and hatch into larvae at rate b multiplied by the density limiting factor L K (1 − ) ∕ . Larvae mature into pupae at rate g L and have natural mortality rate μ L . All pupae mature into adults at rate g P . Adults have natural mortality rate μ A .
The host-only model equations are as follows: Model (1) has two equilibria: the extinction equilibrium, (0, 0, 0), and a nonextinction equilibrium, Since all model parameters are positive, in order for all components of the nonextinction equilibrium to be positive we must have bg To determine conditions for the stability of each equilibrium we look at the eigenvalues of the Jacobian matrix evaluated at that equilibrium point. The characteristic equation for the Jacobian matrix evaluated at the extinction equilibrium has the form: . Under the Routh-Hurwitz criteria for stability, all coefficients must be the same sign. If any coefficient's sign differs, then there exist eigenvalues of opposite signs guaranteeing at least one eigenvalue with positive real part and an unstable equilibrium. Since all parameters are positive, it is clear that a > 0 3 , a > 0 2 , and a > 0 1 . If the criterion (2) holds, so the nonextinction equilibrium is positive, then a < 0 0 , and by the Routh-Hurwitz criteria, the extinction equilibrium will be unstable. If then all coefficients of the characteristic equation are positive. The additional Routh-Hurwitz criteria needed for stability is a a a > 1 2 0 , and it is easy to show this condition holds since all parameters are positive. Therefore the extinction equilibrium is stable when (4) holds and unstable when condition (2) holds.
The characteristic equation for the Jacobian matrix evaluated at the nonextinction equilibrium is of the same form as (3) with a = 1 3 , a = bg g μ μ μ 2 . If condition (2) is met, then all coefficients are positive. The nonextinction equilibrium point is stable if a a a > 1 2 0 , which reduces to ( ) This condition is satisfied since all parameters are positive. Thus, the nonextinction equilibrium in the host-only model is stable for all parameters sets for which (2) holds.

| Host-parasitoid model
In this section we extend the model to include parasitoids, modeling four life-stages: eggs (P E ), larvae (P L ), prepupae/pupae (P P ), and adults (P E ). We note that while we will refer to P P as the pupal class of the parasitoid, it includes both the prepupal and pupal stages. The density of parasitoid eggs, P E , increases as eggs are successfully oviposited into host larvae. We let c be the clutch size, and α be the rate of successful oviposition of eggs by parasitoid adults (P A ) into host larvae. All parasitoid eggs (P E ) are assumed to mature into parasitoid larvae (P L ) at rate g P E .
Parasitoid larvae mature into the pupal class (P P ) at rate g P L and die at rate μ P L . Parasitoid pupae mature into adults at rate g P P and die at rate μ P P . Parasitoid adults have a natural mortality rate of μ P A . When model (1) is coupled with the parasitoid equations we obtain the host-parasitoid model equations as follows: Model (5) yields three equilibria: Extinction: (0, 0, 0, 0, 0, 0, 0) Host only: L P A ( *, *, *, 0, 0, 0, 0) Host and parasitoid: ( ) L P A P P P P *, *, *, *, *, *, * Note that here and in the following sections, for convenience, we use the same notation (i.e., L P A *, *, *, etc.) for each of the equilibria even though L* for one equilibrium in general is not necessarily equal to L* of another equilibrium.
The extinction equilibrium (6) is stable exactly when the extinction equilibrium of the hostonly model is stable. Parasitoids cannot persist in the absence of the host and must also go extinct.
The second equilibrium is the host-only equilibrium (7), where L P *, *, and A* are identical to the positive equilibrium values for the host-only model and densities of all stages of the parasitoid population are zero. When evaluated at equilibrium (7), the Jacobian matrix of model (5) has a block structure where the upper 3 × 3 block is the same as the Jacobian for model (1) evaluated at the nonextinction equilibrium. All eigenvalues were determined to have negative real part provided condition (2) holds (so L*, P*, and A* are positive). The Routh-Hurwitz criteria can also be applied to the lower 4 × 4 block of the Jacobian matrix. Again, a necessary condition for stability is that all coefficients of the characteristic equation have the same sign. All coefficients for the characteristic equation for the lower 4 × 4 block are strictly positive except for the constant coefficient, a 0 , which is negative if the following condition holds: When this condition, also known as the invasion criteria, is met, the host-only equilibrium (7) is unstable and parasitoids are able to grow from low initial numbers when introduced into a stable host population. Biologically, condition (9) is met when the total number of offspring produced per adult parasitoid that survive to adulthood exceeds one, allowing the parasitoid population to persist.
The third equilibrium is the coexistence equilibrium (8), where both the host and parasitoid coexist. The equilibrium values at this steady state are as follows: When the extinction and host-only equilibria are unstable, the coexistence equilibrium (8) may be stable, or the hosts and parasitoids may also enter into a limit cycle; this has been observed to occur in simulations for certain parameter values, and has also been seen in previous models such as the Nicholson and Bailey model (Godfray & Shimada, 1999;Mills & Getz, 1996;Nicholson & Bailey, 1935).

| Host-parasitoid-hyperparasitoid model
In this section, we incorporate two competing hyperparasitoid species to produce the full hostparasitoid-hyperparasitoid model (schematic in Figure 1). One species of hyperparasitoid attacks the larval stage of the parasitoid, and the other species attacks the pupal stage of the parasitoid (recall that P P includes both parasitoid prepupae and pupae). We model two stages (juvenile and adult) for each hyperparasitoid species. H J represents the density of larvae-attacking hyperparasitoid juveniles and H A represents the density of larvae-attacking hyperparasitoid adults. Similarly, Ĥ J and Ĥ A are the densities of the pupae-attacking hyperparasitoid juveniles and adults, respectively.
Larvae-attacking hyperparasitoid adults are assumed to attack parasitoid larvae at rate β L . Larvae-attacking hyperparasitoid juveniles emerge from eggs oviposited into parasitoid larvae. The hyperparasitoid juveniles will mature into adults of the same species at rate g H J , and incur mortality at rate μ H J . Note that this mortality rate incorporates both natural mortality as well as any increased mortality due to nicotine. Larvae-attacking hyperparasitoid adults have mortality rate μ H A . Pupae-attacking hyperparasitoid adults oviposit their eggs into the pupal stage of the parasitoid at rate β P , producing pupae-attacking hyperparasitoid juveniles. These larvae will mature into hyperparasitoid adults at rate g Ĥ J , with mortality rate μ Ĥ J . Pupae-attacking hyperparasitoid adults incur mortality at rate μ Ĥ A .

L P A P P P P H H H H
*, *, *, *, *, *, *, *, *,ˆ*,ˆ* There are no equilibria where parasitoids and/or hyperparasitoids are present at positive levels in the absence of the host, since the host is necessary for these populations to survive. Similarly, there are no equilibria where hyperparasitoids persist in the absence of parasitoids. We will focus our analysis of this model on the stability of the three equilibria containing hyperparasitoids: two exclusion equilibria where only one hyperparasitoid species persists (14 and 15), and the coexistence equilibrium where both hyperparasitoid populations persist (16).

| Invasion criteria
Here we derive criteria for when one species of hyperparasitoid can invade the system with the other hyperparasitoid species at equilibrium. The species that is established is referred to as the "resident" while the species being introduced is referred to as the "invader." If a resident population is at an F I G U R E 1 Schematic for Model (10). Solid lines indicate maturation to an older stage or reproduction.
Dashed lines indicate parasitism, with α, β L , and β P representing the attack rates of primary parasitoid adults on host larvae, larvae-attacking hyperparasitoid adults on parasitoid larvae, and pupae-attacking hyperparasitoid adults on parasitoid pupae, respectively. Dotted lines indicate the production of eggs/juveniles resulting from parasitism. Larvae-attacking hyperparasitoid juveniles (circled in red) suffer increased mortality when nicotine is present equilibrium that is unstable, then when a small number of invaders is introduced their population can grow. This may result in either coexistence or extinction of the resident species (Murdoch et al., 2013), with coexistence occurring when both species can invade the other.
If larvae-attacking hyperparasitoids are the resident species and the system is at equilibrium (14), the criterion for pupae-attacking hyperparasitoids to be able to invade is: is the equilibrium density of parasitoid pupae. If each individual from the invading species can produce more than one offspring that survives to become an adult over its lifetime, then it is possible for the species to invade. Total offspring for the pupae-attacking hyperparasitoid is determined by the product of the successful oviposition rate, β P , equilibrium density of parasitoid pupae, P* P , and the expected lifespan of the hyperparasitoid, Note the total number of offspring produced by the invading species is dependent on the density of their host species at the stage the hyperparasitoid attacks, in this case the pupal stage of the primary parasitoid, P* P . We remark that this density is also dependent on properties of the resident species that attacks an earlier stage of the primary parasitoid.
Substituting P* P into (17) with all parameters for the hyperparasitoid populations equal results in the condition: Therefore pupae-attacking hyperparasitoids can invade a population of larvae-attacking hyperparasitoids if the maturation rate of the primary parasitoid larvae, g P L , is greater than the combined rate at which parasitoid pupae leave the stage, either through maturation (g P P ) or death (μ P P ). We note that condition (18) requires the pupal stage of the primary parasitoid be longer than the larval stage. When nicotine is present it will increase the mortality rate of the larvae-attacking hyperparasitoid juveniles, μ H J . The only place this parameter appears in the invasion criteria (17) is in the numerator of P* P . By increasing (or decreasing) μ H J we increase (or decrease) the left-hand side of the invasion criteria (17). Therefore intensifying the effects of nicotine on the larvae-attacking hyperparasitoid resident makes it easier for the pupae-attacking hyperparasitoid to invade. Increasing the successful oviposition rate of the pupae-attacking hyperparasitoid relative to the larvae-attacking hyperparasitoid will also help to satisfy the invasion criteria.
Next we consider the case where pupae-attacking hyperparasitoids are the resident species, at equilibrium (15). The condition for larvae-attacking hyperparasitoids to be able to invade a resident population of pupae-attacking hyperparasitoids is: Again, we can interpret the left-hand side of condition (19) biologically as the total number of new adults produced by each larvae-attacking hyperparasitoid adult over its lifetime. Note that the effect of increasing μ H J is reversed in the larval-attacker invasion criteria compared to the pupae-attacker invasion criteria. The mortality rate of the larvae-attacking hyperparasitoid juveniles, μ H J , is now in the denominator of the left-hand side of (19) and when the parameter is increased to simulate stronger effects of nicotine it will decrease the left-hand side, making it harder to reach the threshold of one required for invasion, and therefore lowering the ability of the larvae-attacking hyperparasitoids to invade a system where the pupae-attacking hyperparasitoids are established.

| MECHANISMS FOR COEXISTENCE
From Section2.3.1, we know that increasing the death rate of the larvae-attacking hyperparasitoid juveniles (μ H J ) due to nicotine will increase the ability of the pupae-attacking hyperparasitoids to invade a system where larvae-attacking hyperparasitoids are established. Increasing the attack rate of the pupae-attacking hyperparasitoid (β P ) relative to the larvaeattacking hyperparasitoid (β L ) can also provide an advantage for the pupae-attacking hyperparasitoid and help the species to persist. In this section, we further explore the effect of μ H J and β P on the dynamics of model (10), focusing on the hyperparasitoid populations.
Baseline parameter values for hosts, parasitoids, and hyperparasitoids are given in Table 1. Parameters are either estimated from the literature or reasonable values are assumed. Mortality rates were calculated from estimates on the percent of individuals that successfully reached the next stage or expected lifespan. For host pupae and parasitoid eggs, these stages are not subject to parasitism and we make the simplifying assumption that their mortality rates are zero; that is, all individuals mature into host adults or parasitoid larvae respectively. Larvae-attacking hyperparasitoid parameters are taken to be the same as those for pupae-attackers. We note that under the parameter values in Table 1, the positive coexistence equilibrium of the Host-Parasitoid Model (8) is stable. Also, with these parameters and in the absence of nicotine ), larvae-attacking hyperparasitoids exclude pupae-attacking hyperparasitoids in the full Host-Parasitoid-Hyperparasitoid Model (10) as shown in Figure 2. Hosts, primary parasitoids, and larvae-attacking hyperparasitoids all reach a positive steady state while pupae-attacking hyperparasitoids go extinct.
Increasing μ H J relative to μ Ĥ J causes equilibrium (14) (hosts, parasitoids, and larvae-attacking hyperparasitoids) to lose stability. As μ H J is increased from a baseline of 0.016, equilibrium dynamics change from exclusion of pupae-attacking hyperparasitoids, to coexistence of the hyperparasitoid species, to exclusion of the larvae-attacking hyperparasitoids (see Figure 3a). Increasing μ H J further has no impact on the model dynamics after the larvae-attacking hyperparasitoids are extinct.
Differences in the rate of successful oviposition between hyperparasitoid species can also affect the outcome of competition. If pupae-attacking hyperparasitoids have an increased rate of successful oviposition relative to the larvae-attacking hyperparasitoids, for example if β = 0.00002 L and ZIMMERMAN ET AL.
Natural Resource Modeling | 11 of 24 β = 0.00004 P , then the species can coexist on a nicotine-free plant (see Figure 3b). Increasing μ H J from the nicotine-free baseline will continue to increase the total density of the pupae-attacking hyperparasitoid and decrease the total density of the larvae-attacking hyperparasitoid until the latter goes extinct. Exclusion of the larvae-attacking hyperparasitoid occurs at a lower value of μ H J when β P is increased relative to β L . We next look at how changing the oviposition rate of the pupae-attacking hyperparasitoid relative to the larvae-attacking hyperparasitoid affects dynamics, for varying levels of nicotine. Let k be the ratio of the successful attack rates of the pupae-attacking hyperparasitoids and the larvae-attacking hyperparasitoids, so k =  (14) (exclusion of pupae-attacking hyperparasitoids) losing stability for k just below 2, when the oviposition success rate of the pupae-attackers is almost twice that of the larvae-attackers (see Figure 4a). Table 1 and μ = 0.016 HJ . Initial conditions were 10 adult hosts, 10 adult primary parasitoids, and 10 adults of each hyperparasitoid species. Larvae-attacking hyperparasitoids exclude pupae-attacking hyperparasitoids (a) (b) F I G U R E 3 The total population size (juveniles plus adults) of larvae-attacking and pupae-attacking hyperparasitoids at equilibrium is plotted against the bifurcation parameter μ HJ , the mortality rate of larvaeattacking hyperparasitoid juveniles. Other parameters are as in Table 1, with initial conditions as in Figure 2. There is a large range of k for which the hyperparasitoids coexist, either at equilibrium or via stable limit cycles. Increasing k far enough will eventually lead to the exclusion of the larvaeattackers (k 5.5 ≈ ). On a nicotine-producing plant with μ = 0.08 H J , we see similar dynamics when varying the attack rates of the hyperparasitoid species (see Figure 4b). Larvae-attacking parasitoids exclude pupae-attacking parasitoids for low values of k (note there is coexistence here for k = 1). The hyperparasitoids coexist at equilibrium (16) for intermediate values of k (note we no longer see oscillations), and for k high enough larvae-attacking parasitoids are excluded by pupaeattacking parasitoids. This exclusion threshold occurs at a lower value of k (k 3 ≈ ) than in the nicotine-free system. The equilibrium population size for pupae-attacking hyperparasitoids increases with k (and β P ) until larvae-attacking hyperparasitoids go extinct. Continuing to increase k once pupae-attacking hyperparasitoids have excluded larvae-attacking hyperparasitoids has a detrimental effect on the total population size of the pupae-attacking hyperparasitoid due to decreased numbers of primary parasitoids as a result of the higher attack rate.

F I G U R E 4
The total population size ( juveniles plus adults) of larvae-attacking and pupae-attacking hyperparasitoids is plotted against k, the ratio of successful oviposition rate of pupae-attacking hyperparasitoid adults to the successful oviposition rate of larvae-attacking hyperparasitoid adults. Here β = 0.00002 L , other parameters are as in Table 1, and initial conditions are as in Figure 2. (a) μ μ = = 0.016 For μ = 0.08 H J , increasing k (or, equivalently, β P ) always results in either constant or decreased equilibrium population sizes for the larvae-attacker, as we change from the larvaeattacker excluding the pupae-attacker, to coexistence, and finally to the pupae-attacker excluding the larvae-attacker. On a plant with even higher levels of nicotine, so μ = 0.36 H J , both species go extinct if k (or, equivalently, β P ) is too low. As k is increased from zero, the pupaeattacking hyperparasitoids are the first to be able to persist (k 0.38 ≈ ). Continuing to increase k further results in a small region of coexistence (k 0.42 − 0.73 ≈ ), after which the pupaeattacking hyperparasitoid again excludes the larvae-attacking hyperparasitoid. We note the latter region includes the scenario of equal attack rates (k = 1). Figure 5 shows the outcome of competition between the two hyperparasitoid species as the mortality rate due to nicotine, μ H J , and the attack rate of the pupae-attacking hyperparasitoid, β P , vary. If β P is sufficiently low and μ H J is sufficiently high, neither species of hyperparasitoid can persist. For low values of β P and μ H J , larvae-attacking hyperparasitoids exclude pupae-attacking hyperparasitoids. As either parameter is increased, pupae-attacking hyperparasitoids gain a fitness advantage over their competitor, either by incurring a lower mortality rate or increased attack rate.

| EVOLUTION OF PLANT CHOICE
In this section we use an evolutionary game theoretic approach to investigate how plant-specific oviposition rates of competing larvae-attacking and pupae-attacking hyperparasitoid species may evolve in an environment containing both nicotine-free plants (e.g., tomato) and plants that produce nicotine (e.g., tobacco), referred to here as nicotine-producing plants. Although C. congregata adults can fly, we make the simplifying assumptions that each plant type has an associated host and primary parasitoid population that does not utilize the alternative plant type. Both species of hyperparasitoid have access to both types of plants, with overall oviposition rates β L and β P for the larvae-attacking and pupae-attacking hyperparasitoids, respectively. We define the strategy parameter p to be the fraction of attacks by larvae-attacking hyperparasitoids that occur on parasitized hosts on nicotine-free plants, and similarly the strategy parameter q is the fraction of attacks by pupae-attacking hyperparasitoids that occur on parasitized hosts on F I G U R E 5 Persistence outcomes as a function of β P and μ HJ . μ = 0.016 ĤJ , β = 0.00002 L , and other parameters are as in Table 1 with initial conditions as in Figure 2 nicotine-free plants. It follows that p 1 − and q 1 − are the fraction of attacks made by larvae-and pupae-attacking hyperparasitoids, respectively, on parasitized hosts on nicotine-producing plants.
A model schematic is shown in Figure 6 and   where all variables and parameters are defined as in the single plant model and a "+" denotes populations on the nicotine-producing plant. The subset of the two-plant model (20) representing the dynamics of the larvae-attacking hyperparasitoids is given by ⎤ ⎦ , and the "fitness" matrix, which determines net reproductive output, is ( ) The subset of (20) representing the pupae-attacking hyperparasitoids is The strategies p and q are taken to be hereditary traits that evolve according to the differential equations F I G U R E 6 Schematic for two plant model (20). Solid lines indicate maturation to an older stage or reproduction. Dashed lines indicate parasitism, with α, pβ L , p β (1 − ) L , qβ P , and q β (1 − ) P representing the attack rates of primary parasitoid adults on host larvae on each plant, larvae-attacking hyperparasitoid adults on parasitoid larvae on the nicotine-free plant, larvae-attacking hyperparasitoid adults on parasitoid larvae on the nicotine-producing plant, pupae-attacking hyperparasitoid adults on parasitoid pupae on the nicotine-free plant, and pupae-attacking hyperparasitoid adults on parasitoid pupae on the nicotine-producing plant, respectively. Dotted lines indicate the production of eggs/juveniles resulting from parasitism. Larvae-attacking hyperparasitoid juveniles (circled in red) suffer increased mortality on the plant producing nicotine where G and Ĝ are the fitness functions and denote the spectral radius of G and Ĝ , respectively (Vincent & Brown, 2005). The terms σ p 2 and σ q 2 are the variances of the traits in each population, which determine the speed of evolution. Modeling evolution simultaneously with population dynamics, the full system includes the two-plant model (20) with the equations for trait evolution (23) and (24). Equilibria of this model that optimize fitness will indicate evolutionarily stable strategies. We find that if the death rate of larvae-attacking hyperparasitoids is the same on plants with and without nicotine (μ μ = H H J J + ), then the evolutionary equilibrium for the strategies is p q * = * = 0.5, or equal attack rates on both plants for both species. When the larvae-attacking parasitoids are adversely affected by nicotine this is no longer the outcome of evolution. Figure 7 shows how the strategies p and q evolve over time from an initial condition of p q = = 0.5 (both species using both plants equally) for the case where there is an increased death rate due to nicotine for the larvae-attacking hyperparasitoids on the nicotine-producing plant only (μ = 0.08 , and pupae-attacking hyperparasitoids have a higher attack rate than larvae-attacking hyperparasitoids (β = 0.00002 L , β = 0.00004 P ). Recall these parameters result in a stable coexistence equilibrium for the single plant model (10). The evolutionary equilibrium values of the strategies for the two plant model are now p* = 0.778 and q* = 0.146. That is, larvae-attacking hyperparasitoids attack parasitoid larvae on nicotine-free plants 77.8% of the time, and the other 22.2% of attacks are on parasitoid larvae on nicotine-producing plants. Pupaeattacking hyperparasitoids only attack parasitoid pupae on nicotine-free plants 14.6% of the time, and the other 85.4% of attacks are on parasitoid pupae on the nicotine-producing plant.
To compare the corresponding hyperparasitoid population dynamics for two sets of strategies, model (20)   . While the number of larvaeattacking hyperparasitoid adults remains almost constant, there are more larvae on the nicotine-free plant and fewer on the nicotine-producing plant for the equilibrium values of p and q compared to when p q = = 0.5. The adult pupae-attacking hyperparasitoid population is 23% larger for p p = * and q q = * compared to when p q = = 0.5, with fewer larvae on the nicotine-free plant and larger populations on the nicotine-producing plant. Figure 8 shows how strategy equilibria, along with the corresponding steady state population dynamics, change as the death rate of larvae-attacking hyperparasitoid juveniles on the nicotine-producing plant (μ H J + ) increases. Results are shown for the pupae-attacking hyperparasitoid oviposition rate equal to that of the larvae-attacking hyperparasitoids (β = 0.00002 P ) as well as twice as large (β = 0.00004 P ). If larvae-attacking and pupae-attacking hyperparasitoids have equal attack rates (β β = = 0.00002 L P ), then larvae-attackers will exclude pupae-attackers on both plants when either both plants are nicotine-free or the effects of nicotine on the nicotine-producing plant are small (μ μ < 0.07 . As μ H J + increases, fewer larvae-attacking hyperparasitoid juveniles can survive on the nicotine-producing plant and p*, the fraction of attacks on the nicotine-positive plant, increases slightly. Once μ H J + is large enough, both species of hyperparasitoids can coexist on the nicotineproducing plant. Since pupae-attacking hyperparasitoids can never persist on the nicotine-free  Table 1 plant, all of their attacks are on the nicotine-producing plant (q* = 0). As μ H J + increases, the equilibrium number of larvae-attacking hyperparasitoid juveniles on the nicotine-free plant is unchanged, equilibrium population sizes of the pupae-attacking hyperparasitoid on the nicotine-producing plant increase and equilibrium population sizes of the larvae-attacking hyperparasitoid juveniles on the nicotine-producing plant decrease. This results in a decrease in the adult larvae-attacking hyperparasitoid population size and an increase in the percentage of attacks on the nicotine-free plant, p*, as μ H J + increases. Once μ H J + reaches a value of approximately 0.29 the larvae-attacking hyperparasitoids can no longer persist on the nicotine-producing plant, and 100% of their attacks are on the nicotinefree plant (p = 1), resulting in complete spatial segregation of the species.
When β = 0.00002 L and β = 0.00004 P , both hyperparasitoid species coexist in the absence of nicotine. Both species have equal attack rates (p q * = * = 0.5) and population sizes on both Again, equilibrium population sizes on the nicotine-free plant are not affected by μ H J + . As μ H J + increases, p* increases and q* decreases as the juvenile population sizes supported on the nicotine-producing plant decrease for the larvae-attacking hyperparasitoid and increase for the pupae-attacking hyperparasitoid.
Once μ H J + exceeds a value of approximately 0.14, the effects of nicotine are too strong for the larvae-attacking hyperparasitoids to persist on the nicotine-producing plant, and 100% of their attacks are on the nicotine-free plant (p = 1). Here both species still coexist on the nicotine-free plant, but the nicotine-producing plant contains only pupae-attacking hyperparasitoids. Since pupae-attacking hyperparasitoids are not affected by nicotine, p* and q* remain constant for any further increase in μ H J + .

| DISCUSSION
We have modeled a system of hosts, parasitoids, and two competing hyperparasitoids, incorporating trade-offs between the competing hyperparasitoid species that affect the population dynamics of the system and the outcome of competition between the hyperparasitoids. While the larvae-attacking hyperparasitoids have the advantage of attacking an earlier stage of the primary parasitoid, they incur a cost of increased mortality at the juvenile stage in the presence of nicotine. On the other hand, the pupae-attacking hyperparasitoids must wait for primary parasitoids to escape parasitism by larvae-attacking hyperparasitoids and make it past the larval stage to the pupal stage before they can attack, but they do not incur any cost due to nicotine. We find that larvae-attacking hyperparasitoids will exclude pupae-attacking hyperparasitoids in the absence of nicotine, when there are no differences in the species other than the stage of the parasitoid in which they deposit their eggs. Similarly, Briggs found that two parasitoid species attacking different life-stages of their host, eggs, and larvae, could not coexist when larvaeattacking parasitoids could not attack already-parasitized hosts (Briggs, 1993). We find that the two hyperparasitoid species can coexist when the cost to the larvae-attackers, due to nicotine, becomes high enough, or when pupae-attackers have an additional competitive advantage such as an increased oviposition rate. With equal attack rates among hyperparasitoid species, coexistence can be attained if the cost of increased mortality due to nicotine for the larvae-attacking hyperparasitoid outweighs the benefit of attacking the earlier stage of the parasitoid. As μ H J increases, eventually a level is reached where the larvae-attackers can no longer exclude the pupae-attacking hyperparasitoids and coexistence is possible. This threshold value of μ H J is lower when pupae-attacking hyperparasitoids also have an additional advantage, such as increased attack rates on the parasitoid. If μ H J continues to increase, the cost of mortality becomes too great to overcome and eventually the larvae-attacking hyperparasitoid will not be able to persist. The mortality rate of larvae-attacking hyperparasitoid juveniles may depend on the concentrations of nicotine being produced by the plant and subsequently transferred by the host to the primary parasitoid. The impact of a given concentration of nicotine may also vary with species and how well-adapted they are to nicotine. Species with no prior exposure may incur greater costs than species that have had prior exposure and evolved adaptations.
Coexistence of larvae-and pupae-attacking hyperparasitoid species has been observed to occur on nicotine-free plants. For example, both larvae-and pupae-attacking hyperparasitoids have been reared from cocoons from single broods of C. congregata collected from hosts that fed on nicotine-free plants in the field; however, only prepupae-attacking hyperparasitoids have been reared from C. congregata collected from M. sexta on tobacco (Karen Kester, personal communication). The advantage larvae-attacking hyperparasitoids gain from attacking the earlier stage of the parasitoid would prevent coexistence on a nicotine-free plant, given equal attack rates. However, coexistence can occur when the pupae-attacking species has the advantage of a higher successful rate of oviposition (β β > P L ). The attack rates β L and β P are based on the success of the hyperparasitoids both in search and handling time. Ovipositing into the larva of the parasitoid through the cuticle of the host is a relatively rare trait among hyperparasitoids. Possibly, handling times are greater for larvae-attacking hyperparasitoids due to defensive behavior of the parasitized host. The hyperparasitoids can also coexist in the presence of nicotine; however, if the effects of nicotine become too strong, the larvae-attacking hyperparasitoid can no longer persist. As the successful attack rate of the pupae-attacking hyperparasitoids is increased relative to that of the larvae-attacking hyperparasitoids on a nicotine-producing plant, we see that exclusion of the larvae-attackers generally occurs for weaker effects of nicotine (lower values of μ H J ). The scenario where pupae-attacking hyperparasitoids have increased oviposition rates and larvae-attacking hyperparasitoid juveniles suffer increased mortality due to nicotine results in model dynamics consistent with observations: coexistence of both hyperparasitoid species on nicotine-free plants, and the absence of larvae-attacking hyperparasitoids on plants producing high enough levels of nicotine.
We note that our model has been parameterized so larvae-attacking hyperparasitoids can parasitize all internal instars of the larval stage of the primary parasitoid and pupae-attacking hyperparasitoids can parasitize both the prepupal and pupal stages. In reality, the stages that specific hyperparasitoids attack may be narrower. Future work includes extending the model to explicitly account for individual larval instars and the prepupal stage alone.
To further explore the impact of nicotine on hyperparasitoid diversity, we used evolutionary game theory to determine optimal oviposition strategies for both species when faced with the decision of whether to attack parasitoids on a plant with or without nicotine. When overall attack rates are equal between hyperparasitoid species, the optimal strategy for the pupaeattacking hyperparasitoid is to never attack any parasitoids on the nicotine-free plant. However, this limited diversity outcome is a result of the two species being unable to coexist with equal attack rates in the absence of nicotine, with pupae-attacking hyperparasitoids being excluded.
If pupae-attacking hyperparasitoids have an attack rate sufficiently greater than that of larvae-attacking hyperparasitoids, both species can coexist in the absence of nicotine. In the case that nicotine does not affect larvae-attacking hyperparasitoids, both plants are effectively the same and both species make 50% of their ovipositions on each plant. As mortality due to nicotine on the nicotine-producing plant, μ H J + , increases, larvae-attacking hyperparasitoids make a greater percentage of ovipositions on the nicotine-free plant, and pupae-attacking hyperparasitoids increase the percentage of attacks on nicotine-producing plants. Eventually μ H J + is increased high enough that larvae-attacking hyperparasitoids cannot persist on the nicotine-producing plant, and 100% of their ovipositions occur on the nicotine-free plant. At this point, the pupae-attacking hyperparasitoid is attacking parasitoids on both plants, so we observe both species on the nicotine-free plant but only pupae-attacking hyperparasitoids on the nicotine-producing plant.
Our model assumes that each plant has separate populations of hosts and primary parasitoids with equal parameter values. This assumption is reasonable for populations that are nicotineadapted but may break down for populations that are not nicotine-adapted. For example, C. congregata associated with M. sexta on tobacco is unaffected by low concentrations of nicotine in the host diet that otherwise cause mortality in C. congregata associated with M. sexta on tomato (Kester & Barbosa, 1991). When offered tomato and tobacco plants simultaneously in the field, the tobacco-associated population parasitizes hosts on tobacco and tomato equally, whereas the tomato-associated population prefers hosts on tomato (Kester & Barbosa, 1994). Relaxing the assumption that all populations of the primary parasitoid are adapted to nicotine would introduce differences between the nicotine-producing and nicotine-free primary parasitoid populations and therefore alter the optimal distribution of hyperparasitoid attack rates between the two plants.

| CONCLUSION
The dynamics of populations interacting over multiple trophic levels are complex and intertwined. Here we have shown that the result of competition between two hyperparasitoid species at the fourth trophic level can depend upon the presence or absence of a plant allelochemical, such as nicotine, at the primary trophic level. The hyperparasitoid attacking the earlier life-stage of the primary parasitoid is favored in the absence of nicotine, but coexistence of the hyperparasitoids is possible if the species attacking the later life-stage has an additional advantage such as a higher oviposition rate. The presence of nicotine at the primary trophic level also promotes coexistence of the hyperparasitoid species by increasing the mortality rate of the larvae-attacking hyperparasitoids. However, if the effects of nicotine are too strong, larvae-attacking hyperparasitoids will no longer be able to persist.
In a system with two hyperparasitoid species attacking parasitoids on two different plants, each species will distribute attacks between both plants equally if they are both nicotine-free, with coexistence on both plants possible when pupae-attacking hyperparasitoids have a higher attack rate. When one of the plants produces nicotine, however, this distribution is altered. As the levels of nicotine on one of the plants increases, larvae-attacking hyperparasitoids increase the fraction of their ovipositions on the nicotine-free plant, while pupae-attacking hyperparasitoids favor the nicotine-producing plant. When deleterious effects of nicotine on the larvae-attacking hyperparasitoid are high enough, this species attacks the nicotine-free plant exclusively, while the pupae-attacking hyperparasitoid continues to attack parasitoids on both plants. In this case the diversity of hyperparasitoid species is greater on nicotine-free plants compared to plants producing nicotine, consistent with observations.

ACKNOWLEDGMENTS
A portion of this study was completed as part of Mark Zimmerman's masters thesis at Virginia Commonwealth University. This study was supported by a Simons collaboration grant (SR, 426126).