Temporal change in floral availability leads to periods of resource limitation and affects diet specificity in a generalist pollinator

Generalist species are core components of ecological networks and crucial for the maintenance of biodiversity. Generalist species and networks are expected to be more resilient, and therefore understanding the dynamics of specialization and generalization in ecological networks is a key focus in a time of rapid global change. Whilst diet generalization is frequently studied, our understanding of how it changes over time is limited. Here we explore temporal variation in diet specificity in the honeybee (Apis mellifera), using pollen DNA metabarcoding of honey samples, through the foraging season, over two years. We find that, overall, honeybees are generalists that visit a wide range of plants, but there is temporal variation in the degree of specialization. Temporal specialization of honeybee colonies corresponds to periods of resource limitation, identified as a lack of honey stores. Honeybees experience a lack of preferred resources in June when switching from flowering trees in spring to shrubs and herbs in summer. Investigating temporal patterns in specialization can identify periods of resource limitation that may lead to species and network vulnerability. Diet specificity must therefore be explored at different temporal scales in order to fully understand species and network stability in the face of ecological change.

they may experience periods of increased vulnerability as a result.
Understanding the dynamics and structure of ecological networks is therefore a key focus of ecology in a time of rapid global change such as habitat loss, climate change and biological invasions (Burkle & Alarcon, 2011).
Generalization and specialization can be defined and measured in many ways (Devictor et al., 2010).Generalist species are classed as those with a large niche breadth, becoming more specialist as the niche narrows.However, generalization in a species or population may arise through each individual having a similar generalist diet, or through each individual utilizing a small subset of the resources utilized by the whole population (Araújo et al., 2011;Bolnick et al., 2003).Individual specialization is thought to help generalist species switch resources quickly (Szigeti et al., 2019), with species with large variation in diet between individuals thought to be more robust to ecological change (Forsman & Wennersten, 2016).
Establishing the prevalence and temporal consistency of individual specialization therefore helps to further understand how vulnerable generalist species may be in the face of ecological change.
Optimal foraging theory (OFT) can be used to help understand the mechanisms behind specialization.OFT states that the selection of resources by individuals is determined by the overall net energetic gain, which is calculated from the value of the reward against the energy required to locate and extract it (Araújo et al., 2011).It predicts that when preferred resources are abundant, individuals utilize few resources, and when preferred resources become limited, dietary breadth is increased as the diet is supplemented with lower value resources (MacArthur & Pianka, 1966;Svanbäck & Bolnick, 2007).
During periods of low resource availability, individuals must not only increase diet breadth but also become opportunistic as the search time for rewarding resources increases, resulting in increased interindividual variation in diet due to stochastic variation in resource discovery.Individual specialization may allow generalist species to switch resources easily if their preferred resources are not available, but using alternative resources which are less preferable could impact species and result in periods of resource limitation.Temporal patterns in individual specialization can be compared to periods of resource limitation detected through physiological responses (Herring et al., 2011) or food stores (Mattila & Otis, 2006) to assess how they are linked.
Identifying preferred resources requires an estimate of how much or little resources are used compared to their availability.
Generalists are expected to utilize resources in proportion to their abundance (Feinsinger et al., 1981), but resource selection will also be influenced by other traits such as nutritional value (Hicks et al., 2016).If the abundance of resources is quantified, its relationship to relative use may be tested to identify true preferences.If foragers select resources in proportion to their abundance, then individual specialization is driven by preference for the most abundant resources.If foragers select resources disproportionally to their abundance, then other factors may be important.
One of the most well-studied ecological networks are those concerning the interactions between plants and pollinators (Bascompte & Jordano, 2006;Burkle & Alarcon, 2011;Memmott, 1999), with their resilience against ecological change having huge implications for pollination ecosystem processes (Klein et al., 2007).A key component of many plant-pollinator networks globally is the generalist honeybee Apis mellifera (Hung et al., 2018).Honeybees have a highly developed recruitment strategy which enables them to communicate the location and quality of resources to other members of the colony, and as such foraging decisions occur at both an individual and colony level (Seeley, 1995).Honeybees are often referred to as super-organisms as colony-level decisions are stronger than those of the individual (Moritz & Fuchs, 1998), and therefore the foraging decisions of a colony are comparable to those of an individual in other species.Although colonies within the same apiary have been shown to utilize different resources (de Vere et al., 2017), how this degree of intercolony specialization changes through time as a result of fluctuating resources is poorly understood.Identifying how diet specificity fluctuates temporally in honeybees therefore has implications for plant-pollinator networks and can be used to make predictions about generalist species in other ecological networks.
Identifying preferences for resources based on their abundance or other traits may be key to understanding any temporal patterns in intercolony specialization.Given the highly developed recruitment strategy of honeybees (Seeley, 1995), it is expected that there would be a positive association between the abundance of a plant in the landscape and its relative use by honeybees.If plant abundance is not the sole driver of selection, some taxa will be visited more, or less, than expected when considering their availability in the landscape.In addition, trees and shrubs have been found to be an important resource for honeybees (de Vere et al., 2017;Jones et al., 2022), but whether this is simply due to their relative abundance in the landscape compared to herbs is poorly understood.If honeybees are selecting trees over herbs, their relative use will not be in proportion to their availability in the landscape.Furthermore, given that A. mellifera is native to the UK (Carreck, 2008), we may expect it to show a preference for native plants over horticultural plants.If honeybees show a preference for native over horticultural plants, they will proportionally use more native plants than are available in the landscape.
Managed honeybee colonies can require supplemental feeding due to low food stores (Mattila & Otis, 2006) and, consequently, periods of resource limitation have long been described by beekeepers (Seeley, 1995).A prominent food shortage is thought to occur between periods of strong nectar flow in spring and summer, anecdotally known as the "June gap" within the UK (Crane, 1976;Suryanarayana & Singh, 1989).These periods of resource limitation can therefore be used to investigate temporal diet specialization in a generalist pollinator.
Pollen DNA metabarcoding is a powerful tool for characterizing floral visitation by pollinators (reviewed in Bell et al., 2022;Lowe, Jones, Brennan, Creer, & de Vere, 2022).The use of DNA metabarcoding for botanical identification of honey is well established (Hawkins et al., 2015;Milla et al., 2021;Valentini et al., 2010;Wirta et al., 2021), but few have applied this method to explore ecological questions related to honeybee foraging (de Vere et al., 2017;Jones, Brennan, et al., 2021;Lucek et al., 2019).Here, we use DNA metabarcoding to explore temporal patterns in diet specificity in A. mellifera by characterizing honey samples from multiple colonies over two years.Whether generalist species, such as honeybees, show temporal specialization has implications for species and network vulnerability to ecological change.We investigate temporal specialization in a generalist and relate it to periods of resource limitation, specifically asking the following questions: 1. Do honeybees show diet specificity, and does it vary through time?
2. If honeybees show temporal patterns of specificity, are these linked to any periods of resource limitation? 3. Is temporal diet specificity related to floral abundance, or preferences for (i) trees, shrubs or herbs, or (ii) native or non-native plants?

| Study site
The study was conducted at the National Botanic Garden of Wales, UK (51°50′33.4″N, 4°08′49.2″W), a diverse landscape (230 ha) consisting of formal gardens and organic farmland, designated as a National Nature Reserve (Waun Las NNR) (Figure S1).The botanic garden contains over 5000 plant species and varieties from throughout the world, including many horticultural plants grown throughout Western Europe.It is set within an agricultural area in West Wales, UK, with semi-improved grassland comprising the major habitat.Two apiaries are located within the study site, 1 km apart, one within the botanic garden with close access to horticultural plants and the other at the edge of the nature reserve.All colonies had identical management practices, housed within British Standard National hives.Each month colony size was estimated by counting the number of frames of bees in each hive (Table S1).Food availability was determined by recording whether honey stores were present (Table S2).

| Floral surveys
To estimate the availability of floral resources, the outside areas of the botanic garden and nature reserve were surveyed monthly from April to September in 2018 and 2019.The site was split into 287 survey zones, each of 18 ha, and classified into one of three habitat types: broadleaved woodland and hedgerows; horticultural; and grassland (Figure S1).Each of the zones was mapped using qgis version 3.6.1 and r version 4.0.2.For each zone, citizen scientists walked around the perimeter to record all plant species in flower, along with an estimate of the proportion of the survey zone taken up by available flowers of each plant species.Transects were walked through larger survey zones (fields) to capture as much information as possible.An estimate of the total area covered by each plant taxon was determined by multiplying the percentage cover values by the area of the zones it occupied.

| Honey sampling and DNA extraction
During the same time period as floral surveys were undertaken, honey was sampled from six honeybee colonies (three from each apiary) during the last week of each month.Thirty millilitres of honey was sampled from a newly capped comb in each colony using a sterile 50-ml centrifuge tube.Wax was removed using sterile forceps and 10 g of each honey sample was weighed into a new sterile 50-ml centrifuge tube.A modified version of the Qiagen DNeasy 96 Plant Kit was used for DNA extraction.The 10 g of honey was made up to 30 ml with molecular-grade water (Sigma) and placed in a water bath at 65°C for 30 min, shaking each sample at 10 min intervals.
Samples were transferred to 50-ml Nalgene round-bottomed tubes and placed in a high-speed centrifuge (Sorvall RC-5B) at 15,000 rpm for 30 min.The supernatant was discarded, and the pellet resuspended in 400 μl of buffer, made up of 400 μl AP1 from the DNeasy 96 Plant Kit (Qiagen), 80 μl proteinase K (1 mg ml −1 ) (Qiagen) and 1 μl RNase A (Qiagen).Samples were incubated for 60 min at 65°C and then disrupted in a TissueLyser II (Qiagen) at 30 Hz for 4 min.The remaining steps were carried out according to the manufacturer's protocol, excluding the use of the QIAshredder and the second wash stage.The OneStep PCR Inhibitor Removal Kit (Zymo Research) was used to purify the DNA extract, which was then diluted 1 in 10, prior to amplification.
Primers were appended with universal tails to attach custom indices in the second-round PCR.The first PCR used a final volume of 20 μl: 2 μl template DNA, 10 μl of 2× Phusion Hot Start II High-Fidelity Mastermix (New England Biolabs UK), 0.4 μl (2.5 μm) forward and reverse primers, and 7.2 μl of PCR-grade water.A negative control was included within each PCR to test for cross-contamination and reagent contamination.Thermal cycling conditions for the first rbcL PCR were: 98°C for 30 s, 95°C for 2 min; 95°C for 30 s, 50°C for 30 s, 72°C for 40 s (40 cycles); and 72°C for 5 min, 30°C for 10 s.Thermal cycling conditions for the first ITS2 PCR were as follows: 98°C for 30 s, 95°C for 10 min; 95°C for 30 s, 56°C for 30 s, 72°C for 1 min (40 cycles); and 72°C for 10 min.
Agarose gel electrophoresis was used to visualize PCR products and to assess whether amplification was successful.The PCR was carried out three times, and the products of each PCR were pooled.The pooled products were purified using Illumina's 16 S Metagenomic Sequencing Library preparation protocol using a 1:0.6 ratio of product to Agencourt AMPure XP beads (Beckman Coulter).
A second PCR was carried out on the purified product to anneal custom unique and matched i5 and i7 indices to each sample (Ultramer, Integrated DNA Technologies).The second PCR used a final volume of 25 μl consisting of: 5 μl of pooled purified first-round PCR product, 12.5 μl of 2× Phusion Hot Start II High-Fidelity Mastermix (New England Biolabs UK), 1 μl of i5 and i7 Index Primer, and 6.5 μl of PCR-grade water.Thermal cycling conditions for the index PCR were as follows: 98°C for 30 s; 95°C for 30 s, 55°C for 30 s, 72°C for 30 s (eight cycles); and 72°C for 5 min, 4°C for 10 min.
The index PCR product was compared to the cleaned-up product from the first PCR on a 1% agarose gel to confirm amplification of index tags.A second clean-up stage was followed according to the Illumina protocol with a 1:0.8 ratio of product to beads.Products were quantified using a Qubit 4.0 fluorescence spectrophotometer (Thermo Fisher Scientific) and pooled at equal concentrations to create the final library for sequencing.The negative PCR controls from each plate were sequenced along with the products.Sequencing took place on an Illumina MiSeq using a 2 × 300-bp kit.

| Bioinformatic analysis
The Illumina sequence reads were processed using a custom bioinformatic pipeline (Ford & Jones, 2020).Initially, raw sequences were trimmed to remove low-quality regions, paired and merged.Only sequences >450 bp (rbcL) and 350 bp (ITS2) were used in downstream analysis.Identical reads were dereplicated within each sample and clustered at 100% similarity across all samples with singletons (sequence reads occurring once across all samples) removed.In the UK, the Barcode Wales and Barcode UK projects provide 98% coverage of all native flowering plants and conifers using three plant DNA barcode markers, rbcL, matK and ITS2, allowing reliable identification at the species and genus level (de Vere et al., 2012;Jones, Twyford, et al., 2021).Sequences were compared to a custom, curated reference library, containing 5887 plant species, for identification (Jones, 2020;Jones, Brennan, et al., 2021).The reference library represented UK native species (Stace, 2019), naturalized and alien species (Preston et al., 2002), and all horticultural species from the IRIS BG database at the National Botanic Garden of Wales.Specieslevel coverage was 57% for rbcL and 52% for ITS2, and genus-level coverage was 96% for rbcL and 84% for ITS2 (Jones, 2020).

| Assigning plant taxa
Sequences were compared against the reference library using blastn, recording the top 20 BLAST hits and grouping together sequences with identical BLAST results across all 20 hits.The taxonomic identifications of these grouped sequences were then automatically assigned based on the highest bit score.If the top bitscore belonged to a species, the sequence was assigned to that species.If the top bitscore belonged to different species within the same genus, a genus designation was made for that sequence.If the top bitscore matched multiple genera of the same family, then a family designation was made.Sequences returning top bitscores of multiple families within different orders were removed, assuming that these were poor-quality sequences.The identification of the sequences was then manually checked to ensure botanical veracity, relating to the plant's presence within the landscape, whilst considering the discrimination ability of each marker (Jones, Twyford, et al., 2021).
Following the assignment of identifications, a consensus was reached to combine the taxa identified by rbcL and ITS2 at differing taxonomic resolution.A rule-based, objective and conservative approach was used which assigned taxa to the higher taxonomic level (Methods S1).For example, taxa identified to the genus level with one marker and species using the other were assigned to genus for the consensus.The numbers of reads for each consensus taxon within a sample were then summed to combine the results from both markers.The proportion of reads per taxon per sample was used as a measure of relative read abundance.The relationship between the relative read abundance of each matched taxon, identified using both markers within a sample, was tested using Spearman's rank correlation with Holm correction for multiple testing.
Plants identified to species or genus were assigned a native status and form, while those identified to family were not categorized.Native status was assigned as either "native and near native," "naturalized" or "horticultural" based on Stace (2019).The category "native and near native" comprised native species and also taxa identified to genus level that include native species and horticultural varieties which are functionally similar, such as Prunus spp.Naturalized plants were those which have been introduced and become widespread and self-perpetuating in the wild.All remaining non-native plants were classified as horticultural.Taxa were grouped into three form categories: "tree," "shrub" or "herb."Herbs were defined as non-woody species, shrubs were defined as woody species <5 m tall and trees were defined as woody species >5 m, following the Royal Horticultural Society classification (Bricknell, 2010).The plant taxa found were designated with four categories of abundance, according to the total proportion of sequences contributed each month: major (≥10% sequences), secondary (≥1% and <10%), minor (≥0.01%and <1%) and occasional (>0% and <0.01%).

| Statistical analyses
The data obtained through DNA metabarcoding are considered semiquantitative, owing to the stochasticity of detecting rare taxa and potential biases in sampling, along with extraction, amplification and sequencing of DNA (Bell et al., 2022;Lowe, Jones, Brennan, Creer, & de Vere, 2022).Relative read abundance was used for all analyses (Deagle et al., 2019), either using the proportion of taxa as a percentage or, for the models, the number of sequences, controlling for sequencing depth by setting the total number of sequences per sample as an offset (comparable to proportion) (Jones, Brennan, et al., 2021;Lowe, Jones, Brennan, Creer, & de Vere, 2022).To analyse how the floral composition of honey changes through time, the function manyglm within the r package mvabund (Wang et al., 2012) was used to run a generalist linear model, using all of the plant taxa present in each sample.The data best fit a negative binomial distribution due to the strong meanvariance relationship within the data and the high number of zeros in the abundance data occurring from plant taxa absent in a sample (Figure S2).The number of sequence reads for each plant taxon was set as the multivariate response, with the effect of month (measured as the number in the calendar), year and location of hives included as predictor variables.The number of reads per sample was included as an offset to control for differences in the number of sequences between samples (Deagle et al., 2019;Jones, Brennan, et al., 2021).To visualize temporal changes in the composition of the honey, non-metric multidimensional scaling (NMDS) ordination was used based on the proportion of reads returned for each plant taxon.Ordinations were carried out using the metaMDS function in the vegan package in r (Dixon, 2003), using Bray-Curtis dissimilarity indices.Post hoc pairwise comparisons were run using the anova.manyglmfunction (Wang et al., 2012) to identify which foraging seasons (spring, summer, autumn) significantly differed from each other with respect to taxonomic composition.
To visualize interactions between each colony and plant taxon, bipartite networks were constructed with the proportion of sequences in a sample as a measure of relative abundance (Deagle et al., 2019).We used network metrics to quantify the changes in colony-plant interactions through the year.Individual specialization (IS) was measured as the mean proportional similarity (PS i ) between one colony's diet compared to the combined diet of all colonies during the same sampling period (Feinsinger et al., 1981).
IS was calculated using the "PsiCalc" function in the RinSp package, with a Monte Carlo resampling simulation (2000 iterations) to test its significance (Zaccarelli et al., 2013).The value of IS = 1 when the average diet of each colony is directly proportional to the diet of all colonies (meaning the diets are more similar to each other) and decreases toward 0 when each colony has a distinct diet (Bolnick et al., 2002).We calculated the proportional generality (G), a quantitative measure of the mean number of effective plant taxa per colony, as an additional measure of diet breadth, using the "bipartite" package (Dormann & Strauss, 2014).A higher value of G indicates a broader diet breadth whilst a value at or near zero would indicate few resources are used (Blüthgen et al., 2006;Cusser et al., 2019).All graphical and statistical analyses were carried out in r (version 4.0.2).

| Comparison with floral surveys
All plant taxa identified using DNA metabarcoding were matched to those recorded in floral surveys, for comparison of relative abundance in the honey and landscape.Plants which were identified to species using DNA were matched to the same species in the landscape.To account for the differing taxonomic levels at which plants were identified, those identified at the genus level in the honey were matched to all species belonging to that genus in the floral surveys.
Those identified as one or more genera were also matched by grouping each matched genus in the floral survey.Plants which were identified to the tribe or family level in honey samples were omitted from comparison with floral surveys.
The relationship between the proportion of sequence reads in each honey sample and the proportion of area flowering through the season (for plant taxa contributing over 1% of sequence reads per month) was assessed using Spearman's rank correlation with Holm correction for multiple testing.To test whether plant taxa were used more or less than expected by chance given their relative abundance in the landscape, null models were created using the package econullnetr for preference analysis (Vaughan et al., 2018).Included in these models were plants which either contributed more than 1% of sequences or over 1% of total flowering area in any given month.
Changes in the use of plant type (relating to form and status) over time were investigated using two multivariate generalist linear models, both with season and year included as predictor variables (Figures S3 and S4).The response variable was the number of reads returned (categorized by plant form for one model and categorized by plant status for the other).The total number of sequences per sample was retained as an offset in both models.The abundance and richness of plants in the honey (characterized by native status and form) were compared to the floral survey data each month using Fisher's exact tests.These were (i) comparisons between the relative read abundance and the proportion flowering in the landscape each month and (ii) comparisons between the generic richness of honey samples and the generic richness of the landscape each month.

| RE SULTS
DNA metabarcoding of 54 honey samples returned 6,984,378 sequences after stringent quality control (3,996,872 rbcL and 2,987,506 ITS2).The mean number of sequences per sample was 129,340 (SD = 24,488), ranging from 57,495 to 167,861.The rbcL marker detected 99 taxa at the following taxonomic ranks: three families, one tribe, 75 genera and 20 species.ITS2 identified 120 taxa consisting of 93 genera and 27 species.Following the matching of taxa identified by both markers at varying discrimination, there were 143 plant taxa using the rbcL and ITS2 regions combined with 18% of taxa identified to species, 79% to genus, 2% to family and 1% to tribe.Sequence reads for both markers separately and combined can be found in the Supporting Data.There was a strong correlation between the abundance of sequences found within each sample using both rbcL and ITS2 (n = 59) (Spearman correlation coefficient r s = .604,p < .001)(Figure S5).

| Do honeybees show diet specificity, and does it vary through time?
Both month (LR 1,52 = 799.3,p < .001)and year (LR 1,51 = 404.0,p < .001)were found to be good predictors of plant composition within the honey, irrespective of apiary location (LR 1,50 = 412.0,p = .069).The most abundantly used plants were the same in 2018 and 2019 but 25 taxa were unique to 2018 and 36 to 2019 (Table S4).
NMDS showed that samples collected in the same month were most similar to each other, with samples collected later in the year becoming increasingly divergent (Figure 1).Post hoc pairwise comparison revealed that plant composition in the honey differed between all pairs of seasons (spring vs. summer: sum-of-LR = 571.9,p < .001;Honeybees were found to use 143 plant taxa, covering 66 families and 139 genera.However, despite the diversity of taxa used, only 15 taxa were found to contribute more than 1% of total sequence reads across all samples (Table 1).Included in this list are plants used abundantly by honeybees in a particular season and those with long flowering periods, which are used through the year.
Over both years, 10 plant taxa were classified as major (≥10% of reads in at least one month), 30 as secondary taxa (≥1% and <10% of reads in at least one month), 66 as minor taxa (≥0.01%and <1% of reads in at least one month) and 37 as only ever occurring occasionally (>0% and <0.01% of reads in at least one month) (Table S4).
Over the study period, honeybees had access to a total of 1497 plant taxa covering 613 genera, distributed across native habitats consisting of semi-improved grassland, woodland and hedgerows, along with planted horticultural areas and amenity grassland (Lowe, Jones, Brennan, Creer, Christie, & de Vere, 2022).There was an average of 260 genera in flower each month (SD = 55.25), but honeybees only used a minority of available genera (Figure 3), with each colony choosing between 2% and 20% of available genera per month (Table S5).
Patterns in diet specificity were visualized through the structure of seasonal bipartite networks (Figure 2a,b).Interindividual or, as in this case, intercolony specialization (IS) measures the degree of similarity between the diet of a colony and the total diet of all colonies during the same sampling period and equals one when all colonies have the same diet, decreasing toward zero when each colony has its own unique diet.During spring (April and May) IS values were greater than 0.6 (p < .001),which means that each colony's diet was strongly similar to the total dietary niche of all colonies sampled during the same time period (Table 2).This is a result of all colonies abundantly using flowering tree species such as cherries and plums (Prunus spp.), along with sycamore and maples (Acer spp.) and herbs such as dandelion, Taraxacum officinale, and Pulmonaria spp.(Figure 2a).During this same time period, the average number of effective plant taxa per colony (G) was between 5 and 6, except for in May 2018, where diet breadth was higher, and G was 9.85 (Table 2).
In June, the values of IS decreased to 0.42 (p < .001) in 2018 and to 0.56 (p < .001) in 2019, demonstrating that the diets of colonies diverge and, on average, the diet of each colony in comparison to the diet of all colonies together became more distinct than in previous months (Table 2).Whilst spring-flowering taxa such as Acer spp.
and Taraxacum officinale are retained in the diet of all colonies, the introduction of Rubus fruticosus into the diet of some colonies marks

1.17
Note: Per month, taxa contributing to ≥10% of sequences are categorized as major plant taxa, <10% to ≥1% are secondary taxa, <1% to ≥0.01% are minor taxa, and <0.01%are occasional.Those noted as >0.00 were present but at very low values.An extended list of all taxa found in honey samples can be found in Table S4.
a shift in resource use (Figure 2a).2).There was higher intercolony variation in July 2018 compared to July 2019 due to colonies using Cirsium/Centaurea/Hypochaeris spp.and Coreopsis spp.
Intercolony foraging differences increased in August, illustrated by lower IS values in comparison to July (Table 2), as colonies retained the use of resources such as Cirsium/Centaurea/Hypochaeris spp., Filipendula ulmaria, R. fruticosus and Trifolium repens, whilst a subset of colonies began using Impatiens glandulifera (Figure 2b).Diet breadth measured by generality (G) remained the same between July and August in 2018 (G = 4.2), but increased slightly between July (G = 2.5) and August 2019 (G = 3.6).
Honey was not collected in September 2019 due to poor weather conditions.However, during September 2018, intercolony variation remained similar to that in August (IS = 0.76, p < .001),but diet breadth (G) reduced to 2.7, the lowest value seen throughout the season, with the majority of colonies using I. glandulifera most abundantly (Figure 2b).In summary, we see that honeybees  Note: Intercolony specialization (IS) = 1 when the average diet of a colony is directly proportional to the diet of all colonies during the same sampling period and nears 0 when each colony has its own unique diet.The metric G, generality, is a quantitative measure of the average number of effective plant taxa per colony per month.Honey was not sampled in September 2019 due to poor weather.
at some points through the year, whilst at other times resources differ.

| If honeybees show temporal patterns of specificity, are these linked to any periods of resource limitation?
Over the foraging season, honeybees experienced periods of food shortage, represented by a lack of honey stores in the hive (Table S2).
Shortages occurred early in the season in 2018 (April and May-three hives), due to colonies building up after winter.We then see periods  (Figure S6).Preference analysis showed that of 75 plant taxa, which either contributed more than 1% of sequences or over 1% of total flowering area in any given month, 31 were used more than expected given their abundance in the landscape, including the two key plants R. fruticosus and I. glandulifera (Figure S7).
Twenty-five taxa were used less than expected given their abundance in the landscape, including Cirsium/Centaurea/Hypochaeris spp.and Ranunculus/Ficaria spp.The remaining 19 taxa were used as expected (Figure S7).
Honeybee forage preferences for trees, shrubs and herbs varied by month (LR 1,48 = 74.79,p < .001)and year (LR 1,47 = 11.33,p < .007)(Figure 4).In April and May, the majority of honeybees' diet was composed of trees and herbs.Throughout June and July, the use of trees decreased, and honeybees replaced these with shrubs.The switch from trees to shrubs occurred earlier in the season in 2018 than in 2019, with honeybees using mostly shrubs in June 2018 and mostly trees in June 2019.In August and September, the majority of honeybees' diet was composed of herbs and shrubs.Each month, the relative DNA sequence read abundance of each plant form identified in honey samples was significantly different to the relative abundance flowering within the landscape, with honeybees using trees and shrubs more than expected throughout the year (Figure 4).
Although herbs were used less than expected based on their relative abundance in the landscape, these make up a large proportion of honeybees' diet through the year (Figure 4).A similar trend was found when comparing the generic richness of each plant form found in honey compared to its availability in the landscape each month (Figure S8).
Overall, throughout the year, the most abundantly used plant taxa were native and near-native (Figures S9 and S10).However, the relative use of plants in each status category did vary by month (LR 1,48 = 71.24,p < .001)and year (LR 1,47 = 8.21, p = 0.041), with increased use of the naturalized I. glandulifera causing a reduction in the use of native and near-native plants in August (2018 and2019) and September (2018) (Figure S9).There was no significant difference found between the relative abundance of each plant status (native and near-native, horticultural, naturalized) in honey samples and their relative abundance in the landscape in any month, except during August and September 2018 (Figure S9).However, when comparing generic richness of each status category in honey samples to generic richness availability in the landscape each month, honeybees used more native and near-native genera than expected by chance in all months except April 2019 (Figure S10).

| DISCUSS ION
In this study, we identified that a generalist forager can exhibit temporal diet specialization through the year.Periods of specialization were found to correspond to periods of resource shortage, as predicted by OFT.We identified preferences for particular resources which drive temporal specialization and could make generalists more vulnerable to periods of resource limitation.We highlight the importance of considering temporal differences when assessing diet specificity and the necessity of understanding the mechanisms behind periods of resource limitation, in order to ensure species and their networks are robust in a time of rapid global change.
The floral resources selected by honeybees covered a diverse taxonomic breadth, with 66 plant families represented in honey samples, demonstrating diet generalization.Despite this generalization, we found that honeybees showed temporal variation in the degree of specialization, considering both intercolony dietary variation and generality as a measure of diet breadth.The use of a wide range of plants is expected in social bees with long foraging periods, as it allows species to adapt and survive in response to changing availability of resources through the year (Ogilvie & Forrest, 2017).
The prevalence of this short-term specialization and long-term generalization allows honeybees to switch resources quickly and easily in response to changes in floral abundance, making them more robust than long-term specialists (Szigeti et al., 2019).Changes in plant use in response to varying plant abundance have also been recorded in bumblebees, hoverflies and butterflies, widening their forage choice as the density of their preferred resource decreases (Goulson, 1999). | Although there was a high diversity of plants available through the year, we found periods of resource limitation for honeybees, demonstrated by a lack of honey stores.Periods of resource limitation were most prevalent in June, supporting anecdotal evidence from beekeepers that there is a "June gap" in resource availability (Crane, 1976;Suryanarayana & Singh, 1989).The nectar collected by honeybees is processed and stored as honey to be used when foraging is not possible, for example during winter or periods of cold, wet weather.As a result, honey stores increase during periods of high productivity.Under OFT, individuals increase dietary breadth when the availability of preferred resources is low, which simultaneously results in an increase in diet variation (MacArthur & Pianka, 1966).Seasonal food shortages are well documented in agricultural habitats across Europe and North America, following mass flowering of insect-pollinated crops (Couvillon et al., 2014;Jachuła et al., 2021;Timberlake et al., 2019).We reveal a shift in major resource use between spring and summer, with periods of re-  (Seeley, 1995), but the relationship between patch size and recruitment in social species is not fully understood (Goulson, 1999).Flowering trees and hedgerow species are large three-dimensional flowering patches in the landscape, which may provide a denser resource compared to two-dimensional flowering patches such as grasslands (Donkersley, 2019).Whilst honeybees demonstrated a high use of native and near-native plants, this was not more than expected given their relative abundance in the landscape.We did, however, find that honeybees used a greater number of native plants than expected given their generic richness in the landscape.Therefore, we show that whether honeybees prefer native over non-native plants is dependent on the measure used.
The "June gap" in floral resources (Crane, 1976;Suryanarayana & Singh, 1989) has been suggested to be a result of reduced floral availability between preferred resources of spring-flowering trees and summer-flowering shrubs and herbs (Balfour et al., 2018).By relating temporal patterns in intercolony specialization to periods of resource shortage and preferences for trees and shrubs, our work supports this theory.Whilst we found that honeybees are generalist and can switch between resources in response to phenological changes in floral availability, the lack of preferred, abundant resources in June is illustrated by increased intercolony variation in resource use.This period also corresponds to the point where the dominant plant form used by honeybees changes from trees (April and May) to shrubs (July).If honeybees were able to switch between these resources successfully without experiencing resource limitation, colonies would remain productive throughout the floral season.This suggests that whilst there is an abundance of floral resources available through the year, honeybees require successive periods of preferred, abundant resources which are not always provided in this landscape.These periods of resource limitation may also impact wild pollinators if preferred resources are shared with honeybees, by increasing exploitative competition between species (Wignall et al., 2020).Further, it is likely that any gaps in floral resource availability will have a greater impact on species that lack the ability to store resources and those with short flight seasons such as solitary bees (Ogilvie & Forrest, 2017;Timberlake et al., 2019).
A key area of further work is to determine both the quantity and the quality of nectar and pollen resources, through the season, to further characterize periods of resource dearth in both agricultural and horticultural landscapes (Timberlake et al., 2019).Whilst foraging in social insects is influenced by the abundance of resources, it is important to note that communication between individuals, environmental factors, and phenotypic or genotypic differences in individuals and colonies may also play a role in resource selection (Frank & Linsenmair, 2017).
Here, we show that a generalist species can show temporal specialization, which has implications for the vulnerability of species and networks in the face of ecological change.We highlight the importance of considering temporal patterns when investigating diet specificity, as network stability may be over-estimated if temporal variation is not considered.Our identification of periods of resource shortage in a generalist demonstrates that whilst individual specialization may allow species to switch resources easily in a changing environment, there is still a requirement for a phenological overlap in the availability of preferred resources to adequately support species.Therefore, assuming generalist species do not show specificity and that all resources are assumed to be equally rewarding could result in periods of resource limitation being overlooked and vulnerability in the face of ecological change being under-estimated.
summer vs. autumn: sum-of-LR = 231.3,p = .005;autumn vs. spring: sum-of-LR = 729.2,p < .001),indicating a transition in honeybee forage choice between these periods.The plants used by each colony varied through the year, resulting in phenological shifts discernible at the colony and network level (Figure 2a,b).F I G U R E 1 Non-metric multidimensional scaling (NMDS) plots showing plant composition of honey samples in relation to month of collection for 2018 and 2019, with each sample labelled with the plant taxon that was most abundant in terms of proportion of sequence reads

F
Plant composition of honey samples from up to six honeybee colonies (H1-H6), collected in 2018 and 2019.The size of each hive bar (left) is proportional to the number sampled and the size of each plant bar (right) is the relative read abundance of DNA sequences recovered.The width of the connecting ribbon is based on the relative read abundance of plant sequences.Plant taxa which contributed >5% of total sequences in any month are coloured and labelled (see key in panel b).Details of all the interactions, including those not labelled, are available in the Supporting Data.(a) Samples collected from April to June in 2018 and 2019.Missing samples were not collected due to limited stores b) Samples collected from July to September in 2018 and July to August 2019.Samples were not collected in September 2019 due to poor weather.TA B L E 1 List of plant taxa found in honey samples using rbcL and ITS2 DNA barcode markers which are represented by >1% of total sequence reads across all samples Diet breadth, measured by generality (G), remained comparable to previous months, with values of 6.0 and 7.8 for 2018 and 2019 respectively.By July, R. fruticosus became the dominant taxon used by most colonies (Figure 2b), quantified by low diet breadth in comparison to previous months (2018: G = 4.2, 2019: G = 2.5).This resulted in a reduction in intercolony variation, represented by an increase in IS values between June and July (Table are generalists, using a wide range of plants, although only a small number of these are used abundantly.Honeybees show temporal patterns of diet specificity with colonies sharing similar resources F I G U R E 3 (a) Total number of genera in flower in floral surveys from April to September.(b) Proportion of genera found in DNA compared to those in flower.No honey samples were collected in September 2019 due to poor weather TA B L E 2 Network metrics quantifying temporal changes in foraging behaviour

3. 3 |
Is temporal specificity related to floral abundance, or preferences for (i) trees, shrubs or herbs, or (ii) native or non-native plants?When considering plant taxa that contributed over 1% of sequence reads in at least one month (major and secondary taxa, n = 36), no significant relationship was found between the proportion of DNA sequences in the honey and flowering abundance within the study site in any month, with the exception of July 2019 (Spearman correlation coefficient r s = .68,p = .002)and August 2019 (Spearman correlation coefficient r s = .59,p = .04) source limitation corresponding with increases in intercolony variation in diet as predicted by OFT.Our findings are supported byRequier et al. (2015) who found that the diets of honeybee colonies in France were most dissimilar in June, during a period of low food supply between peak flowering of two crop blooms.We therefore demonstrate that individual colony specialization within honeybee diets allows utilization of a wide range of resources, but the lack of preferred, abundant resources drives periods of resource limitation with low productivity.In order to further explore which resources were preferred by honeybees, and consequently driving periods of resource limitation through changes in their availability, we assessed the relative use of plant taxa in relation to their abundance.We found a lack of relationship between the use of floral resources and their relative abundance in the landscape, indicating the influence of other drivers on floral selection.When assessing preferences for types of plants, honeybees were found to use trees more than expected given their abundance in spring and more shrubs than expected in summer.Honeybees use patch-based recruitment strategies to discover floral resources

F
Relative abundance of each plant form in honey (DNA) compared to the relative abundance in flower (floral surveys) each month for 2018 and 2019.Using Fisher's exact test, a significant difference was found each month across both years.2018: April (p < .001),May (p < .001),June (p < .001),July (p < .001),August (p < .001),September (p = .018).2019: April (p < .001),May (p < .001),June (p < .001),July (p < .001),August (p < .001) food shortage during June (three hives in 2018, one in 2019) and in August (one hive in 2018).The periods of food shortage in June corresponded to periods where individual colony specialization was most prevalent (lowest IS values) across both years (Table2).The food shortage in August 2018 corresponded to another increase in individual colony specialization from July, but this difference is less pronounced than seen between May and June in both years.