Nicholson-Bailey Model

The Nicholson-Bailey model is a classic discrete-time, coupled difference-equation model in population ecology that describes the intergenerational dynamics of an insect host (H) and a specialist parasitoid (P) with synchronized, non-overlapping generations (Nicholson & Bailey, 1935). Formulated in 1935 by Australian entomologist Alexander John Nicholson and physicist Victor Albert Bailey, the model serves as the discrete-generation counterpart to the continuous-time Lotka–Volterra predator–prey equations, providing the mathematical foundation for insect population dynamics, parasitoid foraging theory, and classical biological pest control across terrestrial ecosystems.

Unlike true predators—which kill and consume multiple prey items across their lifetime—an insect parasitoid (principally parasitic wasps of the order Hymenoptera, such as Ichneumonidae, Braconidae, and Chalcidoidea, as well as tachinid flies of the order Diptera) has a free-living adult stage in which females search the environment and oviposit one (solitary) or several (gregarious) eggs into, on, or near the immature stages (eggs, larvae, or pupae) of an arthropod host. The developing parasitoid larva feeds on the tissues of the single host individual, eventually killing the host and emerging as a winged adult in the next generation (t + 1). Because host reproduction and parasitoid emergence in temperate climates are synchronized into discrete seasonal generations, Nicholson and Bailey modeled their interaction using coupled first-order difference equations.

Let H_t and P_t denote the population densities of host larvae and searching adult female parasitoids, respectively, in generation t. Nicholson and Bailey (1935) derived their coupled recurrences from four mechanistic assumptions:

1. Unbounded Host Geometric Growth in the Absence of Parasitoids: If no parasitoids are present (P_t = 0), each surviving host produces λ (λ > 1) net offspring in the next generation (H_t+1 = λH_t), without intraspecific resource competition (compare the density-dependent Beverton–Holt and Ricker models).
2. Linear (Unsaturated) Random Search ("Area of Discovery" a): Each adult female parasitoid searches randomly and independently with a constant lifetime searching efficiency a (a > 0, termed the area of discovery by Nicholson), without handling-time saturation (a linear Type I functional response where the total number of random encounters per host across the parasitoid population is E_enc = aP_t).
3. Poisson Distribution of Encounters Across Hosts: Because encounters occur randomly and independently at mean rate μ = aP_t per host, the probability p(k) that an individual host is encountered and parasitized k times (k = 0, 1, 2, ...) follows the Poisson distribution:
   • Poisson Encounter Probability: p(k) = (aP_t)^k exp(−aP_t) / k!
   • Fraction of Hosts Escaping Parasitism (k = 0): f(P_t) = p(0) = exp(−aP_t)
   • Fraction of Hosts Parasitized at Least Once (k ≥ 1): 1 − p(0) = 1 − exp(−aP_t)
4. Parasitoid Offspring Recruitment (Numerical Response): Every host parasitized one or more times (k ≥ 1) dies before reproducing and yields an average of c viable adult female parasitoids in generation t + 1 (c = 1 for a solitary parasitoid with a 1:1 sex ratio, or c > 1 for gregarious species).

Combining these four assumptions yields the canonical Nicholson-Bailey coupled difference equations:

• Host Generation t + 1: H_t+1 = λ H_t exp(−aP_t)
• Parasitoid Generation t + 1: P_t+1 = c H_t [ 1 − exp(−aP_t) ]

Eigenvalue Modulus and Diverging Oscillations

To solve for the non-trivial interior equilibrium (H^∗, P^∗) where both host and parasitoid densities remain constant across generations (H_t+1 = H_t = H^∗ and P_t+1 = P_t = P^∗), set the first equation to H^∗ = λ H^∗ exp(−aP^∗):

• Equilibrium Parasitoid Density (P^∗): P^∗ = ln(λ) / a
• Equilibrium Host Density (H^∗): H^∗ = [ λ ln(λ) ] / [ ac (λ − 1) ]

Because λ exceeds 1, a unique positive coexistence equilibrium always exists. However, to determine whether small perturbations around (H^∗, P^∗) decay back to equilibrium or amplify over time, one evaluates the Jacobian community matrix (J) at (H^∗, P^∗) using the Jury stability criteria for discrete-time systems (which require that both eigenvalues μ₁, μ₂ lie inside the unit circle in the complex plane, |μ| ≤ 1).

For the Nicholson-Bailey Jacobian at (H^∗, P^∗), the determinant is:

• Determinant of Jacobian at Equilibrium: Det(J) = λ · ln(λ) / (λ − 1)

By elementary calculus, for every host reproductive rate λ above 1, the inequality λ ln(λ) > λ − 1 holds universally, meaning Det(J) exceeds 1 for all λ above 1. Because the product of the two complex conjugate eigenvalues equals Det(J) = |μ|² > 1, the eigenvalues always have modulus greater than unity.

Consequently, the interior equilibrium (H^∗, P^∗) of the unmodified Nicholson-Bailey model is unstable: any tiny deviation from equilibrium initiates diverging oscillations of progressively increasing amplitude, culminating in the deterministic extinction of the host (or of the parasitoid during the deep host trough, followed by unchecked geometric outbreak of the surviving hosts). In a classic laboratory test inside multi-cage arenas, Syunro Utida (1957) observed rapid oscillatory crashes in single-patch cultures of the azuki bean weevil (Callosobruchus chinensis) and its braconid parasitoid (Heterospilus prosopidis), or the greenhouse whitefly (Trialeurodes vaporariorum) and its aphelinid parasitoid (Encarsia formosa), unless spatial heterogeneity intervened.

Because host and parasitoid species coexist stably across the biosphere, theoretical ecologists (Michael P. Hassell, Robert M. May, Ilkka Hanski, and William W. Murdoch) identified five biological mechanisms that modify the escape function f(P_t) and shift Det(J) below 1:

1. Robert May's (1978) Negative Binomial Model and the CV² Rule

In a celebrated 1978 paper, Robert M. May replaced the random Poisson assumption with a negative binomial distribution of parasitoid attacks across habitat patches, governed by clumping parameter k (k > 0, where smaller k denotes stronger spatial aggregation):

• Negative Binomial Escape Function: f(P_t) = [ 1 + (aP_t) / k ]^−k

As k → ∞, the negative binomial converges to the unstable Poisson model exp(−aP_t). However, May proved analytically that whenever the clumping parameter satisfies 0 < k < 1—meaning parasitoid attacks are strongly aggregated across patches so that many hosts experience low risk while high-density patches absorb heavy parasitism—the interior equilibrium is globally asymptotically stable for any host reproductive rate λ above 1. In 1990, Michael P. Hassell, Stephen W. Pacala, and Robert M. May generalized this into the CV² > 1 rule: regardless of whether spatial heterogeneity is host-density-dependent (aggregative numerical response) or host-density-independent, a host–parasitoid interaction is stabilized whenever the squared coefficient of variation (CV² = Variance / Mean² = 1/k) of searching parasitoid density across hosts exceeds 1.

2. Spatial Lattice Metapopulations and Spiral Waves (Hassell, Comins, & May, 1991)

Even when every individual habitat patch follows the unstable Nicholson-Bailey equations, Michael P. Hassell, Hugh N. Comins, and Robert M. May (1991) demonstrated on a two-dimensional spatial lattice of coupled patches with limited neighborhood dispersal (d_H, d_P) that the system self-organizes into rotating spiral waves and spatial chaos. Because adjacent patches oscillate out of phase within the metacommunity, regional persistence is maintained without any local patch sitting at a static equilibrium.

1. Nicholson, A. J., & Bailey, V. A. (1935). "The balance of animal populations.—Part I." Proceedings of the Zoological Society of London, 105(3), 551–598. https://doi.org/10.1111/j.1096-3642.1935.tb01680.x
2. Beddington, J. R., Free, C. A., & Lawton, J. H. (1975). "Dynamic complexity in predator-prey models framed in difference equations." Nature, 255(5503), 58–60. https://doi.org/10.1038/255058a0
3. May, R. M. (1978). "Host-parasitoid systems in patchy environments: a phenomenological model." Journal of Animal Ecology, 47(3), 833–844. https://doi.org/10.2307/3674
4. Hassell, M. P., Comins, H. N., & May, R. M. (1991). "Spatial structure and chaos in insect population dynamics." Nature, 353(6341), 255–258. https://doi.org/10.1038/353255a0