Using ITS2 metabarcoding and microscopy to analyse shifts in pollen diets of honey bees and bumble bees along a mass‐flowering crop gradient

Worldwide pollinator declines lead to pollination deficits in crops and wild plants, and managed bees are frequently used to meet the increasing demand for pollination. However, their foraging can be affected by flower availability and colony size. We investigated how mass‐flowering oilseed rape (OSR) can influence the pollen resource use of small and large honey bee (Apis mellifera L.) and bumble bee (Bombus terrestris L.) colonies. Colonies were placed adjacent to strawberry fields along a gradient of OSR availability in the landscapes. We used ITS2 metabarcoding to identify the pollen richness based on ITS2 amplicon sequencing and microscopy for quantification of target pollen. Bumble bees collected pollen from more different plant genera than honey bees. In both species, strawberry pollen collection decreased with high OSR availability but was facilitated by increasing strawberry flower cover. Colony size had no effect. The relationship between next‐generation sequencing‐generated ITS2 amplicon reads and microscopic pollen counts was positive but pollen type‐specific. Bumble bees and, to a lesser degree, honey bees collected pollen from a wide variety of plants. Therefore, in order to support pollinators and associated pollination services, future conservation schemes should sustain and promote pollen plant richness in agricultural landscapes. Both bee species responded to the availability of flower resources in the landscape. Although honey bees collected slightly more strawberry pollen than bumble bees, both can be considered as crop pollinators. Metabarcoding could provide similar quantitative information to microscopy, taking into account the pollen types, but there remains high potential to improve the methodological weaknesses.


| INTRODUC TI ON
Pollinators not only contribute to the yield and quality of many crops, but also deliver pollination services to many wild plants, enhancing seed set (Klatt et al., 2014;Klein et al., 2007;Wietzke et al., 2018). However, pollinators worldwide are at risk due to multiple stressors, such as agricultural intensification, habitat loss, and accompanying reduction in the diversity and abundance of host plant species Vaudo et al., 2015). At the same time, global pollinator-dependent crop production is intensifying and the demand for pollination services is increasing (Aizen & Harder, 2009). When native and domesticated pollinators are rare or absent, farmers are exposed to high economic risks due to reduced pollination rates (Potts et al., 2016).
The use of managed bees, such as honey bees (Apis mellifera L.) and bumble bees (Bombus terrestris L.), can contribute to meeting the pollination demand in crop production. Both species are generalist pollinators that visit a great variety of plant species and are suitable for the pollination of many crops (Kleijn et al., 2015). Crop pollination can be promoted through the use of honey bee pollination services offered by beekeepers or by placing commercially available bumble bee colonies in or next to crops. The use of bumble bees in particular has become a widely used practice for pollination in glasshouses, but is also frequently used in fields (Gosterit & Baskar, 2016;Velthuis & van Doorn, 2006), while keeping honey bees is an established practice in many field crops (Klein et al., 2007;Potts et al., 2010).
Bees collect pollen because it provides proteins, vitamins and minerals for larval development (Thorp, 2000). Bees often forage on a large number of different plant species to meet their dietary requirements (Alaux et al., 2017;Di Pasquale et al., 2013;Leidenfrost et al., 2020) and they can balance nutrient deficits by collecting either greater amounts or diversity of pollen (Danner et al., 2017;Hendriksma & Shafir, 2016). Colony development can be enhanced through increasing intake of higher amounts and/or better pollen quality in terms of pollen diversity or species composition (Brodschneider & Crailsheim, 2010;Kämper et al., 2016;Vaudo et al., 2018). To counteract ongoing bee declines and to sustain vital populations of managed and wild pollinators in agricultural landscapes, it is important to understand the temporal and spatial dynamics of pollen resource exploitation (Bertrand et al., 2019;Kämper et al., 2016), the effects of pollen richness on reproductive success (Hass et al., 2019;Requier et al., 2017;Vaudo et al., 2018), and plantpollinator interactions and networks (Arceo-Gomez et al., 2020;Bell et al., 2017). High flower and pollen constancy of bees may imply high carry-over of the respective pollen on the stigma, and hence could be used as an indicator of pollination efficiency (Gyan & Woodell, 1987;Marzinzig et al., 2018;Montgomery, 2009).
Pollen richness can be investigated using ITS2 (internal transcribed spacer) metabarcoding. In comparison to traditional methods (e.g., microscopy), metabarcoding can provide a higher resolution of taxon richness and has a higher throughput with a predictable cost and time frame (Bell et al., 2019;Keller et al., 2015;Smart et al., 2017). However, metabarcoding is restricted in its applicability in pollen quantification. Several factors may affect the quantification results obtained by next-generation sequencing (NGS)-based amplicon sequencing, such as DNA extraction, PCR (polymerase chain reaction) amplification or barcode copy number (Peel et al., 2019), which can lead to over-or underestimation of actual pollen counts (Baksay et al., 2020;Pornon et al., 2016;Richardson et al., 2015).
The number of metabarcoding reads may be more likely to be positively related to microscopy-based counts for commonly occurring taxa than for rare pollen taxa in the samples (Smart et al., 2017).
Hence, the abundance of metabarcoding amplicons could be used as an estimate of relative abundance but should be applied with respect to the investigated taxa and research questions (Danner et al., 2017;Nürnberger et al., 2019;Smart et al., 2017). Complementary to metabarcoding, microscopic pollen analysis can be used for pollen quantification as fewer processing steps are needed that could bias the pollen counts. However, great expertise is needed to identify all plant taxa present in pollen pellets using microscopy , but target pollen can be counted by nonexperts, especially when pollen grains have a characteristic surface structure (Beug, 2015). However, rare pollen taxa may not be detected with microscopy, which can also depend on the number of pollen grains counted per sample (Lau et al., 2018;Smart et al., 2017). Up to now, only a few studies have compared quantitative results from metabarcoding and traditional microscopy and indicate that outcomes are not necessarily correlated and can depend on plant species and species composition (Bell et al., 2019;Richardson et al., 2015).
Bumble bees build colonies with up to 600 individuals (e.g., B. terrestris) while honey bee colonies (A. mellifera) can achieve colony sizes of up to 80,000 individuals (Amiet & Krebs, 2012). Thus, both bee species can provide many foraging and pollinating individuals. Yet, using large colonies for pollination services does not guarantee a high number of individuals in target crop fields because both honey bees and bumble bees are able to explore landscapes extensively due to their large foraging radii (Steffan-Dewenter & Kuhn, 2003;Westphal et al., 2006). However, bumble bees tend to exploit more diverse resources in the close surroundings of their colonies (Wolf & Moritz, 2008) because they are not able to communicate as effectively as honey bees do by waggle dance (Couvillon, 2012;Dornhaus & Chittka, 2001).
Honey bees in particular but also some bumble bee species, such as B. terrestris, preferably exploit highly rewarding mass-flowering crops, such as oilseed rape (OSR, Brassica napus L.; Rollin et al., 2013).
The availability of mass-flowering resources within the foraging ranges can lead to lower bee densities in minor flowering crops, such as strawberries (Fragaria × ananassa Duch.) (Bänsch, et al., 2020a;Grab et al., 2017) and hence potentially affect pollination services.
Simply increasing the number of managed bee colonies in crop fields is presumably not the optimal solution for pollination management because this would result in higher costs for farmers and high managed bee densities may negatively affect wild bees (Herbertsson et al., 2016;Mallinger et al., 2017). As an alternative, the selection of pollinator species and adaptation of colony size may be more efficient in directing foraging bee pollinators to crops in agricultural landscapes. To our knowledge, research on this possible adaption of colony size in pollinator management is still very rare but can have severe implications for pollination management. A few studies have suggested that foraging can vary between colonies of different sizes (e.g., shorter foraging distance of smaller colonies ;Beekman et al., 2004;Boecking & Kreipe, 2015;Westphal et al., 2006).
In general, we focus on the following questions in the present study: Which pollen resources do managed pollinators use in agricultural landscapes and to what extent? Is the pollen foraging behaviour bee species-specific and influenced by (mass-)flowering resource availability and colony size? Do quantitative analyses of pollen samples based on ITS2 metabarcoding and microscopy yield comparable results?
More specifically, we hypothesized that honey bees and bumble bees will use flowering crops (i.e., strawberry and OSR) as major pollen resources in agricultural landscapes because flowering crops provide ample pollen and can cover large areas. Furthermore, we expected that bumble bees collect pollen from a greater number of plant genera than honey bees and that they collect pollen from resources close to the colonies (i.e., strawberry) in greater proportions and more frequently due to less effective communication.
Considering the bees' preference for mass-flowering crops, we hypothesized that pollen richness and the proportion of strawberry pollen, as a minor flowering target crop, will decrease in the pollen loads of both bee species with increasing OSR availability. Larger colonies were expected to collect greater pollen richness while smaller colonies should collect greater proportions of crop pollen close to the hive. Due to several factors that influence outcomes of pollen quantification based on metabarcoding (Bell et al., 2019;Peel et al., 2019), we hypothesized a positive and pollen type-specific relationship between the number of metabarcoding reads and microscopic pollen counts for the two target crop species, namely strawberry and OSR.
We investigated the pollen foraging preferences of two managed bee species, A. mellifera and B. terrestris, during the coflowering of two economically important crops (strawberry and OSR) by placing large and small honey bee and bumble bee colonies next to strawberry fields in agricultural landscapes that represented a gradient of mass-flowering OSR availability. Making use of the advantages offered by ITS2 metabarcoding, we identified the plant genera that were present in mixed pollen samples collected by the honey bee and bumble bee colonies. Traditional microscopy was applied to quantify the proportions of strawberry and OSR pollen in the samples. Finally, we analysed the pollen type-specific relationships between the number of metabarcoding reads and microscopic pollen counts for the two crop species as only a few studies have compared pollen quantification by metabarcoding and microscopy so far, and these did not achieve clear results (Bell et al., 2019;Richardson et al., 2015).

| Study location
We studied the pollen foraging behaviour of bee colonies adjacent to nine commercial strawberry fields in central Germany in the surrounding regions of the cities Göttingen (southern Lower Saxony) and Kassel (northern Hesse) in 2016 ( Figure S1). Most study fields were managed for public pick-your-own harvesting and usually had three different strawberry varieties to extend the harvesting season. Study field sizes were on average 2.23 ha (±0.92 SD, range 1.01-3.61 ha). The strawberry study fields were surrounded by an agricultural landscape matrix which we mapped using a geographical information system (ESRI ArcGIS, version 10.3.1). We classified the land cover types into strawberry fields, OSR fields, other cropland (mainly nonentomophilous crop fields with adjacent field margins, single trees and country lanes), seminatural habitats (e.g., hedges, fallow land and meadow orchards), forests and urban areas. The landscapes were mapped within a typical honey bee and bumble bee foraging distance of about 2,000-m radius around our study fields (Härtel & Steffan-Dewenter, 2014;Westphal et al., 2006).
Availability of OSR was calculated as the product of OSR land cover (ha) and OSR flower cover (%). OSR flower cover was estimated at each observation date within the nearest OSR field to the strawberry field along a transect of 50 × 4 m. All OSR fields had relatively uniform germination and flowering within our study landscapes. To validate the OSR gradient, we assessed the OSR land cover in autumn 2015 (winter OSR plants can already be identified in autumn).
Strawberry flower cover within the adjacent field was determined by counting the open flowers along 2 m of the row with greatest flower abundance, as this area was probably the most attractive for bees.

| Experimental set-up
We established an experiment with small and large honey bee (Apis mellifera carnica Pollmann 1879) and bumble bee (Bombus terrestris Linnaeus, 1758) colonies adjacent to strawberry fields in order to study pollen foraging in relation to OSR availability and strawberry flowering. One small and one large colony of each species was placed at the edge of each of the nine study fields (over a distance of ≤5 m to the strawberry cultivation). Hence, we studied 36 colonies. Large honey bee colonies had around 20,000 workers at the beginning of the season and one queen. Small honey bee colonies were built as nuclei with around 4,000 workers from additional large colonies. Honey bee colony sizes were estimated following the "Liebefeld Method" by visually estimating the number of adults on the comb surface (Dainat et al., 2020;Imdorf et al., 1987). All small colonies successfully raised their own queen, which emerged a few days after experimental set up in the field. Even immature queens produce queen substance pheromone and stimulate pollen collection in foragers (Boch, 1979;Free et al., 1984). Hence, we do not expect large differences in foraging behaviour due to the queens' age. At the end of the experiment, large honey bee colonies comprised ~30,000-35,000 workers and small colonies ~6,000-8,000 workers. Small and large bumble bee colonies (B. terrestris) were purchased from a German bumble bee breeder (STB Control). The colonies consisted of one queen bee and 40 and 80 workers, respectively. To monitor the development of the bumble bee colonies, we weighed them as the number of individuals is difficult to quantify (e.g., they may hide in the complex structure of their nests; Lefebvre & Pierre, 2006). We monitored the colony weight during the first observation round (small colonies: mean 1,045.33 g ± 45.68 SE and large colonies: mean 1,155.11 g ± 50.59 SE) and third observation round (small colonies: mean 1,533.50 g ± 69.85 SE and large colonies: mean 1,754.75 g ± 114.32 SE). Data collection began on May 6, 2016 with the beginning of the strawberry blossom.

| Pollen sampling and preparation
We collected pollen loads of honey bees and bumble bees in front of their colonies during the strawberry flowering period on, if possible, five observation days per study landscape. Because we could not analyse pollen from each colony at each sampling date (e.g., due to low colony activity and hence pollen material below 0.05 g), the number of samples per colony type can differ (small honey bee colonies n = 34, large honey bee colonies n = 40, small bumble bee colonies n = 38, large bumble bee colonies n = 37). We set a threshold of 0.05 g to have enough pollen material for the metabarcoding process (0.015 g), repetitions in case something went wrong (e.g., contamination) and microscopy. The study period lasted from May 6 to June 6, 2016, depending on the microclimatic conditions within the study landscapes.
Pollen loads from honey bees were collected using commercial pollen traps installed in front of the colony. The traps guide the bees through a 5-mm grid that removes pollen loads from the hind tibia. Traps were installed in front of each colony for 30 min on each observation day.
Pollen loads from bumble bees were collected by capturing, if possible, five individual bumble bees in front of their colonies with an insect net and placing them into marking cages, respectively. Pollen was removed from the hind tibia with tweezers and stored in 1.5-ml reaction tubes. Bumble bees were released after this procedure. To account for foraging preferences of bees for either pollen or nectar resources at certain times of the day we varied sampling times across landscapes systematically at each visit (i.e., visiting each field equally in morning hours and afternoon hours between 9 a.m. and 5 p.m., respectively).
Pollen sampling was only conducted on days with low wind speed, no rain and a minimum temperature of 14°C.
We pooled the pollen loads of each observation date by colony and homogenized them in 70% ethanol (ratio 1:4 pollen: 70% ethanol).
From this mixture, we prepared 1-ml aliquots in 1.5-ml reaction tubes for microscopic and molecular pollen analysis. The tubes were centrifuged for 10 min at 15,400 g. We then removed and discarded the supernatant ethanol. Afterwards, samples were dried with open lids in a fume cupboard with an air throughput of 1,000 m 3 /hr for 72 hr.
To study the pollen collection constancy (homogeneity in pollen samples) of in-field foragers in strawberry fields, we collected honey bee (n = 37) and bumble bee individuals (n = 36) visiting strawberry flowers. Standardized transect walks were conducted at two locations within strawberry fields (2 × 50 m in 30 min at the edge and in the centre of the field) for the two most common strawberry varieties (Sonata and Honeoye) on five observation dates in seven out of the nine landscapes in 2015. Transect times were varied systematically between 9 a.m. and 5 p.m. (i.e., visiting each field equally in morning hours and afternoon hours, respectively). The collected bees were washed with 70% ethanol individually in 1.5-ml reaction tubes in order to remove the pollen from the bees' bodies.
Subsequently, the bees were removed, and the pollen-ethanol mixture was centrifuged for 10 min at 15,400 g to remove and discard the residual ethanol and dried as described in the paragraph above.

| Pollen analysis
We analysed the richness and amount of target pollen grains in pollen loads using two methodologies: ITS2 metabarcoding and microscopy.
To assess the relationship and variability between the quantitative outcomes for different pollen types, we compared the number of ITS2 amplicon reads and pollen grain counts for Fragaria and Brassica pollen types in the samples. Both methods have the same taxonomic resolution in our study. We assume that ITS2 amplicon reads from Fragaria and Brassica indicate strawberry and OSR, as they were the most common flowering plant species belonging to these genera in our study landscapes during the study period. Fragaria vesca flowers at the same time as Fragaria × ananassa but flowers in much lower abundancies in our study landscapes, and mainly in woody habitat structures. Other congeneric species of Brassica type flower most likely only later in the year (e.g., in flower strips).
We used metabarcoding of ITS2-region PCR amplicons to quantify the richness of pollen collected by small and large honey bee and bumble bee colonies. In the present study, we used the advantages of metabarcoding techniques (e.g., high efficiency and resolution; Baksay et al., 2020;Bell et al., 2019) for qualitative high-throughput identification of PCR-amplified ITS2 sequences from plant genera present in pollen loads collected from honey bees and bumble bees to study pollen richness. Finally, the standard protocol was followed until the DNA was eluted with 50 µl of elution buffer. DNA concentration and quality were measured using a Nanodrop (Thermo Scientific).
Final extension was performed at 72°C for 5 min.
All reactions were checked for successful amplifications and contaminations by gel electrophoresis (1.5% agarose gels stained with ethidium bromide, 120 V for 30 min). Triplicate PCR products were pooled per sample and purified using the QIAquick PCR purification Kit. Then, 500 ng of each PCR product was used for library preparation using NEBNext Ultra II DNA Library Prep Kit for Illumina according to the manufacturer's protocol (New England Biolabs). Paired-end sequencing (2 × 150 bp) was performed on a NextSeq500 platform using a Mid-output flowcell (150 cycles).
To increase the accuracy of assignment of amplicon sequencing reads to plant-specific ITS2 sequences, we extracted all ITS2 sequences from a global eukaryota database (Förster, 2015) that have previously been described for plants occurring in Lower Saxony, Germany (Garve, 2004(Garve, , 2007. The resulting subset was made nonredundant by clustering identical entries with vsearch (version 2.9.1; Rognes et al., 2016) and then used to create a magicBLAST database (version 1.4; Boratyn et al., 2019). After blasting the ITS2 amplicon reads against this database, all paired reads that both aligned to a database entry (plant ITS2 sequence) with at least 50 bp each and had a similarity >98% were kept.
For each matching read, we calculated an alignment score by multiplying the alignment length with the alignment identity; the scores for the forward and reverse read were summed to get the final score for each read-pair. Read-pairs that matched several entries were ordered by this score. Only the top scoring match (plant species) per read was counted. As some plant species have very similar ITS2 sequences and we therefore cannot unambiguously distinguish them on a species level, we decided to use amplicon readbased identification at the genus level only. If there were multiple scoring matches with an identical score, we decided on the match with higher reliability based either on personal observations of flowering plants in the field or otherwise a distribution atlas of plants in Lower Saxony (Garve, 2007). Only in the case of Hirschfeldia incana L., we decided to reject the first match because several other matches for this read indicated the genus Brassica (e.g., Brassica napus L., Brassica rapa L.). In comparison with other studies (e.g., Danner et al., 2017;Nürnberger et al., 2019) that calculated the relative abundance of plant taxa based on metabarcoding sequences, we based our analysis on the presence/absence of ITS2 sequences of certain plant genera only, as the quantitative output (e.g., number of reads) can be biased (Sickel et al., 2015). Pollen richness represents the number of plant genera occurring in one pollen sample.
Strawberry and OSR pollen grains were quantified by microscopy for colony samples and the in-field foragers. One aliquot per sample was diluted with distilled water (ratio 1:4 pollen: 70% ethanol). One drop of the pollen-water mixture was applied to a microscopic slide together with one drop of Kayser's gelatine stained with fuchsine and fixed with a cover slide. We counted 200 pollen grains at 400× magnification on each slide (Lau et al., 2018). For this, we randomly selected one or when necessary more rows on the slides until we reached a number of 200 pollen grains and categorized the pollen into strawberry pollen, OSR pollen and others according to a selfmade reference collection and a determination key (Beug, 2015).

| DATA ANALYSIS
All statistical analyses were performed with the software r version 3.6.2 (R Core Team, 2016). Continuous explanatory variables (i.e., OSR availability and strawberry flower cover) were scaled to a mean of zero and a standard deviation of one to improve convergence of the models using the scale function (R Core Team, 2016). We found only weak if any correlation between fixed effects (−0.3 > r < 0.3; Hinkle et al., 2003). Data were visualized using the package ggplot2 (Wickham, 2017) and mixed model fit was visualized using the package effects (function allEffects, Fox & Weisberg, 2019).

| Effects of bee species, colony size and massflowering resource availability on pollen richness
To analyse the effects of bee species, colony size and flowering resource availability on pollen richness we fitted generalized linear mixed effects models using the glmmTMB package (Brooks et al., 2017). In a first step, we fitted a global model with pollen richness as the response variable and bee species (honey bee/bumble bee), colony size (small/large), OSR availability, and strawberry flower cover and all two-way interactions as explanatory variables. Colony ID nested in study landscape was included as a random effect. Two global models containing all explanatory variables were fitted with Poisson and negative binomial error distributions, respectively, and compared using the second-order Akaike Information Criterion (AICc) (Burnham & Anderson, 2004). The smaller the AICc, the better is the fit of the model.
We decided to use the negative binomial model because the AICc was lowest and over-dispersion was detected in the Poisson model. We applied the multimodel inference approach on our global model using the function dredge (package MuMIn; Burnham & Anderson, 2002; Barton, 2018), which creates a list of candidate models with all possible model combinations. To avoid overfitting, we limited the number of parameters in each candidate model to three (Crawley, 2007). The appropriateness of model assumptions was assessed by plotting residuals vs. fitted values. We ranked the models by AICc. All models within delta AICc (dAICc) <2 from the best fitting model were considered to have substantial empirical support and are reported together with the null model (Table 1). The relative importance of each explanatory variable was assessed using the sum of Akaike weights (Σw i ) over all candidate models which included the respective variables (function importance, package MuMIn; Barton, 2018). We considered all explanatory variables in the best fitting models (dAICc < 2) if Σw i > 0.2 to explain the effects on our response variables. We applied post hoc comparisons using the function emmeans to test for differences between bee species and colony sizes with alpha of 0.05 (Lenth, 2017).
To analyse how pollen composition in pollen loads of small and large honey bee and bumble bee colonies differ, we performed presence/ absence-based nonmetric multidimensional scaling (NMDS, function nmds, package vegan; Oksanen et al., 2019) using Bray Curtis dissimilarity (Clarke et al., 2006). Differences between colony types were tested with an analysis of variance using distance matrices including study landscape as a strata variable (function adonis2, package vegan ; Oksanen et al., 2019).

| Pollen collection
In order to describe the pollen resource utilization of honey bee and bumble bee colonies, we created heat maps displaying the number of ITS2 amplicon reads obtained from metabarcoding. Using only the presence of amplicon reads for certain plant genera, we extracted the five most common plant genera in pollen samples of small and large colonies summarized over the whole study period (i.e., frequency). These data are merely described and not statistically analysed.

| Colony level
Microscopic pollen counts were used to determine the proportion of strawberry and OSR pollen. The effects of bee species, colony size, OSR availability and strawberry flower cover on the proportions of collected strawberry pollen were analysed with generalized linear mixed effects models (glmmTMB, Brooks et al., 2017). We fitted global models with the proportion of strawberry as the response using the cbind function. Explanatory variables were bee species, colony size, OSR availability, strawberry flower cover and all two-way interactions. Because a binomial model resulted in over-dispersion, the global model was fitted using a betabinomial error distribution.
We included colony ID nested in study landscape as random effects to account for our nested study design. We then followed the multimodel inference approach as described in the previous paragraph, and again, allowed only three parameters in candidate models.

| In-field foragers
To analyse the proportion of strawberry pollen in the pollen loads of honey bees and bumble bees (B. terrestris) in the field, we fitted generalized linear mixed effects models with a beta-binomial error TA B L E 1 Summary of best fitting candidate models (dAICc < 2) and null models for the response variables (a) pollen richness and (b) strawberry pollen collection Note: Candidate models are a subset from the global models with a maximum of three parameters. Global models for both pollen richness and strawberry pollen collection include OSR availability, bee species, colony size, strawberry flower cover and all two-way interactions as explanatory variables. The response variable for pollen richness is the number of plant genera in pollen samples based on metabarcoding and the response variable of strawberry pollen collection is the proportion of strawberry pollen found in colony pollen samples counted with the microscope. Model estimates for models indicated with † are shown in Table S2. Explanatory variables: OSR = OSR availability, species = bee species (honey bee/bumble bee), size = colony size (small/large), straw_fc = strawberry flower cover. a dAICc = 0.004. distribution (glmmTMB, Brooks et al., 2017). The proportion of strawberry pollen was used as the response variable using the function cbind and bee species as the explanatory variable. The location of bee collection in the field (edge/centre) was nested in the study landscape as a random effect. We applied post-hoc comparisons to test for differences between bee species with alpha of 0.05 (function emmeans, Lenth, 2017).

| Relationship between quantitative outcomes of metabarcoding and microscopy
To test for pollen type-specific associations between quantitative outcomes of ITS2 metabarcoding and microscopy, we fitted linear regression models for log ITS2 amplicon reads + 1 versus log microscopic pollen grain counts + 1 for each crop type (i.e., strawberry and OSR) separately.

| Effects of bee species, colony size and massflowering resource availability on pollen richness
The effects of bee species, colony size, OSR availability and strawberry flower cover on the richness of pollen, collected by honey bee and bumble bee colonies, were explained by several models with empirical support (dAICc < 2) (Table 1a). Based on the relative importance of each explanatory variable, assessed using the sum of Akaike weights (Σw i ), we found that greatest importance was assigned to the effect of bee species, indicated by a high Σw i of 1, followed by the effect of size and OSR availability (Σw i = 0.45 for both, Table S1). Bee species and OSR availability were included in the best fitting model (dAICc = 0, Table 1a). Pollen richness revealed by metabarcoding was 4.9 times higher in bumble bee colony samples compared to honey bees ( Figure 1a) and increased with increasing OSR availability (Figure 1b; Table S2a). Colony size also affected pollen richness, but to a smaller extent. In general, large colonies collected about 20% more different pollen genera than small colonies ( Figure S3a). Strawberry flower cover (Σw i = 0.2) correlated negatively with pollen richness, independently of bee species and colony size. However, this effect is only of minor importance as indicated by the low sum of Akaike weight and low effect size (see model estimates in Table S2a, Figure S3b). The interaction between species and size was included within the best fitting models, but the sum of Akaike weight was quite low (= 0.11) and hence not considered to have a substantial effect on our response variable. We calculated model-averaged coefficients which support our results ( Figure S4a).

| Pollen community composition
The taxonomic composition of the pollen samples originating from small and large colonies was very similar within bee species but differed significantly between bee species (R 2 = .46, p = .005, stress value = 0.18; Figure S5).

| Pollen collection
ITS2 sequences of Fragaria and Brassica were consistently identified in pollen loads of returning honey bees continuously during the study period (Fragaria sequences were found in 57 samples and Brassica sequences in 60 samples from 77 pollen samples in total, Table 3).
Other pollen resources were typically restricted to a shorter time period (e.g., Salix) or differed between study landscapes. Sequences of Brassica were also identified on many observation dates in bumble bee pollen loads (59 out of 75 pollen samples), while a reduced frequency was observed for the collection of Fragaria pollen (43 out of 75 pollen samples). Compared with honey bees, bumble bee colonies collected pollen from more diverse resources ( Figure S6a,b)

| Colony level
The ITS2 sequences of Fragaria were among the five most frequent genera in the pollen samples of all small colonies and of the large honey bee colonies. However, the ITS2 sequences of Fragaria were not recorded within the most frequent genera collected by the large bumble bee colonies (Table 3).
Based on microscopic quantification, strawberry pollen grains amounted on average to 26.30% of the pollen collected by honey bee colonies and 18.58% of the pollen collected by bumble bee colonies, while the collection of OSR pollen was below 8% for both bee species and both colony sizes (Table 4).
We found two models with empirical support explaining the effects of bee species, colony size, OSR availability and strawberry flower cover on the proportion of collected strawberry pollen (Table 1b). In the best fitting model (dAICc = 0), strawberry flower cover, bee species and OSR availability were included. Based on assessment of the relative importance of each explanatory variable, strawberry flower cover was identified as the most important predictor variable (Σw i = 1), followed by OSR availability (Σw i = 0.94) and bee species (Σw i = 0.56; Table S3). Strawberry pollen collection increased with increasing strawberry flower cover (Figure 2a). Increasing OSR availability decreased the proportion of collected strawberry pollen independently of bee species or colony size ( Figure 2b). We found a higher proportion of strawberry pollen in honey bee samples by 69.61% compared to bumble bee samples (Figure 2c). Colony size and interactions were not included in the best fitting model (dAICc < 2). The results are supported by the model-averaged coefficients ( Figure S4b).

| In-field foragers
Honey bee foragers in strawberry fields collected a 1.27 times greater proportion of strawberry pollen compared to B. terrestris foragers (p < .001, Figure 3).

| Relationship between quantitative outcomes of metabarcoding and microscopy
We found positive relationships between the number of ITS2 amplicon reads and microscopic pollen counts which differed for the F I G U R E 1 Effect of bee species (a, Σw i = 1) and oilseed rape (OSR) availability (b, Σw i = 0.45) on pollen richness (i.e., number of plant genera). Pollen richness in pollen loads is higher in colony samples collected from bumble bees than from honey bees. Different letters indicate significant differences obtained from post hoc Tukey tests (significance level of .05). Predicted values and 95% confidence interval (in black) from mixed effect models are shown. Furthermore, pollen richness increases with high OSR availability (b). The regression line was obtained from mixed model estimates (model R1, see Table S2a) and the 95% confidence region is shown (in grey). Note that pollen richness is shown on a log scale in both plots [Colour figure can be viewed at wileyonlinelibrary.com] two pollen types (see Figure 4 for details on intercept and slope), being stronger for strawberry pollen (R 2 = .69) than for OSR (R 2 = .15). The positive intercepts for both plant genera (i.e., when microscopic pollen counts = 0) indicate that ITS2 metabarcoding was able to detect pollen when microscopy failed. When no pollen grains were detected by microscopy, the average number of ITS2 reads for strawberry was 3.2 (95% confidence interval [CI] 1.5-5.9) and for OSR 113.4 (95% CI 62.4-203.4).

| Pollen richness
Honey bee pollen demand was met continuously by crop plants (i.e., OSR and strawberry). However, the majority of analysed pollen in this study (>70% in many samples) was collected from noncrop plants. Bumble bees (Bombus terrestris) exploited a greater richness of plant genera in agricultural landscapes compared to honey bees and hence can contribute to the pollination of a greater variety of plant species. Although both bee species are generalists and able to collect pollen rewards and nectar from many plant species (Amiet & Krebs, 2012), honey bees are known to focus on only a few species (de Vere et al., 2017). This is probably due to their ability to communicate the most profitable resources in the landscape using the waggle dance (Couvillon, 2012), and the fact that individual foragers alternate pollen and nectar resources only to a limited extent during foraging trips (Keller et al., 2005). Additionally, pollen composition differed greatly between bee species. Large amounts of pollen came from TA B L E 3 Five most common pollen resources of small and large honey bee and bumble bee colonies across all landscapes and observations dates

Monotropa 29
Fragaria 26 Cornus 24 Note: The study period lasted from May 6 to June 6, 2016. Frequency is the number of samples in which the ITS2 sequence was detected and is based on presence/absence data of plant genera in pollen samples.

TA B L E 4
Mean percentages (±SE) of strawberry and oilseed rape pollen in pollen loads of small and large honey bee and bumble bee colonies woody and herbaceous plant genera which can be found in, for example, hedgerows and field groves and have been highlighted by other studies as important sources for bee nutrition and colony growth (de Vere et al., 2017;Kämper et al., 2016;Requier et al., 2015).
Interestingly, we found evidence that high landscape-wide OSR availability increased the pollen richness collected by both bee species, supporting the findings of Requier et al. (2015). Bees may have focused on a greater diversity of pollen rather than on quantity (Leonhardt & Blüthgen, 2012). In contrast, increasing local strawberry flower availability appears to reduce the collected pollen richness, presumably because bees focused on the resources next to their colonies.

F I G U R E 2
Effects of (a) strawberry flower cover (Σw i = 1), (b) OSR availability (Σw i = 0.94) and (c) bee species (Σw i = 0.56) on the proportion of strawberry in pollen loads (n = 157). High strawberry flower cover increased the proportion of strawberry pollen independently of bee species and colony size (a). High OSR availability decreased the proportion of strawberry pollen loads of both species and colony sizes (b). The regression lines were obtained from mixed model estimates (model P1, see Table 2b) and the 95% confidence region is shown (in grey). Honey bees collected greater proportions of strawberry pollen compared to bumble bees (c). Different letters indicate significant differences obtained from post hoc Tukey tests (

| Collection of strawberry pollen
Both bee species collected a fairly high amount of strawberry pollen. However, with high increasing OSR availability, they collected less strawberry pollen. As honey bees can communicate the location of most profitable resources using the waggle dance (Couvillon, 2012), a shift to mass-flowering or other high-reward flower patches is likely. However, we found only limited pollen foraging on OSR, which is in accordance with Garbuzov et al. (2015), but in contrast to Danner et al. (2017). Those contrasting results F I G U R E 3 Effect of bee species on the proportion of strawberry pollen in pollen loads from the hind tibia of in-field foragers. Strawberry pollen had a greater share in pollen loads of honey bees (n = 37) than in bumble bees (n = 36  (Rollin et al., 2013). At the same time, however, they are known to favour foraging close to their colonies to reduce energy costs (Lihoreau et al., 2010;Seeley, 1995). At the time of high OSR flowering, we identified ITS2 sequences of several other genera (e.g., Salix, Prunus and Taraxacum) in the pollen loads using amplicon metabarcoding. Hence, bees have collected pollen on coflowering plants that may be more attractive pollen resources than OSR. However, they may have used OSR as a nectar source.
In addition, local strawberry flower availability in the adjacent strawberry field benefitted the strawberry pollen foraging of both bee species. Thus, high strawberry flower cover will facilitate the pollen collection of the respective crop but may result in pollinator dilution in the field.
In contrast to Boecking and Kreipe (2015), we did not find that small colonies collected greater proportions of strawberry pollen in pollen samples. However, based on the frequency of pollen in the samples, strawberry was collected by small colonies more frequently than by large colonies. Based on these descriptive data, we found first indications that the foraging behaviour of small and large bumble bee colonies can differ. In addition, honey bee colonies collected strawberry pollen more frequently than bumble bee colonies. Additional data from in-field honey bee foragers showed a higher flower constancy for honey bee individuals than bumble bee individuals (Rollin et al., 2013). A high flower and pollen constancy is likely to be linked to pollination success and even higher seed set (Gyan & Woodell, 1987;Montgomery, 2009 (Marzinzig et al., 2018), it may depend even more on the functional diversity of traits (Hoehn et al., 2008).
The proportion of fertile achenes per fruit, which can be linked to higher fruit weight (Klatt et al., 2014), will be higher if several bee species, not just the honey bee, visit the flower (Chagnon et al. 1993).

| Relationships between quantitative outcomes of metabarcoding and microscopy
In general, we observed pollen type-specific, Marzinzig et al., 2018) while a more specific assessment of rare pollen species would probably need a pollen grain count of 500 (Lau et al., 2018). As our correlation analysis was conducted with strawberry and OSR pollen, which are major pollen resources in our study landscapes, we are confident that both pollen species are well represented in our samples and that our data sets provide a sound basis for the analysis. In comparison to microscopy, ITS2 metabarcoding is more advantageous in that it achieves a high taxon richness, allows for a higher throughput with a predictable cost and time frame, and does not need specific expert knowledge in palynology (Bell et al., 2019;Keller et al., 2015).
New developments in microscopic pollen detection using deep learning techniques (Gallardo-Caballero et al., 2019) or in full-length amplicon or genome sequencing with, for example, nanopore sequencing techniques (Lang et al. 2019;Leidenfrost et al., 2020;Peel et al., 2019) could improve the weaknesses of both approaches (e.g., time expenditure in microscopy or quantification accuracy in molecular methods). However, studies are needed to compare and evaluate the accuracy of those new developments.

| CON CLUS IONS
We demonstrate here that honey bee and bumble bee (Bombus terrestris) colonies differ substantially in their pollen resource use in agricultural landscapes. Bumble bees collected pollen from a much larger variety of plant genera compared to honey bees. Thus, conservation schemes should consider bees' foraging preferences by taking diverse plant communities into account to promote pollinators and associated pollination services for wild and crop plants. Annual flowering crops and in particular floral resources in permanent landscape elements, such as hedges, are important in fulfilling the foraging requirements of bees. Both honey bee and bumble bee foragers adapted their foraging behaviour to the availability of mass-flowering resources, which could affect the provisioning of pollination services to minor flowering crops. Honey bees carried slightly more strawberry pollen and less diverse pollen loads than bumble bees, but the consequences for pollination services need to be studied in more detail. If bee densities are low, farmers can use managed bee colonies for crop pollination.
However, we would instead recommend designing pollinator-friendly agricultural landscapes that provide species-rich flower resources for wild and managed pollinators, which in turn can provide pollination services to crops and wild plant species.
ITS2 metabarcoding is a suitable method to study the richness of bee-specific pollen diet using mixed pollen samples of unknown plant communities. However, associations between quantitative outcomes of microscopic pollen grain counts and ITS2 amplicon reads were pollen type-specific and weak, and large proportions of variation were not explained. Our results can contribute to ongoing discussions that apply and test different methods to quantify pollen grain counts (Baksay et al., 2020;Pornon et al., 2016). Given the growing interest in both microscopic (Gallardo-Caballero et al., 2019) and molecular (Baksay et al., 2020;Leidenfrost et al., 2020) pollen analyses for pollen identification and quantification, our study highlights that the methods should be chosen carefully and in a targeted manner.

ACK N OWLED G EM ENTS
We thank the contributing strawberry farmers for allowing us to keep our honey bees and bumble bees on their fields. We also thank Susanne Jahn and Brigitte Jünemann for their valuable help