An energetics-based honeybee nectar-foraging model used to assess the potential for landscape-level pesticide exposure dilution

Landscape and resources

Foraging takes place in a landscape consisting of multiple resource patches. Each resource patch is assumed to be internally homogeneous and to contain a single resource. Resource patches may be fields with mass-flowering crops providing nectar or semi-natural habitats characterised by a dominant nectar-providing species.

Resources have a dynamic density of nectar-providing flowers F (m⁻²). In absence of nectar consumption, F remains F₀. The density of open flowers is the product of plant density and the number of open flowers per plant. Resources are characterised by the average amount of nectar g (mg) a honeybee can obtain from a flower, and by the typical sugar content of its nectar (g g⁻¹).

Resource patch selection

Being a central-place forager, a honeybee foraging trip comprises the travel back and forth between hive and resource patch, and the searching for and handling of flowers in the patch. Within the field, floral resources are collected assuming a type II functional response (Holling, 1959), so the number of flowers with nectar that are visited per time unit (s) amounts to (1)${\displaystyle f=\frac{aF}{1+ahF}}$ with a the attack rate (m² s⁻¹) and h the handling time per flower (s). The rate of nectar collection (mg s^-1) is given by fg. The time t_L it takes to collect nectar to full capacity is (2)${\displaystyle t_L=\frac{\gamma }{fg}}$ where the capacity is given by γ (mg). Note that unlike e.g., Schmid-Hempel, Kacelnik & Houston (1985) we assume that foragers always collect a full load.

Flight time and flight costs are proportional to the distance from hive to the resource patch. The duration of a foraging trip equals the sum of travel times and the time spent at the patch: (3)${\displaystyle t_{trip}=2\frac{D}{v}+\frac{\gamma }{fg}}$ with D being the distance from the hive to the field (m), and v the flight velocity (m s⁻¹).

The energy expenditure EE (J) ignoring basic metabolism can be calculated as the sum of the travel costs EE_travel (travel time × energetic costs per time unit) and the costs while loading nectar in the field EE_field. The latter term is made up by flight costs while searching for nectar flowers and costs while sitting on the flowers extracting the nectar. The latter is ignored, as it is approximately an order of magnitude smaller than flight costs (a value of 0.0042 J s⁻¹ was applied in Schmid-Hempel, Kacelnik & Houston (1985)). Energy expenditure at the field thus becomes: (4)${\displaystyle {EE}_{field}=\left(t_L-\frac{\gamma }{g}h\right)e_F=\frac{\gamma e_F}{gaF}.}$

The equation is simplified by substituting Eqs. (2) and (1), thus eliminating the handling time. By assuming that handling the flower doesn’t take energy, h does not play a role in the energetic balance. In Eq. (4) average flight cost e_F (J s^-1) is used, as during foraging in the field the individual state changes gradually from unloaded to loaded. The value of e_F is obtained from loaded e_F,L and unloaded e_F,U flight costs (Seeley, 1986). The total travel costs are given by: (5)${\displaystyle {EE}_{travel}=\frac{D}{v}\left(e_{F,U}+e_{F,L}\right)=2\frac{D}{v}{e}_F.}$ If the same route is followed back and forth, the energy costs can be averaged.

The total energy expenditure EE_total (J) for a foraging bout sums to: (6)${\displaystyle {EE}_{total}=\left(t_L-\frac{\gamma }{g}h\right)e_F+2\frac{D}{v}{e}_F=\left(t_L-\frac{\gamma }{g}h+2\frac{D}{v}\right)e_F=\left(\frac{\gamma }{gaF}+2\frac{D}{v}\right)e_F.}$

The yield of a trip in terms of energy, energy intake EI (J), depends on the energy content of the collected nectar of resource type R, e_R (J mg⁻¹): (7)${\displaystyle EI=\gamma e_R.}$ With t_trip, EE and EI, the basic ingredients for a “decision-making process” for a colony are specified, and costs (both in energy and time) and yields for specific foraging locations in a landscape with multiple fields can be compared. In theory, there are different currencies that might be optimized (Stephens & Krebs, 1987). For honeybees, considerable effort has been invested in deciding whether the relevant currency is the net rate of energy delivery ((gain − cost)/time) or the net energetic efficiency ((gain − cost)/cost). Experimental data, e.g., Seeley (1994) but see De Vries & Biesmeijer (2002), and models fit to experimental data (Schmid-Hempel, Kacelnik & Houston, 1985), as well as several of the other foraging models discussed by Becher et al. (2013) indicate that net efficiency is the most appropriate currency. We therefore assume that net energetic efficiency NEE (8)${\displaystyle NEE\mathrm{=}\frac{EI-{EE}_{total}}{{EE}_{total}}}$ defines the attractiveness of each patch. In the SO version, the patch with maximum NEE is selected as the single resource patch. In the RL version, the recruitment of foragers for a resource patch is explicitly modelled, following Camazine & Sneyd (1991) and Seeley, Camazine & Sneyd (1991) and extending their model in a somewhat similar way as done by Cox & Myerscough (2003). In the RL-model, NEE is translated into the probability of foragers abandoning resources and recruiting new foragers among the follower-bees present at the dance floor (Supplemental Information 2).

Resource detection

Resource patches are detected primarily by scouts. Scouts make up around 10% of the unemployed foragers: novice foragers and experienced foragers that have recently abandoned a depleted resource patch (Seeley, 1995). Detection of an attractive resource patch may depend on patch distance and size (Dauber et al., 2010), as well as on the number of active scouts. Encounter rates may be obtained from movement models or from dispersal functions summarizing movement simulations (Heinz et al., 2007) as σ_i values: the probability that a single scout encounters patch i. The probability for patch i being detected by at least 1 of s scouts amounts to (9)${\displaystyle P_i=1-{\left(1-{\sigma }_i\right)}^s.}$ When P is set to 1 for each resource patch–as we do in the studied cases–the hive population is assumed to have perfect knowledge of its environment. This setting is appropriate when dealing with a landscape with abundant resources in the immediate surroundings of the hive (Seeley, 1995).

Resource acquisition

The colony’s hourly acquisition of a selected resource depends on the absolute number of foragers n dedicated to this resource. For the SO-model, n equals the total number of foragers active in a time step. For the RL-model, n is a fraction of this total number, as multiple patches may be exploited simultaneously. Foragers may make several trips to the same patch, depending on trip duration t_trip and the time between trips t_UD spent on unloading and dancing in the hive. Thus, the number of foraging trips per hour an individual forager can make, amounts to b = Δt∕(t_trip + t_UD) with Δt representing the time step in seconds.

Not the full load of nectar γ collected in the patch will arrive at the hive as some of the sugar will be consumed on the way back. The amount of energy used during the return flight, E_c, is approximately $e_L\frac{D}{v}$. The amount of nectar arriving at the hive thus decreases to $\gamma -\frac{E_c}{e_R}$ per trip and when summed over all foragers exploiting the patch during this time step becomes: (10)${\displaystyle n\left(\gamma -\frac{E_c}{e_R}\right)b.}$

Resource dynamics

By calculating resource patch selection per hour we may account for: (1) any diurnal pattern in the anthesis of flowers; (2) depletion of floral resources occurring in small resource patches; and (3) diurnal patterns in the number of active foragers in the hive. For simplicity’s sake, we assume a binary pattern in anthesis (all flowers either open or closed).

Resource dynamics may also result from nectar consumption by competing species or foragers from other colonies. A constant background density of competitors, Z, may be defined of species that exploit resources in a similar way as honeybees (same functional response). The main reason to incorporate competition is that it allows us to explore honeybee foraging in a setting of quickly depleting resources.

Within a day there may be also renewal of nectar resources. With a non-zero renewal rate, r, a fraction of the pool of visited (empty) flowers will each time step turn into nectar-providing flowers again.

We deal with resource depletion by subtracting after each time step the visited number of open flowers from the initially present number. Hourly dynamics in open flower density F are given by: (11)${\displaystyle F_{t+1}=F_t-\frac{n\gamma }{gA}b-\frac{a{F}_t}{1+ah{F}_t}Z+r\left(F_0-F_t\right)}$ with A representing field size (m²). A fixed number of flowers γ∕g is visited for a full load. Lowered flower density may change the NEE-based ranking of resource patches in the next time step.

Dilution

Exposure dilution is caused by foraging on other resources than a treated resource when these other resources are not or to a lesser extent contaminated. Dilution is defined relative to a reference concentration in treated resources, e.g., the sugar-based concentration X_P in the nectar of a treated field with resource R: (16)${\displaystyle X_P=C_P\frac{e_{SUGAR}}{e_R}}$ where C_P refers to the concentration (on a wet-weight base) in nectar on the treated field, often referred to as the predicted environmental concentration (PEC) resulting from a certain application rate of the chemical. It is multiplied by e_SUGAR∕e_R to obtain the sugar-based concentration. Dilution factors follow from actual X_m or X_e as: (17)${\displaystyle {\phi }_m=X_m∕X_P\mathrm{and}{\phi }_e=X_e∕X_P.}$ With a patch-specific C we can represent variability in the applied dose for the treated crop, deal with a substance that is applied on different crops with different doses, and include patches representing off-crop or off-field habitats that are exposed through spray drift only.

For an agricultural landscape containing fields with a target crop as the only resources, Eqs. (14) and (15) can be further simplified. With a single resource e_R,i equals e_R and when application of a pesticide on the target crop leads to a constant concentration C_P in the nectar, C_i equals C_P. Dilution obtained applying Eq. (14) then simplifies to: (18)${\displaystyle {\phi }_e=\frac{\sum _{i=1}^L{n}_i\gamma b_i{C}_P}{\sum _{i=1}^L{n}_i\frac{\left(\gamma e_R-E_{c,i}\right)}{e_{SUGAR}}{b}_i}\frac{1}{C_P\frac{e_{SUGAR}}{e_R}}=\frac{\sum _{i=1}^L{{\delta }_{i}n}_i\gamma b_i}{\sum _{i=1}^L{n}_i\left(\gamma e_R-E_{c,i}\right)b_i}}$ with δ_i = 1in case the field selected in time step i is sprayed and δ_i = 0 if it is unsprayed. Similarly, for the case with metabolization and no enrichment (Eq. (15)) we obtain: (19)${\displaystyle {\phi }_m=\frac{\sum _{i=1}^L{{\delta }_{i}n}_i\left(\gamma -\frac{E_{c,i}}{e_R}\right)b_i{C}_P}{\sum _{i=1}^L{n}_i\frac{\left(\gamma e_R-E_{c,i}\right)}{e_{SUGAR}}{b}_i}\frac{1}{C_P\frac{e_{SUGAR}}{e_R}}=\frac{\sum _{i=1}^L{{\delta }_{i}n}_i\left(\gamma e_R-E_{c,i}\right)b_i}{\sum _{i=1}^L{n}_i\left(\gamma e_R-E_{c,i}\right)b_i}.}$

Eqs. (18) and (19) show that for this simplified case, dilution factors do not depend on C_P but simply reflect how many unsprayed fields are selected. When all selected fields are treated φ_m will be one, while φ_e may even exceed one. Whether selected resource patches were treated or not depends on the application scenarios. In a probabilistic approach, with p representing the probability for a field of being treated, dilution will on average approximate p.

90^th percentiles

In an exposure assessment the model will be applied for all possible locations of an apiary within the area of use of a substance. Following recent guidance (EFSA, 2013), all sites at edges of target crop fields are considered to be potential locations and the target crop patch adjacent to the hive location is always assumed to be sprayed or treated. From application of the model at each site, we construct cumulative distributions of exposure or dilution from which statistical descriptors can be obtained, e.g., the 90th percentile that quantifies the dilution obtained at 90% of all sites, or the exposure that is exceeded at 10% of the sites.

Energetics-based model

From the energetics-based model (Eqs. (6) and (7)) we derived thresholds for the exploitation of resources depending on distance and resource characteristics, as done by Cresswell, Osborne & Goulson (2000). For resources at equal distance we derived a simple equation to compare their attractiveness. For landscapes consisting of only large-scale crop fields, we derived rules of thumb for the selection of foraging locations.

Landscape representation

When it is assumed that nectar is collected from all resources in the landscape proportional to their attractiveness (EFSA, 2013), exposure dilution is always to be expected in landscape mosaics with multiple resource patches. In the RL-model dilution is more likely than in the SO-model. To explore whether landscape-level exposure mitigation can be effective, we tested for both models three scenarios, all based on the idea that by providing alternative floral patches foraging effort can be diverted away from exposed patches. Construction of flower strips and pollinator-friendly management of semi-natural habitats are such measurements that may contribute to the persistence of pollinator populations in the agricultural landscape (Garibaldi et al., 2014; Haaland, Naisbit & Bersier, 2011; Wratten et al., 2012). Also, the presence of different (early and late) mass-flowering crops has been suggested to enhance pollinator density (Riedinger et al., 2014).

We parameterized the model for two flower species (Table 2), oilseed rape (Brassica napus, OSR) as a common attractive mass-flowering crop, and white clover (Trifolium repens) being an important floral resource, common along roadsides and field margins (Sponsler & Johnson, 2015). Patches with these resources provided the building blocks of the landscapes tested in the scenarios. In all scenarios OSR represented the target crop. From the model it followed that for any of the scenarios to be effective the alternative resource had to be equally or more attractive (larger NEE at equal distance) than the OSR field. For the scenario study, using a hypothetical landscape, we selected for OSR a conservative and for clover an optimistic estimate of the relevant coefficients (Table 2 and Supplemental Information) and checked whether this condition was met (see ‘Results’). In real landscapes clearly more resource types will be present, in a range of densities. To understand the impact of some of this variability we varied flower densities and sugar content over a range of values.

We defined landscape structure from geographic data for the northern part of Flevoland (The Netherlands). Geo-referenced data for crops were obtained from the land-use database LGN6 (grid-based, 25 m resolution) (Hazeu et al., 2011). Data on the presence of off-field habitats (road side verges, ditch sides, shores, other semi-natural elements) were obtained from the vector-based dataset TOP10NL (PDOK, 2015). We assumed that fields of LGN6 category “other crops” represented OSR. For the “Alternative Fields” scenario, fields of another (randomly chosen) LGN6 category represented the alternative mass-flowering crop (Fig. 1). For the “Off-field Habitats” scenario, off-field habitats were used as specified in the spatial data set. All sites (cells in the 25 m grid) that bordered target crop fields were considered as potential beehive locations (Fig. S1). To avoid excessive computing times the model was finally run on a random selection of 10% of these sites (identical for all scenarios and runs). Tested landscapes contained 1,207 fields with OSR (7,057 ha). The alternative crop was present on 454 fields (3,375 ha). Off-field habitat patches were small and numerous: 58,137 patches comprising 3,410 ha.

We assumed all resources to have open flowers during the whole day, a foraging day of 10 h, no competitors being present (Z = 0), no renewal of resources (r = 0), 100 active foragers and coefficients as in Tables 1 and 2. By definition, sites were next to a large resource patch. Therefore, considering only resources within 2 km distance (D_max = 2000) and assuming perfect knowledge of their presence (P = 1 for each resource patch) seemed reasonable. For simplicity’s sake we assumed that the hypothetical pesticide was metabolized with the consumed sugar on the way back to the hive, so without dilution the dilution factor would simply be one. The nearby crop field was assumed to be always treated with the pesticide.

Exposure mitigation

(i) Alternative Fields: dilution can result from the presence of attractive but untreated fields, either with the target crop or with alternative attractive (mass-flowering) crops. For crop fields other than the one nearest to the hive, treatment occurred following a probability p. The average dilution factor for a site was then calculated from 100 different random realizations of the series of treated fields. The analysis was done for a range of values for p (0–1, with an interval of 0.1). A variant of this mechanism would be the presence of fields with another, attractive, crop type. We tested this by defining another common crop in Flevoland as a hypothetical alternative mass-flowering crop, largely identical to oilseed rape (Fig. 1). The energetic attractiveness of this crop was manipulated to range from less to more attractive than OSR, by adjusting its sugar content (0.8, 0.9, 1, 1.1 and 1.2 times OSR sugar content). Thus, the alternative crop could also represent another OSR variety.

(ii) Flower Strips: dilution can result from the presence of highly attractive flower strips around target crop fields. This was simulated by adding to each target crop field up to 4 flower strips at its edges, each strip of length 1/4 of field perimeter. Strips contained white clover (Trifolium repens, Table 2) in densities comparable to clover fields. Presence of each strip was randomly set, with a probability p (0–1, with an interval of 0.1). The area of a flower strip was set to a prescribed width w multiplied by strip length. Width values of 1, 2, 5 and 10 m were tested. The area contained in the strips was subtracted from the crop field area. Thus, strips were strictly in-field off-crop habitats.

(iii) Off-field Habitats: dilution can result from the presence of off-field semi-natural habitats when these are managed in an appropriate way. We tested this using geo-data for fields (as above) with OSR and for semi-natural elements (Fig. 1). We assumed a single common flower species (white clover) to be representative for all off-field habitats. As the presence and size of these habitats were fixed and defined by the geo-data, we tested their impact on dilution for a range of resource quality values. Open flower density determined quality and was set to low, medium or high value (0.5, 1 or 1.5 times clover field density). The three quality levels were tested for a range of values for p, here representing the probability that an off-field element was considered a nectar resource patch.