Predicting changes in bee assemblages following state transitions in North American drylands

Drylands worldwide are experiencing ecosystem state transitions: the expansion of some ecosystem types at the expense of others. Bees in drylands are particularly abundant and diverse, with potential for large compositional differences and seasonal turnover across ecotones. To better understand how future ecosystem state transitions may influence bees, we compared bee assemblages and their seasonality among three dryland ecosystem types of the southwestern U.S. (Plains grassland, Chihuahuan Desert grassland, and Chihuahuan Desert shrubland). Using passive funnel traps, we caught bees during two-week intervals from March through October during 2002 – 2014. The resulting dataset included 302 bee species and >70,500 individuals. Bee abundance, composition, and diversity differed among ecosystems, indicating the potential for future ecosystem state transitions to alter bee assemblage composition in drylands. We also found strong seasonal turnover in bee species, suggesting that bee phenological shifts may accompany ecosystem state transitions. Common rather than rare species drove the observed trends, and both specialist and generalist bee species were indicators of each ecosystem type or month; these species could be informative sentinels of community-wide responses to future shifts. Our work suggests that predicting the consequences of global change for bee assemblages will require accounting for both within-year and among-ecosystem variation.


32
Drylands worldwide are experiencing ecosystem state transitions: the expansion of some 33 ecosystem types at the expense of others 1,2 . These transitions include encroachment of C3 shrubland into (Table 1). Transitions among these ecosystem types are predicted to occur under climate change, with replace Plains grassland [45][46][47] . In our study, the two Chihuahuan Desert sites were separated by ~2 km; assemblage-level patterns were driven by common or rare species, we ran all abundance, composition, 141 and diversity analyses (described below) on the full dataset, on a dataset with singleton bee species 142 (those caught only on a single transect, in a single month) removed, and finally on a subset of the dataset 143 containing only the bee species that were present in >5% of the samples.

144
Overview. Analyses addressed our two key questions within one set of statistical models 145 (described below). First, (1) How much do bee assemblage abundance, composition, and diversity differ 146 among major southwestern U.S. ecosystem types? was determined by the statistical significance and 147 magnitude of the effect of ecosystem type in our models. We also compared the effect size of ecosystem 148 type against the effect size of month of sampling to estimate the relative importance of inter-ecosystem 149 versus seasonal variability. Then, to address (2) Do dryland ecosystem types differ in their degree of 150 seasonal variation in bee abundance, composition, or diversity? we evaluated whether the interaction 151 between ecosystem type and month of sampling was statistically significant, indicating that ecosystems 152 differed in the seasonality of bee abundance, composition, or diversity. 7 sampling, and the random effect of transect, which was nested within ecosystem type to account for the 156 repeated measures design, using perMANOVA (version 1.0.3) with 9999 permutations of residuals under 157 a reduced model. We additionally examined whether ecosystem types or months differed in bee 158 assemblage dispersion using permDISP in Primer 53 . We visualized assemblage composition with non-159 metric multidimensional scaling analysis (NMDS) implemented with 500 restarts in Primer. For each 160 ecosystem type, we assessed bee species turnover among months, as well as the rate of community 161 change, using the <codyn> package in R 54 . Finally, to identify which taxa contributed most to bee

165
Bee diversity and abundance. For bee diversity, we calculated the Shannon diversity index (H'), 166 species richness, and evenness (Pielou's J) using the <vegan> package in R 56 . We then used linear 167 mixed effects models to examine the influences of ecosystem type, sampling month, and their interaction 168 (fixed effects), as well as transect identity (random effect nested within ecosystem type), on these three 169 responses, as well as on total bee abundance (function lmer, <lme4> package in R) 57 . When there was a 170 significant ecosystem type x sampling month interaction, we tested a priori contrasts for pairs of months 171 within each ecosystem type and for pairs of ecosystem types within each month using Tukey-Kramer 172 multiple comparisons in the <emmeans> package in R 58 .

175
The dataset 176 We captured a total of 70,951 individuals representing 302 species during the 13 years of 177 monthly trapping (see Supplementary Table S1 for a full species list). Species were distributed across 6 families and 56 genera (Supplementary Table S1 and Fig. S2). Our dataset was dominated by a small number of abundant species and contained a large number of rare species (Supplementary Fig. S3). The 180 most commonly collected species were Lasioglossum semicaeruleum (36% of all collected specimens), 8 and Eucera lycii (3%). Amongst the collected species, 30% were singletons, and 58% were found in <5% 183 of all samples.

184
Bee assemblage composition: temporal variation surpassed differences among dryland ecosystem types 185 Variation among ecosystems. All ecosystems significantly diverged in bee assemblage 186 composition, and this pattern was present during all months (Table 2

192
Indicators of variation among ecosystems. We identified 43 bee species as ecosystem indicators 193 according to their Dufrene-Legendre (DL) indicator species values (Table 3)

203
Temporal variation. The month of sample collection explained an order of magnitude more 204 variation in bee assemblage composition than did ecosystem type (Table 2, Fig. 3). Generally, 205 assemblages diverged between the early and late months of the year and converged during the middle of 206 the summer. Across ecosystems, the pair of months most divergent in bee composition was March versus 9 Months additionally differed from one another in the magnitude of assemblage dispersion, a 210 metric that captures the degree of beta-diversity across both sites and transects (Table 2, Fig. 3). The 211 strongest differences in beta-diversity were between March or June, which had the smallest multivariate 212 dispersions (mean ± s.e., March: 21.0 ± 1.5, June: 20.4 ± 0.8), against October, which had the largest 213 average dispersion across ecosystems (29.8 ± 1.8).

214
Bee abundance and diversity: temporal variation exceeded variation among dryland ecosystems another in total bee abundance (Table 4). Bee abundance was on average 43% lower in the Chihuahuan

217
Desert shrubland relative to the two grassland sites from March through July (Fig. 4a). However, 218 abundances within the ecosystems converged in August, and abundance differences disappeared in 219 September and October (Fig. 4a), as indicated by a significant interaction between ecosystem type and 220 month of collection (Table 4: Ecosystem x Month, P < 0.0001).

221
Diversity. Ecosystems also diverged in bee diversity as measured by the Shannon index and sampling, the Chihuahuan Desert shrubland ecosystem had the highest bee Shannon diversity and 225 evenness, with these diversity metrics 5% (Shannon diversity) and 2% (evenness) higher than in the 226 Plains grassland. In turn, Plains grassland diversity metrics were 16% (Shannon diversity) and 12% 227 (evenness) higher than the Chihuahuan Desert grassland. In contrast, on average across months, the 228 ecosystems did not significantly differ in bee species richness (Table 4, Table 4: Ecosystem x Month -Shannon diversity: P < 0.0001, richness: P = 0.0137, evenness: P 231 < 0.0001), indicating that dryland ecosystem types differed in their degree of seasonal variation in bee 232 diversity (Question 2). Specifically, Shannon diversity was greater in the Chihuahuan Desert shrubland 233 than in the Desert grassland in all months except for March; differences in Shannon diversity were largest 234 in May and September, when Shannon diversity respectively was 38% and 33% higher in the The largest difference occurred in October, in which Shannon diversity was 31% higher in the 238 Chihuahuan Desert shrubland than Plains grassland. However, this trend was reversed in both March and September, with greater Shannon diversity in the Plains relative to Chihuahuan Desert grassland in 242 all of these months (Fig. 4b).  Abundance. Like species composition, total bee abundance also varied seasonally across the 261 three ecosystem types (Table 4), and ecosystem types exhibited differing trends in total abundance over 262 the course of the season (Fig. 4a). In the Chihuahuan Desert grassland, bee abundance increased from 11 by a ~50% decrease in abundance between April and May (df = 84, t = 11.12, P < 0.0001) and a 66% 266 increase in abundance between May and June (df = 84, t = -6.99, P < 0.0001). Between July and August, 267 while bee abundance decreased ~30% within both the Chihuahuan Desert grassland (df = 84, t = 4.62, P 268 = 0.0004) and Plains grassland ecosystems (df = 84, t = 5.93, P < 0.0001), it increased by 40% within the 269 Chihuahuan Desert shrubland (df = 84, t = -3.86, P = 0.0053) (Fig. 4a). Across ecosystem types, bee 270 abundances were generally lower in September and October relative to all other months (Fig. 4a).

271
Diversity. Within each ecosystem, most months had similar levels of species richness, with some 272 exceptions (Fig. 4c). Notably, there was a sharp decline in richness between August and October across 273 all three ecosystems (Fig. 4c). During this period, richness declined by 70% within the Chihuahuan Desert 274 grassland (df = 84, t = 11.18, P < 0.0001) and by 60% within both the Plains grassland (df = 84, t = 10.92, Shannon diversity and evenness diverged among ecosystems (Fig. 4b,d). Patterns in total abundance,

277
Shannon diversity, richness, and evenness were all largely driven by common rather than rare species 278 (see Supplementary Fig. S5).

286
In contrast, certain bee taxa were only characteristic of a given month within one or two  (Supplementary Tables S3-S6). Finally, 295 25 species were indicators of a particular month in one or two ecosystem types, and were then indicators 296 of a different month, often the following one, in the other ecosystem(s) (Supplementary Table S3).

299
We found large variation in bee assemblages and their seasonality among three dryland 300 ecosystem types of the southwestern U.S. These results indicate the potential for future ecosystem state 301 transitions to alter bee assemblage composition in drylands. Overall, ecosystem types in our study had 302 similar levels of bee species richness but differed from one another in species evenness and composition.

303
These results imply that state transitions could alter the presence/absence and relative abundances of 304 bee species in our system, bringing about substantial assemblage reordering.

305
Our data suggest that the most probable state transitions in our southwestern U.S. drylands -306 shrub encroachment into Desert grassland and Desert grassland encroachment into Plains grassland 59,60 307 -could substantially reshape bee communities. For the grass-to-shrub transition, our results suggest that, 308 on average over the season, total bee abundance will decrease while richness, Shannon diversity, and 309 Pielou's evenness will increase. In contrast, our findings predict that richness, Shannon diversity, and 310 evenness will decrease while total abundance will remain relatively unchanged if Desert grassland 311 replaces Plains grassland. The simultaneous occurrence of these state transitions could therefore 312 substantially alter the distribution of bees and their ecosystem services across the landscape.

398
In our dataset, common rather than rare bee species drove the trends in abundance, diversity, 399 and composition that we observed over both space and time. Considering these species' ecologies may

442
and thus susceptibility to phenological shifts in response to increasing temperature.

443
We also identified bee species that were characteristic of particular times of year across 444 ecosystem types. Monitoring these species could enable the detection of broad, cross-site phenological 445 shifts that may occur in the future. For instance, Osmia, Anthophora, and Diadasia species may be 446 monitored to consider shifting pre-monsoon bee phenology, and Perdita species may be used to study 447 shifts in emergence timed with monsoon rains. Future publications using these data will investigate inter-448 annual bee assemblage differences and relationships with climate variables.

449
Finally, we identified bee species that were characteristic of particular months only in specific 450 ecosystems. Among these, certain specialist bee species may be candidates for detecting important 451 phenological shifts within ecosystems, identifying phenological differences among ecosystems, and