Insect biomass is not a consistent proxy for biodiversity metrics in wild bees

Recent studies have reported on dramatic cases of aerial insect population declines by focusing on the measure of the total biomass of caught insects. However, there is currently no consensus about how biomass patterns among sites and habitats might consistently capture the subtleties of changes in aerial insect community structure. Here, we investigated the relationship between the total biomass of wild bees collected using pan traps in urban, agricultural, and semi-natural habitats on one hand, and a spectrum of biodiversity metrics on the other hand, particularly species richness (SR), alpha diversity, functional diversity (FD) and three different forms of phylogenetic diversity (PD). Our results indicate that although biomass is significantly and highly correlated with the abundance of wild bees, it is generally significantly but only moderately and non-linearly correlated to the various facets of wild bee diversity among habitats. By contrast, we also found that all three measures of PD used are consistent across habitats, suggesting that a taxonomic hierarchy based on Linnaean classification could be used as a proxy for the measurement of PD in wild bees, particularly in other well-studied areas such as Western Europe where a multi-gene molecular phylogeny is unavailable as yet. Collectively, our results illustrate the clear limitations of biodiversity monitoring through measures of trapped insects biomass. We advocate for more robust measures of biodiversity trends in wild bees, requiring both standardized surveys, and the identification of caught specimens down to the species level to capture the subtleties of species, traits-based and phylogeny-based community changes across habitats or time. Scaling out this approach is an essential prerequisite for more global conservation planning tailored to the ecological requirements of the targeted insect species.


Introduction
A couple of long-term field studies have recently reported on dramatic declines of aerial insect abundance, both in temperate (Hallmann et al., 2017) and in tropical ecosystems (Lister and Garcia, 2018), raising the alarm about the impacts of such large-scale losses of biodiversity and the associated ecosystem services underlying human well-being (Pimm et al., 1995;Daily, 1997). These studies made headlines around the world, echoing fears of ecosystem collapse, at a time when both the general public and the mass media are increasingly aware of the Anthropocene's central role in both climate change and the Earth's ongoing sixth mass extinction (Ellis et al., 2010;Dirzo et al., 2014). This research was quickly followed by a review on the worldwide decline of the insect fauna, in which Sánchez-Bayo and Wyckhuys (2019a) hypothesized that humanity may witness the extinction of 40% of the world's insect species over the next few decades, and by another report on a 35% decrease in global butterfly populations over the past 20 years (Wepprich et al., 2019).
Long-term studies aiming to identify the underlying causes of mass insect defaunation, or "insect Armageddon" as it was also coined, are de facto rare and biased towards NW Europe and North America (see also Conrad and Woiwod, 2007;Martins et al., 2013;Taron and Ries, 2015). Investigating trends of species richness and abundance globally in insects is particularly methodologically challenging, and there is merit in leading such research with a focus on insects, a hyper diverse group of organisms often less seen as targets of conservation compared to vertebrates (Leandro and Jay-Robert, 2019). These efforts are urgent given the need to understand the key roles of these organisms for the proper functioning of life on Earth. Yet, parallel to this, some important methodological limitations of these recent studies were also pointed out, including the neglected impact of changing landscape structure around the study sites on the local biodiversity, the lack of standardized sampling protocols including sampling effort, and overall the far-fetched conclusions reached by the authors, despite the available evidence (see e.g. Willig et al., 2019;Lister and Garcia, 2019). The review by Sánchez-Bayo and Wyckhuys (2019a) in particular stirred a vivid wave of criticisms, with international teams of researchers questioning the validity of their conclusions after highlighting similar methodological limitations, including the non-adoption of literature search standards (Komonen et al., 2019;Mupepele et al., 2019;Simmons et al., 2019;Thomas et al., 2019; but see also Sánchez-Bayo and Wyckhuys, 2019b).
All these studies and their criticisms illustrate the lack of existing knowledge on insect population trends, and the paucity of long-term and large-scale monitoring programmes (Macgregor et al., 2019), particularly for insect groups that are key for ecosystem functioning and food production, such as bees. For example, the first IUCN Red List of European bees (Nieto et al., 2014) has concluded that over half of Europe's ca. 2000 species of bees are "data deficient", meaning that too little or no information is available on the abundance and distribution of these species or of the nature of the threats they face, despite their pivotal importance for our ecosystems, our pollinator-dependent crops, and ultimately both our food security and our well-being (IPBES, 2016;Potts et al., 2016). Despite all the attention that bees have received ever since awareness increased around the "Colony Collapse Disorder" affecting the western honey bee (Apis mellifera) (see e.g. VanEngelsdorp et al., 2009), and despite the resulting increase of awareness on the sheer diversity (over 20,000 species at the worldwide scale to date, see Ascher and Pickering, 2020) and threats to wild bees (Biesmeijer et al., 2006;Nieto et al., 2014), our "forgotten pollinators" (Buchmann and Nabhan, 1996) still and typically suffer from a combination of multiple shortfalls: (i) a Linnean shortfall, given the discrepancy between formally described species and the number of species that actually exist, even in well-surveyed regions such as Europe and N-America, (ii) a Raunkiaeran shortfall associated to the relatively limited knowledge and lack of compilation of species life-history traits, and finally (iii) an Eltonian shortfall stemming from a relatively shallow understanding of the biotic interactions and ecological requirements of most species (see Hortal et al., 2015).
Addressing the urgent need for rigorous demographic data on the populations of multiple species of wild bees to identify the drivers of biodiversity change also requires deploying systematic monitoring programmes with standardised protocols and survey timing. Mainstreaming the use of meaningful, modern community ecology metrics and proxies such as functional and phylogenetic diversity, and structure that adequately capture trends in the state of biodiversity (McGill et al., 2006(McGill et al., , 2015, while also correlating highly with the provision of ecological processes and services (e.g. Cadotte et al., 2008Cadotte et al., , 2009Cadotte et al., , 2012Martins et al., 2015) also requires coordinated efforts and mainstreaming recent methodological advances. For bees, there is currently no consensus about how bee biomass patterns among sites and contrasting habitats might reliably capture the subtleties of bee community structure. This gap in our knowledge is particularly important in a context of (i) standardized passive sampling and monitoring techniques readily available (see e.g. Westphal et al., 2008;O'Connor et al., 2019;Prendergast et al., 2020) that can potentially be deployed across a wide range of sites and habitats by non-professional volunteer citizens, (ii) the fact that the conservation of the biodiversity of bees as a target insect group to promote the ecosystem services they provide is increasingly on the political agenda, and (iii) a more global need for scientificallyinformed conservation actions to safeguard biodiversity.
In this study, we attempt to contribute to the evaluation of the relationship between the total biomass of wild bees and a spectrum of biodiversity metrics, particularly species richness (SR), alpha diversity, functional diversity (FD) and three different forms of phylogenetic diversity (PD). We systematically surveyed wild bees with passive pan traps in urban, agricultural and semi-natural habitats in Belgium during their peak activity in Spring (April-May-June in 2015-2016-2017-2018), as this period of the year coincides with the highest levels of abundance and species richness of these insects in NW Europe. Specifically, using a standardized sampling effort, we examined the relevance of using the biomass of trapped bees measured simply by weight as a proxy for species richness and other biodiversity metrics, both generally and for each habitat category individually. We then discuss how these results should contribute to informing on the sampling design of a pan-European monitoring scheme as part of the European Pollinator Initiative, to help the EU and its Member States set strategic objectives, develop actions to address the decline of pollinators, and ultimately contribute to global conservation efforts.

Bee surveys
Bees were surveyed at 44 sites in Brussels-Capital Region and Wallonia in Belgium. Each site was surveyed in April, May and June of the same year, and all sites were surveyed between 2015 and 2018. Bees were sampled at each site using a combination of standardized white, yellow, and blue UV-reflecting pan traps (diameter = 16 cm, height = 5 cm), altogether nine pan traps per site for a sampling period of seven consecutive hours, between 09:00 and 16:00 CET. Pan traps are one of the most popular sampling techniques for bees, yielding hundreds to thousands of specimens when used in the field, and they therefore represent an adequate sampling technique to investigate the patterns of insect biomass and how they relate to other diversity metrics. The pan traps were filled with 300 ml of tap water with a few drops of the same colourless, odour-free ionic surfactant, and were placed on the ground, uncovered by surrounding vegetation. The collection dates in April, May and June at each of the 44 sites surveyed in 2015-2018 across urban (n = 15 sites), agricultural (n = 15 sites) and semi-natural (n = 14 sites) habitats in Belgium are listed in Table S1. Each site was sampled once a month.
All specimens were stored in 70 • ethanol after each collection date, and they were pooled for each site*collection date according to the colour of the pan trap they were collected in. The preparation of the specimens for their identification in the lab involved a first wash in a FLOUREON 600 ml (42,000HZ) ultrasonic glass cleaner machine filled with tap water and a few drops of the same colourless, odour-free ionic surfactant. This step was followed by an immersion of the specimens in tap water at room temperature to wash away the remaining soap on their body. The specimens were first dried on highly absorbent paper towels, and then systematically blow dried in tea strainers using a hand hairdryer for 1-2 min to restore the fluffy aspect of their hairiness prior to pinning. These steps were essential to facilitate the identification process.
All specimens collected were pinned, labelled, identified down to the species level and digitized except the specimens belonging to the subgenus Micrandrena (genus Andrena, family Andrenidae) as their reliable identification is known to be problematic due to their small body size and micro-morphological resemblance. This makes the identification of these bees still the subject of debates among the few experts for this group (e.g. Schwenninger, 2009;Dardón et al., 2014). We also excluded honey bees (Apis mellifera) for the analyses presented below, because (i) their presence effectively reflects the presence of apiaries, and (ii) we are interested in describing the biodiversity patterns of wild bees (i.e., non-Apis bees in Europe).

Computation of ecological/behavioural traits and biomass
We followed the approach adopted by Normandin et al. (2017) by compiling ecological, life history and morphological traits for all wild bee species recorded in this study. We classified the bee species according to six functional categorical traits chosen to take into account the diverse attributes of wild bee ecology that are known to influence their link to their environment and their functional role as pollinators. The variables that concerned our species were: pollen specialization (oligolectic, polylectic), pollen transportation (accidental, leg and body, corbiculae, legs only, underside, crop), tongue length (short, long), seasonal activity (spring, summer, spring and summer, all year around), sociality (solitary, communal, cleptoparasite, primitively eusocial, social parasite, eusocial), nesting behaviour (underground nest builder, underground nest parasite, renter of pre-existing cavities above ground, nest parasite of a social species, cleptoparasite of a nest established in a pre-existing cavity below ground). Along with these six functional categorical traits, we also computed a quantitative functional trait for each species, namely the inter-tegular distances (ITD; in mm) on five female specimens of each species with a digital caliper (precision = 0.01 mm). The ITD measures the shortest distance between the tegulae, i.e. the small scale-like sclerites covering the base of the forewing in bees and many other insects; this trait is generally regarded as a good proxy for body size (Cane, 1987) and bee foraging range (Greenleaf et al., 2007;Zurbuchen et al., 2010 and references therein).
We then used a recent allometric model developed by Kendall et al. (2019) which takes into account biogeography, phylogenetic relatedness, and intra-specific variation to use the above-mentioned ITDs to predict bee dry body weight (in mg) for each species. This model was used to calculate the total biomass of the insects caught and it is available in the pollimetry package (https://github.com/liamkendall/p ollimetry) in RStudio (version 1.1.456) (RStudio Team, 2016) for R (R Core Team, 2018).

Computation of biodiversity metrics
The data were first aggregated (i.e., all pan trapped bees were pooled monthly for each site) to retain a single line of species occurrence records for each month (i.e., April vs. May vs. June). This database was used to examine patterns of abundance, species richness, alpha diversity, functional diversity and three forms of phylogenetic diversity among months and habitats using a nested ANOVA through mixed effects models fit by maximizing the restricted log-likelihood (REML) in the nlme package (version 3.1-137) (Pinheiro et al., 2018). We first used the "Month" variable as a fixed effect and the "Habitat" variable as random effect to test the impact of the timing of the surveys on the abundance and species richness of wild bees, and we then used the "Habitat" variable as a fixed effect and the "Month" variable as random effect to test the impact of the habitat type on the abundance and species richness of wild bees. Post-hoc tests in the form of Tukey contrasts (multiple comparisons of means) were performed using the multicomp package (Hothorn et al., 2008).
We computed the Simpson's diversity index D (or α-diversity), a composite quantitative measure that reflects dominance-rarity pattern of the community (i.e. by habitats), taking the relative abundances (evenness) of different species into account (Borcard et al., 2018).
The mixed matrix of qualitative and quantitative traits for each species was then converted into a Gower distance matrix with the gowdis function in the FD package (version 1.0-12) (Laliberté and Legendre, 2010;Laliberté and Shipley, 2014). This Gower distance matrix was used to prepare a functional hierarchical cluster, or functional dendrogram, with the hclust function in the stats package of R (R Core Team, 2018). The functional diversity (FD) index, a measure that quantifies the diversity of functional (or behavioural/ecological) traits within species assemblages, was then computed for each community (i.e., study site) with the BAT package (version 1.6.0.) (Cardoso et al., 2018) by summing up the branch lengths on this functional dendrogram connecting species present in each community (Schleuter et al., 2010).
To explore the patterns of phylogenetic diversity (PD) among sites and habitats, we computed the PD index as described for FD above with the BAT package (version 1.6.0.) (Cardoso et al., 2018), using three distinct phylogenetic classifications: a phylogeny based on the hierarchical Linnaean classification, a phylogeny based on sequences of the COI barcode gene, and a multi-gene phylogeny.
The phylogeny based on the hierarchical Linnaean classification follows the higher-rank classification proposed by Danforth et al. (2006), and we followed Michener (2007) for the taxonomy from family to subfamily, tribe, genus, subgenus, and species. The species names used here are in agreement with the nomenclature used in the IUCN Red List of European Bees (Nieto et al., 2014) and the species checklist of EU bees (Michez et al., 2019). The R-package ape (version 5.0) (Paradis et al., 2004) was used to build a polytomous, ultra-metric tree comprising all species recorded in situ. Branch lengths were calculated for the whole tree by setting the p-parameter to 1, following Hoiss et al. (2012).
The two molecular phylogenies were constructed using (i) just the mitochondrial gene cytochrome oxidase I (COI), and (ii) the mitochondrial gene cytochrome oxidase I (COI), the nuclear genes Elongation Factor-1alpha (EF1α) and Long-Wavelength opsin (LW Rh). These are common molecular markers used for the identification of bee species and their phylogenetic relationships (Mardulyn and Cameron, 1999;Danforth et al., 2006;Magnacca and Brown, 2012;Hedtke et al., 2013;Schmidt et al., 2015). We extracted available COI, EF1α and LW Rh sequences from GenBank platform for all Belgian bee species, i.e. a total of 355 barcodes, 141 EF1α and 163 LW Rh sequences were retrieved (see Table S2 for all gene accession details). We then used a selection of five Crabronidae species as outgroups (Pison chilense (SPINOLA, 1851), Philanthus triangulum (F.), Bembix troglodytes (HANDLIRSCH, 1893), Sphecius speciosus (DRURY, 1773) and Crabro scutellatus (SCHEVEN, 1781)), since this family of solitary wasps is recognized as sister group of the bee clade .
We used two separated molecular datasets to reconstruct a COI phylogenetic tree and a multi-gene phylogenetic tree based on the sequences of the available COI, EF1α and LW Rh genes. All sequences of each molecular datasets were aligned using ClustalX v.2.1. (Larkin et al., 2007), with defaults parameters and pairwise deletion for gap treatments. The obtained alignment contained 1481 nucleotide characters for COI tree and 5295 nucleotide characters for multi-gene tree. We used jModelTest 2.1.10. (Darriba et al., 2012) to explore best nucleotide substitution model on our both aligned DNA datasets. Generalized Time-Reversible model with invariable sites and gamma model of rate heterogeneity (GTR + I + Γ) (Tavaré, 1986) was selected as best nucleotide substitution model for our aligned COI sequences, and also for each partition in our aligned multi-gene sequences. The maximum likelihood (ML) method was conducted in RAxML v.7.7.1. (Stamatakis et al., 2008;Stamatakis, 2014) on the CIPRES Science Gateway (Miller et al., 2010) (http://embnet.vital-it.ch/raxml-bb/). One hundred rapid bootstrap inferences were executed and followed by a thorough ML search from both aligned datasets. Bipartition information from both best-known ML trees (in Newick format) was used to process phylogenetic community structure analysis.
The results of the nested ANOVA's on patterns of abundance and species richness is presented in the manuscript, while the nested ANOVA's on alpha diversity, functional diversity and three forms of phylogenetic diversity are presented in the supplementary material.

Computation of species, functional and phylogenetic turnover metrics
To examine in greater detail, the species, functional and phylogenetic structure of wild bee communities among habitats, we aggregated the dataset to retain a single line of species occurrence records for each habitat (urban vs. agricultural vs. semi-natural), and we used the BAT package (version 1.6.0.) (Cardoso et al., 2018) to explore the species, functional and phylogenetic turnover.
We first ran an analysis of similarities (ANOSIM) using the average Bray-Curtis distances among samples and 10,000 permutations to test statistically whether there is a significant difference in community composition among habitat categories with the vegan package (version 2.0-5) (Oksanen et al., 2019).
We then computed the Sørensen index of beta diversity, βsør, a measure of total dissimilarity, and its partition into its two major components: (1) species replacement among habitats, i.e., turnover: βsim; and (2) species loss/gain among habitats, i.e., nestedness: βsne; following the formula βsør = βsim + βsne (Baselga, 2010;Legendre, 2014). The relative values of βsim and βsne allow to quantitatively evaluate the relative importance of species turnover versus nestedness among habitats. The same procedure was used for the analysis of the species, functional and phylogenetic turnover; for the functional and the phylogenetic turnover analyses we used the functional dendrogram and the multi-gene phylogeny described above, respectively. A Venn diagram illustrating the number of bee species recorded in each habitat category was also created using the VennDiagram package (version 1.6.9) (Chen, 2014).

Pairwise comparisons of biomass and biodiversity metrics
The pairwise comparisons of biomass and biodiversity metrics described above were visualized with the ggpairs function in the GGally package (version 1.4.0.) (Emerson et al., 2012;Schloerke et al., 2018); a loess regression is fitted to the observed data with 95% confidence band intervals around the fit.
The scatter plot produced also allows the computation Pearson's correlation coefficient for each pair of variables, both irrespective of the habitat types, and for each individual habitat type (i.e., urban, agricultural and semi-natural).
All the analyses described in the methods section were performed in RStudio (version 1.1.456) (RStudio Team, 2016) for R (R Core Team, 2018).

Patterns of wild bee biodiversity and community structure
In this study, we recorded a total of 2490 specimens belonging to 125 species of wild bees (excluding Apis mellifera): 67 species/792 specimens in urban habitats, 64 species/874 specimens in agricultural habitats, and 47 species/824 specimens in semi-natural habitats. A total of 21 species of wild bees are shared among all habitats, i.e. 31.3-44.7% of the total number of species found in sites from each habitat category. Urban habitats hosted the highest number of unique species (n = 18 species), followed by semi-natural (n = 15 species) and agricultural sites (n = 14 species) (Fig. S2).
The top 5 species of urban habitats were all generalist, groundnesting species in the families Andrenidae (genus Andrena) and Halictidae (genus Lasioglossum), except Nomada fabriciana (family Apidae), the cuckoo bee associated to Andrena chrysosceles and Andrena bicolor (species #10 and #6 on the rank of urban habitats, respectively, see Fig. S1). Likewise, the communities of agricultural sites were dominated by species in the families Andrenidae (genus Andrena) and Halictidae (genus Lasioglossum), with Andrena haemorrhoa (family Andrenidae) as the dominant species found in apple orchards and market gardening micro-farms. Semi-natural habitats (here, sandy grasslands) were associated (i) to high abundances of ground-nesting solitary bees, particularly Andrena vaga (family Andrenidae), a spring specialist species on willows (Salix spp., Salicaceae), but also (ii) to high abundances of Dasypoda hirtipes (family Melittidae), an early summer specialist species on Cichorioideae (Asteraceae), and Halictus sexcinctus (family Halictidae), a rare generalist species restricted to sandy habitats (grasslands, quarries) of Southern Belgium that shows a social polymorphism in which different populations can exhibit solitary, communal, or more advanced social structure (Richards, 2001;Pauly, 2019).
The number of shared/unique species between and among habitat types are illustrated in the Venn diagram presented in Fig. S2, and the ANOSIM analysis conducted by pooling the monthly data for each site showed a significant difference in community composition among habitat categories (ANOSIM, R = 0.1924, p < 0.001). The analysis of the overall differentiation in community structure among habitat categories illustrated in Table 1 shows that the turnover component (βsim) is primarily responsible for the total species/functional/phylogenetic dissimilarity among habitats (βsør) ( Table 1).
Our results clearly indicate that the sampling period has a significant impact on the levels of wild bee abundance and species richness, with the month of April associated with significantly higher number of specimens and species (Fig. S3). Overall, the three habitat categories were associated with equal number of specimens collected in pan traps (Nested ANOVA through mixed effects model fit by maximizing the restricted log-likelihood (REML) with "Habitat" as a fixed effect, and "Month" as random effect: p = 0.1735) (Fig. S3). By contrast, we found significantly different levels of species richness associated with these similar levels of specimen abundance: semi-natural habitats being associated to a significantly lower species richness compared to both urban (Post-hoc test of Nested ANOVA, p = 0.0379) and agricultural (Post-hoc test of Nested ANOVA, p = 0.0353) habitats. No significant difference in species richness was observed between urban and agricultural habitats across months (Post-hoc test of Nested ANOVA, p = 0.9995) (Fig. S3).
Likewise, our results indicate that there was no significant effect of the habitat type on the alpha diversity on a monthly basis (Nested ANOVA through mixed effects model fit by maximizing the restricted log-likelihood (REML) with "Habitat" as a fixed effect and "Month" as random effect: p = 0.095), but that the sampling period has a significant impact on the levels of wild bee functional and phylogenetic diversity (Nested ANOVA's through mixed effects model fit by maximizing the restricted log-likelihood (REML) with "Month" as a fixed effect and "Habitat" as random effect). The month of April is systematically associated with significantly higher diversity for both FD and PD (Post-hoc tests of Nested ANOVA's, all p < 1e-04), whereas May and June display equal levels of functional and phylogenetic diversity across habitat types (Fig. S4).

Biomass correlation with biodiversity metrics
Total insect biomass of wild bees (in log form, measured as the community-level total dry body weight of all specimens) was highly significantly and positively correlated with the abundance of insects (measured as the number of specimens collected) (Pearson's correlation coefficient = 0.966, p < 0.001) (Fig. 1) and showed a good consistency among habitat categories (Pearson's correlation coefficients range: 0.911-0.997) (Figs. 1 & S7). However, insect biomass was significantly but only moderately correlated to species richness (Pearson's correlation coefficient = 0.521, p < 0.001), functional diversity (FD) (Pearson's correlation coefficient = 0.490, p < 0.001) and the three forms of phylogenetic diversity (PD) tested here (Pearson's correlation coefficient PD taxo = 0.539; Pearson's correlation coefficient PD COI = 0.535; Pearson's correlation coefficient PD multi = 0.493, all with p < 0.001) (Figs. 1 & S7). Biomass and insect abundance were invariably significantly, but negatively, correlated to Simpson's D, suggesting that an increase in the biomass of insects collected reflects an increase of abundance of insects and, to some extent only, of species richness. This, in turn, translates into a decreased species diversity in the communities, presumably because of an increasing numerical imbalance (or decreased evenness) among species, since Simpson's D index of diversity takes into account both the number of species present, as well as the relative abundance of each species.
Perhaps more importantly, we provide evidence that these correlations between total insect biomass and species/functional/phylogenetic biodiversity metrics are not only moderate, but they also vary non- Fig. 1. Scatter plot matrix illustrating the pairwise relationships between biomass and other biodiversity metrics of wild bees collected once a month in April, May and June with UV-painted pan traps in urban (URB., n = 15 sites), agricultural (AGR., n = 15 sites) and semi-natural (SNAT., n = 14 sites) habitats in Belgium. The scatter plot also shows Pearson's correlation coefficient for each pair of variables, both irrespective of the habitat types (« Cor », in black) and for each individual habitat type (i.e., urban, agricultural and semi-natural). The diagonal of the scatter matrix represents the density distribution and the right-hand box plots illustrate the variation of each continuous variable, also coloured by habitat type. The results indicate that although biomass (in log form, measured as the community-level total dry body weight of all specimens) is highly correlated with the abundance of insects (measured as the number of specimens collected), it is generally significantly but weakly correlated to the various facets of wild bee diversity (Pearson's correlation coefficient ~ 0.5). Furthermore, the results indicate that these correlations vary by almost a factor of 2, depending on the habitat type considered. The results also show that the measure of phylogenetic diversity (PD) is consistent among the different trees used (i.e., taxonomic hierarchy based on Linnaean classification vs. COI tree vs. multi-gene tree) (Pearson's correlation coefficient ~0.95), suggesting that a taxonomic hierarchy based on Linnaean classification could be used as a proxy for the measurement of PD in wild bees in NW Europe, and perhaps also in other areas where a considerable share of known species still await sequencing of DNA barcodes and other nuclear genes. See Methods section for details, and Fig. S4 for the results of significance tests on the Pearson's correlation coefficients.
linearly and by almost a factor of 2, depending on the habitat type considered (Figs. 1 & S7). We observed a trend in the measure of biomass being systematically almost twice as weakly correlated with species/functional/phylogenetic biodiversity metrics in semi-natural habitats compared to agricultural or urban habitats (Figs. 1 & S7).
Last, our results also show that all three measures of PD are highly significantly and consistently correlated (all Pearson's correlation coefficients > 0.95 and associated p < 0.001), and that our measure of FD is significantly and positively correlated to our measures of phylogenetic diversity (PD) (Pearson's correlation coefficient with PD taxo = 0.936; Pearson's correlation coefficient with PD COI = 0.876; Pearson's correlation coefficient with PD multi = 0.959, all with p < 0.001), with a significant congruence among habitat categories (Figs. 1 & S7). These results suggest that the combination of behavioural/ecological traits used in this study for the computation of FD is highly redundant with measures of PD, perhaps because (i) of the high number of character states used, and (ii) the overall tendency of wild bees to exhibit relatively high levels of phylogenetic conservatism in their fundamental behavioural/ecological traits (Michener, 2007), including in Europe (Michez et al., 2019).

Community turnover and habitat types
Our results suggest that the turnover among these habitat categories is responsible for the non-linear changes in the correlation between insect biomass and species/functional/phylogenetic biodiversity metrics: the unpredictable species replacement between habitats and the aforementioned unpredictable trends in biomass imply that there is an equal likelihood of turnover between small and large species of bees. This complexity and unpredictability are also illustrated by Shortall et al. (2009) who analysed almost 30 years of trap data from four sites in the United Kingdom, representing one of the few long-term, dedicated insect-monitoring programs in the world. They found significant declines overall in large insects, particularly large flies, but no significant trends for other groups, such as aphids, with location-specific trends. Similar complex and noisy trends were observed in the UK by Bell et al. (2015) who investigated life-history traits dependent changes in aphid phenology under climate change. The extent to which the results observed for wild bee diversity in different habitats in Belgium apply to other insect taxa in the same and other region remains to be investigated more thoroughly, although other reports indicate that bees and other insects are often poor surrogates for each other (Yong et al., 2020;Oberprieler et al., 2020).
Our findings are also in agreement with previous reports on the role of urban green spaces to support high levels of wild bee diversity (see Hernandez et al., 2009;Banaszak-Cibicka and Zmihorski, 2012;Normandin et al., 2017;Theodorou et al., 2020), including anthropophilic species less regularly found outside of urban settings (Fig. S1). This phenomenon is presumably the result of the wide variety of plant genera and families found in urban private and community gardens, as well as in parks, botanical gardens and urban farms such as in our study sites (listed in Table S1) (MacIvor et al., 2014;Pardee and Philpott, 2014;Langellotto et al., 2018;Potter et al., 2019;Egerer et al., 2020).
The wild bee diversity presented in this study is the product of homogenous collection techniques. All species were collected using the same methodology, in Spring and using pan traps. How species are sampled is another key factor in determining biodiversity metrics. Here we show a strong correlation between biomass and abundance for bees, but all specimens were collected using pan-traps and we would not necessarily expect the same relationships given net sampling, particular with inexperienced samplers who are likely to sample a greater proportion of larger species compared to pan traps (Westphal et al., 2008). It is important to note though that we do not compare collection techniques in this study, but we rely on data collected only using passive pan traps. Therefore, the abundance metrics collected, whilst comparable across sites, are likely to be biased towards certain genera and not represent an accurate measure of total wild bee abundance (Portman and Cariveau, 2020). To best approximate true abundances of species, we would need a combination of traps, transect walks and netting (Portman and Cariveau, 2020;Prendergast et al., 2020).

Biomass vs. other biodiversity metrics
In this study, we provide the first attempt to evaluate the relationship between the total biomass of insects caught in "passive" traps and several important biodiversity metrics for wild bees. Our results, based on standardized field surveys in contrasting habitat types in Spring, a period corresponding to the highest levels of abundance and species richness of wild bees in NW Europe, show that biomass and abundance will not consistently allow researchers to estimate meaningful measures of wild bee species richness, or functional and phylogenetic diversity. Furthermore, our results provide evidence that similar levels of specimen abundance can be associated with highly contrasting and significantly different levels of species, functional and phylogenetic turnover.
Our results are in agreement with the recent study by Seibold et al. (2019) who showed that the magnitude of insect decline in forest and grassland sites across Germany is dependent on the biodiversity metric used. Furthermore, these authors also showed that percentage losses in abundance and biomass did not correspond to equivalent losses in species richness. Besides, as with our results, the relationship between these metrics is generally expected to be influenced by the surrounding landscape context (e.g. Hopfenmüller et al., 2014;Kratschmer et al., 2019;Ekroos et al., 2020).
Monitoring losses of total biomass or total abundance of insects can be appealing to non-experts, and can be seen as a useful approach: however, a loss of biomass may be more consequential than local declines in species diversity, as common insect species contribute the most to ecosystem services, such as pollination (see e.g. Kleijn et al., 2015). Additionally, at the ecosystem level, it is often temporally and economically unfeasible to monitor all insects (Montgomery et al., 2020) despite the need for more in-depth studies and a "one health" approach as we increasingly realize that -and how much -human and insect well-being intertwine . However, in wellknown taxa, there are a number of examples whereby biomass is not informative of other population and biodiversity measures. For example, Hallmann et al. (2020) found stronger declines in smaller, rarer ground beetle species, suggesting that declines in biomass might be restricted to some groups of organisms only at the community level. Likewise, Homburg et al. (2019) showed that significant declines in species richness and phylogenetic diversity are unrelated to biomass loss in carabid beetles, and Saint-Germain et al. (2007) showed that the strength of the correlation between biomass and abundance of coleopterans varied significantly between the chosen taxa and sampling method.
Recent studies also suggest that even measures of species richness alone are unlikely to be completely instructive of biodiversity changes. Turnover in species and the functional traits associated to those species over geographic and temporal gradients are generally considered to be more informative of the health and functioning of an ecosystem (Hillebrand et al., 2018). Species identity is crucial to understand the dynamics of a system, and in particular how insect-provided services may be threatened . For example, Heleno et al. (2009) found that insect abundance did not decrease with the increase of alien plant species, but biomass did, because large insects where replaced by smaller insects. Therefore, these studies and our results presented here are calling to consider biomass, abundance, species richness, species identity, functional diversity and phylogenetic diversity as necessary metrics that should be used as complementary tools in ecosystem management, and not individually.
An additionally important conclusion from our results on biodiversity metrics is the apparent interchangeability between the different phylogenetic trees used (i.e., taxonomic hierarchy based on Linnaean classification vs. COI tree vs. multi-gene tree) (Pearson's correlation coefficient ~0.95), suggesting that the performance of a taxonomic hierarchy based on Linnaean classification could be used as a proxy for the measurement of PD in wild bees. This is a result of pivotal importance, reflecting the fact that the present-day Linnaean classification of bees is largely consistent with multi-gene phylogenies, owing to (i) the consensus in the scientific community on the higher-rank classification of wild bees ever since Michener (2007), and (ii) the fact that the contemporary Linnaean classification of bees also evolved in recent decades to incorporate the results of advances in the systematics of wild bees based on molecular markers Schmidt et al., 2015). Therefore, the results we provide here demonstrate the operability of using trees based on the Linnaean Classification for the analysis of patterns in PD for bees, and they represent an impetus to analyse broad phylogenetic trends in regions with a well-known but as yet not sequenced wild bee fauna. This includes NW Europe and perhaps also other areas where a considerable share of known species still await sequencing of DNA barcodes and other nuclear genes.

Bee diversity metrics and pollination
Here, we illustrate the importance of measuring biodiversity metrics associated to species identity for wild bees. Biomass may be essential for understanding some ecosystem services, such as changes in the prevalence of insects as a food source (Morse, 1971), but because the functional role of many insect species is poorly understood, several other metrics are urgently needed (see e.g. Dorchin et al., 2018). Pollination services provided by wild bees are linked strongly to species identity and functional role, with only certain species important for crop pollination (Kleijn et al., 2015), but important crop pollinating species are likely to exhibit turnover along a geographic gradient . Wild bees are also vital for wildflower pollination, and these relationships can be highly specialized (Rasmussen et al., 2020). Even though overall hypotheses predict a global decline for wild bees, the diversity, distributions and population trends of many wild bees remain understudied globally (IPBES, 2016; Bartomeus et al., 2019;Ascher and Marshall, 2020), and there are also examples of areas that show a positive trend in wild bee diversity (Herrera, 2019). Consequently, measuring the biomass of wild bees alone cannot capture the subtleties of these patterns of community structure that are key both for pollination services and for conservation.

Conclusion
The ideal situation to measure insect decline requires continued monitoring of communities in a variety of habitats to provide long-term time-series data which can provide robust trends, which in turn can be directly linked to proposed drivers of decline (Didham et al., 2020). Our results suggest that these long-term series must also include relevant diversity metrics and cannot depend on one alone. Therefore, biomass has been rightly criticized as a singular measure of insect decline, particularly as the mainstream media has shown a propensity to emphasize the magnitude of these results, to the detriment of more nuanced studies (Saunders et al., 2020a). A greater responsibility is needed by scientists to present accurately the relevance of their methods and results to insect decline (Saunders et al., 2020b). Habel et al. (2019) argued that for a general consensus of decline of European insects more studies are needed that agree on long-term negative trends of diversity, abundance and biomass, and incongruence or absence of one or more of these measures makes generalizations difficult. We show that this is a feasible goal for wild bees and our results emphasize that the sampling design of a proposed pan-European monitoring scheme as part of the European Pollinator Initiative needs to ensure that all relevant metrics are measured. This is vital to help the EU and its Member States set strategic objectives, develop actions to address the decline of pollinators, and ultimately contribute to global conservation efforts. The effort required to formulate detailed baselines is immense (Karlsson et al., 2020), but meaningful trends in species functional role and turnover is vital to properly understand long-term trends in wild bees and other insects.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. to the DNF agents B. Allard, B. Vandoren, C. Gillot, C. Pierlot, F. Gruselin, P. Fery, P. Havart, P. Protin, P. Verté and S. Terweduwe for their dedicated contribution to the fieldwork. NJV and his team received financial support for this study from Bruxelles Environnement (BE/ IBGE) in the framework of the Brussels Bee Atlas (WildBnB), from the Walloon Region through a research grant delivered by the Direction générale opérationnelle de l'Agriculture, des Ressources naturelles et de l'Environnement (DGO3) for the "Modèle permaculturel" project on biodiversity in micro-farms, as well as from the FNRS/FWO joint programme "EOS -Excellence Of Science" for the project "CliPS: Climate change and its impact on Pollination Services (project 30947854)" and the LIFE 11NAT/BE/001060 LIFE Herbages project. NJV supervised the sampling, designed the study, led the writing of the manuscript and the statistical analyses; all authors contributed either to the collection, analysis or interpretation of data, and all authors reviewed, edited and approved the final version of the manuscript prior to submission. None of the funding sources had any role in the study design, in the collection, analysis and interpretation of data, in the writing of the report and in the decision to submit the article for publication.