Complex urban environments provide Apis mellifera with a richer plant forage than suburban and more rural landscapes

Abstract Growth in the global development of cities, and increasing public interest in beekeeping, has led to increase in the numbers of urban apiaries. Towns and cities can provide an excellent diet for managed bees, with a diverse range of nectar and pollen available throughout a long flowering season, and are often more ecologically diverse than the surrounding rural environments. Accessible urban honeybee hives are a valuable research resource to gain insights into the diet and ecology of wild pollinators in urban settings. We used DNA metabarcoding of the rbcL and ITS2 gene regions to characterize the pollen community in Apis mellifera honey, inferring the floral diet, from 14 hives across an urban gradient around Greater Manchester, UK. We found that the proportion of urban land around a hive is significantly associated with an increase in the diversity of plants foraged and that invasive and non‐native plants appear to play a critical role in the sustenance of urban bees, alongside native plant species. The proportion of improved grassland, typical of suburban lawns and livestock farms, is significantly associated with decreases in the diversity of plant pollen found in honey samples. These findings are relevant to urban landscape developers motivated to encourage biodiversity and bee persistence, in line with global bio‐food security agendas.

widespread pesticide use are among the most powerful drivers of pollinator decline globally (Gill et al., 2012;Goulson et al., 2015;Jachuła et al., 2022). In the case of honey bees, rates of infection and colony collapse have been attributed to several interacting factors, of which loss of plant forage diversity and abundance has been identified to be an important cause (Branchiccela et al., 2019;Requier et al., 2018;Smith et al., 2013;Vanbergen & The Insect Pollinators Initiative, 2013).
The main sources of nutrition for honey bees are the floral rewards including nectar, which is the primary source of carbohydrate, and pollen, which is the main source of protein, both of which are collected by worker bees (Brodschneider & Crailsheim, 2010).
Nutritional requirements of individuals vary by role in the colony, as foragers and nurse bees require different nutrition, and overall foraging intensity is modulated at the colony level (Altaye et al., 2010;Seeley, 1986). Efficient colony maintenance and brood rearing require not only a sufficient quantity of pollen and associated nutrients, but there are also notable benefits associated with a diverse pollen diet. Colonies reliant upon a single monofloral crop, such as those often found in agricultural habitats, experience a brief glut of pollen. However, they may struggle for sufficient nutrition at other times and are particularly susceptible to a failed crop or inclement weather (Dolezal et al., 2019;Topitzhofer et al., 2019;Vanbergen & The Insect Pollinators Initiative, 2013). A diverse diet of plants with flowering times spread throughout the season offers security against these risks, allowing increased temporal stability of nutrient availability (Avni et al., 2014). It also has a direct nutritional benefit compared to a mono-floral diet, with a diverse diet better able to meet differential nutritional requirements of the different roles within the colony (Paoli et al., 2014). In addition, a poly-floral diet can increase the immunocompetence of bees and is indirectly associated with an increased number of female offspring and reduced disease and pesticide susceptibility of the colony (Alaux et al., 2010;Centrella et al., 2020). While diet diversity is important for health, it is not the only factor, and total pollen and nectar are also critical.
A sufficient biomass of collected pollen is clearly fundamentally important for growth and development, particularly in the spring when stored pollen stocks are low and foraging levels can be inconsistent (DeGrandi-Hoffman et al., 2016).
Urban intensification and expansion may provide a relatively novel opportunity for wild bees and other pollinators (Ayers & Rehan, 2021;Hall et al., 2017;Turo & Gardiner, 2019;Wenzel et al., 2020). Further increases in the development of urban environments are projected to continue (Chen et al., 2020;Gao & O'Neill, 2020;Seto et al., 2012), and understanding pollinator ecology and behavior in response to changing habitat is necessary.
There has been a recent increase in the number of urban beekeepers, with many utilizing gardens and rooftops for their hives (Lorenz & Stark, 2015). Reduced colony mortality, fewer parasitic invasions, and increased colony longevity and reproductive output are all characteristics reported for bee colonies in more urban, compared to rural, environments (Baldock et al., 2019). These benefits are largely attributed to the availability of floral resources and the lower concentration of pesticides (Botías et al., 2017;Samuelson et al., 2018). Agricultural landscapes are often capable of providing comparable or larger quantities of pollen than urban areas, leading to high food accumulation in a hive (Sponsler & Johnson, 2015).
However, the diversity of diet is likely to be lower than that collected when foragers are able to access urban landscapes. Urban areas are known to support diverse populations of wild native bees (Baldock et al., 2015;Casanelles-Abella et al., 2022) alongside honey bees, which like many other bee species, are generalist foragers able to take advantage of the floral diversity available in urban areas.
Studying accessible, managed honey bee colonies as models for their wild counterparts is therefore a powerful tool to better understand urban pollinator ecology (Giannini et al., 2015;Lowenstein et al., 2019).
Urban environments, in general, can be considered to be rich in plant diversity, including a mixture of native, those widespread and not introduced by human activity, neophyte, those not native but in the wild and naturalized, and non-native species, which may include garden plants, recent non-naturalized escapees, and contemporarily invasive species (Aronson et al., 2017;Baldock et al., 2015;Gaertner et al., 2017;Grimm et al., 2008). This richness is due to the spatial heterogeneity of the areas, which provide niches for opportunistic seedlings, and also the presence of a broad range of cultivated plants in private gardens, parks, allotments, urban food production, and across green infrastructure (Frankie et al., 2005;Garbuzov & Ratnieks, 2014;Knapp et al., 2012;Matteson & Langellotto, 2009).
The degree of urbanization and habitat fragmentation can greatly alter the availability and diversity of floral resources for pollinators (Levé et al., 2019;McKinney, 2002). Even generalist species, such as honey bees, exhibit selectivity in the plant species visited depending upon the needs of the colony at specific times and the availability of resources (Hawkins et al., 2015;Lowenstein et al., 2019;Nottebrock et al., 2017;Requier et al., 2015;Ruedenauer et al., 2020;Salisbury et al., 2015;Vaudo et al., 2015). Foraging distances vary depending on the level of landscape complexity surrounding hives and have been reported to often be shorter in complex and urban or suburban landscapes when foraging for pollen, but the pattern does not continue when foraging for nectar (Garbuzov et al., 2015;Steffan-Dewenter & Kuhn, 2003). Longer foraging flights of over 9.5 km are known to occur, with foraging strategy theorized to be linked to patch size and quality (Beekman & Ratnieks, 2000).
Plant taxa contributing to the forage of a hive are generally characterized by identification of pollen sourced from hive pollen traps, isolated from honey, or through physically tracking foraging bees (Carvell et al., 2007;Dimou et al., 2006;Valentini et al., 2010).
Metabarcoding of DNA; species identification through the analysis of complex, mixed community DNA (Deiner et al., 2017;Hebert et al., 2003;Statnikov et al., 2013), has benefits over methods based on morphology and was popularized for bee forage analysis as DNA barcoding became more prevalent in plants (Dunning & Savolainen, 2010;Kress et al., 2005;Newmaster et al., 2006). While plant species identified from pollen loads give a direct measure of the plants visited by bees when collecting pollen, information from honey-extracted plant DNA can be used to describe plants visited for both pollen and nectar collection over a longer period Hawkins et al., 2015;Louveaux et al., 1978). Some foraging is known to target pollen only and may therefore be missed when honey-based sampling is used (Synge, 1947). A number of gene regions (e.g. rbcL, trnL and ITS2) have been identified for use as metabarcodes in plants and a multi-gene region metabarcoding approach has been recommended to increase the discriminatory power and broaden the range of species detection, as specific gene regions show biases in detection range and level of plant taxon identification Burgess et al., 2011;Hollingsworth et al., 2011;Kress et al., 2009).
In this study, we used DNA metabarcoding of honey samples to determine the diet diversity of honey bees in and around Greater Manchester, UK, and asked how the surrounding landscape composition might influence the diversity and composition of the community of plants visited.
Our aims were to use the data from DNA metabarcoding along with GIS analysis of land cover to address three questions: 1. Does land cover diversity surrounding an apiary predict honey bee diet diversity? 2. Does the proportion of a specific type of land cover surrounding an apiary predict honey bee diet diversity?
3. Does the proportion of urban land cover surrounding an apiary predict the proportions of native, non-native, and neophyte plants in the honey bee diet? 2 | MATERIAL S AND ME THODS

| Sample collection
Honey samples were sourced from 14 Apis mellifera (European honey bee) hives; 12 from Greater Manchester, one from Warrington in Cheshire, and one from Rossendale in Lancashire, across an urban gradient ( Figure 1). Extraction of honey varied by the apiary, but the majority of samples were processed with a standard honey extraction method whereby cells are uncapped, the honey is removed from cells by centrifugation in a tangential extractor, and finally filtered to remove large particulates in the honey. Generally, samples were taken from pools of extracted honey collected from multiple frames, but in some instances, comb, chunk, or unfiltered honey were sampled (Table 1). Comb-honey is honey that has not been removed from the cells, chunk-honey is a blend of extracted honey and combhoney, and unfiltered honey follows the traditional extraction methodology with the omission of the final filtering stage. A single sample was procured from each apiary in 2014/15, although the extraction date is unknown. Honey was sourced from small, independent apiaries that often produce a single harvest per year, with harvest typically occurring in late summer.

| Landscape analysis
The hive locations were identified by postcode of apiaries, at the re- Locations of apiaries were mapped using QGIS v.2.14.0 (QGIS Geographic Information System, 2020) based on the geographic coordinates (Table 1). To adequately describe land cover across a range of spatial scales most typically used by foraging bees (Garbuzov et al., 2015;Sponsler & Johnson, 2015), the landscape surrounding each apiary was characterized using buffers with radii of 500, 1000, 2500, and 5000 m. The majority of foraging flights for both pollen and nectar occur within a radius of ≤5000 m, although less foraging behavior has been known to occur further afield when necessary (Beekman & Ratnieks, 2000;Couvillon et al., 2015). The GB 25m raster land cover data set for the study area was obtained that identifies 21 landscape classes (Rowland et al., 2017). The proportion of each of the 21 land cover class at the different spatial scales was determined using LECOS (Jung, 2013), a QGIS plugin for calculating patch-based landscape metrics (Tables 1 and 2 and Figure 1). To test for associations between proportions of land cover types in the buffer zones, a Spearman's correlation coefficient was calculated for each pair using the R "psych" package (Revelle, 2017).

| Total DNA extraction
Total DNA was extracted from 40 g honey using a modified protocol for the DNeasy Plant Mini Extraction Kit (Qiagen) described in Hawkins et al. (2015). Each honey sample was homogenized by stirring with a sterile stirrer. Four subsamples of 10 g were diluted in TA B L E 1 Names, abbreviated identifiers, sample type, latitude, longitude, and proportion of urban land cover in the 500-m, 1000-m, 2500m-and 5000-m buffers surrounding the postcode of the apiary calculated by GIS analysis of the 14 hives used in the study.

| DNA amplification and sequencing
DNA was amplified using two sets of PCR primers, one amplifying the chloroplastic rbcL gene (Hollingsworth et al., 2009) and one amplifying a plant-specific variant of the internal transcribed spacer region 2 (ITS2) of the nuclear ribosomal region (Chen et al., 2010).
Herein, they are referred to as rbcL and ITS2 pLant (ITS2p), respectively (Table 3). Control amplifications using negative control DNA extractions as the DNA template produced no visible bands on 2% electrophoresis gels and were not progressed further. The PCR products were prepared for sequencing using a two-stage PCR protocol Illumina HiSeq 2500 as a 2 × 300 bp rapid run using a V2 flowcell.
Raw DNA sequence data are available from the N.C.B.I. sequence read archive under accessions SRR16143791 to SRR16143818.
Any ASVs assigned a taxonomy outside the plant kingdom were removed. Any ASV lacking species-level resolution was searched against the N.C.B.I. nucleotide database using megablasts for highly similar sequences (Camacho et al., 2009). Species-level taxonomy was assigned where the best match was achieved against a voucher specimen of a single species. If multiple species in the same genus were equally probable, then genus-level taxonomy was assigned. If multiple genera were equally probable, then family-level taxonomy was assigned, and so on. In the data of each gene region, ASVs with identical taxonomic assignments were collapsed together and ASV counts combined.
Low-frequency incidences of collapsed ASVs were removed from individual samples; where the percentage of reads associated with an ASV in a single sample was <0.03% of total reads associated with that ASV. Low-frequency ASVs within each sample (<1% of total sample reads) were also subsequently removed (Taberlet et al., 2018). Species missing from these sources were checked for availability in online, UK-based garden centers and assigned UK plausibility accordingly. Genus-level assignments were filtered based on the plausibility of the genus being present in the UK based on the same records as the species-level data. Implausible taxa were removed, and abundance matrices from the two gene regions were combined into a single, unweighted presence-absence matrix to maximize the detection range. Only ASVs achieving genus-level assignment or better were retained. Species-level assignments were assigned a category based upon their status in the UK, as either native, nonnative, or neophyte.

| Plant community composition
Per-apiary, mean rbcL reads passing quality trimming and denoising were 120,773.6 (SD: 59,371. (Impatiens glandulifera) and two unresolved genera. Shannon diversity ranged from a minimum of 1.1 to a maximum of 3.14 ( Table 4).
Taxonomic resolution differed between the two gene regions with entire genera unable to be resolved to species in both. Due to the different specificities of each gene region, only Hydrangea and Juglans were unable to be classified to species level in either region, and as such the combination of data from both regions gives broad spectrum taxonomic detection. Unspecified Hydrangea spp.
and Juglans spp. reads accounted for 0.97% and 0.37% of raw reads, respectively.
Considering the ASVs with species-level assignment, native, non-native, and neophyte species were each found in most, but not all hives, and the number of species in each category was indepen- The most common plant species was Impatiens glandulifera, a neophyte, which was found in honey from every apiary. Also very common were Olea europaea and Rubus armeniacus, both nonnatives, and Trifolium repens, a native, as well as Impatiens spp. and Rubus spp. which did not achieve species-level resolution. These common ASVs were all found in honey from >50% apiaries. Other trees and shrubs detected included those frequently found in towns and cities such as Quercus spp., Tilia spp., Malus spp., Buddleja spp.,  F I G U R E 2 Shannon diversity of plant taxa within a hive is significantly negatively correlated with Shannon diversity of land cover in the surrounding 5000 m, a distance likely to cover the majority of foraging flights by bees. Each data point represents a single hive. Pearson's correlation coefficient is presented along with the linear regression line. The shaded area represents the 95% confidence interval.

F I G U R E 3
Proportion of urban land in the 5000-m buffer surrounding a hive (a distance likely to cover the majority of foraging flights by bees) is significantly positively correlated with the Shannon diversity of the plant taxa detected in the honey of the hive. Pearson's correlation coefficient is presented along with the linear regression line. Each data point represents a single hive. The p-value has been adjusted to control for false discovery. The shaded area represents the 95% confidence interval.
of very little vegetation cover and also includes areas such as docks, car parks, and industrial estates. Improved grassland includes high production grassland, characterized by a lack of winter senescence and is sometimes heavily grazed.

| Influence of landscape components on honey bee diet
We found plant taxa richness to be significantly negatively correlated with diversity of land cover surrounding apiaries. Habitat heterogeneity theory tells us that a larger, more heterogeneous environment provides a greater number and wider variety of available habitats or niches and is therefore likely able to support a more diverse flora and fauna (Kallimanis et al., 2008) in apparent contradiction to this finding. In terms of diversity of plant forage available, highly heterogeneous urban landscapes can, in some instances, host more diverse plant communities than landscapes consisting of more diverse but homogenous, land cover types.
In this, they can offer an attractive refuge for a diverse community of bees and other pollinators (Daniels et al., 2020;Hall et al., 2017;Hülsmann et al., 2015;Kowarik, 2011;Lowenstein et al., 2019;Somme et al., 2016;Theodorou et al., 2020). This is further supported by our findings as we recorded a significant positive association between the proportion of urban land cover surrounding apiaries and the diversity of plants detected in the honey. Manchester city center has a disproportionately large number of high-density residential properties (62.6% of all housing in the city) (Baker et al., 2018), and we see that in other comparable cities, the presence of many smaller gardens, cultivated or left wild, provide a diverse forage for bees (Gaston et al., 2005;Lowenstein & Minor, 2016). The post-industrial cityscape also contains many brownfield sites described as being characteristically long-term derelict, vacant, and/or contaminated (Dixon et al., 2010), as well as verges, canal towpaths, and other unmanaged areas. Unmanaged areas, urban meadows, and private gardens are very often occupied by native "weed" species, many of which are highly prized sources of pollen and nectar (Sponsler & Johnson, 2015;Turo & Gardiner, 2019;Weaver, 1965). These species often provide their floral rewards either very early or very late in the bee foraging season, providing high value nutrition when forage availability might otherwise be low (Hicks et al., 2016).
Garden escapees, alongside wild opportunistic seedlings in urban areas, have the added advantage of being unlikely to be treated with pesticides common in agricultural, horticultural, and floricultural trades (Goulson et al., 2018;Lentola et al., 2017). Studies in other urban areas have shown higher plant diversity in the private gardens of diverse areas with both ornamental and weed species contributing to the complexity, and there is ample scope for high-density vertical planting (green walls, planters on balconies of high-rise buildings) alongside relatively little use of ornamental F I G U R E 4 Analysis of the relationship between the proportion of improved grassland (IG) in the 5000-m buffer surrounding a hive (a distance likely to cover the majority of foraging flights by bees) and the Shannon diversity of the plant taxa detected in the honey showed a significant negative correlation. Pearson's correlation coefficient is presented along with the linear regression line. The p-value has been adjusted to control for false discovery rate. Each data point represents a single hive.
lawns in highly urbanized municipal planting schemes and domestic gardens (Aronson et al., 2017;Knapp et al., 2012;Lowenstein & Minor, 2016). Several studies now show that urban and suburban environments appear to support a greater diversity of pollen in the diet than that provided by other surrounding land cover types (Lucek et al., 2019;Richardson et al., 2021). Our study adds to and supports this body of work, expanding it to include other sites.
Every hive in our study had multiple types of land cover in the surrounding buffers, with larger buffers more likely to capture a diverse range of land cover types. As such, even those hives in which the smallest (500 m) buffer was dominated by a single land cover type, the proportion of that dominant type was reduced at larger buffers.
In the cases of four hives with very high proportions of urban land cover in the 500-m buffer, the proportion of urban land decreased as the buffer size around the hive increased. Foragers from these hives, therefore, have ample scope to access other types of land cover on foraging flights, most commonly suburban areas, the secondary land cover type surrounding these hives. Previous studies have demonstrated that urban bees in the UK can access ample floral resources at close proximity and forage mostly within a smaller range (500 m to 1.2 km) than the maximum distances recorded by foraging bees (~12 km); however, this does appear to be seasonal and not universal to all urban areas (Beekman & Ratnieks, 2000;Garbuzov et al., 2015;Sponsler & Johnson, 2015;Steffan-Dewenter & Kuhn, 2003).
Honey sampling is recognized to provide detection of plants over a broad temporal range and here enabled the detection of a wide range of plants Hawkins et al., 2015;Louveaux et al., 1978). It is reasonable to argue that our analysis describes forage collected in highly urbanized areas at some point of the foraging year, but we are unable to describe the foraging range of bees from these hives conclusively. The pooling of honey from multiple combs during a traditional extraction, albeit from within an apiary, will very likely increase pollen diversity in a given sample due to honey being made at different stages of the foraging season. While our samples were a mixture of extraction types, removal of non-traditional (single frame) samples did not yield any notable differences in our results, and as such all samples were included in analyses. The patterns we report in the data may well be specific to samples collected at a particular time of the foraging season, as forage diversity is known to vary throughout the year (Requier et al., 2015). (very likely Impatiens glandulifera, but species resolution not possible for this ASV), and Rhododendron spp., all of which are known to be rich sources of pollen (Hicks et al., 2016). Clearly absent from this rural hive are some of the woody species that make up a large component of the diet of other hives. Notable by their absence are the Oleaceae, Fabaceae, Fagaceae, and Brassicaceae, which although not uniformly present in every other hive, are all common families across the data set. The reliance of bees from this rural hive on a very small number of taxa is of concern. The loss of flower resources due to farming intensification is recognized as an important driver in pollinator declines (Potts et al., 2010). Furthermore, the importance of the introduction or restoration of flower-rich habitats in improved grasslands in order to enhance biodiversity for pollinators has also been established previously (Orford et al., 2016).

| Plant metabarcoding technical considerations
We found that the combination of laboratory and bioinformatics methods employed produced many false positive results at the taxon assignment stage, in common with many metabarcoding studies (Ficetola et al., 2016;Porter & Hajibabaei, 2018;Zinger et al., 2019). For example, the data described plant species in our samples that were unlikely to be growing in the region. A database of plausible taxa in the ecosystem is therefore invaluable for quality control and should be generated and evaluated with the highest possible level of stringency. In the British Isles, we have an extremely well-characterized and barcoded native flora (Ratnasingham & Herbert, 2007;Stace, 2010), but the range of plants visited by bees and other pollinators is often much more diverse. A generalist forager such as a honey bee will forage exotic cultivars in gardens, and some invasive species are known to provide the majority of nectar and pollen to a hive when available (Donkersley et al., 2017), as we found in the present study. Where research is focused upon urban and suburban ecosystems, in particular, it is important to adjust the criteria by which plausible plant taxa are filtered, as a simple filter which only passes native, or naturalized, plant species will not suffice. Initial exploratory analysis of the data revealed that implausible taxonomies were much more likely to be assigned to ASVs present at extremely low relative abundances. This is potentially due to errors inserted at a low frequency into amplicons during PCR amplification, DNA sequencing or as other artifacts of the data analysis method.
To handle these low-frequency, implausible taxa, we removed any low-frequency ASVs during data processing (Taberlet et al., 2018).
After this blunt-edged, but highly effective, data processing step, relatively few implausible ASVs remained, and any that remained were also subsequently removed by our filter against the database of plausible taxa.

| CON CLUS ION
The

ACK N OWLED G M ENTS
We are very grateful to all the beekeepers who kindly collected and donated honey samples from their hives, as well as the Manchester

CO N FLI C T O F I NTE R E S T
The authors declare no conflicts of interest.