DOI: 10.6086/d1x10s
title: Diversity and turnover of wild bee and ornamental plant assemblages in commercial plant nurseries

description: In human-modified landscapes, such as agricultural areas, understanding
 how habitat characteristics influence the diversity and composition of
 beneficial organisms is critical both to conservation efforts and to
 modeling ecosystem services. Assessing turnover, or changes in species
 identity across sites or through time, is crucial to determine how shifts
 in community composition relate to changes in ecosystem services. For
 pollinators like wild bees, factors affecting their diversity have been
 well studied, but variables influencing temporal turnover, particularly
 across seasons within a year, remain relatively poorly understood. To
 investigate how local and landscape characteristics correlate with bee
 diversity and turnover across seasons, we recorded wild bee and flowering
 ornamental plant assemblages at 13 plant nurseries in southern California
 between spring and autumn over two years. Nurseries cultivate a broad
 diversity of flowering plant species that differ widely across sites and
 seasons, providing a field setting to test for correlations between
 turnover and diversity of plants and bees. As expected, we documented
 strong seasonal trends in wild bee diversity and composition. We found
 that local habitat factors, such as increased cultivation of native
 plants, were positively associated with bee diversity in sweep netting
 collections, whereas we detected moderate influences of landscape-level
 factors like proportion of surrounding natural area in passive trap
 collections. We also detected a moderate positive correlation between
 overall gains in plant species and gains in bee species at nurseries
 across consecutive seasons. Our results have important considerations for
 the conservation of wild bees in landscapes dominated by ornamental
 plants, and highlight the utility of plant nurseries for investigating
 hypotheses related to diversity and turnover in plant-pollinator systems.
description: Study sites and bee collection We collected
 wild (non-<em>Apis</em>) bees at 13 wholesale or retail plant
 nurseries (each &gt; 2 ha in area) in southern California, USA (Table
 S1, Fig. S1). We sampled each nursery from March 2016 – May 2018 once per
 season, except winter: once each in spring (March–May), summer
 (June–August), and autumn (September–November) (Table S2). We used two
 methods to collect bees. First, we deployed blue vane traps (SpringStar,
 Woodinville, WA) suspended from 1.2-m hooks in established sampling plots
 (see below) at a density of roughly one plot per hectare of nursery
 property (2-12 plots per nursery). We collected vane traps 72 h after
 deployment. Second, a single observer used a sweep net to collect bees
 visiting flowering plants in nurseries for 30-minute periods on each day
 we deployed and collected vane traps, spending no longer than 5 minutes
 collecting from a given plant species. Non-native honey bees,
 <em>Apis mellifera</em>, were abundant on many plants at these
 nurseries and were intentionally not sweep netted, but were sometimes
 collected in the vane traps. We identified all wild bees to species (or
 morphotaxon, if they could not be confidently assigned to species) using
 keys (Michener 2007, Ascher and Pickering 2020), reference collections,
 and assistance from bee specialists at the University of California
 Riverside Entomology Research Museum. Local and
 landscape features We surveyed flowering plants that
 were currently blooming, regardless of floral abundance, within a 15-m
 radius (~700 m<sup>2</sup>) plot around each vane trap during
 each nursery visit. We designated each plant species as native or
 non-native to California referencing the CalFlora database (calflora.org).
 We further classified five nurseries as “native nurseries” because they
 (1) advertised themselves as cultivating mostly plants native to
 California, and (2) on average, over two-thirds of plant taxa documented
 in plots were California-native. We classified the eight other nurseries
 as “conventional nurseries”. We estimated floral cover as percent cover
 within plots by dividing the plot into quadrants and assigning a rank to
 each quadrant using an ordinal scale from 0 to 4: ‘0’ = no flowers of any
 species, ‘1’ = few flowers with sparse cover, ‘2’ = more than ‘1’ but
 covering less than 50% of quadrant, 3 = &gt; 50% cover but including
 large patches with no flowers, and ‘4’ = near 100% cover, then averaging
 these values for the plot. Using QGIS (QGIS Development Team 2020) and the
 2016 National Land Cover Database (mrlc.gov), we calculated the proportion
 of natural area inside 500-m and 1-km buffers around each nursery’s
 perimeter. We treated land cover classes 41–43, 52, 71, 90, and 95 as
 “natural” landscape, excluding classes denoting water and developed,
 agricultural, or barren land (11, 21–24, 31, 81, and 92). We also recorded
 from local weather stations the daily high temperature for each of the
 four calendar days that vane traps were deployed, which we then
 averaged. Statistical analyses We
 conducted all statistical analyses in R version 3.3.3 (R Core Team 2021).
 All means are reported ± SEM.<em> </em>The dataset represented
 58 sampling events for blue vane traps (of which 21 were at native
 nurseries) and 80 sampling events for sweep netting (32 at native
 nurseries). For seven of these sweep netting sampling events, no wild bees
 were detected. We excluded all <em>Apis mellifera</em> (633
 specimens) from the vane trap dataset so as to focus on patterns in wild
 bee assemblages. Comparisons of analyses with the full dataset and the
 dataset excluding honey bees revealed that their removal had minimal to no
 influence on the outcomes of our statistical models (Table S3).
 We performed individual-based rarefaction
 using <em>vegan</em> (Oksanen et al. 2018) and
 <em>iNEXT</em> (Hsieh et al. 2016) to examine sampling
 completeness of wild bees across sites, nursery types, and collection
 methods. For each sweep netting and vane trap sampling event, we
 calculated bee Shannon diversity and evenness. For each survey plot, we
 calculated the richness of blooming plants and the proportion of richness
 represented by native plant species, then averaged across plots within a
 given sampling event. In statistical models, we substituted proportional
 native plant richness with the residuals of this variable obtained from a
 regression against nursery type, as these predictors were collinear
 (variance inflation factor &gt; 3) (Graham 2003). The same was also
 done with average daily high temperature in blue vane trap models, which
 was collinear with season. We constructed four separate linear mixed
 models (LMMs) using <em>lme4</em> (Bates et al. 2015) for each
 of bee Shannon diversity and evenness and each of two collection methods.
 We included as fixed effects season and nursery type (and their
 interaction), blooming plant richness, proportional native plant richness,
 percent floral cover, daily high temperature, nursery area, and
 proportional surrounding natural area at 1 km. Models in which natural
 area was quantified at 500 m were not significantly better than those at 1
 km (all ΔAIC &lt; 2.0), thus we only report models using the 1-km
 buffer data. Nursery, study year, and date of sampling served as random
 effects in LMMs. For LMMs, we confirmed lack of multicollinearity (VIF
 &lt; 3) using ‘vif’ (<em>car</em>) (Fox and Weisberg
 2011). Using the function ‘emmeans’ in the package <em>emmeans
 </em>(Lenth 2019), we conducted post-hoc Tukey’s tests to compare
 treatment means for statistically significant factors and corrected for
 multiple comparisons. We next examined variation in the
 composition of plant and wild bee assemblages in nurseries. For bees, we
 conducted these analyses only on vane trap samples, as sweep net samples
 were not large enough to sufficiently represent assemblage composition. We
 tested for differences in the multivariate centroids of samples using
 permutational multivariate analyses of variance (permANOVA) on
 dissimilarity matrices using ‘adonis’ (<em>vegan</em>). For
 bees, we used Bray-Curtis dissimilarity as it accounts for the relative
 abundances of species, whereas for plants, we used Jaccard dissimilarity
 as plant species were recorded in plots as present or absent. We included
 the aforementioned fixed effects and stratified permutations (10,000) by
 nursery. We also examined beta diversity by analyzing multivariate
 homogeneity in dispersion of plant and wild bee assemblages around their
 centroids using ‘betadisper’ (<em>vegan</em>) between nursery
 types, seasons, and years. We performed Mantel tests using ‘mantel’
 (<em>vegan</em>) to explore correlations between dissimilarity
 matrices of plant and bee (vane trap samples) assemblages across nurseries
 separately in spring, summer, and autumn. We used
 function ‘RAC_change’ in the package <em>codyn</em> (Hallett
 et al. 2016) to calculate relative changes in plant and bee richness
 through time (Kazenel et al. 2020), pooling bees from both collection
 methods. Temporal species turnover can be partitioned into gains and
 losses (Faleiro et al. 2018). Function ‘RAC_change’ quantifies changes in
 richness between a pair of time periods as the number of species gained or
 lost divided by the total unique species present in both time periods. We
 calculated these changes for each nursery between consecutive seasons in
 which samplings occurred (spring to summer, summer to autumn, and autumn
 to the following spring, over two years), totaling six inter-season
 comparisons. To determine if turnover in bee species was related to that
 of plants, we constructed LMMs with changes in bee richness (as gains or
 losses) as dependent variables, and the corresponding plant turnover
 metric, nursery type, category of inter-season comparison, and their
 interaction as fixed effects. Finally, to determine if
 any bee species in our sweep netting or vane trap collections were
 significantly associated with local or landscape factors, we conducted a
 multi-level pattern indicator species analysis using function ‘multipatt’
 (<em>indicspecies</em>) (De Cáceres and Legendre 2009). This
 analysis identifies species whose occurrence is associated with groups of
 sites categorized by a factor. We conducted analyses for season, nursery
 type, proportional native plant richness (which we dichotomized as
 &gt; or ≤ 50% native), floral cover (&gt; or ≤ the average percent
 cover of all plots), and proportional natural area within 1 km of the
 nursery (&gt; or ≤ 50% natural area).
description: Descriptions of all sheets and variables in file are given in the
 first sheet ("metadata").
description: # Diversity and turnover of wild bee and ornamental plant assemblages in
 commercial plant nurseries | sheet
 name | variable name | type of variable | possible values | description of
 variable | | :--------------- | :------------ | :--------------------- |
 :-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- | :-------------------------------------------------------------------------------------------- | | bees | year | factor | 2016, 2017, 2018 | calendar year in which the bee was collected | | season | factor | spring, summer, autumn | season in which the  bee was collected | | | nursery.code | factor | various | the nursery at which the bee was collected | | | nursery.nativity | factor | conventional, native | the type of nursery | | | trap.method | factor | bvt, netting | whether the bee was captured via blue vane traps ("bvt") or sweep netting | | | species.code | factor | various | unique code denoting the species of bee | | | remove | factor | \[blank], yes | whether to remove the specimen ("yes") prior to running statistical models on the data; these 11 bees were from blue vane traps that were unintentionally disturbed (knocked over) while deployed, and are thus inappropriate to include in quantitative models | | | plants | year | factor | 2016, 2017, 2018 | calendar year in which the plant was recorded | | season | factor | spring, summer, autumn | season in which the plant was recorded | | | nursery.code | factor | various | the nursery at which the plant was recorded | | | plot | factor | 1 to 12 | the sampling plot at the nursery | | | taxon | factor | various | the plant taxon | | | net\_samplings | year | factor | 2016, 2017, 2018 | calendar year in which the bee was collected | | season | factor | spring, summer, autumn | season in which the  bee was collected | | | date.collected | factor | various | calendar date on which bee was collected | | | nursery.code | factor | various | the nursery at which the bee was collected | | | nursery.nativity | factor | conventional, native | the type of nursery | | | temp | numeric | various | temperature (°F) | | | plant.richness | numeric | various | average plot-level blooming plant richness for that sampling event | | | floral.cover | numeric | 0 to 100 | average plot-level % floral cover for that sampling event | | | nat.area.1km | numeric | 0 to 1 | proportion of natural area surrounding nursery at 1-km buffer | | | nursery.area | numeric | various | area of nursery, in hectares (ha) | | | nat.resids | numeric | various | average plot-level proportion of blooming plant richness represented by native taxa (residuals) | | | shannon.index | numeric | various | Shannon diversity index (H') for sweep netted bees | | | evenness | numeric | various | Pielou's evenness index (J') for sweep netted bees | | | bvt\_samplings | year | factor | 2016, 2017, 2018 | calendar year in which the bee was collected | | season | factor | spring, summer, autumn | season in which the  bee was collected | | | nursery.code | factor | various | the nursery at which the bee was collected | | | nursery.nativity | factor | conventional, native | the type of nursery | | | retrieve.date | factor | various | calendar date on which vane trap was collected | | | plant.richness | numeric | various | average plot-level blooming plant richness for that sampling event | | | floral.cover | numeric | 0 to 100 | average plot-level % floral cover for that sampling event | | | nat.area.1km | numeric | 0 to 1 | proportion of natural area surrounding nursery at 1-km buffer | | | nursery.area | numeric | various | area of nursery, in hectares (ha) | | | nat.resids | numeric | various | averate plot-level proportion of blooming plant richness represented by native taxa (residuals) | | | temp.resids | numeric | various | temperature (residuals) | | | shannon.div | numeric | various | Shannon diversity index (H') for vane trapped bees | | | evenness | numeric | various | Pielou's evenness index (J') for vane trapped bees | | | turnover | time.numeric | factor | 1 to 7 | earlier of the two seasons being comparied (1=spring 2016, 2=summer 2016, \[…] 7=spring 2018) | | time.numeric2 | factor | 1 to 7 | the latter of two seasons being compared | | | nursery.code | factor | various | the nursery at which seasons are being compared | | | gains\_bees | numeric | various | proportional gains in bee species between time periods | | | losses\_bees | numeric | various | proportional losses in bee species between time periods | | | year | factor | 2016, 2017 | calendar year of earlier time period | | | season | factor | spring, summer, autumn | season of earlier time period | | | nursery.nativity | factor | conventional, native | the type of nursery | | | gains\_plants | numeric | various | proportional gains in blooming plant species between time periods | | | losses\_plants | numeric | various | proportional losses in blooming plant species between time periods | | ## Code/Software This is an optional, freeform section for describing any code in your submission and the software used to run it. Describe any scripts, code, or notebooks (e.g., R, Python, Mathematica, MatLab) as well as the software versions (including loaded packages) that you used to run those files. If your repository contains more than one file whose relationship to other scripts is not obvious, provide information about the workflow that you used to run those scripts and notebooks.