Increasing synchrony opposes stabilizing effects of species richness on terrestrial communities

Ecological theory has predicted that species richness should stabilize communities, with mechanisms including species synchrony and population variability determining the net impacts. While these theories have been supported empirically, results can be sensitive to taxonomic bias as studies are often focussed on plants. Trophic differences between consumers and primary producers can lead to varying stabilizing effects of species richness. Here, we compared the impact of species richness on community variability in four taxonomic groups: terrestrial birds, mammals, invertebrates and plants.


| INTRODUC TI ON
Community variability, the inverse measure of community stability, is often measured using aggregated abundance or biomass of all species through time (Olivier et al., 2020;. Understanding drivers of community variability is crucial for predicting ecosystem function provision (Ives & Carpenter, 2007;Renard & Tilman, 2019). Theories predict that community variability should be reduced in species-rich communities relative to species-poor communities (Doak et al., 1998;Ives & Carpenter, 2007;Tilman, 1999). Two central mechanisms have been hypothesized as the mechanistic basis of the relationship between community variability and species richness: (1) species synchrony and (2) population variability (Downing et al., 2014;Thibaut & Connolly, 2013).
Population variability represents the average variability of the abundance of each species within a community and is expected to positively correlate with community variability (Thibaut & Connolly, 2013). Theoretical studies indicate that species richness should have a positive effect on population variability within competitive communities (Lehman & Tilman, 2000;Tilman, 1999).
Species synchrony, which captures the extent of correlations in the abundance of different species within communities over time, is known to increase community variability through weakened compensatory dynamics (Thibaut & Connolly, 2013). Species synchrony is expected to exhibit a negative relationship with species richness as species have variable responses to environmental conditions, and interspecific competition can be stronger in species-rich communities (Loreau & de Mazancourt, 2013;Tredennick et al., 2017).
Most studies examining relationships between species richness and community variability have focussed on simplified ecosystems with single trophic levels, typically primary producers such as algae and vascular plants (Campbell et al., 2011;Gross et al., 2014;Houlahan et al., 2018;Jiang & Pu, 2009;Xu et al., 2021). Whether findings from these studies apply to consumers are unclear, but biological differences between consumers and primary producers have been suggested to alter the impacts of species richness on ecosystem properties (Duffy et al., 2007;Ives et al., 2005), with support from meta-analyses finding that effects of species richness on ecosystem functions varied across trophic levels Griffin et al., 2013;Hooper et al., 2012;Lefcheck et al., 2015;O'Connor et al., 2017;Soliveres et al., 2016;Srivastava et al., 2009).
While primary producers rely on common resources such as light, consumers have more dissimilar requirements among species (Bengtsson et al., 1994;Duffy et al., 2007;Ives et al., 2005). This should lead to more differential responses to environmental conditions and, therefore a more negative species richness-synchrony relationship in consumers (Loreau & de Mazancourt, 2013). Alternatively, more dissimilar resource requirements can reduce interspecific competition, leading to a less negative species richness-synchrony relationship in consumers (Loreau & de Mazancourt, 2013). However, the latter is less likely, as differential responses to environmental conditions are known to be more important than interspecific competition in regulating synchrony within natural communities (Mutshinda et al., 2009;Thibaut et al., 2012;Tredennick et al., 2017).
Reduced interspecific competition among consumers can result in a weaker increase in population variability along species richness gradients (Lehman & Tilman, 2000;Loreau & de Mazancourt, 2013;Tilman, 1999). Overall, net species richness effects on community variability should be more stabilizing in consumer communities, as reduced population variability and synchrony decrease community variability.
Meta-analyses have provided support that species richness stabilizes communities through reducing synchrony within communities, while effects of species richness on population variability can strongly vary across studies (Campbell et al., 2011;Gross et al., 2014;Houlahan et al., 2018;Jiang & Pu, 2009;Xu et al., 2021).
These datasets, however, are strongly biased towards or only involve primary producer communities, so the conclusions may not apply to consumer communities. Pinpointing the mechanisms that generate results variability across studies is hindered by a lack of studies focussed on other organisms or bias in study design across ecosystems. For example, in a recent meta-analysis (Xu et al., 2021), only two multitrophic studies were available to evaluate the generality of species richness-synchrony relationships. For other relationships, the effects of trophic complexity were confounded by ecosystem type (aquatic/terrestrial) and investigational approach (observational/experimental).
With the increasing availability of time series data globally (Dornelas et al., 2018), general patterns between species richness and community variability at large spatial scales and across multiple taxa can be better identified. We conducted a cross-taxon analysis on the open-access BioTIME database, including time series data on birds, mammals, invertebrates and plants across the globe (Dornelas et al., 2018) to examine drivers of community variability for these taxa. A recent study using BioTIME detected weakly negative relationships between species richness and community variability but did not examine taxonomic differences in detail (Dallas & Kramer, 2022). Here we used the BioTIME dataset to assess whether the community variability and species richness relationship of primary producers (i.e. plants) can be generalized to different consumer (i.e. animal) taxa while eliminating confounding factors that typically exist in currently available meta-analyses (e.g. Xu et al., 2021).
Based on most previous theoretical (Ives & Carpenter, 2007;Lehman & Tilman, 2000;Loreau & de Mazancourt, 2013;Thibaut & Connolly, 2013) and empirical findings (Xu et al., 2021), we hypothesized that species richness would reduce synchrony and community variability but increase population variability of different taxa. As consumers should have more dissimilar resource requirements compared with plants (Duffy et al., 2007;Ives et al., 2005), we expect species richness to have less positive impacts on population variability in animal communities due to weakened interspecific competition (Loreau & de Mazancourt, 2013).
Effects on synchrony should be more negative, driven by more differential responses to environmental conditions (Loreau & de Mazancourt, 2013). Ultimately this should lead to a more stabilizing effect of species richness on the community variability of consumers.

| ME THODS
To examine the drivers of community variability of different taxa, we utilized the open-access BioTIME dataset. The dataset includes species abundance data from marine, freshwater and terrestrial realms based on consistent methodology within each time series (Dornelas et al., 2018). We analysed time series that recorded count or density data and spanned ≥5 years to ensure temporal variability could be reliably quantified. Hereafter we collectively term these data as abundance data. Additionally, we only included taxa with ≥10 studies for all analyses. This resulted in only terrestrial time series being analysed, as there were insufficient studies for marine and freshwater primary producer communities, making comparisons with consumers from those environments difficult. An additional advantage to analysing terrestrial time series is that the environmental history of sites can be obtained through global land cover datasets. This allows for the controlling of effects of longterm directional changes in community composition generated by other factors such as land-use change, which is recommended for studies examining the effects of biotic variables on community variability Valerio et al., 2022). We only included communities that had <10% changes in agricultural, urban and forest land within the study period at the landscape scale, defined here as 96 km 2 grids the sites belong to within a global hexagonal system (Barnes & Sahr, 2017).
This definition follows a previous study showing the impacts of deforestation on terrestrial communities based on BioTIME at this scale (Daskalova et al., 2020). To minimize temporal mismatches with land cover data (see below), we only included data between 1992 and 2015. Furthermore, we excluded data if the grain size of the communities were larger than the buffer (96 km 2 ) used to extract agricultural, urban and forest land cover. Our protocol yielded data for four taxonomic groups (birds, mammals, terrestrial invertebrates and plants). Taxa classifications followed the BioTIME database.
We defined sampling strategies as methods to collect community data in each survey. If two communities had the same sampling strategy but were sampled in different years or had different sampling efforts, we still considered them as having the same strategies. For example, in a butterfly study with six sites, sampled through walking a 700 m transect, the start of the sampling period differed across sites, and each year had a different sampling effort (Bonebrake et al., 2016). We examined the data and pooled some studies according to this definition (Table S1).
Communities within the same study were sampled in different years in many studies. Therefore pooling data across communities/ locations to form aggregated time series would introduce temporal variability in community metrics due to differences in sampling years. We then quantified temporal variability for each community, defined as the smallest sampling units based on sampling descriptions of each study, such that we avoided pooling data across studies)-29 out of 50 studies contained annually surveyed data. By contrast, others were sampled irregularly with a mean sampling interval of 1.5 years (see Table S1 and Figure S1 for a list of studies included, their methodology and locations).

| AB IOTIC VARIAB LE S
We obtained land-use and climatic data at each location based on available spatial coordinates to control for environmental variation between sites to avoid these variables confounding effects of species richness, synchrony and population variability on community variability. Global land cover data between 1992 and 2015 at 300 m resolution were obtained from ESA CCI to obtain the extent of agricultural, urban and forest cover for each location in different survey years (https://www.esa-landc over-cci.org/?q=node/175). As some spatial coordinates were approximate, we followed Daskalova et al. (2020)'s methodology. We used a 96 km 2 hexagonal grid system to obtain the agricultural and urban land cover proportion for each community, although we did not combine community data within the same hexagon. We then calculated the proportion of natural habitats as 1 − Proportion Agriculture − Proportion Urban in each sampling year, before averaging them for each community. Some datasets provided a single spatial coordinate for communities from a large region (≥96 km 2 ), and we obtained more accurate coordinates for each sample by examining associated publications or sources (e.g. official websites of LTER sites). If not possible, community data aggregated across an area larger than 96 km 2 were excluded. We used the 0.5° resolution climatic data from Jiang et al. (2017) since it provides inter-and intra-annual climatic variability based on climatic data from 1901 to 2012. We obtained the annual mean of temperature and precipitation, as well as their intra-and inter-annual variations at each location.
Because of the coarse resolution, the land-use and climatic variables included should be considered as reflecting the landscape and regional factors but not local environmental conditions. Communities within the same hexagonal grid or sharing the exact coordinates would have identical regional abiotic conditions. Since these community metrics were quantified at the smallest sampling units, any signal from landscape or regional abiotic conditions would indicate their effects in regulating community variability or its associated mechanisms at more local scales.

| COMMUNIT Y ME TRIC S
For communities in each study, we first pooled different samples within each year. Thus, all metrics captured inter-annual but not intra-annual signals. As some communities had varying sampling efforts between years, each year's total count or species density was divided by the number of surveys within each sampling year before calculating any community variables. The number of surveys was quantified as the number of unique sampling dates each year. In all studies, only one taxonomic group was studied, so community metrics were never averaged across taxonomic groups. For each study, we obtained all community variables at the smallest sampling levels to minimize pooling across sampling units with different survey years. Note that some communities had different values for community metrics but identical values for abiotic drivers, as some studies provided only one set of spatial coordinates for all communities.
We used the coefficient of variation (CV) to quantify temporal variability, an inverse measure of stability, at community and population levels. We chose to use CV as it facilitates comparisons with previous studies. Increasing CV indicates higher variability, which ranges from 0 to ∞. Community variability was quantified as the CV of annual community abundance, while population variability was quantified as averaging CVs of the annual abundance of each species. For synchrony, we used Gross's synchrony since it is mathematically independent of species richness, and thus, the results would not be affected by statistical artefacts Gross et al., 2014). Gross's synchrony ranges from −1 to 1, with increasing values indicating more synchronous communities. Gross's synchrony is undefined when only one species is detected throughout the time series or if at least one species has the same abundance throughout all years. Therefore 272 communities (~4%) were excluded from the analyses. Common species weighed more heavily than rare species in synchrony and population variability metrics, as we weighted species using their total abundance per sample across sampling years (Table S2).
We also obtained average species richness across the years of each community. Sampling effort (i.e. number of sampling dates) could vary across years in the same time series; thus we constructed rarefaction curves to obtain a standardized measure of species richness in each year before averaging. Specifically, we used the species data from the smallest sampling unit of each year to build samplebased rarefaction curves. We obtained rarefied species richness for each year based on the lowest sampling effort along the time series.
The rarefied species richness of each year was then averaged and used in subsequent analyses. Matuschek et al., 2017). Nevertheless, this also meant that abiotic and methodological variables may have inflated type-I and type-II errors (Harrison et al., 2018). Therefore we only focussed on interpreting the results of community variables. No random slopes were correlated with random intercepts to facilitate model convergence (Matuschek et al., 2017). For all individual models, we started with the maximal model, meaning that all random intercepts and random slope of community metrics were included. In cases where there were convergence or overfitting (as indicated by singular fit) issues, we removed random effects with zero variances and refitted the model again (Barr et al., 2013). If the issue persisted, we removed the weakest random effect sequentially until no issues emerged.
Individual models' final random effect structures are presented in Table S4. If the directed separation test indicated a missing relationship between community synchrony and population variability, we added a correlated error among them, as the causality of this relationship was unexpected from our model. To assess richness's net effects on community variability, we built a linear mixed model with richness, abiotic and methodological variables as predictors for each taxon. The selection of random effect followed the protocols above.  (Legendre & Legendre, 2012). Individual linear mixed models were built using R package lme4 (Bates et al., 2015), while piecewiseSEM (Lefcheck, 2016) was used to build structural equation models. We adjusted p-values from the models using small-sample (CR2) corrections to address heteroskedasticity using clubSandwich (Pustejovsky, 2022). iNEXT was used to build rarefaction curves (Hsieh et al., 2016).

| RE SULTS
In all structural equation models, species synchrony and population variability increased community variability, with both predictors showing the highest standardized path coefficient and, therefore relative importance in regulating community variability By contrast, a weaker negative relationship with marginal significance (p = .081, Figure 5d) was detected in plants. The magnitude of standardized path coefficients was consistently higher in richnesssynchrony relationships than in richness-population variability relationships, although the differences were minor in all taxa (0%-11%,

| DISCUSS ION
By comparing community variability patterns across taxa, our results demonstrate that species richness can drive destabilizing and stabilizing mechanisms, as it increases synchrony but reduces population variability. Other studies have also found that species richnesssynchrony and species richness-population variability relationships have opposite effects on community variability (Xu et al., 2021). In many cases, the net stabilizing effects of richness are driven by reducing synchrony rather than population variability (Xu et al., 2021).
Thus, our results provide empirical support for how stabilizing effects of species richness (if any) can be facilitated by reducing population variability, which can be countered by increasing synchrony, and potentially explain the weak net impacts of species richness on natural communities (Houlahan et al., 2018). We also identified inconsistencies among taxa in the importance of synchrony relative to Note that identical colour between panels represent different studies. Clusters and data points from the same study have the same colour. All icons are from public domain license and are available at Phylopic 2.0 (https://beta.phylo pic.org/). population variability, as reflected in their differences in standardized path coefficients. Particularly, the importance of synchrony was higher or comparable to population variability in all animal communities but not plant communities. These results indicate that trophic differences can alter how community variability is regulated (Ives et al., 2005;Thébault & Loreau, 2005). Our use of terrestrial and observational data generates stronger inferences about the effects of trophic differences, as in previous studies, their results could be confounded by investigational approaches and ecosystem types (Xu et al., 2021).  (Downing et al., 2014). Biases towards temperate regions in our dataset can also reduce interspecific competition's importance (Paquette & Hargreaves, 2021), leading to more positive richness-synchrony relationships. The results can also be explained by temporal scale, as increasing length of time series shifts temporal dynamics from synchronous to asynchronous (Downing et al., 2008). As synchrony tends to increase community variability, our results suggest the richness-synchrony relationship is unlikely to be a general mechanism stabilizing communities-instead, it can destabilize some communities. Different colours within each panel present different studies. Note that identical colour between panels represent different studies. Clusters and data points from the same study have the same colour. All icons are from a public domain license and are available at Phylopic 2.0 (https://beta.phylo pic.org/). with theoretical models (Lehman & Tilman, 2000;Loreau & de Mazancourt, 2013;Tilman, 1999) and experimental results (Xu et al., 2021). One explanation is that findings from most experimental communities or theoretical models are based on single trophic levels, while in natural communities, species are also experiencing trophic interactions, which can stabilize communities under field conditions (Brose et al., 2006;Jiang & Pu, 2009;McCann, 2000). Particularly, if species richness increases weak interactions in food webs, population variability should be reduced (Downing et al., 2014). The results do not support the hypothesis that terrestrial environments should be characterized by positive species richness-population variability relationships (Xu et al., 2021). The hypothesis is based on terrestrial environments lacking strong trophic interactions that destabilize communities (Shurin et al., 2002;Shurin et al., 2006), such that the stabilizing impacts of species richness through increasing weak trophic interactions would be less pronounced. However, recent studies have shown that trophic interaction strength strongly varies across terrestrial environments (Zvereva & Kozlov, 2022), suggesting the hypothesis may only apply to specific systems.
The net effect of species richness on community variability in each taxon depends on the relative importance of synchrony, population variability and effects of species richness unmediated through these two components. Theoretically, the net effect of species richness should also depend on whether species richness exhibits a stronger effect on synchrony or population variability (Xu et al., 2021). However, we found that the magnitude of standardized path coefficients of richness-synchrony and richness-population variability relationship is comparable in all taxa; thus, they are not the primary driver explaining the variable net impacts across taxa. Note that identical colours between panels represent different studies. Clusters and data points from the same study have the same colour. All icons are from a public domain license and are available at Phylopic 2.0 (https://beta.phylo pic.org/). variability models (Eschenbrenner & Thébault, 2022;Thibaut & Connolly, 2013).
Nevertheless, we note that the importance of synchrony was more similar or even higher than population variability in animal communities compared with plant communities. This finding is inconsistent with recent studies suggesting synchrony is less critical than population variability in animals compared with plants (Danet et al., 2021;Eschenbrenner & Thébault, 2022). Empirical studies on the relative importance of synchrony and population variability across taxa are limited, with mixed results on the importance of synchrony relative to population variability across trophic levels (Eschenbrenner & Thébault, 2022;Lamy et al., 2020;Olivier et al., 2020;Siqueira et al., 2022). The pattern shown here is unlikely explained by animal communities having multitrophic structures, as this should promote the importance of population variability over synchrony (Eschenbrenner & Thébault, 2022). Additionally, our invertebrate dataset is dominated by herbivores such as butterflies and grasshoppers, therefore, single trophic communities. The differential sample size is also unlikely to be the primary driver of the results, as invertebrates (N = 494) have a similar sample size to plants (N = 720), yet the importance of synchrony relative to population variability in invertebrates is more similar to birds (N = 3389) and mammals (N = 2160). Thus, we believe biological differences between animals and plants, such as more dissimilar resource requirements, can best explain our results (Bengtsson et al., 1994;Duffy et al., 2007;Ives et al., 2005).
The primary aim of our study was to evaluate potential differences in community variability patterns across taxonomic groups rather than a comprehensive test on the effects of diversity for different taxa. The latter requires analysing additional metrics that capture other dimensions of diversity (Olivier et al.,  All axes are at logarithm scale. The black lines indicate predictions based on fixed effects only, with solid, dash-dotted and dotted lines representing significant (p < .05), marginal (.05 ≤ p≤ .1) and insignificant results (p > .1), respectively. Coloured lines are effects of species richness in each lowest level cluster (Birds: Study; Mammals: Subblock; Invertebrates: Block; Plants: Subblock). Note that identical colours between panels represent different studies. Blocks, subblocks and data points from the same study have the same colour. All icons are from a public domain license and are available at Phylopic 2.0 (https://beta.phylo pic.org/). Valerio et al., 2022) and conducting analyses at multitrophic community levels, which is not possible with the BioTIME dataset because most time series are focussed on single taxonomic groups (Jiang & Pu, 2009;Xu et al., 2021). Particularly, asynchronous dynamics between resources and consumers (Downing et al., 2008) can be important in stabilizing natural communities, yet taxon-specific analyses may fail to capture these processes.
However, analysing stabilizing effects of species richness for different organisms in specific environments can still be useful in assessing generalities of different hypotheses related to community variability, particularly when used for comparisons with previous studies, which have typically focussed on species richness for individual taxonomic groups (Houlahan et al., 2018;Xu et al., 2021).
Additionally, taxon-specific analyses can be important in revealing community variability patterns for managing specific taxonomic groups, such as those that offer important ecosystem functions (van Klink et al., 2019).
Despite being suggested as a statistical inevitability (Doak et al., 1998), our results only detected an overall stabilizing effect of species richness on community variability in plants only. Multiple mechanisms exist for species richness to stabilize or destabilize communities (Downing et al., 2014;Thibaut & Connolly, 2013), and the varying importance of different mechanisms across taxonomic groups further complicates the net effects of species richness.
Thus, our multi-taxa analysis questions the generality of the relationship between species richness and community variability, given that empirical support is mostly from datasets based on or biased towards primary producers (Houlahan et al., 2018;Xu et al., 2021).
It also suggests that high richness, especially of consumer taxa, is unlikely to provide strong stabilizing effects in changing environments. Understanding how environmental conditions drive various stabilizing and destabilizing mechanisms in different organisms, especially synchrony and population variability, may provide more robust insights into community variability. This can help explain or predict the impacts of environmental changes on community variability of different taxa, which has been found to show inconsistent responses in other multi-taxa studies (e.g. Blüthgen et al., 2016, Olivier et al., 2020.

ACK N O WLE D G E M ENTS
We thank all data contributors and curators of the BioTIME database. TPNT was supported by a General Research Fund (17115119) award from the Hong Kong Research Grants Council to TCB and LCP, and the University of Toronto Scarborough Postdoctoral Fellowship Programme.

CO N FLI C T O F I NTER E S T S TATEM ENT
The authors declare no conflict of interest.

PE E R R E V I E W
The peer review history for this article is available at https:// www.webof scien ce.com/api/gatew ay/wos/peer-revie w/10.1111/ ddi.13696.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available in DataDryad, at http://doi.org/10.5061/dryad.9p8cz 8wk6. These data were derived from the following resources available in the public domain: BioTIME dataset (https://bioti me.st-andre ws.ac.uk/), ESA CCI land cover (https://www.esa-landc over-cci.org/?q=node/175) and Jiang et al. (2017)