1 Introduction

The complexity of altruistic societies and their inherent inequalities among individuals have puzzled evolutionary biologists ever since Darwin (Ratnieks and Wenseleers 2008). In particular, disparate theories have been proposed to explain how reproductive cooperation arises and is maintained in some animal societies (Nowak et al. 2010). Perhaps the most remarkable example of reproductive cooperation is the division of labor in eusocial insects. In eusocial species of Hymenoptera (bees, ants, and wasps), not all individuals get to reproduce. For instance, in honeybees, one queen (or a small group of queens) is responsible for laying eggs, while up to millions of workers collect provisions for the colony and cooperatively raise and protect the royal offspring (Bourke 1988).

The strength of reproductive cooperation in eusocial societies varies largely across taxa. For instance, in some ant species, workers have vestigial ovaries and therefore are fully sterile (Holman et al. 2010). In many other hymenopterans, however, workers retain functional ovaries and can lay unfertilized eggs that develop into males, but usually do not do so due to regulatory mechanisms in queenright colonies (Bourke 1988). Queens produce pheromones to signal dominance, which (i) inhibit the development of ovaries in workers and (ii) stimulate workers to attack reproductive workers and destroy (by eating) their laid eggs, a behavior known as “worker policing.” Although the physiological mechanisms underlying the inhibition of workers’ reproduction are well-studied, the contribution of worker policing to hive population dynamics and persistency has received much less attention. Curiously, even for honeybees, Apis mellifera, the best-studied eusocial organism, this is not an exception.

Interestingly, despite these control mechanisms, some workers still lay eggs, a phenomenon called anarchy in the beehive (Oldroyd et al. 1994). The underlying causes of anarchistic behavior in workers, and particularly its consequences for the stability of the beehive, are not fully understood so far (Bourke 1988). If anarchic workers lay eggs in the presence of a functional queen (i.e., queenright colonies), the abundance of sexually active males, known as drones, in the colony tends to increase over time because unfertilized eggs necessarily develop into drones. Normally, about one in a thousand drones is worker-derived (Visscher 1989). However, a gradual increment in the proportion of drones derived from eggs laid by anarchic workers may lead to drone overpopulation, which can directly affect colony’s fate. The loss of the natural balance in the drone population is energetically costly for the colony because they do not gather neither nectar nor pollen but still consume resources. Therefore, an overpopulation of drones may result in the collapse of the hive sustainability (Winston 1991).

Anarchistic behavior in workers has an additional cost to the colony associated with the fact that anarchic bees are less productive (Dampney et al. 2004). Therefore, both by reducing the drone population and favoring “good,” diligent workers, worker policing seems to be a beneficial regulatory mechanism for the long-term colony persistence (Ratnieks and Visscher 1989). However, this patrolling behavior adds to a large list of worker’s duties and, if done too often, could also lead to a reduction in hive productivity because workers are collecting fewer resources. Nonetheless, it is expected that the benefits of worker policing outweigh its costs. Because this patrolling behavior is triggered by queen’s pheromones, this chemical signaling has to reach the whole colony to inhibit anarchistic behavior (Visscher 1989). For example, confining the queen to a specific area of the hive favors workers to lay eggs (Dampney et al. 2004; Visscher 1998). Thus, we hypothesize that the rate of anarchic workers increases in bigger colonies as queen’s pheromones may not be as efficient to control the whole hive as it is in smaller colonies.

Here, we developed a dynamical mathematical model to understand the dynamics of anarchy in workers and its consequences to the hive persistence. The model considers two different costs of the anarchistic behavior to the hive productivity: (i) the reduction in energetic gain associated with increments in the population of drones and (ii) the reduction in normal workers’ productivity due to the time allocation to policing. Mathematical models have previously been used to gain insights into the consequences of behavioral changes within the beehive (Khoury et al. 2011, 2013), but to the best of our knowledge, none has addressed the effects of anarchistic behavior to the beehive stability and long-term persistence. Altogether, our results suggest that worker policing is overall beneficial for the beehive; however, when anarchy rate exceeds a certain threshold, policing is not able anymore to compensate the costs. There are certain circumstances where worker policing could even be detrimental to the hive. Our model suggests that the balance between worker policing and frequency of anarchic workers’ behavior can provide new insights into the maintenance of eusociality in honeybees.

2 Materials and methods

2.1 Colony dynamics

We propose a mathematical model for the dynamics of the worker bee population including anarchic ones. We define two groups in this population: worker bees that are not affected by anarchism W (“normal workers”) and anarchic workers A (“anarchist workers”). We assume that colony increase (i.e., new offspring) is proportional to the resources acquired by the hive. This assumption is common and standard in population ecology (Cain et al. 2011; Loreau 2010; Rockwood 2006). The per capita resource rate captured by a worker bee is denoted by β. When there are worker-laid eggs in the colony, worker policing takes place (Ratnieks and Visscher 1989; Visscher 1998), reducing productivity because workers allocate time to policing duties. We model this reduction by multiplying β by (1 − ρ), where ρ represents the fraction of time that the bees spend performing worker policing. It has been observed that anarchic bees do not contribute to worker policing, but their baseline resource production is less than in an unaffected bee (Dampney et al. 2004). We denote the per capita resource production rate of an anarchist bee by δ. Adding the contribution from the unaffected bees to the affected bees, we find the production rate or productivity of the colony to be P = β(1 − ρ) W + δA. We then introduce the population death rate of workers as – μNW – μγAW / (1 + ρθW). In the term μNW, N is the total population (N = W + A), and μ represents the intensity of the population autoregulation processes (Loreau 2010). Thus, μNW is the death rate due to population autoregulation. The term μγAW / (1 + ρθW) is the increase in the death rate due to anarchistic behavior, where γ is the per capita cost of anarchy and θ is the per capita worker policing efficiency. This term can be understood as follows: anarchistic behavior breaks the normal sex ratio by generating additional costs for the colony, i.e., increment in the population of drones (Visscher 1998). We model this extra cost to be proportional to the number of anarchist bees γA. On the other hand, this extra cost is reduced by the worker policing done by the normal worker population W (Dampney et al. 2004). The reduction due to worker policing is modeled by dividing the cost γA by the factor (1 + ρθW). By looking at the previous expression, we notice that if there is no worker policing ρ = 0, then there is no reduction in the cost produced by the anarchist. Finally, we assume that this resulting cost increases the death rate proportionally to itself. The previous assumption is also common for ecological models (Cain et al. 2011; Loreau 2010). Finally, it has been hypothesized that failure or weakening in the queen’s control via pheromones (e.g., too big colonies) led to anarchistic behavior (Bourke 1988), because controlling large colonies is harder than small ones; thus, we modeled the transition rate from worker bee to an anarchist bee to be proportional to the colony size αNW, where α is the spontaneous recruitment rate. In other words, αNW is the number of workers that become anarchist per unit of time, thus leaving the population W and entering the population A. Equation (2) describes the population dynamics of anarchic bees A, and the term μNA is the death rate due to population autoregulation. In a similar way as in equation (1), the term μγAA / (1 + ρθW) is the increase in the death rate due to the costs of anarchistic behavior. In summary, our model consists of the following equations:

$$\frac{dW}{dt}=\beta \left(1-\rho \right)W+\delta A\hbox{-} \mu NW\hbox{-} \frac{\mu \upgamma AW}{1+\rho \theta W}-\alpha NW$$
(1)
$$\frac{dA}{dt}=\alpha NW-\mu NA-\frac{\mu \upgamma {A}^2}{1+\rho \theta W}$$
(2)
$$N=W+A$$
(3)

We assume that new offspring emerge as normal workers (W) and that population autoregulation processes and costs of anarchism affect equally anarchic and non-anarchic bees.

2.2 Analysis of the system

To find the number of stationary points and their stability, we plot the nullclines of the system. The W-nullcline is the set of points in the space (W, A) where the equation:

$$0=\beta \left(1-\rho \right)W+\delta A-\mu NW-\frac{\mu \upgamma AW}{1+\rho \theta W}-\alpha NW$$
(4)

is satisfied. Correspondingly, the A-nullcline is the set of points where the equation:

$$0=\alpha NW-\mu NA-\frac{\mu \upgamma {A}^2}{1+\rho \theta W}$$
(5)

is satisfied. We can find the stationary state, i.e., dW / dt = 0 and dA / dt = 0, by looking where this two nullclines intersect. We then performed a stability analysis of the stationary points. A stationary point is stable if the system tends to return to this stable equilibrium even if there is a perturbation that deviates the populations from this equilibrium. In addition to this, we perform numerical simulations of the model using the function odeint of the scipy library from Python v3 for different test values of ρ and θ. The values of all other parameters are shown in Table I.First, we analyzed the population growth dynamics as the colony grows. We studied the stationary total hive population N and productivity P for different values of ρ to find its optimum value. This is done by taking the stationary value of N = W + A and P = β(1 − ρ)W + δA in simulations where t >> 1. The simulation was stopped when there was no observable change in this dynamical variable. We perform a simulation for different values of ρ starting from ρ = 0 and taking increments of ∆ρ = 0.0002 until a value of ρ = 0.02 was reached. We repeated this process for low γ = 1, medium γ = 2, and high γ = 5 per capita generation costs of anarchy. We also plot the total population cost of anarchy versus anarchy frequency (α). Finally, to find the most effective pattern of worker behavior in the presence of anarchy, we analyzed the optimum policing ρ versus anarchy per capita frequency α.

Table I Parameter values used in the model and the description of how these values were chosen
Full size table

3 Results

Based on the nullclines of the system (Figure 1), there is only one nontrivial stationary point. The stability analysis reveals that this nontrivial stationary point is always stable. Stability means that the hive will return to its stationary state after a small external perturbation as a consequence of its internal dynamics. Thus, we can conclude that anarchy in the beehive promoted by hive size cannot drive the colony to collapse.

Figure 1.
figure 1

a Nullclines. Continuous line: W-nullcline, dotted line: A-nullcline. b Population dynamics in the hive. Parameter values γ = 1, ρ = 0.0054, all other parameters as shown in Table I.

Full size image

According to the dynamics described by the model, the anarchistic behavior emerges just after the colony reaches a certain size (Figure 1b). Then, when the hive population stabilizes, the subpopulation of anarchic bees also stabilizes.

When the per capita costs from anarchy in the beehive are low (γ = 1), then worker policing has not a real benefit for the colony. This pattern can be observed in Figure 2a where the increment in hive population and productivity is marginal and even can decrease by making policing. By comparing the maximum of the curves on Figures 2a, b, and c with its value when ρ = 0, it becomes clear that as the per capita cost of anarchy increases (subfigure c), more beneficial are the effects of worker policing.

Figure 2.
figure 2

Stationary population and productivity for different values of policing fraction rate ρ. a When side costs of anarchy are low γ = 1, b side costs of anarchy are moderate γ = 2, and c side costs of anarchy are high γ = 5.

Full size image

Interestingly, maximum productivity (red line in Figure 2) is obtained with a lower worker policing rate than the rate needed to obtain the maximum of the population (blue line in Figure 2). These two different maxima are related to two possible optimum behaviors. If the colony prioritizes the viability of the next generation (blue line), then they have to perform policing with high frequency to decrease death rate, thus obtaining the largest possible population size. On the other hand, if the colony adjusts the worker policing to obtain maximum productivity, it means that bees prioritize the fertility of the queen because the productivity is directly related to offspring production rate. Thus, workers should perform less policing in this latter scenario.

We have related the increase of death rate due to anarchy behavior μγA / (1 + ρθW) with the cost of anarchy syndrome. From Figure 3a, we observe that there is an inflection in the resultant anarchy cost when the per capita anarchy frequency γ is increased. The inflection is caused because the intensity of policing is limited (Figure 3b). Once the maximum possible policing is reached, the resulting cost of anarchy cannot be further mitigated by the workers. Furthermore, if the beneficial effects of worker policing are low θ = 0.001, then the optimum policing ρ is reduced as per capita anarchy frequency γ increases, as shown in Figure 3c. The shifting to the left of the optimal worker policy means that if policing is not efficient, it is more beneficial to allocate more time to other duties.

Figure 3.
figure 3

a Net cost of anarchy = μγA / (1 + ρθW) at steady state as anarchy frequency α increases. b Optimal worker policing fraction rate ρ as a function of anarchy frequency α when worker policing is efficient at reducing the cost of anarchy θ = 0.003. c Optimal worker policing fraction rate ρ as a function of anarchy frequency α when worker policing is inefficient at reducing the cost of anarchy θ = 0.001.

Full size image

4 Discussion

The causes and consequences of anarchy in bees have challenged biologists for decades, and recently, we have started to gain insights into the underlying mechanisms driving this peculiar break in eusociality in bees (Barron et al. 2001; Beekman and Oldroyd 2003, 2008; Dampney et al. 2004; Oldroyd and Ratnieks 2000). Here, we used this accumulated biological knowledge about anarchy in bees to construct a dynamic mathematical model to unravel how anarchy and worker policing interact to determine hive productivity and stability.

Our model reveals that even if anarchy is costly for the colony, it cannot result in a collapse of the hive. Furthermore, this phenomenon would be hard to emerge in small-sized colonies because queen’s pheromones avoid workers laying eggs. In turn, worker policing is overall beneficial for the beehive as it tends to increase population size and/or fertility (Figure 3). This result emphasizes the importance of the queen’s regulatory mechanisms via pheromones to signal dominance and inhibit worker reproduction and promote worker policing (Hoover et al. 2005a, b; Wossler et al. 2001).

However, when excessive time is allocated to worker policing, this behavior can be detrimental to hive productivity. Therefore, there is an optimal effort to be invested in policing, which depends on the intrinsic costs associated with the anarchistic behavior, particularly the reduction in productivity due to larger drone populations and lower worker efficiency by anarchic bees. By comparing the panels in Figure 3, it is clear that the maximum productivity benefits are obtained with a higher policing frequency (ρ) for a higher intrinsic cost of anarchy γ. Figure 3a shows that if anarchistic behavior has low intrinsic costs for the colony (e.g., in small beehives), then policing has little or no benefits. On the other hand, if the intrinsic cost of anarchy is high, then worker policing has a big impact on mitigating this cost (see Figure 3c). Current literature reports empirical observations of different frequencies of anarchic behavior in beehives. For example, Barron et al. (2001) explain that worker policing, in conjunction with the pheromonal systems of the hive, maintains the reproductive division when the anarchy behavior is rare. But there are anarchic colonies where worker reproduction is more common. In particular, a standard beekeeping practice is to confine the queen to one area of the hive using a queen excluder. Colonies with male brood in regions of the hive off-limits to the queen are therefore believed to be anarchic.

We quantitatively analyzed the consequences of this natural variability in anarchism syndrome prevalence by analyzing the dependence of the effective cost of anarchistic behavior on the anarchy frequency. Our model suggests that worker policing prevents the anarchy cost to grow uncontrolled and dominate the hive, but there is an inflection point at a critical frequency of anarchy (Figure 3a). After this point, the costs associated with anarchistic behavior grow faster because worker policing reaches its maximum possible value (see Figure 3b). Therefore, even if anarchy behavior increases, policing intensity does not change. Furthermore, if the efficiency of worker policing in reducing the costs of anarchy is low, a decrease in the intensity of policing by workers is energetically a more convenient strategy (see Figure 3c). Conversely, when anarchy frequency is high, reduction in productivity of workers due to policing behavior overcomes the benefits of this patrolling behavior.

An important emerging property of our model is that the maximum population size and maximum productivity (and, thus, fertility) are achieved at different rates of policing. This difference implies that the balance between the intensity of worker policing and the frequency of anarchy can uncover whether the colony optimizes reproductive success in the next generation (i.e., fertility) or the life expectancy of the next generation (i.e., viability) (Sober 2001). Understanding when fertility or viability is optimized by beehives is an important next step to elucidate dynamics in eusocial societies.

Our mechanistic model offers new insights into how different levels of cooperation and conflict affect productivity and stability in eusocial animals. However, the functioning and dynamics of eusocial societies are inherently complex, and we suggest that future theoretical studies should also incorporate other direct and indirect effects known as relevant drivers of beehives’ productivity (e.g., interspecific competition, parasites, kin structure).