# How wild bees find a way in European cities: Pollen metabarcoding unravels multiple feeding strategies and their effects on distribution patterns in four wild bee species

Joan Casanelles-abella, Stefanie Müller, Alexander Keller, Cristiana Aleixo, Marta Alós Orti, François Chiron, Nicolas Deguines, Tiit Hallikma, Lauri Laanisto, Pedro Pinho, et al.

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Joan Casanelles-abella, Stefanie Müller, Alexander Keller, Cristiana Aleixo, Marta Alós Orti, et al. How wild bees find a way in European cities: Pollen metabarcoding unravels multiple feeding strategies and their effects on distribution patterns in four wild bee species. Journal of Applied Ecology, 2022, UN Decade on Ecosystem Restoration, 59 (2), pp.457-470. 10.1111/1365-2664.14063. hal-03385649v2

HAL Id: hal-03385649

https://hal.science/hal-03385649v2

Submitted on 19 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
Journal of Applied Ecology

MR JOAN CASANELLES ABELLA (Orcid ID: 0000-0003-1924-9298)

CRISTIANA ALEIXO (Orcid ID: 0000-0002-3293-476X)

MS MARTA ALÓS ORTI (Orcid ID: 0000-0002-5722-6583)

MR LOÏC PELLISSIER (Orcid ID: 0000-0002-2289-8259)

DR MARCO MORETTI (Orcid ID: 0000-0002-5845-3198)

Article type : Research Article

Editor : Margaret Stanley

Corresponding Author Email ID: joan.casanelles@wsl.ch

**Title:** How wild bees find a way in European cities: pollen metabarcoding unravels multiple feeding strategies and their effects on distribution patterns in four wild bee species

**Authors:**

Joan Casanelles-Abella 1,2 *

Stefanie Müller 1,3

Alexander Keller 4

Cristiana Aleixo 5

Marta Alós Orti 6

François Chiron 7

Nicolas Deguines 7,8

Tiit Hallikma 6

Lauri Laanisto 6

Pedro Pinho 5

Roeland Samson 9

Piotr Tryjanowski 10

Anskje Van Mensel 9

Loïc Pellissier 2,11 +

This article has been accepted for publication and undergone full peer review but has not been
through the copyediting, typesetting, pagination and proofreading process, which may lead to
differences between this version and the Version of Record. Please cite this article as doi:
[10.1111/1365-2664.14063](https://doi.org/10.1111/1365-2664.14063)

This article is protected by copyright. All rights reserved
1 Biodiversity and Conservation Biology, Swiss Federal Institute for Forest, Snow and Landscape Research WSL, Birmensdorf, Switzerland

2 Institute of Terrestrial Ecosystems, ETH Zurich, Zurich, Switzerland

3 Department of Evolutionary Biology and Environmental Studies, University of Zurich, Switzerland

4 Organismic and Cellular Interactions, Biocenter, Faculty of Biology, Ludwig-Maximilians-Universität München, Martinsried, Germany

5 Centre for Ecology, Evolution and Environmental Changes (cE3c), Faculdade de Ciências, Universidade de Lisboa, Lisboa, Portugal

6 Institute of Agricultural and Environmental Sciences, Estonian University of Life Sciences, Tartu, Estonia

7 Université Paris-Saclay, CNRS, AgroParisTech, Ecologie Systématique Evolution, Orsay, France

8 Laboratoire Ecologie et Biologie des Interactions, Equipe Ecologie Evolution Symbiose, Université de Poitiers, UMR CNRS 7267, France

9 Laboratory of Environmental and Urban Ecology, Department of Bioscience Engineering, University of Antwerp, Antwerp, Belgium

10 Poznan University of Life Sciences, Department of Zoology, Poznań, Poland

11 Land Change Science, Swiss Federal Institute for Forest, Snow and Landscape Research WSL, Birmensdorf, Switzerland

<footnote>* Corresponding author</footnote>

<footnote>+ Joint senior position</footnote>
ABSTRACT

1. Urban ecosystems can sustain populations of wild bees, partly because of their rich native and exotic floral resources. A better understanding of the urban bee diet, particularly at the larval stage, is necessary to understand biotic interactions and feeding behaviour in urban ecosystems, and to promote bees by improving the management of urban floral resources.

2. We investigated the larval diet and distribution patterns of four solitary wild bee species with different diet specialization (i.e. *Chelostoma florisomne*, *Osmia bicornis*, *Osmia cornuta* and *Hylaeus communis*) along urban intensity gradients in five European cities (Antwerp, Paris, Poznan, Tartu and Zurich) using two complementary analyses. Specifically, using trap-nests and pollen metabarcoding techniques, we characterized the species' larval diet, assessed diet consistency across cities, and modelled the distribution of wild bees using species distribution models (SDMs).

3. Our results demonstrate that urban wild bees display different successful strategies to exploit existing urban floral resources: not only broad generalism (i.e. *H. communis*), but also intermediate generalism, with some degree of diet conservatism at the plant family or genus level (i.e. *O. cornuta* and *O. bicornis*), or even strict specialization on widely available urban pollen hosts (i.e. *C. florisomne*). Furthermore, we detected important diet variation in *H. communis*, with a switch from a herbaceous pollen diet to a tree pollen diet with increasing urban intensity.

4. Species distribution modelling indicated that wild bee distribution ranges inside urban ecosystems ultimately depend on their degree of specialization, and that broader diets result in less sensitivity to urban intensity.

5. *Policy implications*. Satisfying larval dietary requirements is critical to preserving and enhancing wild bee distributions within urban gradients. For high to intermediate levels of feeding specialization, we found considerable consistency in the preferred plant families or genera across the studied cities, which could be generalized to other cities where these bees occur. Identifying larval floral preferences (e.g. using pollen metabarcoding) could be helpful for identifying key plant taxa and traits for bee survival and for improving strategies to develop bee-friendly cities.

**KEYWORDS:** cavity-nesting bees, feeding behaviour, remote sensing, trap-nests, urbanization, urban biodiversity, species distribution models

1. INTRODUCTION

Wild bees are responsible for major ecosystem functions and make many contributions relevant to people, including pollination and maintenance of ecosystem stability, and they represent social and cultural values
(e.g. Potts et al., 2016). Over the last decades, bee populations have dramatically decreased (Zattara & Aizen, 2021). Multiple causes have been identified (Goulson et al., 2015), and the loss of floral resources is among the most important (Goulson et al., 2015). Wild bees critically depend on large amounts of nectar and pollen to survive and reproduce during their life cycle, and most species display some degree of fidelity to specific plant taxa (Goulson, 1999; Vanderplanck et al., 2014). Consequently, diet specialization and diet preference are two key traits determining the sensitivity to land-use changes (e.g. urbanization; Dharmarajan et al., 2021) and the distribution patterns of wild bees (e.g. Fournier et al., 2020).

Urban ecosystems can harbour large and diverse bee communities, helping to preserve and promote wild bee diversity. Although urbanization has major negative impacts on biodiversity (Theodorou et al., 2020), a significant number of wild bee species can thrive in cities. Documentation of wild bees in urban ecosystems has frequently indicated diverse wild bee communities (Baldock et al., 2015; Casanelles-Abella, Chauvier, et al., 2021), although this ultimately depends on each species' traits and its response to urbanization. Urban ecosystems are warmer (Roth et al. 1989), have higher landscape heterogeneity (Turrine & Knopp, 2015) and are generally less polluted by pesticides (Scheyer et al. 2007) than intensive agricultural areas. Moreover, while intensified agricultural systems have impoverished floral resources, in cities these resources might be maintained, thanks to social investment, high availability of woody species (e.g. street trees) in highly urbanized areas, and the presence of flower-rich habitats (Somme et al., 2016; Tew et al., 2021). In both public and private urban greenspaces, there are important efforts to establish and maintain flowering plant assemblages, with each phase reflecting the preferences and needs of the specific owners and managers (Harrison & Winfree, 2015). Therefore, there is a major opportunity to promote wild bee fitness and reproduction by increasing and improving wild bee habitats.

Urban ecosystems can induce dietary changes in species, due to their distinct availability of various food resources. In urban ecosystems, natural food resources are complemented with anthropogenic food resources (Faeth et al., 2005), whose accessibility is modulated by each species' diet specialization. Urban floral resources are especially diverse in cities, due to gardening and horticultural activities, with many native and exotic species planted for different purposes. Some of these species provide additional sources of food for pollinators within their range of foraging preference, phenology and trait matching (Garbuzov & Ratnieks, 2014; Harrison & Winfree, 2015). Therefore, generalist wild bee species with a broad dietary range might be better able to exploit the existing urban resources, access and forage on a greater variety of patches, and consequently be more widely distributed.

Knowledge on bee diet preferences could reveal which plant species are important for their survival and reproduction and could be translated into important decisions concerning the planning and management of floral resources, e.g. what species to plant. Plant identification with DNA metabarcoding techniques has increased in diet studies and provides new knowledge about the feeding preferences of animals, which can
help us to understand their distribution along environmental gradients (Pitteloud et al., 2021). So far, diet preferences have mostly been assessed indirectly through observations of adult bee plant visitation (e.g. Marquardt et al., 2021) or through the morphological identification of pollen grains (Haider et al., 2013; Sedivy et al., 2011). Nonetheless, specific sampling methodologies, such as trap-nests that allow standardized sampling (Staab et al., 2018), combined with metabarcoding techniques promise to be a powerful tool to characterize and study larval bee diets (Bell et al., 2016; Keller et al., 2015). Trap-nests target larval pollen and thus can better describe bee diet preferences than measurements of adult visitation, while pollen metabarcoding techniques can identify a larger number of taxa with a higher taxonomic resolution than pollen morphological identification; this is particularly useful in urban ecosystems with unique and rich plant pools. Metabarcoding techniques reduce the need for taxonomic expertise associated with pollen morphological identification, thus broadening its application across multiple sites.

Here, we investigated the larval diet and distribution of four widespread wild bee species of urban ecosystems, representing a gradient of decreasing diet specialization (*Chelostoma florisomne*, *Osmia bicornis*, *Osmia cornuta* and *Hylaeus communis*), along urban intensity gradients in five European cities (Antwerp, Paris, Poznan, Tartu and Zurich). In particular, we asked the following questions: (1) What is the taxonomic and trait-based composition of the bee diets in different urban areas? (2) How consistent are the bee diets across urban areas? (3) How does diet specialization influence the bee distribution in urban ecosystems? We expected that specialized bee species (i.e. *C. florisomne*) would have strong preferences for specific plant taxa and thus a highly consistent diet across urban areas and within urban gradients. Conversely, we predicted that more generalist bees (i.e. *O. cornuta*, *O. bicornis* and *H. communis*) would have a more flexible diet and be capable of switching to alternative floral resources, including exotic taxa, and thus have a higher turnover in the diet composition (at the plant family, genus and species level) and a less consistent diet. Finally, we hypothesized that bee species with greater diet specialization would have low flexibility in terms of switching their diet to other plant taxa and thus would be more sensitive to urban intensity.

# 2. MATERIALS AND METHODS

## 2.1 Cities and study sites

We investigated wild bee diets in five cities in Europe: Antwerp (Belgium), Greater Paris (France, hereinafter referred as Paris), Poznan (Poland), Tartu (Estonia), and Zurich (Switzerland), covering a large part of the climatic variability in mainland Europe. Site selection followed Casanelles-Abella, Frey et al. (2021). Overall, we selected sites from the urban green areas mapped and defined in the pan-European Urban Atlas (EEA, 2012). We used an orthogonal gradient of patch area and connectivity. In particular, we calculated connectivity using the proximity index (PI), which considers the area and the distance to all nearby patches with a favourable habitat, within a given search radius. We considered as favourable habitat the land cover classes urban green areas, forest and low urban density with < 30% impervious
surface. We set the search radius to 5 km from each focal patch, the maximum value possible with the available cartography. This resulted in the final selection of 80 sites: 32 in Zurich and 12 in each of the remaining four cities (Fig. S1, Table S1). We maintained a minimum distance of 500 m between selected sites (except for two sites in Zurich selected based on their position in the patch and connectivity gradient, which were separated by 260 m).

## 2.2 Bee sampling

At each site we installed trap-nests in trees, and in three cases (one in Paris, one in Tartu and one in Zurich) in other vertical structures (e.g. lamp post). We constructed trap-nests with reeds and cardboard tubes (Fig. S1). Our sampling trap-nests consisted of a standardized wood box with three plastic pipes 15 cm in diameter and 20 cm long. We assembled the first two pipes using 200–300 reeds from *Phragmites australis* (Cav.) Trin. and 5–10 bamboo reeds with diameters of 1–10 mm and a length of 20 cm to cover all requirements of the cavity-nesting bee community. We assembled the third pipe only using cardboard tubes of 7.5 mm diameter specific for *Osmia* spp. (WAB Mauerbienenzucht, Konstanz, Germany). We installed trap-nests at 2.5–3.5 m height with direct sunlight and SE or SW exposition, and kept them in the field from January until October 2018. In October, we collected the trap-nests and stored them at ca. 5°C until February 2019, and then transferred them to a new room at ambient temperature to recreate spring-like conditions. Bees hatched and were identified to the species level from February to June 2019.

## 2.3 Study organisms

We collected pollen from the nests of four solitary bee species, *Chelostoma florisomne* (Linnaeus, 1758), *Hylaeus communis* (Nylander, 1852), *Osmia bicornis* (Linnaeus, 1758), and *Osmia cornuta* (Latreille, 1805). These species encompass a gradient of diet specialization (i.e. the number of different plant families exploited as resources, from the oligolectic *C. florisomne* to the highly polylectic *H. communis*), differ in phenology, and are common species in urban areas in Europe. In our study, each species was present in at least three of the studied cities. For more details about the ecology of these four wild bee species, see Text S1.

## 2.4 Pollen identification

We extracted a total of 464 pollen samples (Table S1) from undeveloped cells (i.e. where the larva had died) in nests where at least one adult had emerged and thus taxonomic identification of the bees was possible. Specifically, for *C. florisomne* we used 121 samples distributed in 3 cities and 18 sites (2 in Antwerp, 1 in Tartu and 15 in Zurich), for *O. cornuta* we used 66 samples distributed in 3 cities and 20 sites (6 in Antwerp, 5 in Paris and 9 in Zurich), for *O. bicornis* we used 176 samples distributed in 5 cities and 37 sites (3 in Antwerp, 10 in Paris, 8 in Poznan, 1 in Tartu and 15 in Zurich), and for *H. communis* we used 101 samples distributed in 5 cities and 33 sites (4 in Antwerp, 6 in Paris, 6 in Poznan, 9 in Tartu and 8 in Zurich). DNA metabarcoding (isolation, amplification and sequencing) was performed by AllGenetics laboratories (AllGenetics & Biology SL, A Coruña, Spain). We followed the method described by Sickel et
al. (2015) and Vierna et al. (2017) to produce a pooled amplicon library on the ITS2 genomic region for the Illumina platform (Illumina, 2017). See Text S3 for details on the laboratory procedure of pollen metabarcoding. Bioinformatics followed mainly the procedure described in Campos et. al. (2021) with minor modifications: We used VSEARCH v2.14.2 (Rognes, Flouri, Nichols, Quince, & Mahé, 2016) to join paired ends of forward and reverse reads. We also used VSEARCH to remove reads shorter than 200 bp, complete quality filtering (EE < 1 Edgar & Flyvbjerg, 2015) and de-novo chimera filtering, and define amplicon sequence variants (ASVs), as previously done successfully for pollens (e.g. Wilson et al., 2021). The ITS2 rDNA reads were first directly mapped with VSEARCH global alignments and an identity cut-off threshold of 97% against a floral ITS2 reference database generated with the BCdatabaser (Keller et al., 2020), which consisted of plants recorded within the study regions. For the remaining unclassified reads, we first used global alignments against a global reference database (Ankenbrand 2015; Keller et al., 2015). For reads that were still unclassified, we used SINTAX (Edgar, 2016) to assign taxonomic levels as deep as possible but a maximum of genus level with the same global reference database. In total, 82% of species recorded at the sites were present in the local database, 83% of species in the global database (direct classification). Further, 92% of genera were covered by the global database for hierarchical classification. Please note that the global database contains 112.115 unique species, 11.321 genera and 710 families in total, with a very high likelihood of coverage for species and genera of any interest for anthropogenic use, including exotic garden species.

## 2.5 Environmental variables

We assembled variables that were potential drivers of bee diets and distributions and that represented different aspects of the urban environmental gradients. Specifically, we focused on proxies of stress (particularly thermal stress), amount of habitat and resource availability at different spatial scales. We inferred resource availability at the local scale by performing floristic inventories on standardized plots, as explained in Casanelles-Abella, Frey et al. (2021) and in Supplementary Text S4. Furthermore, we collected information on two functional plant traits sensitive to bee-plant interactions, i.e. growth form (Tables S2–S3) and blossom type (Tables S2–S3) using information available in Casanelles-Abella, Frey et al. (2021). See Text S4 and Tables S2–S3 for additional information on the definition of the traits.

We used local and landscape connectivity metrics, local land cover metrics and landscape remote sensing based indices to infer thermal stress and the amount of available habitat, particularly regarding resource availability. As connectivity metrics we used patch size and the proximity index. We obtained the local land cover metrics by mapping grasslands, artificial surfaces, bare land, coniferous trees and deciduous broad-leaved trees and then calculating their proportions at different spatial scales (i.e. 8, 16 and 32 m) from the focal trap-nest (see Text S5 for additional details). Finally, we used remote sensing based indices on land surface temperature, impervious surfaces, soil, water and vegetation at different spatial scales (50, 100, 200, 400, 800, 1600 m). Specifically, we used land surface temperature (LST), the urban index (UI), the colour index (CI), the normalized difference water index (NDWI) and the normalized difference vegetation
index (NDVI), which can be used to characterize existing vegetation and urban infrastructure. In addition, we performed a principal component analysis (PCA) on the explanatory variables to define new meaningful underlying variables while reducing the dimensionality of the data set (see section 2.6). See Text S6 for details on the calculation of the remote sensing based indices and Fig. S2 for the distribution of values of each predictor in each city.

## 2.6 Statistical Analysis

We conducted all analyses in R version 4.0.2 (R Core Team, 2019) with RStudio version 1.4.1106 (RStudio Team, 2020).

### Species diet analysis

We performed taxonomic and trait-based metrics on the bee diets at the city and site level. Specifically, we computed the proportion of different plant taxa at the family, genus and species level (Table S4). Furthermore, we calculated the species, genus and family richness and the Shannon diversity index. Concerning trait-based responses, we calculated the proportion of the different categories of the three studied traits (Table S5). For each of the four studied bees, we performed Pearson correlations to investigate the relationships between the taxonomic and the trait-based diet metrics with the proxies of urban intensity, habitat amount and resource availability. We assessed these relationships: (1) for each single city and (2) for all the cities combined.

We calculated the pairwise correlations between cities for each bee species to study diet consistency. Specifically, we first assembled binary trophic interaction matrices between the four bee species and the plant taxa at the family, genus and species level and then calculated the Pearson correlations of the binary trophic interaction matrices between pairs of cities for each bee and plant level. However, the trophic interaction matrix for a given city, and thus the pairwise correlations between cities, is influenced by the available plant pool. To avoid effects of plant composition, we first created a list with the plant pool occurring in each city at the family, genus and species level. We used the plant species sampled within a 100-m buffer by Casanelles-Abella, Frey et al. (2021) and complemented with the plant species recorded in GBIF at each city for the period 2000–2018 (Fig. S3). If a plant family, genus or species was missing in one of the plant pools of a pair of cities, we removed the interaction when performing the correlations. Moreover, we computed a Chi-squared ($\chi^2$) test on the family and trait composition between cities’ plant species pools (Text S7, Table S6).

### Species distribution of urban gradients

We studied bee distribution patterns with species distribution models (SDMs). We assembled occurrence matrices indicating the occurrence of each bee species in the different sites. From the candidate environmental predictors, we evaluated the statistical relevance of each predictor using the predictive power ($D^2$; Table S7) and then manually picked three predictors that had correlations <0.7 to avoid
collinearity (Fig. S4) for each bee separately. We used an ensemble of two common modelling techniques to account for model uncertainty and specificity (Buisson et al., 2010). Specifically, we used two regression-based models, i.e. generalized linear models (GLMs) and generalized additive models (GAMs), and two tree-based models, i.e. gradient boosting machines (GBMs) and random forests (RFs), that show a higher complexity in their fitting procedures than GLMs and GAMs. We used city as a fixed factor to account for the nested structure of the data with a binomial probability distribution. We parameterized each modelling technique in the following way: we calibrated GLMs with first-order polynomials, GAMs with a spline smoothing term of intermediate complexity (k=4), RFs with a node size of 5 (nodesize = 5) and 1000 trees, and GBMs with an interaction depth of 1, a shrinkage of 0.001 and 1000 trees. We ran the models using the R packages *mgcv* version 1.8-30, *RandomForest* version 4.6-14 and *gbm* version 2.1.5. We randomly split the species records of the four bees into two sets containing 80% of the data for model calibration and 20% of the data for model evaluation. We repeated the procedure five times. We assessed model performance with the True Skill Statistic (TSS; Allouche, Tsoar, & Kadmon, 2006). TSS evaluates model skill in distinguishing absences from presences. The predictive performance of the different models was deemed acceptable when TSS > 0.4, following a commonly used minimum threshold (Thuiller et al., 2019). Thus, we discarded models with TSS values lower than 0.4. We used the selected models of each studied bee species to predict the probability of occurrence over the environmental space of the studied cities.

# 3. RESULTS

## 3.1 Species diet analysis

A total of 41 plant families, 93 genera and 135 species were identified from the nests of the four bee species (Fig. 1, Tables S8–S9). Over half of the species were native (55%), there were more herbs (42%) than trees (34%), and dish-bowl blossoms were more common (56%) (Figs S6–S8, Table S9). The number of plant species per bee nest was similar among bee species (Table S9). The total number of collected plant taxa varied greatly among bee species, reflecting their differences in diet specialization: 1 family and 4 species in *C. florisomne*, 12 families and 33 species in *O. cornuta*, 18 families and 51 species in *O. bicornis*, and 32 families and 81 species in *H. communis*. (Fig. 1, Table S9). At the city level, we found dominance patterns in pollen abundance for some bees (Fig. 1). In *O. bicornis*, the most abundant species in pollen were *Quercus robur* (Antwerp, 70%), *Acer pseudoplatanus* (Paris, 64%; Poznan 44%; and Zurich 33%). In *H. communis*, *Styphnolobium japonicum* was the most abundant species in pollen but only in Paris (52%) and Poznan (32%), with the vast majority of species representing 1-14% of the pollen abundance. Interestingly, in *C. florisomne*, the most abundant *Ranunculus* spp. in pollen changed between cities (*R. acris* in Antwerp, *R. repens* in Tartu and *R. bulbosus* in Zurich). Finally, in *O. cornuta*, no species made up more than 37%, being the most abundant ones *A. pseudoplatanus* (Antwerp, 37%; Paris, 21%; Zurich, 24%) and *Prunus lusitanica* (Paris, 33%). In addition, very few nests were dominated by a single plant species and mostly in *C. florisomne* (Fig. S6).
We found different levels of diet conservatism across cities at the plant family and plant genus level, according to the bee specialization degree and taxonomic resolution of the plant taxa (Fig. 2). At the family and genus level, diet consistency was high for *C. florisomne* and declined with broader feeding niches, particularly at the genus level (Fig. 2a, Table S10). Conversely, we found major variation at the plant species level, which was particularly prominent for the broad generalist *H. communis*, which switched from herbaceous pollen to tree pollen with increasing urbanization (Fig. 3b).

The extent to which bee diet taxonomic and trait-based composition were conserved also varies according to the degree of specialization. *Chelostoma florisomne* had the most conserved diet, composed exclusively of native *Ranunculus* spp. (Figs 1, 2a, Tables S8–S9). *Osmia cornuta* primarily collected the pollen of native tree and shrub species, mainly with dish-bowl or brush type blossoms, from the families Sapindaceae, Salicaceae and Rosaceae (Figs 1, S7–S10, Tables S4, S5, S8–S11). Nevertheless, in Paris and Zurich we also found a considerable proportion of Ranunculaceae (Fig. 1b). Both *C. florisomne* and *O. cornuta* taxonomic and trait-based metrics showed no or little variation along urban intensity gradients (Fig. S11, Table S12). In *O. bicornis*, native tree species with dish-bowl or brush blossoms from the families Sapindaceae and Fagaceae represented a large part of the diet (Figs 1, S7–S10, Tables S4, S5, S8–S11), but we found some variation across cities concerning the remaining species in the diet (Figs. 1b, S10).

*Hylaeus communis* had the most diverse and variable diet. The Fabaceae family represented 34% in Paris and 42% in Poznan of the species found in the larval diet (Fig. 1b) and a minor part in the remaining cities (Fig. 1b). Furthermore, exotic species were more frequent for *H. communis* than for the other three bee species (Figs 3a, S7, Tables S5, S9). In addition, we found family richness, species richness and pollen diversity to be positively correlated with NDVI for *H. communis* for all cities except Poznan (Fig. S11, Table S12). Finally, in Paris, an important part of the diet was trees with flag type blossoms (Fig. 3b), in part due to the contribution of *Styphnolobium japonicum* (Fig 1a, Table S8), which became more dominant in the diet with increasing urban intensity (e.g. decreasing NDVI and increasing UI and CI at different scales; Fig. 3, Table S12).

## 3.2 Species distributions along urban gradients

The PCA conducted on the explanatory variables returned two main axes that explained 38% and 12% of the variation, respectively. The first axis was composed of remotely sensed variables, with larger values of the PC axis indicating less vegetation (i.e. higher UI, CI and LST, and lower NDVI; Fig. S12) independently of the landscape scale considered. The second axis was mostly composed of local land cover variables and metrics representing the available floral resources. Specifically, larger values on the PC axis indicated larger proportions of grasslands and lower proportions of deciduous trees and artificial surfaces, independently of the local scale considered (Fig. S12).
We found two distinct distribution patterns of the four bee species along urban intensity gradients. The first type of response was composed of *C. florisomne*, *O. cornuta* and *O. bicornis*. The probability of occurrence of these three bee species decreased rapidly with increasing urban intensity at the landscape scale (Figs 4a–b, S13) and increased with higher proportions of grasslands at the local scale (Fig. 4c–d). Strikingly, the probability of occurrence peaked at low proportions of broad-leaved trees (Figs 4c–d, S13), even though the diets of *O. cornuta* and *O. bicornis* were composed largely of tree pollen. By contrast, the probability of occurrence of *H. communis* remained constant for larger values of urban intensity (Figs 4a–b, S13). Moreover, the probability of occurrence was positively affected by the amount of deciduous broad-leaved trees at local scales (Figs 4c–d, S13).

4. DISCUSSION

With biodiversity loss increasing at an unprecedented speed (Leclère et al., 2020) and urbanization expanding worldwide, there is a pressing need to better understand the niche requirements of different taxa in order to satisfy them in urban ecosystems through targeted planning and management (Aronson et al., 2017). Here, using trap-nests and pollen metabarcoding techniques, we studied the larval diet and distribution patterns of four wild bee species in five European cities.

The larval diet of the studied oligolectic and intermediate polylectic bees was conserved in terms of plant family composition, representing successful strategies in urban ecosystems as an alternative to generalism. As hypothesized, diet consistency decreased with decreasing wild bee degree of specialization and was relatively high at the plant family level and to lower extend at the plant genus level. Moreover, for *C. florisomne*, *O. cornuta* and *O. bicornis*, bee diet at the plant family level in the studied cities was consistent with prior studies in non-urban ecosystems (e.g. Haider et al., 2014). Diet conservatism regarding plant family composition has also been observed in other polylectic bee species such as *Bombus* spp. (Wood et al., 2021) and *Osmia* spp. (Haider et al., 2014; Vaudo et al., 2020), indicating that diet conservatism might be widespread despite a lack of data. Finally, diet consistency at the plant species level, with the exception of highly specialized *C. florisomne*, was very low, possibly because the specific trait values (e.g. pollen properties imposing cognitive or physiological restrictions; Vaudo et al., 2020) driving plant-pollinator interactions in less specialized wild bees might be relatively consistent within plant families or genera.

Across different animal taxa, broad feeding niches (i.e. low diet conservatism) have been identified as a key characteristics to pass the urban ecological filtering, and hence to thrive in cities (Fournier et al., 2020), with some species even undergoing rapid phenotypic changes to broaden their diet (e.g. Eggenberger et al., 2019). While widely distributed species inside cities are typically generalists, especially when highly mobile and with specific nesting and sociality modes (i.e. cavity nesting or social), intermediate specialization might be an alternative advantageous strategy in urban ecosystems. Focusing on wild bees,
species whose preferred plant families are selected and facilitated by stakeholders (in Central European cities, e.g. Rosaceae, Fagaceae, Salicaceae; Ossola et al., 2020) could have very diversified resources, due to the co-occurrence of spontaneous and cultivated native and exotic plant species. Here, we show that this is the case for *O. cornuta* and *O. bicornis*. Consequently, wild bees could still switch pollen hosts at the species level and better exploit existing resources while avoiding the costs associated with broad niches, as some pollen types are indigestible or toxic (Eckhardt et al., 2014; Praz et al., 2008). On the other hand, strict specialization, as in *C. florisomne*, can still be a successful strategy in urban ecosystems when pollen hosts are highly facilitated and widespread in urban ecosystems. Nonetheless, strict specialists are vulnerable to the partial or total loss of pollen hosts, due to e.g. urban sprawl across the habitat of the bee or pollen host species, changes in social investment when gardening or managing greenspaces, or the arrival of pests.

The degree of specialization in bees is associated with distribution patterns along urban gradients. Our results suggest that increasing specialization leads to a higher sensitivity towards increased urban intensity. Diet specialization determines the possibility of occupying new patches, and greater specialization represents a strong limitation when the nutritional requirements are not met (e.g. with agricultural intensification; Peters, Keller, & Leonhardt, 2021). Rarity in bumblebees has been associated with narrower feeding niches (Goulson, 1999), and this is also the case for other wild bees (e.g. Deguines et al., 2016). Some types of urban greenspaces consistently contain high bee diversity, including several specialists (Baldock, 2020; Salisbury et al., 2015). Still, only a handful of usually broad generalists can colonize various types of urban greenspaces. Hence, these species are widespread within urban ecosystems (Casanelles-Abella, Chauvier et al., 2021; Fournier et al., 2020), even though not only the degree of specialization but also other correlated functional traits are involved in the response to urbanization (e.g. stress tolerance, dispersal; Harrison & Winfree, 2015). To promote a larger number of bees in previously unoccupied areas of the city, wild bee habitat must be strategically increased, including enhancing the availability of high floral resource diversity. A first step to achieve this is to map the existing floral resources within the different urban land covers, making use of ongoing inventories (Ossola et al., 2020) or sampling schemes from research projects (Baldock et al., 2019; Casanelles-Abella, Frey et al., 2021). These products could be combined and compared with wild bee diet preferences to detect where and what kind of planning actions can be taken. For instance, strategically increasing plant diversity rather than overall quantity with targeted taxonomic and trait groups is a measure that has been successful in other ecosystem types (e.g. agroecosystems; Sutter et al., 2017), although other wild bee requirements (e.g. nesting mode) must also be satisfied to successfully promote them (Requier & Leonhardt, 2020).

## 4.1 Limitations and prospects
Pollen metabarcoding only yields relative abundance data and can be subject to PCR or taxon dependent biases as discussed in several studies. (e.g. Bell et al., 2016; Keller et al., 2015) The semi-quantitative abundances obtained in metabarcoding analyses are however still useful to identify the relative proportions of taxa and differentiate dominant, common and rare contributions to a mixed pollen sample by showing correlations (even though not perfect) between read and grain numbers (Keller et al., 2015). It is further agreed that pollen metabarcoding can identify to deep taxonomic levels, is able to detect rare taxa and is well comparable between studies (Bell et al., 2016). In addition, the prevalence of certain plant families and species across cities (Fig. 1), and the existing information on *C. florisomne*, *O. cornuta* and *O. bicornis* (Haider et al., 2013, 2014; Sedivy, Praz, Müller, Widmer, & Dorn, 2008) showing patterns similar to those observed here support our findings. Combining molecular and morphological approaches would provide more robust relative abundance estimates as suggested by Sickel et al. (2015). Further, if taking dry weight of the collected pollen samples is an option (e.g. Lihoreau et al., 2012), these relative data could be converted to absolute abundances per sample. Furthermore, we could only use morphological plant traits to investigate diet composition and shifts along urban gradients. However, bee choices when selecting larval food are likely to be influenced by the nutritional properties of the pollen, such as the content of protein, sugar or other essential nutrients (Vanderplanck et al., 2014; Vaudo et al., 2020). Evidence seems to indicate important differences in the nutritional value of the available floral resources between urban land cover types (Tew et al., 2021) and possibly between ornamental plant species (Garbuzov & Ratnieks, 2014). Incorporating nutritional traits in future dietary studies will better elucidate the mechanisms behind diet composition and trophic niche shifts.

Finally, our results open future research directions. First, further studies could consider whether there are signals of feeding behavioural changes due to anthropogenic factors (e.g. plant or management preferences by stakeholders), as seen in other taxa (e.g. Cucherousset et al., 2012), and if these are reflected in the variation of other functional traits (e.g. tongue length; Eggenberger et al., 2019). To achieve this, it would be critical to study and compare bee diets in both urban and non-urban ecosystems. Second, assessing the effects of intraspecific diet variability on fitness (e.g. number and sex-ratio of offspring, adult and larval survival) would provide important mechanistic insights for anticipating how plant management choices might affect urban bee communities.

## 4.2 Policy implications

The uncovered taxonomic and functional diet preferences can support the planning and management of urban greenspaces to promote wild bees, particularly in wild bee species where the diet patterns are consistent across cities. Typical common urban weeds (e.g. *Taraxacum officinale*, *Bellis perennis*, *Trifolium pratensis*, *T. repens*) that are important floral resources in urban areas (Baldock et al., 2019; Kanduth et al., 2021; Larson et al., 2014) contributed little or nothing to the larval diet of the studied bees. Conversely, both native and exotic woody species proved to be a widely used floral resource. For example, *O. cornuta* and *O. bicornis* collected large amounts of tree pollen from different plant families, which might enable them to exploit secondary pollen hosts by mixing high quality pollens with less digestible ones (e.g.
*Ranunculus* spp., Asteraceae; Eckhardt et al., 2014; Praz et al., 2008). Thus, the maintenance of different vegetation and urban habitat types (e.g. meadows, street trees, shrublands) is of major importance for preserving bees. In addition, our results indicate that the occurrence of specific plant taxa (e.g. at the family or genus level) or trait values is more important than the origin of the plants (Harrison & Winfree, 2015). For example, non-native plants have been observed to retain the blooming time of their original region (Godoy et al., 2009), and thus plant species from the Northern Hemisphere (e.g. planted urban trees or shrubs such *Acer* spp., *Salix* spp., *Crataegus* spp and, *Quercus* spp.) may bloom in synchrony with natives species, providing additional resources for wild bees. Furthermore, urban trees have been shown to be an important resource for several wild bee species (Somme et al., 2016). For instance, in our study, *H. communis* increasingly collected more pollen from ornamental trees (mainly *Styphnolobium japonicum*) with decreasing amounts of greenspace. Because the tree distribution in cities is mostly driven by anthropogenic factors, it represents an important point of action for greening strategies, specifically in densely urbanized city areas with limited herbaceous vegetation. Overall, floral preferences obtained from pollen metabarcoding, particularly when combined with existing information on available floral resources, could help to improve current strategies for developing bee-friendly cities (e.g. in the EU, see Wilk et al., 2020). In particular, characterizing bee diets could inform planning, management and decision making (e.g. what species, genera or families to plant), involving stakeholders, policymakers, nurseries and plant centres, for urban greenspaces to preserve and further promote wild bees in urban ecosystems.

## ACKNOWLEDGEMENTS

This study was supported by the Swiss National Science Foundation (project 31BD30_172467) within the ERA-Net BiodivERsA project "BioVEINS: Connectivity of green and blue infrastructures: living veins for biodiverse and healthy cities" (H2020 BiodivERsA32015104). We thank AllGenetics and especially A. Vizcaíno for conducting the pollen metabarcoding in the middle of a pandemic. We thank the local authorities of Zurich – particularly Grün Stadt Zurich (GSZ), Antwerp, Paris, Poznan and Tartu for supporting this study. Moreover, we acknowledge L. Roquer, J. Bosch and A. Rodrigo for their input on the trap-nest and pollen collection. We thank E. Eggenberg, K. Kilchhofer, P. Bischof, R. Veenstra, Ł. Myczko, A. Zanetta, D. Frey and L. Villarroya-Villalba for helping collect and process the bee data. We deeply thank B. Peters, C. Bello and C. Waldock for their input on data analysis, A. Müller for his guidance on bee biology and nutrition, C. Pitteloud for her help preparing the metabarcoding analyses, M. T. J. Johnson for his input in urban eco-evolutionary dynamics, and J. Bascompte for providing valuable comments on the manuscript.

## CONFLICT OF INTEREST

The authors declare no conflict of interest.

## AUTHORS’ CONTRIBUTIONS
JCA, SM, MM, and LP conceived the ideas. MM, FC, LL, RS, PP and PT developed the site selection methodology. JCA collected the bee data. AK performed the bioinformatic analyses on the pollen. JCA, SM, ND, TH and MAO collected the plant data. MAO, TH and LL provided the local land cover metrics. CA and PP provided the remote sensing data. JCA and SM analysed the data. JCA led the manuscript writing, with editing and supervision by MM and LP. All the authors reviewed and revised the manuscript and gave final approval for publication.

## DATA AVAILABILITY STATEMENT

Data on diet composition supporting this paper available via EnviDat v.0.1 repository
https://doi.org/10.16904/envidat.249 (Casanelles Abella, J., et al., 2021)

## REFERENCES

Allouche, O., Tsoar, A., & Kadmon, R. (2006). Assessing the accuracy of species distribution models: Prevalence, kappa and the true skill statistic (TSS). *Journal of Applied Ecology*, *43*(6), 1223–1232.
doi: 10.1111/j.1365-2664.2006.01214.x

Ankenbrand, M. J., Keller, A., Wolf, M., Schultz, J., & Förster, F. (2015). ITS2 database V: Twice as much. *Molecular Biology and Evolution*, *32*(11), 3030–3032. doi: 10.1093/molbev/msv174

Aronson, M. F., Lepczyk, C. A., Evans, K. L., Goddard, M. A., Lerman, S. B., Maclvor, J. S., ... Vargo, T. (2017). Biodiversity in the city: key challenges for urban green space management. *Frontiers in Ecology and the Environment*, *15*(4), 189–196. doi: 10.1002/fee.1480

Baldock, K. C. R. (2020). Opportunities and threats for pollinator conservation in global towns and cities. *Current Opinion in Insect Science*, *38*, 63–71. doi: 10.1016/j.cois.2020.01.006

Baldock, K. C. R., Goddard, M. A., Hicks, D. M., Kunin, E., Mitschunas, N., Osgathorpe, L. M., ... Memmott, J. (2015). Where is the UK's pollinator biodiversity? The importance of urban areas for flower-visiting insects. *Proceedings of the Royal Society Biological Sciences*, *282*, 20142849. doi: 10.1098

Baldock, K. C. R., Goddard, M. A., Hicks, D. M., Kunin, W. E., Mitschunas, N., Morse, H., ... Memmott, J. (2019). A systems approach reveals urban pollinator hotspots and conservation opportunities. *Nature Ecology and Evolution*, *3*(3), 363–373. doi: 10.1038/s41559-018-0769-y

Bell, K. L., De Vere, N., Keller, A., Richardson, R. T., Gous, A., Burgess, K. S., & Brosi, B. J. (2016). Pollen DNA barcoding: Current applications and future prospects. *Genome*, *59*(9), 629–640. doi: 10.1139/gen-2015-0200

Buisson, L., Thuiller, W., Casajus, N., Lek, S., & Grenouillet, G. (2010). Uncertainty in ensemble forecasting of species distribution. *Global Change Biology*, *16*(4), 1145–1157. doi: 10.1111/j.1365-2486.2009.02000.x
Campos, M.G., Anjos, O., Chica, M., Campoy, P., Nozkova, J., Almaraz-Abarca, ... Carreck, N.L. (2021)
Standard methods for pollen research. *Journal of Apicultural Research*, 60(4), 1–109
https://doi.org/10.1080/00218839.2021.1948240

Casanelles-Abella, J., Chauvier, Y., Zellweger, F., Villiger, P., Frey, D., Ginzler, C., ... Pellissier, L. (2021).
Applying predictive models to study the ecological properties of urban ecosystems: A case study in
Zürich, Switzerland. *Landscape and Urban Planning*, 214, 104137. doi:
10.1016/j.landurbplan.2021.104137

Casanelles-Abella, J., Frey, D., Müller, S., Aleixo, C., Ortí, M. A., Deguines, N., ... Moretti, M. (2021). A dataset of the flowering plants (Angiosperma) in urban green areas in five European cities. *Data in Brief*, accepted, 107243. doi: 10.1016/j.dib.2021.107243

Casanelles Abella, J., Moretti, M., Pellissier, L., Keller, A., Müller, S., Samson, R., Pinho, P., Chiron, F.,
Laanisto, L., Tryjanowski, P., Villarroya Villalba, L. (2021). Larval food composition of four wild bee
species in five European cities. EnviDat. doi:10.16904/envidat.249.

Cucherousset, J., Boulêtreau, S., Azémar, F., Compin, A., Guillaume, M., & Santoul, F. (2012). "Freshwater killer whales": beaching behavior of an alien fish to hunt land birds. *PLoS ONE*, 7(12), 1–6. doi: 10.1371/journal.pone.0050840

Deguines, N., Julliard, R., de Flores, M., & Fontaine, C. (2016). Functional homogenization of flower visitor communities with urbanization. *Ecology and Evolution*, 6(7), 1967–1976. doi: 10.1002/ece3.2009

Dharmarajan, G., Gupta, P., Vishnudas, C. K., & Robin, V. V. (2021). Anthropogenic disturbance favours generalist over specialist parasites in bird communities: Implications for risk of disease emergence. *Ecology Letters*, 24(9), 1859-1868. doi: 10.1111/ele.13818

Eckhardt, M., Haider, M., Dorn, S., & Müller, A. (2014). Pollen mixing in pollen generalist solitary bees: A possible strategy to complement or mitigate unfavourable pollen properties? *Journal of Animal Ecology*, 83(3), 588–597. doi: 10.1111/1365-2656.12168

Edgar, R. (2016a). SINTAX: a simple non-Bayesian taxonomy classifier for 16S and ITS sequences.
*BioRxiv*, 074161. doi: 10.1101/074161

Edgar, R. (2016b). UCHIME2: improved chimera prediction for amplicon sequencing. *BioRxiv*, 074252. doi: 10.1101/074252

Edgar, R. C., & Flyvbjerg, H. (2015). Error filtering, pair assembly and error correction for next-generation sequencing reads. *Bioinformatics*, 31(21), 3476–3482. doi: 10.1093/bioinformatics/btv401

EEA. (2012). Urban Atlas. Retrieved February 16, 2017, from http://land.copernicus.eu/local/urban-atlas/urban-atlas-2012
Eggenberger, H., Frey, D., Pellissier, L., Ghazoul, J., Fontana, S., & Moretti, M. (2019). Urban bumblebees are smaller and more phenotypically diverse than their rural counterparts. *Journal of Animal Ecology*, 88(10), 1522–1533. doi: 10.1111/1365-2656.13051

Faeth, S. H., Warren, P. S., Shochat, E., & Marussich, W. A. (2005). Trophic dynamics in urban communities. *BioScience*, 55(5), 399–407. doi: 10.1641/0006-3568(2005)055[0399:TDIUC]2.0.CO;2

Fournier, B., Frey, D., & Moretti, M. (2020). The origin of urban communities: From the regional species pool to community assemblages in city. *Journal of Biogeography*, 1–15. doi: 10.1111/jbi.13772

Garbuzov, M., & Ratnieks, F. L. W. (2014). Quantifying variation among garden plants in attractiveness to bees and other flower-visiting insects. *Functional Ecology*, 28(2), 364–374. doi: 10.1111/1365-2435.12178

Godoy, O., Richardson, D. M., Valladares, F., & Castro-Díez, P. (2009). Flowering phenology of invasive alien plant species compared with native species in three Mediterranean-type ecosystems. *Annals of Botany*, 103(3), 485–494. doi: 10.1093/aob/mcn232

Goulson, D. (1999). Foraging strategies of insects for gathering nectar and pollen, and implications for plant ecology and evolution. *Perspectives in Plant Ecology, Evolution and Systematics*, 2(2), 185–209. doi: 10.1078/1433-8319-00070

Goulson, D., Nicholls, E., Botías, C., & Rotheray, E. L. (2015). Bee declines driven by combined stress from parasites, pesticides, and lack of flowers. *Science*, 347(6229). doi: 10.1126/science.1255957

Haider, M., Dorn, S., & Müller, A. (2013). Intra- and interpopulational variation in the ability of a solitary bee species to develop on non-host pollen: Implications for host range expansion. *Functional Ecology*, 27(1), 255–263. doi: 10.1111/1365-2435.12021

Haider, M., Dorn, S., Sedivy, C., & Müller, A. (2014). Phylogeny and floral hosts of a predominantly pollen generalist group of mason bees (Megachilidae: Osmiini). *Biological Journal of the Linnean Society*, 111(1), 78–91. doi: 10.1111/bij.12186

Harrison, T., & Winfree, R. (2015). Urban drivers of plant-pollinator interactions. *Functional Ecology*, 29(7), 879–888. doi: 10.1111/1365-2435.12486

Kanduth, L., Chartier, M., Schönenberger, J., & Dellinger, A. S. (2021). Red and white clover provide food resources for honeybees and wild bees in urban environments. *Nordic Journal of Botany*, 39(3), 1–11. doi: 10.1111/njb.03005

Keller, A., Danner, N., Grimmer, G., Ankenbrand, M., von der Ohe, K., von der Ohe, W., ... Steffan-Dewenter, I. (2015). Evaluating multiplexed next-generation sequencing as a method in palynology for mixed pollen samples. *Plant Biology*, 17(2), 558–566. doi: 10.1111/plb.12251

Keller, Alexander, Hohlfeld, S., Kolter, A., Schultz, J., Schultz, J., Gemeinholzer, B., ... Ankenbrand, M. J.
(2020). BCdatabaser: On-the-fly reference database creation for (meta-)barcoding. *Bioinformatics*, 36(8), 2630-2631. doi: 10.1093/bioinformatics/btz960

Larson, J. L., Kesheimer, A. J., & Potter, D. A. (2014). Pollinator assemblages on dandelions and white clover in urban and suburban lawns. *Journal of Insect Conservation*, 18(5), 863-873. doi: 10.1007/s10841-014-9694-9

Leclère, D., Obersteiner, M., Barrett, M., Butchart, S. H. M., Chaudhary, A., De Palma, A., ... Young, L. (2020). Bending the curve of terrestrial biodiversity needs an integrated strategy. *Nature*, 585(7826), 551-556. doi: 10.1038/s41586-020-2705-y

Lihoreau, M., Raine, N.E., Reynolds, A.M., Stelzer, R.J., Lim, K.S., Smith, A.D., Osborne, J.L. & Chittka, L. (2012) Radar Tracking and Motion-Sensitive Cameras on Flowers Reveal the Development of Pollinator Multi-Destination Routes over Large Spatial Scales. *PLoS Biology*, 10(9), e1001392.

Marquardt, M., Kienbaum, L., Kretschmer, L. A., Penell, A., Schweikert, K., Ruttensperger, U., & Rosenkranz, P. (2021). Evaluation of the importance of ornamental plants for pollinators in urban and suburban areas in Stuttgart, Germany. *Urban Ecosystems*, 24(4), 811-825. doi: 10.1007/s11252-020-01085-0

Ossola, A., Hoeppner, M. J., Burley, H. M., Gallagher, R. V., Beaumont, L. J., & Leishman, M. R. (2020). The Global Urban Tree Inventory: A database of the diverse tree flora that inhabits the world's cities. *Global Ecology and Biogeography*, 29(11), 1907-1914. doi: 10.1111/geb.13169

Peters, B., Keller, A., & Leonhardt, S. (2021). Diets maintained in a changing world: Does land-use intensification alter wild bee communities by selecting for flexible generalists? *Ecology and Evolution*, in review.

Pitteloud, C., Walser, J. C., Descombes, P., Novaes de Santana, C., Rasmann, S., & Pellissier, L. (2021). The structure of plant-herbivore interaction networks varies along elevational gradients in the European Alps. *Journal of Biogeography*, 48(2), 465-476. doi: 10.1111/jbi.14014

Potts, S. G., Imperatriz-Fonseca, V., Ngo, H. T., Aizen, M. A., Biesmeijer, J. C., Breeze, T. D., ... Vanbergen, A. J. (2016). Safeguarding pollinators and their values to human well-being. *Nature*, 540(7632), 220-229. doi: 10.1038/nature20588

Praz, C. J., Müller, A., & Dorn, S. (2008). Specialized bees fail to develop on non-host pollen: Do plants chemically protect their pollen? *Ecology*, 89(3), 795-804. doi: 10.1890/07-0751.1

Requier, F., & Leonhardt, S. D. (2020). Beyond flowers: including non-floral resources in bee conservation schemes. *Journal of Insect Conservation*, 24(1), 5-16. doi: 10.1007/s10841-019-00206-1

Rognes, T., Flouri, T., Nichols, B., Quince, C., & Mahé, F. (2016). VSEARCH: a versatile open source tool for metagenomics. *PeerJ*, 24(4), e2584. doi: 10.7717/peerj.2584
Roth, M., Oke, T.R. & Emery, W.J. (1989) Satellite-derived urban heat islands from three coastal cities and the utilization of such data in urban climatology. *International Journal of Remote Sensing, 10*(11), 1699-1720. doi: 10.1080/01431168908904002

Salisbury, A., Armitage, J., Bostock, H., Perry, J., Tatchell, M., & Thompson, K. (2015). Enhancing gardens as habitats for flower-visiting aerial insects (pollinators): Should we plant native or exotic species? *Journal of Applied Ecology*, 52(5), 1156-1164. doi: 10.1111/1365-2664.12499

Scheyer, A., Morville, S., Mirabel, P. & Millet, M. (2007) Pesticides analysed in rainwater in Alsace region (Eastern France): Comparison between urban and rural sites. *Atmospheric Environment*, 41(34), 7241-7252. doi: 10.1016/j.atmosenv.2007.05.025

Sedivy, C., Müller, A., & Dorn, S. (2011). Closely related pollen generalist bees differ in their ability to develop on the same pollen diet: Evidence for physiological adaptations to digest pollen. *Functional Ecology*, 25(3), 718-725. doi: 10.1111/j.1365-2435.2010.01828.x

Sedivy, C., Praz, C. J., Müller, A., Widmer, A., & Dorn, S. (2008). Patterns of host-plant choice in bees of the genus Chelostoma: The constraint hypothesis of host-range evolution in bees. *Evolution*, 62(10), 2487-2507. doi: 10.1111/j.1558-5646.2008.00465.x

Sickel, W., Ankenbrand, M. J., Grimmer, G., Holzschuh, A., Härtel, S., Lanzen, J., ... Keller, A. (2015). Increased efficiency in identifying mixed pollen samples by meta-barcoding with a dual-indexing approach. *BMC Ecology*, 15(1), 1-9. doi: 10.1186/s12898-015-0051-y

Somme, L., Moquet, L., Quinet, M., Vanderplanck, M., Michez, D., Lognay, G., & Jacquemart, A. L. (2016). Food in a row: urban trees offer valuable floral resources to pollinating insects. *Urban Ecosystems*, 19(3), 1149-1161. doi: 10.1007/s11252-016-0555-z

Staab, M., Pufal, G., Tscharntke, T., & Klein, A.-M. M. (2018). Trap nests for bees and wasps to analyse trophic interactions in changing environments—A systematic overview and user guide. *Methods in Ecology and Evolution*, 9(11), 2226-2239. doi: 10.1111/2041-210X.13070

Sutter, L., Jeanneret, P., Bartual, A. M., Bocci, G., & Albrecht, M. (2017). Enhancing plant diversity in agricultural landscapes promotes both rare bees and dominant crop-pollinating bees through complementary increase in key floral resources. *Journal of Applied Ecology*, 54(6), 1856-1864. doi: 10.1111/1365-2664.12907

Tew, N. E., Memmott, J., Vaughan, I. P., Bird, S., Stone, G. N., Potts, S. G., & Baldock, K. C. R. (2021). Quantifying nectar production by flowering plants in urban and rural landscapes. *Journal of Ecology*, 109(4), 1747-1757. doi: 10.1111/1365-2745.13598

Theodorou, P., Radzevičiūtė, R., Lentendu, G., Kahnt, B., Husemann, M., Bleidorn, C., ... Paxton, R. J. (2020). Urban areas as hotspots for bees and pollination but not a panacea for all insects. *Nature Communications*, 11(1), 1-13. doi: 10.1038/s41467-020-14496-6
Thuiller, W., Guéguen, M., Renaud, J., Karger, D. N., & Zimmermann, N. E. (2019). Uncertainty in ensembles of global biodiversity scenarios. *Nature Communications*, **10**(1), 1446. doi: 10.1038/s41467-019-09519-w

Turrini, T. & Knop, E. (2015) A landscape ecology approach identifies important drivers of urban biodiversity. *Global Change Biology*, **21**(4), 1652–1667.

Vanderplanck, M., Moerman, R., Rasmont, P., Lognay, G., Wathelet, B., Wattiez, R., & Michez, D. (2014). How does pollen chemistry impact development and feeding behaviour of polylectic bees? *PLoS ONE*, **9**(1), 1–9. doi: 10.1371/journal.pone.0086209

Vaudo, A. D., Tooker, J. F., Patch, H. M., Biddinger, D. J., Coccia, M., Crone, M. K., ... Grozinger, C. M. (2020). Pollen Protein: Lipid Macronutrient Ratios May Guide Broad Patterns of Bee Species Floral Preferences. *Insects*, **11**(2), 132. doi: 10.3390/insects11020132

Vierna, J., Doña, J., Vizcaíno, A., Serrano, D., Jovani, R., & Chain, F. (2017). PCR cycles above routine numbers do not compromise high-throughput DNA barcoding results. *Genome*, **60**(10), 868–873. doi: 10.1139/gen-2017-0081

Wilk, B., Rebollo, V., & Hanania, S. (2020). *A Guide for Pollinator-friendly Cities: Environment A guide for pollinator-friendly cities: How can spatial planners and land-use managers create favourable urban environments for pollinators?* Guidance prepared by ICLEI Europe for the European Commission

Wilson, R. S., Keller, A., Shapcott, A., Leonhardt, S. D., Sickel, W., Hardwick, J. L., ... Wallace, H. M. (2021). Many small rather than few large sources identified in long-term bee pollen diets in agroecosystems. *Agriculture, Ecosystems and Environment*, **310**(June 2020), 107296. doi: 10.1016/j.agee.2020.107296

Wood, T. J., Ghisbain, G., Rasmont, P., Kleijn, D., Raemakers, I., Praz, C., ... Michez, D. (2021). Global patterns in bumble bee pollen collection show phylogenetic conservation of diet. *Journal of Animal Ecology*, **310**, 1365–2656.13553. doi: 10.1111/1365-2656.13553

Zattara, E. E., & Aizen, M. A. (2021). Worldwide occurrence records suggest a global decline in bee species richness. *One Earth*, **4**(1), 114–123. doi: 10.1016/j.oneear.2020.12.005
FIGURES

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>
**FIGURE 1.** Bee larval diet composition in the studied cities. (a) For each bee species, the collected plant species in each city where the bee species was recorded (three cities for *Chelostoma florisomne* and *Osmia cornuta*, five cities for *Osmia bicornis* and *Hylaeus communis*) are shown. The size of the circle represents the mean relative abundance of plant species contributing to pollen samples per city and bee species. (b) For each bee species, the proportion in the pollen of the different collected plant families in the studied cities is shown (mean relative abundance of plant species contributing to pollen samples per city and bee species). Only families with a proportion in pollen ≥0.01 are plotted, whereas the remaining ones are represented in the category "Other families". Note that the proportion in pollen for *O. bicornis* in Antwerp and Tartu was obtained using only four and one samples, respectively. For each bee and city, we provide the number of sites where pollen samples were taken, and the total number of samples. Information on the computation of the phylogenetic tree can be found in Text S2 and Fig. S5. Data supporting (a) can be found in Table S8. Ast, Asteraceae; Api, Apiaceae; Fab, Fabaceae; Fag, Fagaceae; Sal, Salicaceae; Br, Brassicaceae; Ma, Malvaceae; Sap, Sapindaceae; Ra, Ranunculaceae; Pa, Papaveraceae.
<figure><img src="image_6.png" /><figcaption>**FIGURE 2.** Pairwise correlations of the larval diet composition among cities. For each of the four studied bee species, the city pairwise correlations of the collected plant taxa are shown at the family (a), genus (b), and species (c) level. The colour of the dots indicates the value of the correlation, with lower and higher values in orange and blue, respectively. Note that the correlation values are expressed as absolute values. Note also that the pollen for *O. bicornis* in Antwerp and Tartu was obtained using only four and one samples, respectively. Data supporting Fig. 2 can be found in Table S10.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**FIGURE 3.** Trait-based larval diet composition in *Hylaeus communis*. (a–c) Composition of the diet according to the origin status (a), growth form (b) and blossom class (c) of the plant species in the larval pollen. (d–g) first-order GLMs of the proportion of the different plant trait levels in relation to the first (d–f) and second (g–i) PC axes for the origin status (d and f), growth form (e and g), and (f and i) blossom class. Grey shaded bands indicate 95% confidence intervals. Higher PC1 values indicate less vegetation (lower normalized difference vegetation index) and more artificial surfaces (urban index, land surface temperature, colour index). Higher PC2 values indicate a larger proportion of grasslands and lower proportion of deciduous trees at local scales. PC1 explained 38% of the variation and PC2 explained 12% of the variation. See Figs S7–S9 for the trait-based composition and change along urban gradients of the other three bee species.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>FIGURE 4. Bee distribution along urban gradients. (a and c) Loess smoothing of the mean predicted probability of occurrence of the four bee species in relation to the (a) first PCA axis (PC1) and (c) second PCA axis, performed on the explanatory variables, representing 38% and 12% of the variation, respectively. The mean predicted probability of occurrence results from the predicted probabilities of occurrence of the models with TSS > 0.4. Bands represent 95% confidence intervals. (b and d) Variation in the explanatory variables contributing the most to PC1 (b) and PC2 (d). (e) Larger values of PC1 correspond to higher values of impervious surfaces (urban index, UI), bare land (colour index, CI) and land surface temperature (LST) and lower vegetation cover (normalized difference vegetation index, NDVI) at different landscape scales (i.e. 100 and 400 m). (f) Lower values of PC2 correspond to higher proportions of deciduous trees and lower proportions of grasslands and floral resources (plant species richness) at local scales. Other scales and variables have been omitted here for simplicity (see Fig. S13). See also Fig. S12 for more details on the PCA.</figcaption></figure>