A model of infection in honeybee colonies with social immunity

Abstract

Honeybees (Apis mellifera) play a significant role in the pollination of various food crops and plants. In the past decades, honeybee management has been challenged with increased pathogen and environmental pressure associating with increased beekeeping costs, having a marked economic impact on the beekeeping industry. Pathogens have been identified as a contributing cause of colony losses. Evidence suggested a possible route of pathogen transmission among bees via oral-oral contacts through trophallaxis. Here we propose a model that describes the transmission of an infection within a colony when bee members engage in the trophallactic activity to distribute nectar. In addition, we examine two important features of social immunity, defined as collective disease defenses organized by honeybee society. First, our model considers the social segregation of worker bees. The segregation limits foragers, which are highly exposed to pathogens during foraging outside the nest, from interacting with bees residing in the inner parts of the nest. Second, our model includes a hygienic response, by which healthy nurse bees exterminate infected bees to mitigate horizontal transmission of the infection to other bee members. We propose that the social segregation forms the first line of defense in reducing the uptake of pathogens into the colony. If the first line of defense fails, the hygienic behavior provides a second mechanism in preventing disease spread. Our study identifies the rate of egg-laying as a critical factor in maintaining the colony's health against an infection. We propose that winter conditions which cease or reduce the egg-laying activity combined with an infection in early spring can compromise the social immunity defenses and potentially cause colony losses.

Fig 1. Model diagram.

Upper: Core model. Lower: Extended model. The diagram illustrates class (horizontal arrows) and state (vertical arrows) transitions of honeybees. The queen produces new brood (B). As bees are older, they change their classes from brood (B), nurses (N), nectar-receivers (R₀ and R₁), to foragers (F₀ and F₁). Bees can change into an infection state (iB, iN, iR₀, iR₁) with a certain probability if they receive nectar from infected ones. Nectar-receivers and foragers can also change their nectar-loaded state between unloaded (subscript 0) and loaded (subscript 1). φ represents ‘death’. See the main text for a more detailed description of each transition. The full model is composed of both core and extended parts.

Fig 2. Role of the hygienic behavior (k_rem).

A: Healthy brood and nurses in a colony with the hygienic behavior are maintained against an infection. The simulation is based on Eqs 1–4 with model parameters listed in Table 1 with l₀ = 120 and k_rem = 2.5 × 10⁻³. We treat p_t1 ⋅ k_RN ⋅ iR₁ as a parameter, which is equal to 5 × 10⁻⁴. The infection is introduced (p_t1 ⋅ k_RN ⋅ iR₁ = 5 × 10⁻⁴) at day 20 (vertical purple line). B: All brood and nurses in a colony without the hygienic behavior become infected. The simulation setting is the same as in panel A, except that k_rem = 0. C: The steady-state number of healthy nurses is plotted as k_rem is varied while other parameters are fixed. Bold lines represent a stable steady-state and a thin line represents an unstable steady-state. D: The 2-parameter bifurcation diagram is plotted as both k_rem and p_t0 are varied. The diagram is divided into three regions depending on the state of the colony (healthy, bistable, and infected).

Fig 3. Role of the egg-laying rate (l₀).

A: All brood and nurses in a small-size colony (e.g., colony with a low egg-laying rate) become infected even with the hygienic behavior response. The simulation setting is the same as in Fig 2A, except that l₀ = 70. The infection is introduced (p_t1 ⋅ k_RN ⋅ iR₁ = 5 × 10⁻⁴) at day 20 (vertical purple line). B: The two-parameter bifurcation diagram is plotted between l₀ and p_t0 while other parameters are fixed. The red letters a and b mark the critical values of l₀, below which the colonies are always infected when p_t0 is 0.3 and 0.5, respectively. C: Simulation of three colonies with the same parameter set from Fig 2A, but the simulation begins with 2400 initial healthy brood and different initial numbers of healthy nurses. The infection is introduced (p_t1 ⋅ k_RN ⋅ iR₁ = 5 × 10⁻⁴) at day 0. D: Simulation setting is the same as that of panel C, except that the initial number of healthy brood is 600.

Fig 4. Role of the probability of nurses being infected from infected brood (p_t,rem) and the death rate of infected bees (k_d).

The two-parameter bifurcation diagrams are plotted between A: k_rem and p_t,rem, B: l₀ and p_t,rem, C: k_rem and k_d, and D: l₀ and k_d. The diagrams are divided into regions corresponding to the state of the colony (healthy, bistable, and infected).

Fig 5. Simulations of the full model.

The full model is simulated with parameters listed in Table 1 with l₀ = 120 (panels A and B) and l₀ = 70 (panels C and D). One infected forager returning from nectar-collecting (iF₁) is introduced at day 20 (iF₁(20) = 1) (vertical purple line).

Fig 6. Simulations of the full model under seasonal effects.

The full model is simulated with parameters listed in Table 1 under seasonal effects. During winter (grey stripes), which lasts five months a year, we set l₀ = 0, k_r = 0, and n_N, n_R, and n_F are increased 4-folds. A and B: l₀ during non-winter months = 300. C and D: l₀ during non-winter months = 1000. One infected forager returning from nectar-collecting (iF₁) is introduced at day 360 (vertical purple line).

Fig 7. Stochastic simulations of the full model under seasonal effects.

Stochastic simulations of two independent colonies under an identical condition are shown in upper and lower panels, respectively. Seasonal effects are implemented, as described in Fig 6 and the main text. Both colonies are simulated with the same parameter set (l₀ during non-winter months = 1000 and other parameters from Table 1). One infected forager returning from nectar-collecting (iF₁) is introduced at day 360 (vertical purple line).

Fig 8. Percentages of survived colonies from stochastic simulations.

Percentages of survived colonies are calculated from 1000 repeats of stochastic simulations. Model parameters are from Table 1, except those that are listed in each figure panel. Seasonal effects are implemented as described in the main text. A: k_r during non-winter months (black line) or p_surv (red line) is varied. B and C: l₀ during non-winter months is varied. D: The timing of an infection (t₀) is varied in relation to the onset of spring (i.e., iF₁(t₀) = 1).