Meta-analysis reveals widespread negative associations between species richness and ecological uniqueness

Domain: biology
Published: 2026-03-26
Editorial outcome: major_revision

Abstract
Species richness and ecological uniqueness (defined as local contributions of sites to beta diversity of a region) constitute two crucial components of biodiversity. Understanding the relationship between these two components is critical for effective conservation planning, yet global patterns and underlying drivers remain largely controversial. Here, we conduct a comprehensive global meta-analysis to investigate the patterns of the species richness-ecological uniqueness relationship across taxa groups and to evaluate the relative importance of four ecological hypotheses in driving this relationship (i.e. species pool, dispersal limitation, environmental filtering and spatial grain size). We find that negative richness-uniqueness relationships are prevalent across different taxa groups (e.g. terrestrial plants, freshwater macroinvertebrates, birds, and fishes), and such negative association is robust to different data types including presence absence and abundance data. Notably, relationships based on presence-absence data are primarily predicted by species pool attributes (i.e. the proportion of rare species and gamma diversity), whereas abundance-based relationships are more strongly associated with dispersal limitation. In summary, these findings reveal consistent global patterns and mechanistic underpinnings of the richness-uniqueness relationship, offering key insights for biodiversity conservation. We recommend that conservation strategies should prioritize both species-rich and ecologically unique communities to maximize biodiversity protection globally.

Introduction
Understanding spatial patterns of biodiversity at local and regional scales and their applications to biodiversity conservation is a central topic in ecological research1,2. Prioritization in the conservation planning have traditionally been addressed to species-rich regions3,4. However, those areas with high ecological uniqueness—quantified as their local contributions to beta diversity (LCBD) in a region—has recently received increasing attention as a key indicator of biodiversity5. Studies have revealed the importance of ecological uniqueness in conservation planning because sites with high uniqueness may often harbor functionally unique and/or valuable endemic species that contribute to regional biodiversity despite low species richness6,7. Moreover, sites exhibiting high ecological uniqueness may also contribute disproportionately to regional biodiversity by enhancing trait divergence and niche differentiation of specific species assemblages8. As a result, many conservation efforts have simultaneously emphasized species richness and ecological uniqueness in order to better set the biodiversity hotspots9-11. Meanwhile, the relationship between species richness and ecological uniqueness (hereafter RUR) has emerged as an important yet unresolved issue in ecology and conservation biology9-13. Despite perceived importance, the patterns of this relationship and underlying mechanisms have yet to be synthesized at the global scale.

RUR has important direct implications for biodiversity conservation5,13. If species-poor sites tend to exhibit higher ecological uniqueness of species assemblages10,13,14, conserving only species-rich sites may be insufficient. This is because low-richness sites are more likely to contain high proportions of rare or endemic species which may increase variation in community composition (beta diversity) and regional richness—two important considerations in the conservation planning10,11. In contrast, if species-rich sites are also ecologically unique, conservation efforts targeting on these areas may simultaneously capture both species richness and ecological uniqueness13,15. Despite the potential utility of RUR in conservation planning, empirical studies have shown that RUR can be negative10,16,17, positive13,18 or non-significant13. Moreover, the form of RUR may vary (e.g. linear or quadratic10,12,14), and these variations likely carry important implications for biodiversity conservation strategies across distinct regions and taxa groups. However, there remains no consensus on the general patterns of RUR, and the factors underpinning its context-dependency are not well understood.

Four ecological hypotheses have been proposed to explain variations in RUR. First, the species pool hypothesis posits that attributes such as regional species pool size, variation in local species richness (i.e. breadth of richness) and the proportion of rare species may influence RUR through altering the proportion of shared species13,14,19. For instance, if larger species pool size and larger breadth of richness are mainly driven by larger proportion of shared species (lower species turnover), species-rich communities tend to be dominated by common species and show high compositional similarity with others. Therefore, RUR is more likely to be negative due to low uniqueness exhibited in species-rich communities10,16. Conversely, if increased species pool size and breadth of richness reflect enhanced species turnover, this might result in a more positive RUR where species-rich communities are more likely to harbor unique species combinations13,15. In addition, the proportion of rare species has also been suggested to determine whether RUR is negative or positive12,13. An increase in the proportion of rare species is expected to reduce compositional overlap among communities, with species-rich communities more likely to harbor distinct assemblages and become more unique12,16. Therefore, studies have suggested that RUR may shift from negative to positive with increasing proportion of rare species in regional species pool13,15.

Second, the dispersal limitation hypothesis is expected to affect RUR by regulating species sorting in different communities, as weak dispersal capacity of species or strong habitat isolation may restrict them from tracking their preferred habitat conditions20,21. For instance, studies reported that communities with strong dispersal limitation exhibited higher dissimilar compositions than those with low dispersal limitation, resulting in higher ecological uniqueness20. While existing studies have documented the significant influence of dispersal limitation on ecological uniqueness and species richness20,22, current understanding of how dispersal limitation shape the RUR remains critically underexplored.

Third, the environmental filtering hypothesis argues that environmental conditions, especially climatic conditions could influence RUR through shaping species range sizes, richness patterns or species turnover23,24. For example, studies have suggested that climatically distinct regions may select for species adapted to these conditions, which in turn promote small range sizes in those localities16. This may lead to more distinct species assemblage which could shape patterns of both richness and ecological uniqueness. Although accumulating evidence has demonstrated that environmental filtering (e.g., temperature, precipitation) emerge as key determinants of variations in both ecological uniqueness and species richness16,25,26, how RUR varies across environmental gradients remains poorly understood.

Finally, the spatial grain hypothesis suggests that sampling unit size may shape RUR through modulating habitat heterogeneity, community dissimilarity and/or species rarity27. For instance, larger sampling units are more likely to encompass a wider range of microhabitats, increasing habitat heterogeneity and potentially supporting higher species richness22,28. However, this increase in richness may not correspond to increased compositional uniqueness if additional species are already shared across different communities. Therefore, larger grain size may also potentially homogenize species composition across sites and reducing between-site dissimilarity due to compositional averaging29. Moreover, expanding grain size may alter species rarity as larger grain size are more likely to detect rare species that occur in few communities, thereby changing the distributions of uniqueness and richness14,30. Despite its potential importance, the influence of spatial grain size on RUR remains poorly understood and warrants further systematic investigation.

While the aforementioned hypotheses emphasize the roles of different ecological drivers in shaping RUR, they are not mutually exclusive and may jointly determine RUR13,16,25. However, as different hypotheses predict distinct (even opposite) RURs12-14, identifying the main drivers of the RUR across multiple taxa groups is essential for understanding and predicting RUR but remains unresolved at the global scale. Here, we assembled a global dataset combining published and private sources to assess RUR patterns and their main drivers across multiple taxonomic groups (e.g. birds, terrestrial plants, fishes and freshwater macroinvertebrates) spanning broad latitudinal gradients from the tropics to near polar regions (Fig. 1). We calculated RUR using both presence absence (PA) and abundance (AB) data and examined the influence of four ecological hypotheses on RUR across taxa groups. Our goals are to (1) synthesize global patterns in RUR, (2) identify the key factors that shape RUR, and (3) quantify the relative importance of the four hypotheses in explaining variation in RUR. Our meta-analysis provides a comprehensive framework to better understand the mechanisms underlying biodiversity patterns and offers practical insights for conservation prioritization across diverse taxa groups.

Results
The global patterns of RURs

Our meta-analysis revealed a widespread prevalence of negative RURs across the globe (Fig. 2). Meta regression analyses confirmed significant overall negative effect sizes of RURs for both presence absence (PA) and abundance (AB) data (PA: overall effect size of RUR: -0.488, 95% confidence interval (CI): -0.583 – -0.392, p < 0.001; AB: overall effect size of RUR: -0.362, 95% CI: -0.458 – -0.265, p < 0.001; Fig. 2; Supplementary Table 1). We also found that the overall effect sizes of RUR were consistently and significantly negative across all three data sources (i.e. primary studies, data repositories and compiled datasets), and no significant differences were found among these data sources either for PA or AB data types (Supplementary Fig. 1; Supplementary Table 2 and 3).

Taxon-specific analyses demonstrated significantly negative effect sizes of RURs across 11 taxonomic groups in both PA and AB datasets (Supplementary Fig. 2; Supplementary Table 4 and 5), including birds, terrestrial plants, fishes, terrestrial insects, freshwater macroinvertebrates, algae, arthropods, invertebrates, zooplankton, reptile and lichens. Several other groups exhibited negative but non-significant RURs (e.g. aquatic plants, bats, mammals, bivalves and seaweed in PA data; Phytoplankton, amphibians, bats and nematodes in AB data). Only a few groups showed weakly positive effect sizes of RURs, such as aquatic plants (overall effect size: 0.074, 95% CI: -0.136 – 0.284) in AB data (Supplementary Fig. 2; Supplementary Table 4 and 5).

Drivers of RUR

The mixed-effect meta-regression models revealed several key drivers of RURs across different taxa groups at the global scale, which were grouped into four hypotheses (i.e. species pool, dispersal limitation, spatial grain size, and environmental filtering). Under the species pool attributes, the effect sizes of RUR (Fisher’s Z statistics) consistently declined with increasing gamma diversity (β = -0.197, Q = 6.81, degree of freedom (df) = 1, p = 0.009; Fig. 3A; Supplementary Table 6), although the relationship was not significant in AB data (β = -0.075, Q = 1.27, df = 1, p = 0.26; Fig. 3B; Supplementary Table 6). Taxon-specific analyses supported these findings, particularly for birds and fishes (Supplementary Table 7). By contrast, the effect sizes of RUR increased significantly with breadth of richness especially in PA data (β = 0.254, Q = 6.12, df = 1, p = 0.013; Fig. 3C; Supplementary Table 6), with similar positive effect observed for terrestrial plants in PA data (Supplementary Table 8). In AB data, however, the effect sizes of RUR did not show significant association with breadth of richness (β = -0.068, Q = 0.26, df = 1, p = 0.61; Fig. 3D; Supplementary Table 6 and 8). Similarly, RUR effect sizes increased significantly with the proportion of rare species in PA data overall (β = 0.400, Q = 8.03, df = 1, p = 0.005; Fig. 3E; Supplementary Table 6), with significant positive association found for birds (β = 0.778, Q = 6.29, df = 1, p = 0.012, PA data; Supplementary Table 9). However, this positive association showed marginal significant in general AB data (β = 0.618, Q = 3.39, df = 1, p = 0.066; Fig. 3F; Supplementary Table 6), and no significant association was found in any taxon-specific AB models (Supplementary Table 9). Most of above results remained robust after controlling for other factors in the multivariate meta-regression analyses, except for breadth of richness with significant negative association in PA data (Supplementary Table 10A and 10B).

Dispersal limitation was significantly negatively associated with the effect sizes of AB-based RUR (β = -0.554, Q = 9.74, df = 1, p = 0.002; Fig. 3H; Supplementary Table 6), particularly for freshwater macroinvertebrates and terrestrial insects (Supplementary Table 11). This pattern remained robust after accounting for other variables in the AB-based full multivariate models (Supplementary Table 10A and 10B). Spatial grain size significantly negatively predicted the effect sizes of RUR (PA data: β = -0.056, Q = 5.82, df = 1, p = 0.016; AB data: β = -0.052, Q = 5.06, df = 1, p = 0.025; Fig. 4A and 4B; Supplementary Table 6), although the negative associations were generally weak in the full model analyses (Supplementary Table 10A and 10B). However, significantly negative associations were observed for birds based on taxon-specific analyses (Supplementary Table 12).

Environmental filtering, captured by the first two principle components (i.e. PC1 and PC2) of climatic variables, generally showed weak associations with the effect sizes of RUR in single-factor analyses (PA data: β = -0.006, Q = 0.23, df = 1, p = 0.63; β = 0.024, Q = 1.60, df = 1, p = 0.21; AB data: β = -0.014, Q = 2.12, df = 1, p = 0.15; β = -0.017, Q = 1.08, df = 1, p = 0.30; Fig. 4C-F; Supplementary Table 6). However, in the full model analyses, PC1 of climatic variables had a significantly positive association with RUR effect sizes in PA data (Fig. 5A; Supplementary Table 10A and 10B). In taxon-specific analyses, PC1 of climatic variables significantly negatively predicted RUR effect sizes for terrestrial insects (β = -0.100, Q = 8.74, df = 1, p = 0.003, AB data; Supplementary Table 13), while PC2 of climatic variables significantly negatively predicted RUR effect sizes for algae (β = -0.044, Q = 4.65, df = 1, p = 0.031, AB data; Supplementary Table 14) and significantly positively predicted RUR effect sizes for fishes (β = 0.102, Q = 3.93, df = 1, p = 0.047, AB data; Supplementary Table 14).

In terms of variance explained for PA data analysis, the full model analysis suggested that 33% of the variation was explained by seven factors (marginal R squared of multiple meta-regression model: 0.33; Fig. 5A; Supplementary Table 10A). Species pool had the largest explanatory power on RUR effect sizes, where gamma diversity (34.5%) and the proportion of rare species (30.6%) were the dominant contributors, followed by breadth of richness (16.4%), while environmental filtering, spatial grain size and dispersal limitation contributed marginally (< 10%; Fig. 5C; Supplementary Table 10B). For AB data, the model explained 25% of the variation (marginal R squared: 0.25; Fig. 5B; Supplementary Table 10A), with dispersal limitation (21.1%), the proportion of rare species (17.7%), gamma diversity (16.8%), and PC1 of climatic variables (16.6%) as the major contributors (Fig. 5D; Supplementary Table 10B). We found higher contribution of spatial grain size in AB model (15.2%) than in PA model (5.9%; Supplementary Table 10B).

Tests of publication bias

We assessed potential publication bias for full dataset (including primary studies, data repositories and compiled datasets) and different data sources. Funnel plots and Egger’s regression tests showed no significant evidence of publication bias (full dataset: PA data: z = -0.13, p = 0.89; AB data: z = 0.14, p = 0.89; Supplementary Fig. 3; primary studies: PA data: z = -0.55, p = 0.58; AB data: z = -1.43, p = 0.15; Supplementary Fig. 4; data repositories: PA data: z = 0.12, p = 0.91; AB data: z = -1.29, p = 0.20; Supplementary Fig. 5; compiled datasets: PA data: z = -1.86, p = 0.06; AB data: z = -0.96, p = 0.34; Supplementary Fig. 6). Rosenthal’s fail-safe values further confirmed low publication bias (PA data of full dataset: 3124252; AB data: 1794882; PA data of primary studies: 97459; AB data: 28703; PA data of data repositories: 1577588; AB data: 1020582; PA data of compiled dataset: 39593; AB data: 25474), all of which far exceeded the recommended threshold for robust and significant meta-analytic results (PA data of the full dataset: 3450; AB data: 3085; PA data of primary studies: 610; AB data: 385; PA data of data repositories: 2500; AB data: 2310; PA data of compiled dataset: 360; AB data: 410).

Methods
Literature search and data collection

Three distinct approaches were employed here to allow for a broad and representative collection of studies or datasets to be gathered (See Supplementary Note 1: Compiled datasets and Supplementary Data 1 for details):

Peer-reviewed primary studies: First, we searched for peer-reviewed journal articles using the ISI Web of Science using predefined keyword combinations to obtain an extensive collection of primary studies that focused on ecological uniqueness globally. The search was conducted using the ISI Web of Science Core Collection database, accessed via the institutional subscription of the University of Chinese Academy of Sciences. The structure of keyword combinations was defined as:

TS = (biodiversity OR diversity)
AND
TS = ((ecological uniqueness) OR (taxonomic uniqueness) OR (composition* uniqueness) OR (community uniqueness) OR (local contribution* to beta divers*) OR (LCBD) OR (taxonomic distinctness) OR (unique species composition*) OR (distinct assemblage*) OR (distinct communit*))
AND
TS = ((taxonomic beta divers*) OR (taxonomic composition*) OR (species composition*) OR (species combination*) OR (species assemblage*) OR (community composition*) OR (community assemblage*) OR (community dissimilarit*) OR (community variation*) OR (composition* dissimilarit*))
AND
AB = ((ecological uniqueness) OR (taxonomic uniqueness) OR (composition* uniqueness) OR (community uniqueness) OR (local contribution* to beta divers*) OR (LCBD) OR (taxonomic distinctness))

These used keywords were typical terms in empirical studies related to ecological uniqueness and the last search was done on December 1, 2024. This process was conducted in adherence to the Preferred Reporting Items for Systematic Reviews and Meta-analyses (PRISMA) 2020 statement57. Additionally, as a supplementary strategy, we also performed a citation tracking search and investigated all studies that cited Ref 5, on December 1, 2024. This is because this seminal paper provided both the conceptual foundation and terminology of ecological uniqueness (i.e. the metric of local contributions to beta diversity), and which has become a cornerstone in this research domain. Collectively, 1763 records were obtained from these combined search methods (i.e. records of keyword-based searching: 876; records of seminal literature citation tracking: 887) and 1542 records were retained after removing duplicated records (Supplementary Fig. 8).

We established several criteria to exclude several primary studies (Supplementary Fig. 8), these are: (1) We excluded experimental studies because their species combinations are typically artificially manipulated or intentionally set in each community or study unit (Reason 1, 5 studies; Supplementary Fig. 8). (2) We excluded studies that failed to provide taxonomic information at the species level. For instance, several studies limited their taxonomic resolution to the family, subfamily, OTU or other levels, which did not reflect species level richness characteristics (Reason 2, 44 studies; Supplementary Fig. 8). (3) We excluded studies that did not report the distributions of both species richness and ecological uniqueness values, or studies that only report related statistics (e.g. regression coefficients, slope t-tests or model estimates) without providing richness or uniqueness distributions. Although it is methodologically feasible to convert such statistics into correlation coefficients58, this exclusion is necessary because previous empirical evidence demonstrates that the RUR may exhibit either linear or nonlinear (e.g. curve-linear) relationships10,14,25. Without detailed distributions of species richness and ecological uniqueness, we could not reliably characterize the association form of RUR. To empirically test the form of the relationship for each study, it was essential for us to have access to the data distributions. Moreover, relying solely on converted summary statistics would have presupposed a linear association and prevented us from verifying this key assumption. Such incomplete reporting would fundamentally limit our capacity to distinguish between alternative ecological models while introducing substantial uncertainty in the synthesis (Reason 3, 186 studies; Supplementary Fig. 8). (4) We excluded studies that exhibited methodological variation in data acquisition approaches. Specifically, studies combining multiple data sources with inconsistent sampling protocols (e.g., those integrating community data from both literature-compiled records and original field surveys within the same analysis) were excluded to ensure methodological consistency. Furthermore, we also excluded studies lacking standardized sampling protocols, particularly those with inconsistent sampling areas across sampling units. This is because studies have revealed that both species composition and species richness are sensitive to the area of sampling unit30,59. This exclusion criteria was implemented to minimize confounding effects and enhance comparability across studies (Reason 4, 8 studies; Supplementary Fig. 8). (5) We excluded studies that calculated ecological uniqueness metrics based on other analytical units incompatible with our spatial site-based framework. Specifically, studies employing non-spatial compositional units (e.g., temporal observations at different time points, individual organisms from distinct populations, or nodes at food webs) to obtain ecological uniqueness were excluded because these study units fail to capture site-based distributions of ecological uniqueness5. This exclusion criterion ensures conceptual consistency in defining ecological uniqueness as a spatial-community property emerging from β-diversity patterns (Reason 5, 39 studies; Supplementary Fig. 8). (6) We excluded studies that solely reported species richness of specific taxonomic groups (e.g., threatened species, invasive/exotic species, or other preselected categories). This exclusion criterion was implemented because our primary focus lies in understanding RUR at the community or sampling site level (Reason 6, 7 studies; Supplementary Fig. 8). (7) We excluded methodological studies that primarily focused on synthetic/simulated data to quantify species richness or ecological uniqueness, because our analytical framework specifically investigates patterns of ecological uniqueness and species richness derived from empirical observations (Reason 7, 4 studies; Supplementary Fig. 8). (8) We also excluded studies sharing identical datasets (i.e., multiple publications analyzing the same datasets) (Reason 8, 1 study; Supplementary Fig. 8). Using these eight exclusion criteria, 56 studies with 164 unique observations were retained (Supplementary Fig. 8 and Supplementary Data 1) to evaluate RUR across taxa groups globally.

Moreover, it should be noted that both ecological uniqueness and species richness can be derived from raw species composition data (PA or AB data matrices5). During our literature screening, we identified 44 studies that, while not explicitly addressing ecological uniqueness as their primary focus, contained suitable raw datasets. These spatial datasets satisfied essential quality controls including: (1) being observational rather than manipulated experimental studies, (2) employing consistent sampling protocols across surveyed locations or sites, (3) achieving species-level taxonomic resolution, and (4) containing non-overlapping datasets from distinct research efforts. We extracted and processed these datasets to calculate ecological uniqueness, species richness and RUR. The acquisition of these 44 datasets yielded 124 additional observations. Thus, we retained 100 qualified studies following the PRISMA-guided screening process, which comprised 56 case studies explicitly investigating ecological uniqueness and 44 studies containing compatible observational datasets (See Supplementary Data 1 for details). All retained studies were systematically incorporated into the final datasets, collectively yielding 288 independent observations. A PRISMA diagram summarizing the workflow of literature search is included in Supplementary Fig. 8.

Public data repositories: Second, we compiled several additional raw datasets across different taxonomic groups from data repository, because ecological uniqueness and species richness were both computable from raw species composition datasets5. We separately used the following keywords to search for biological community composition datasets in Dryad digital repository (http://datadryad.org/): (taxonomic composition) OR (species composition) OR (species assemblage) OR (community composition) OR (community assemblage) OR (species abundance) OR (presence absence) (Supplementary Fig. 8). The last search was done on December 15, 2024. These search strings yielded a total of 14693 datasets, 9044 datasets were retained after removing duplicated records across these search strings. We initially screened the titles, abstracts, related references to exclude datasets with unrelated topics. We found that 314 records provided raw datasets of species compositions (Supplementary Fig. 8). After that, we also established five criteria to exclude some of these datasets, in order to ensure the quality of the datasets.

(1) We excluded datasets related to experimental studies because the species combinations are artificially manipulated in each sampling unit (Reason 1, 21 datasets; Supplementary Fig. 8). (2) We excluded datasets that failed to provide taxonomic information at the species level. For instance, several datasets reported richness of taxonomic levels of other types (e.g., genus, family, OTU or ASV), which did not meet the requirement of species-level richness in the analysis. Additionally, we disqualified datasets containing species with missing information (e.g., species that could not be identified or lacked pertinent data within the dataset), as these omissions made it difficult for us calculate species richness (Reason 2, 5 datasets; Supplementary Fig. 8). (3) We excluded datasets that exhibited methodological heterogeneity in data acquisition approaches. For example, some datasets obtained part of their community data from literature sources while other parts were based on field surveys. In addition, we also found that sampling sites differed in their sampling area in certain datasets. Such studies were not included in our synthesis. Meanwhile, datasets lacking details of sampling protocols were also excluded (Reason 3, 30 datasets; Supplementary Fig. 8). (4) We excluded datasets based on other sampling unit incompatible with our spatial-community framework (e.g., temporal observations, individual organisms from distinct populations, or food web nodes) (Reason 4, 4 datasets; Supplementary Fig. 8). (5) We also excluded studies sharing identical datasets (i.e., the same dataset used in primary studies or multiple data publications containing the same datasets) (Reason 5, 11 datasets; Supplementary Fig. 8). Following the established exclusion criteria, we retrieved 243 raw datasets from the Dryad repository, encompassing 20 distinct taxonomic groups (e.g., terrestrial plants, fish, birds, terrestrial insects, invertebrates; Supplementary Table 15). We then calculated ecological uniqueness and species richness based on these datasets. This process generated 1415 unique records of observations which were subsequently incorporated into the final dataset. We have also provided the PRISMA diagram summarizing the workflow of dataset search and the list of these datasets (See Supplementary Fig. 8 and Supplementary Data 1 for details).

Following the PRISMA guidelines, two authors independently conducted the literature and dataset screening based on the above exclusion criteria. The two authors examined 3 primary studies and 4 datasets that were only selected by one of the authors. Following a joint re-assessment of the studies and datasets with initial screening discrepancies between the two authors, a full consensus was reached to exclude all of them. We have listed all unselected literature and datasets with the reasons for exclusion (See Supplementary Data 1 for details).

Compiled datasets: Third, we compiled a series of additional raw datasets from our own data. These datasets were mainly acquired through: (1) direct correspondence with principal investigators maintaining unpublished datasets, or (2) using the existing data sources including: GrassPlot60,61, sPlot62 and CESTES63 by ourselves. We included datasets that provided raw data of species abundance (AB) or presence absence (PA). Some of the datasets (i.e. 5 datasets) provided only PA data rather than AB data (See Supplementary Note 1: Compiled datasets). All compiled datasets also underwent rigorous screening against the predefined five-criteria framework originally developed for Dryad dataset selection. In total, 108 datasets with 215 unique observations were assembled and included in the synthesis, 94 datasets of which were published in Ref 43. These datasets include a wide range of taxonomic groups including: Algae (19 datasets), amphibians (2 datasets), aquatic plants (3 datasets), arthropods (9 datasets), bats (2 datasets), birds (10 datasets), fishes (12 datasets), freshwater macroinvertebrates (12 datasets), phytoplankton (11 datasets), terrestrial insects (9 datasets) and terrestrial plants (60 datasets), providing a solid data foundation for the synthesis of RUR (See Supplementary Note 1: Compiled datasets for details).

In summary, our final dataset incorporated a total of 1918 data records from 451 primary studies or datasets (Fig. 1). The final dataset included three types of data sources: Peer-reviewed primary studies (100 studies, 288 data records), public data repositories (243 datasets, 1415 data records) and compiled datasets (108 datasets, 215 data records; Fig. 1). However, for consistency in description, we uniformly referred to each dataset or primary study as a "study" in the subsequent text. For each study, we extracted its geographic coordinates (latitude and longitude), taxonomic group, number of study units (the number of sampling or observation sites). We also recorded the type of data used to calculate ecological uniqueness (i.e., PA or AB data) because abundance weight was proved to be a key factor determining spatial patterns of ecological uniqueness5,16. Additional variables were collected within each study to represent the influence of different hypotheses (See section: Quantification of different hypotheses for details). In general, these studies cover a wide range of taxonomic groups (e.g. terrestrial plants, birds, fishes and freshwater macroinvertebrates) and we divided them into 20 groups (see Supplementary Table 15 for details). The final dataset included 1076 records of PA data and 842 records of AB data.

The calculation of ecological uniqueness

We applied two alternative approaches to calculate ecological uniqueness. Approach 1 was based on total variance in the community dataset5, which the sum of squares from raw species composition matrices (PA or AB data) were calculated and converted into total beta diversity. We then decomposed the total beta diversity into site-specific contributions (i.e. local contributions to beta diversity; defined as Unique_R). Approach 2 involved partitioning a pairwise dissimilarity matrix (e.g., Jaccard or Sørensen) among all study sites, to assess each site’s contribution to overall heterogeneity5,19. This yielded two additional uniqueness indices (Unique_J and Unique_S), based on the decomposition of the Jaccard and Sørensen dissimilarity indices19. The calculation of Unique_R was based on the beta.div function in the “adespatial” package, while the calculations of Unique_J and Unique_S were primarily based on the beta.div.comp and LCBD.comp functions from the same package5,19. All three indices of ecological uniqueness were highly correlated within the same case observations (Pearson’s r > 0.75 for all pairwise comparisons, See Supporting data). Therefore, we prioritized Approach 1 in the subsequent RUR analyses, as it is the most commonly used in studies of ecological uniqueness7,10,12,13,16.

Effect size of RUR

Primary studies have documented different forms of RUR. Through systematic review of previous studies that have reported on the relationship between species richness and ecological uniqueness, we identified linear (first-order) and quadratic (second-order) models as they are predominant forms of RUR characterization (See Appendix 1, which is available on Figshare at: https://doi.org/10.6084/m9.figshare.29270048). To determine the optimal form of RUR, we systematically compared these two formulations across all studies. For each RUR observation across studies, we fitted both linear and quadratic models as follows:

Linear model:
Uniqueness ~ β + β *log_richness + ε (1)
0 1

Quadratic model:
Uniqueness ~ β + β * log_richness + β * (log_richness)2 + ε (2)
0 1 2

Where Uniqueness is ecological uniqueness; log_richness is log-transformed species richness to account for its wide range of values across study units in some cases, and ε is the residual error. We evaluated model performance through Akaike Information Criterion (AIC) comparisons. For studies with sample sizes (i.e. number of study units) less than 30, we additionally calculated corrected AIC (AICc) values to account for small-sample bias. We performed statistical model comparisons between linear and quadratic forms using the anova function in R, which performs an F-test to evaluate whether the additional terms in the more complex (quadratic) model significantly improve the fit. Models were considered statistically equivalent when the resulting p-value exceeded the 0.05 significance threshold. Consequently, our analysis revealed that 32.1% of records (615 out of 1918) exhibited significantly better fit with quadratic model than linear model (p < 0.05; See Appendix 3, which is available on Figshare at: https://doi.org/10.6084/m9.figshare.29270048), while 67.9% (1303 records) showed no statistical difference between linear and quadratic formulations (i.e. p > 0.05). For these cases, the differences in AIC or AICc between linear and quadratic models were also less than 2 units, further supporting their statistical equivalence64 (See Appendix 3, which is available on Figshare at: https://doi.org/10.6084/m9.figshare.29270048). Therefore, we considered the linear model as the most parsimonious representation of RUR patterns for these 1303 observations. For the 615 records where the quadratic model showed a significantly better fit, we examined if these cases were biased toward specific taxonomic groups. We found that quadratic relationships represented 36.1% of PA records and 27.0% of AB records across all taxa groups, with no substantial bias observed across taxonomic groups with more than 10 records (See Supplementary Table 15).

Given the difficulty of deriving standardized effect sizes for non-linear (quadratic) relationships of RUR, and to maintain methodological consistency and cross-study comparability, we restricted subsequent analyses to the 1303 records best represented by linear models. These included of 688 PA data records and 615 AB data records. However, we also performed a qualitative analysis of the nonlinear cases of RUR, better represented by quadratic forms. First, we fitted quadratic models to describe the relationships between ecological uniqueness and species richness in these cases. Then, we calculated the axis of symmetry for each quadratic model and quantified the number of data points distributed on each side of the axis. Based on these patterns, the relationships were categorized as follows: (1) For models with R² greater than 0.1, if all data points were distributed entirely on the left or right side of the axis of symmetry, the relationship is classified as either a "negative association" or "positive association" (i.e., ecological uniqueness decreases or increases with increasing species richness; Supplementary Fig. 9A-D). (2) For models with R² greater than 0.1, if more than 80% of the data points were located on one side of the axis of symmetry, the relationship is classified as a "general negative association" or " general positive association", indicating a predominant but not strictly monotonic trend (Supplementary Fig. 9E-H). We also tested alternative thresholds (e.g., 70% and 90%) and found the results to be similar, and 80% was accordingly selected as the distribution threshold. (3) When fewer than 80% of the data points were located on one side of the axis of symmetry, the relationship was classified as "no general association", indicating that both species-rich and species-poor sites could exhibit high ecological uniqueness. (4) If the model R² is less than 0.1, regardless of the distribution of points relative to the axis of symmetry, the relationship is also classified as "no general association". The details of the distributions of these nonlinear RUR can be found in Appendix 2 (which is available on Figshare at: https://doi.org/10.6084/m9.figshare.29270048).

We then used Pearson’s product–moment correlation coefficients (r) to quantify those linearly modelled relationship between species richness and ecological uniqueness (i.e. RUR). We implemented the following calculation strategy: (1) In several cases, the original species composition data were not available, and only the reported values of ecological uniqueness and species richness were accessible. Upon review, we determined that these studies predominantly employed Approach 1 to calculate ecological uniqueness. Therefore, we directly extracted these values and used them to derive the corresponding RUR. (2) For studies providing both diversity distributions (i.e. species richness and ecological uniqueness) and species composition datasets, we preferentially recalculated ecological uniqueness using Approach 1 described above to ensure methodological consistency and then incorporated them into RUR analyses; (3) In cases where studies contained composition data without reporting distributions of ecological uniqueness and species richness, we implemented the same method as (2). Some studies incorporated more than one estimate of RUR in different locations, biological groups or different seasons. Therefore, these r values were included as separate individual observations within each study and then referred to as an estimate. While correlation does not imply causation, the correlation coefficient provides a simple measure of association between species richness and ecological uniqueness, and thus has been used across studies53,65. The Pearson’s correlation coefficients (r) from each estimate were then normalized using Fisher’s Z transformation as an effect size of RUR:

1 1 + r
Fisher's Z = ln( ) (3)
2 1 - r

Note that a negative effect in Eqn (3) suggests a negative RUR while a positive effect indicates a positive RUR. We weighted each effect size based on the sample size (i.e., the number of sites for which of ecological uniqueness and species richness were calculated) using the formula:

1
v = (4)
i N s - 3

Where N is the sample size (i.e., number of sampling units) of each estimate within each study. This indicates that studies had larger sample sizes would be given higher weights.

Quantifications of different hypotheses

In addition to investigating the global patterns of RUR, we are also interested in exploring the role of four hypotheses with a combination of abiotic and biotic variables (i.e. species pool hypothesis, dispersal limitation hypothesis, environmental filtering hypothesis, and spatial grain hypothesis) in determining RUR. Specifically, to quantify the attributes of species pool, we selected gamma diversity (i.e. total number of species within each study) to represent the size of regional species pool. We then calculated coefficient of variation of species richness among different study units for each observation to represent the breadth of richness, with larger coefficient of variation indicating wider breadth of species richness within each observation. We also investigated the effect of the proportion of rare species within each study, which was defined through dual criteria corresponding to different data types (i.e. PA and AB). For PA data, species detected at less than 10% of sampling units within individual observation were classified as rare. The ratio of these low-frequency species to total number of species was quantified as the proportion of rare species under PA framework. For AB data, we implemented a composite threshold to define the proportion of rare species, and species were designated as rare only if they simultaneously satisfied two conditions: (1) Occurrence in less than 10% of sampling locations, and (2) abundance falling below the 10th percentile rank within the full species abundance distribution. These three metrics provide a multi-dimensional quantification of the species pool hypothesis. Gamma diversity directly measures the size of the regional species pool, representing the total source of potential colonist species available for community assembly34. Breadth of richness specifically reflects the variability of richness within given region - a higher value indicates greater richness variation across sampling units, thereby illustrating how the species pool maintains different community configurations across varying local habitats14. The proportion of rare species reveals the pool's internal structure, with a higher proportion suggesting lower number of shared species or higher species turnover among communities13.

To assess the potential influence of dispersal limitation on RUR, we employed the method of distance-based Moran’s Eigenvector Maps (dbMEM) to model spatial structure across study sites. This is mainly because they effectively detect spatial autocorrelation patterns characteristic of constrained species movement. Dispersal limitation produces a distance-decay relationship in community similarity, where nearby sites tend to share more species (positive spatial autocorrelation). dbMEM variables capture these underlying spatial structures at different scales, from broad-scale gradients to fine-scale patches. This method has been used recently to estimate potential dispersal limitation in various taxa13,66-69. Using georeferenced coordinates from each study (if available), we first constructed a truncated spatial connectivity matrix by retaining pairwise distances below a threshold corresponding to the maximum spanning tree length. Eigen decomposition of this matrix generated orthogonal spatial eigenvectors (i.e. dbMEM variables), with positive eigenvalues representing meaningful multi-scale spatial patterns70-72. We then implemented beta regression models73, employing logit-link functions, to quantify the explanatory power of these spatial eigenvectors, given the bounded nature of ecological uniqueness values (i.e. range from 0–1). The full model structure was specified as:

k
logit(uniqueness) ~ β + ∑β * (dbMEM) + ε
0 i i
i = 1

Where dbMEM are orthogonal spatial eigenvectors derived from truncated distance matrices, each capturing distinct spatial scales of autocorrelation71, k denotes retained eigenvectors with positive eigenvalues, β is the intercept term, β are regression coefficients quantifying the direction and magnitude of spatial eigenvector i's effect on ecological uniqueness in each study or datasets, ε is the residual error. The marginal 𝑅-squared (i.e. Nagelkerke’s pseudo-𝑅2) from these models was interpreted as the potential explanatory capacity of dispersal limitation, representing the proportion of uniqueness variance attributable to spatially structured processes66,72,74. When these variables explain a significant portion of the variation in community composition, the most parsimonious interpretation is that dispersal limitation plays a dominant role in shaping communities.

To assess the effect of environmental filtering, we extracted 19 bioclimatic variables (Bio1-Bio19) from the WorldClim database75 with 30 arc seconds spatial resolution using each study's mean geographic coordinates. This is mainly because of their direct quantification of the critical climatic conditions that determine species' survival and distribution. The species-energy hypothesis posits that climate strongly influences primary productivity, or the total energy available for partitioning within most ecological communities, thereby affecting species richness or species composition via its impact on the number and size of organisms that can be supported76-79. These variables such as the minimum temperature of the coldest month (Bio6) and precipitation of the driest quarter (Bio17) capture the critical extremes of temperature and water availability that often limit species persistence.

Principal Component Analysis (PCA) was applied to mitigate multicollinearity among these variables. The first two principal components (PC1 and PC2) were retained as synthetic climatic predictors in subsequent analyses (see Supplementary Table 16 for details). We first conducted a single-factor analysis and found that PC1 showed strong negative associations with variables related to mean temperature (e.g., Bio1, Bio5, Bio6, Bio10, Bio11) and annual as well as seasonal precipitation (e.g., Bio12, Bio13, Bio16, Bio19), while displaying positive associations with temperature seasonality (Bio4) and annual temperature range (Bio7) (Supplementary Table 16A). This suggests that PC1 captures a transition from warm, wet, and thermally stable environments to cooler, drier regions characterized by greater temperature fluctuations in general. By contrast, the PC2 represents a shift from warm, seasonally dry climates with high temperature variability to cooler, wetter environments with more evenly distributed precipitation throughout the year (Supplementary Table 16A). We then conducted a full model analysis based on a filtered dataset due to differences in data availability (Supplementary Table 16B). In general, the environmental gradients captured by PC1 and PC2 were similar with the single-factor analysis, i.e. PC1 suggested a gradient from cool and relatively dry climates to warm and wet climates while PC2 from warm and seasonal dry climates to cooler and wetter climates, although the gradient direction of the single-factor and full model analyses was different in some cases (e.g. PA data; Supplementary Table 16A and 16B). Finally, we re-performed PCA on the bioclimatic variables when investigating the effect of environmental filtering on RUR based on single-factor taxon-specific analyses (details of these PCA results can be found in Supplementary Table 16C-J).

As for spatial grain size, it was defined as the sampling area of each study unit within each study if reported. The size of spatial grain not only directly determines the number of species observed and the spatial turnover in species composition, but also potentially alters the influence of other ecological processes such as biotic interactions and environmental heterogeneity27. Larger grains tend to integrate greater habitat heterogeneity and buffer fine-scale stochasticity, while smaller grains better capture local biotic and environmental filtering effects. It should be noted that some studies reported sampling area of smaller sampling units. We converted these values into the corresponding area of the study unit used in the study. For example, plots or transects were divided into smaller quadrats in some surveys, and the area of each quadrat was reported in these studies. In such cases, we summed the areas of the smaller quadrats to obtain the total area of each plot or transect because species richness or ecological uniqueness was calculated at the plot or transect level in these studies. Additionally, for certain taxonomic groups (e.g., algae or freshwater macroinvertebrates), some studies reported sampling unit size based on other metrics (e.g. volume or transect length51). We did not include these cases in the current analyses due to the limited availability of such data. We applied a logarithmic transformation to account for the large variation in the magnitude of spatial grain size. These factors were chosen as they may represent the main potential sources of influence on RUR as proposed by previous studies12-14,16, 20.

Statistical analysis

To assess the overall mean effect size across data types (PA and AB data), we implemented an inverse-variance-weighted hierarchical mixed-effects meta-regression model incorporating a multi-level random effects structure. The model specified two independent random effects to account for potential sources of heterogeneity: (1) nested random effects of repeated individual RUR estimates within each study and (2) an independent random effect for taxonomic group to address between-taxa heterogeneity. We also conducted subgroup analyses to examine RURs across different taxonomic groups. For subgroup analyses by data type, we included both the nested structure and taxonomic group as random effects. For subgroup analyses by taxonomic group, we retained only the nested random effect (i.e., estimates nested within studies), as taxonomic group was the grouping variable being analyzed. The Fisher's Z transformation was used as effect size of RUR, while data type and taxonomic groups were used as moderators in the models respectively. The overall effect size was considered to be significant if the 95% confidence intervals did not overlap zero. We calculated I2 across all effect sizes to estimate total heterogeneity80. We partitioned the total I2 at different levels of random effects set above, to quantify the amount of variation in effect size that could be attributed to differences in between-study, within study or taxonomic groups81.

Prior to the analysis of driving factors, we performed a null model analysis to test whether the observed RURs significantly deviated from those of random communities. The null model randomizes species distributions under specific ecological constraints. For presence absence data, the null model analysis explicitly maintains the incidence frequency of each species, the observed species richness of each site and the total number of species occurrences in the study region. By preserving these fundamental constraints, the model randomizes only the co-distribution of species among sites, thereby simulating communities where composition is driven solely by stochastic colonization under the observed geographic and environmental boundaries. This approach allows for a direct test of whether non-random patterns in species co-occurrence exist beyond what is expected by chance, given the underlying ecological setting82.

For abundance data type, the null model strictly holds two key properties of the original abundance matrix constant: 1. Site-level total abundance: The sum of all individuals in each site/community remains unchanged. 2. Species-level total abundance: The sum of all individuals of each species across all sites/communities remains unchanged. This preserves the overall abundance distribution of the species pool. Meanwhile, the null model randomizes the internal structure of the community matrix. Specifically, it reshuffles the individuals among the sites which indicates that the specific association between a given species and a given site (i.e., how many individuals of a species are found in a particular site) is randomized83,84.

Subsequently, we constructed a null distribution of the RUR by creating 999 randomized null communities. By comparing the observed patterns against the distribution of patterns from these randomized null communities, we can test whether the observed RUR patterns differ significantly from what would be expected by chance. The statistical significance of the observed metric was assessed by calculating a two-tailed p-value. This was done by determining the proportion of null values that were as extreme as, or more extreme than, the observed value. A statistically significant deviation from the null distribution was inferred if the observed value fell outside the range defined by these 2.5th and 97.5th percentiles, which corresponds to a p-value of less than 0.05. This indicates that the observed pattern is unlikely to have arisen by random chance, it provides evidence that other ecological processes—such as species interactions, environmental filtering, or dispersal limitation—may actively shape the observed patterns.

Consequently, the null model analysis revealed that a total of 645 records (290 records of presence absence data and 355 records of abundance data) differed significantly from random communities (See Supporting data for details). Therefore, our subsequent analyses examining the predictions of different hypotheses on RUR and their relative importance were conducted exclusively on this subset of records. For these 645 records, we examined the potential influence of seven predictors belong to four hypotheses including: (1) species pool (i.e. gamma diversity, breadth of richness and proportion of rare species); (2) dispersal limitation, (3) environmental filtering (i.e. PC1 and PC2 of climatic variables) and (4) spatial grain size on the effect sizes of RUR across different taxa groups. We first conducted single-factor analyses for each predictor:

Effect ~ β + β * P + δ + δ + δ + ε (5)
RUR 0 1 Study Study/Estimate Taxa

Where Effect is the effect size of RUR; β is the coefficient to be estimated; P is one of the RUR 1 above predictors; δ and δ are the between-study and within-study random effects, Study Study/Estimate δ is the between-taxonomic group random effect, 𝜀 is the sampling error. In specifying the Taxa random effect associated with different taxa groups, we included only those taxa groups for which the number of effect sizes was greater than or equal to 5 (See Supplementary Table 17 for details). This threshold was applied to ensure that each group contributed enough information to reliably estimate the group-level variance component, and to avoid overparameterization arising from sparse data in infrequently represented groups85. We also separately conducted single-factor analyses for taxa groups with over 5 effect sizes (Supplementary Table 17). However, for taxa groups exhibiting low heterogeneity in effect sizes80 (I² < 0.25), we did not perform single-factor analysis. The nested structure of estimate within study was retained whenever applicable. However, the nested structure may not be applicable to some taxa groups because most studies included only one estimate. Therefore, we specified random effects at the study level in the analyses of these groups.

We further examined relative importance of these hypotheses by conducting full-model analyses incorporating different predictors:

Effect ~ β + β * GD + β * BS + β * Pr + β * DL + β * SG +
RUR 0 s1 s2 s3 s4 s5
β * Clim1 + β * Clim2 + δ + δ + δ + ε (6)
s6 s7 Study Study/Estimate Taxa

Where GD is gamma diversity; BS is breadth of richness (i.e. the coefficient of variation of species richness); Pr is the proportion of rare species; DL is dispersal limitation; SG is spatial grain size; Clim1 is the PC1 of climatic variables; Clim2 is the PC2 of climatic variables. All predictors were standardized by subtracting the average value and dividing by the standard deviation. β to β are the standardized coefficients to be estimated within full model analyses. s1 s7 We also excluded those taxa groups with fewer than 5 effect sizes when constructing the full model analyses incorporating all factors across different taxa groups. It should be noted that only four major taxa groups (terrestrial plants, freshwater macroinvertebrates, birds and terrestrial insects) had sufficient number of effect sizes and were included in the full model analysis incorporating all predictors. In addition, we did not perform taxon-specific full model analysis because the sample size was considered inadequate, with the largest group comprising only 44 observations that had complete records for all seven explanatory variables. Multicollinearity of variables was tested before full model analyses using variance inflation factor (VIF) statistics to ensure that all variables with VIF < 5 (Supplementary Table 18). We estimated the relative importance of each predictor in the full model analyses as the ratio of each standardized regression coefficient to the sum of all absolute standardized coefficients in the full model86.

We calculated Q -statistic to test for significance of factors in meta-regressions. The Q - M E M statistic is a Wald-type test of model coefficients, and a significant Q -statistic indicates that the M factors contribute to the heterogeneity in effect sizes. We also calculated marginal R squared for single-factor models and Rmultivariate models to show how much heterogeneity is explained by these factors87. All parameters in the meta-analytical models were estimated using restricted maximum likelihood (REML) estimation to ensure unbiased variance component calculations in nested structures88. We explored potential publication bias in the current meta-analysis using funnel plots, Egger’s regression test for funnel plot asymmetry89, and Rosenthal’s Fail-Safe N test90. The “trim and fill” method was used to both identify and adjust for funnel plot asymmetry from publication bias91. All analyses were performed using R with the “metafor” package88,92. The visualization of results in the meta-regressions was created using the “orchaRd” package93.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "I have reviewed the MS entitled Meta-analysis reveals widespread negative associations between species richness and ecological uniqueness across the globe.\nThe MS is generally well written and easy to follow. However, it lacks some important methodological details and coverage of the relevant literature.\nGeneral comments\nOn one hand, I applaud the exercise of evaluating the richness-uniqueness relationship (RUR) on many datasets across ecosystems and taxonomic groups. The main result that the RUR is generally negative is an important one. One the other hand, I believe that the paper makes a poor job of testing the many different processes that could affect the direction and strength of the RUR. Instead, the paper focusses on species replacement (LCBD_repl) as the core process. In particular, the paper fails to revisit the results of Dansereau, Legendre and Poisot (2022), who extensively tested the properties of the LCBD metric on a large dataset of bird communities.\nDansereau, G., Legendre, P. and Poisot, T., 2022. Evaluating ecological uniqueness over broad spatial extents using species distribution modelling. Oikos, 2022(5), p.e09063.\nIt is clear from Dansereau et al.’s results that the strength and direction of RUR depends on the richness of the species pool, the breadth of the richness gradient across site, the spatial extent of the sampling, as well as the proportion of rare species in the pool (Dansereau et al., 2022). Furthermore, their results clearly show that the canonical RUR is not linear (rather strongly concave), and therefore Pearson’s correlation is perhaps not a proper approximation of its strength and direction.\nSpecific comments\nIt is not clear if the meta regression aims at testing the “strength” OR “direction” of the RUR as both terms are used in the introduction (line 39 and line 71 of the introduction) as well as in the results. There may be different processes affecting the strength or the direction. These differences are not well explained (but see Dansereau et al., 2022).\nJudging from Figure 4, there would be no association between Rate of species replacement and RUR if only strength was tested (absolute Fisher Z). But then, Fisher Z is a metric mostly related to strength (Pearson’s r) rather than direction. This is a problem.\nIt is also not clear why should we expect changes in the strength or direction of the RUR across taxa and ecosystems. Are factors such as the size of the species pool, the breadth of the richness gradient, the spatial extent, or the proportion of rare species, expected to vary across taxa and ecosystems, thus explaining some of the differences? Probably. These are not detailed.\nMore fundamentally, the authors partitioned the LCBD metric (beta diversity) into its components of species replacement (LCBD_repl) and richness difference (LCBD_rich).\nThey then tested if the strength and direction of the RUR is affected by LCBD_repl. First, it is mentioned in the MS that the effect of LCBD_rich was also evaluated (line 437), but I could not find these results. Second, since both the slope and Pearson’s correlation of the RUR depend on LCBD values, and that LCBD_repl also depends on LCBD values (LCBD_tot = LCBD_repl + LCBD_rich), there is a high risk of spurious association here. A full simulation should be conducted to avoid any artifactual findings.\nThe above concerns are major, and therefore I cannot recommend the paper for publications.", "reviewer_2": "MAJOR CONCERNS\nIn its current state, this is not a high-quality meta-analysis and does something I’ve never seen before which is only use snowballing methods to find studies for the synthesis (e.g., screening and including only studies that cite a key paper, here Legendre & Ce Caceres 2013). Typically, this approach is meant to supplement an exhaustive keyword search of multiple databases -- here they skipped the exhaustive search and went straight to screening and including only studies citing Legendre & Ce Caceres (2013). They also supplement sources but do not justify why these sources were chosen or whether they are representative of what’s available. How many studies were missed with these unusual decisions? There is also no information on who screened the studies, or whether there were any practices to validate screening decisions. I suggest the authors have a look at PRISMA guidelines to develop practices that generate high quality syntheses.\nThey also apply a custom weighting that will have unknown consequences for the meta-analysis and publication bias tests because the variance of Fisher’s Z is meant to strictly model sampling error among effect sizes (and the variability among these weights later, in turn, impact the random-effects estimated by the meta-analysis model), but here they’ve confounded this weight and its purpose with aspects of spatial extent without justification. The meta-analysis model should not be modified like this, and if you believe spatial extent impacts effect sizes then treat it as a continuous predictor in the model, and then explore pooled effects as marginal means corrected by spatial extent. Also please report your full meta-analysis model; it’s unclear what random effects were included, what predictors are combined to achieve synthesis, or whether these were implemented correctly with the metafor R package. Results are also incomplete as they emphasize P values but do not report the statistic of these P values. On Line 131, you also cannot do a meta-analysis of absolute values of effect sizes because this moves the null from zero to an unknown null value -- making significant tests meaningless.\nA lot of this work (some of it is quite amazing and remarkable) is simply waste because odd decisions were made that make it impossible to assess the robustness of what was achieved.\nMINOR NOTES\n1) Lines 37-38: do you mean “have yet to be synthesized”\n2) k instead of N should be used to describe number of effect sizes (N is reserved for number of data points used to estimate Z)\n3) report the degrees of freedom for all Q tests\n(Remarks on code availability)\nI could not review the code as it was bundled as a .rar file which my institution does not permit opening or extracting files from.", "reviewer_3": "In this article, the authors carried out a global meta-analysis to investigate the relationships between species richness and ecological uniqueness (RUR), defined as the local contribution to beta diversity (LCBD), and to study the role of species replacement and environmental covariates in determining the direction and strength of RUR across taxa and ecosystems.\nThe methodology of the meta-analysis is well explained and summarized in Figure 2.\nThe resulting dataset is very interesting for studying patterns of community composition across taxa and ecosystems.\nThe authors refer to the local beta diversity contribution (LCBD) as defined by Legendre and De Cáceres (2015).\nIt can be calculated using the beta.div function in the R adespatial package, which I will do to illustrate my issues below.\nI'm not sure that the questions and hypotheses (as illustrated in Figure 1) that are addressed here are sufficiently relevant, and that the way of analyzing and interpreting RUR patterns can support the authors' conclusions.\nMore specifically, my main concern is that the authors have not taken into account the influence of dispersion limitation on RUR.\nThis influence is only mentioned late in the Discussion (L227): \"Moreover, habitat isolation can also have non-negligible effect on the RURs since it limits the capacity of species to disperse between local communities, which are highly correlated with species replacement\".\nIt is known that dispersal limitation leads to a decrease in beta diversity, e.g., as a function of distance (Condit et al. 2002), which in turn should influence LCBD variation in assemblages.\nHere's a specific example:\nConsider a neutral case such as the Hubbell's (2001) mainland-island metacommunity model. The following R code simulates many neutral communities connected to a common mainly, with varying migration rates.\n#devtools::install_github(\"frmunoz/ecolottery/pkg\", build_vignettes = TRUE)\nrequire(ecolottery)\n# Building a neutral mainland assemblage, as in Hubbell (2001)\nJ <- 500;\npool <- coalesc(100*J, theta = theta)$com\n# Simulate a number of neutral communities of same size J=500, for varying migration rate from the mainland, between 0.002 and 1, with 1/500 increment\nm <- 1/500\ncomm1 <- list()\nfor(i in 1:500)\ncomm1[[i]] <- coalesc(J, m*i, pool=pool)\n# Let us calculate the richness of simulated communities\nS1 <- unlist(lapply(comm1, functon(x) nlevels(as.factor(x$com$sp))))\nplot(S1,(1:500)/500, ylab=\"m\")\n# We note here that richness increases with decreasing dispersal limitation, as expected\n# See Reviewer #3 Review Attachment #1\n# Making a community x species table\ntab1 <- matrix(0,nrow=500,ncol=Spool)\ncolnames(tab1) <- levels(as.factor(pool$sp))\nfor(i in 1:500) {\nsad.loc1 <- table(comm1[[i]]$com$sp)\ntab1[i,names(sad.loc1)] <- sad.loc1\n}\n# Calculating the Local Contributions to Beta Diversity of the simulated communities\ncontrib1 <- adespatial::beta.div(tab1)$LCBD\nplot(contrib1, S1,log=\"xy\")\n# See Reviewer #3 Review Attachment #2\ncor.test(contrib1, S1, method=\"spearman\")\n# Spearman's rank correlation rho\n#\n#data: contrib1 and S1\n#S = 39535164, p-value < 2.2e-16\n#alternative hypothesis: true rho is not equal to 0\n#sample estimates:\n# rho\n#-0.8976955\nConsequently, there may be a strongly negative RUR (shown in Reviewer #3 Review Attachment #2) for purely neutral communities that differ only in immigration intensity.\nTo provide a convincing analysis of the drivers of RUR in taxa and ecosystems, the study should thus acknowledge the influence of dispersal limitation and determine/discuss the extent to which dispersal limitation and environmental variation may influence RUR.\nReferences\nCondit, Richard, Nigel Pitman, Egbert G. Leigh Jr., Jérôme Chave, John Terborgh, Robin B. Foster, Percy Núñez, et al.\n« Beta-Diversity in Tropical Forest Trees ». Science 295 (janvier 2002): 666‑69.\nHubbell, S. P. The Unified Neutral Theory of Biodiversity and Biogeography. 1 vol. Princeton and Oxford: Princeton University Press, 2001.\nLegendre, P. & De Cáceres, M. Beta diversity as the variance of community data: dissimilarity coefficients and partitioning. Ecol. Lett. 16, 951-963 (2013)."}, "author_response": "Response to reviewers\nReviewer #1 (Remarks to the Author):\nI have reviewed the MS entitled Meta-analysis reveals widespread negative\nassociations between species richness and ecological uniqueness across the globe.\nThe MS is generally well written and easy to follow. However, it lacks some important\nmethodological details and coverage of the relevant literature.\nGeneral comments\nOn one hand, I applaud the exercise of evaluating the richness-uniqueness relationship\n(RUR) on many datasets across ecosystems and taxonomic groups. The main result that\nthe RUR is generally negative is an important one. One the other hand, I believe that\nthe paper makes a poor job of testing the many different processes that could affect the\ndirection and strength of the RUR. Instead, the paper focusses on species replacement\n(LCBD_repl) as the core process. In particular, the paper fails to revisit the results of\nDansereau, Legendre and Poisot (2022), who extensively tested the properties of the\nLCBD metric on a large dataset of bird communities.\nDansereau, G., Legendre, P. and Poisot, T., 2022. Evaluating ecological uniqueness\nover broad spatial extents using species distribution modelling. Oikos, 2022(5),\np.e09063.\nResponse: We sincerely appreciate your insightful comments and suggestions.\nFollowing your suggestions, we have included more detailed methodological\ndescriptions, such as the definition of RUR forms and the selection of explanatory\nvariables (See Methods and Supplementary Text 5). We have also incorporated\nadditional key literatures (including the study by Dansereau et al. 2022) by re-\nconducting the literature search and screening, to strengthen the contextual\nfoundation of our work.\nSpecifically, we re-conducted the literature search and screening by searching\nfor peer-reviewed journal articles using the ISI Web of Science with\npredefined keyword combinations to obtain an extensive collection of studies\nthat focused on ecological uniqueness globally (See Supplementary Figure S1\nand Supplementary Text 1 for details). Following the PRISMA-guided screening\nprocess, we finally retained 100 qualified studies, which comprised 56 case studies\nexplicitly investigating ecological uniqueness and 44 studies containing\ncompatible observational datasets from literature studies. A PRISMA diagram\nsummarizing the workflow of literature search is included in Supplementary Figure\nS1.\nFigure S1 PRISMA diagram showing the process of locating Peer-reviewed literature studies\nincluded in this synthesis. N: number of studies.\nRegarding the processes that could affect the variation of RUR, we have tested\nthe effects of four distinct hypotheses: (1) species pool, (2) dispersal limitation,\n(3) environmental filtering and (4) spatial grain size and their relative\nimportance on RUR in the current manuscript (Lines 52-100). We selected\nseven predictors related to these hypotheses in the analyses (See Section:\nQuantifications of different hypotheses in Supplementary Text 5 for details).\nOverall, we found that PA (presence-absence) -based RURs were primarily affected\nby species pool attributes (i.e. the proportion of rare species and gamma diversity),\nwhereas AB (Abundance) -based RURs were more strongly influenced by dispersal\nlimitation (See Figures 3-5 in the main text for details).\nReferences:\nDansereau, G., Legendre, P. & Poisot, T. Evaluating ecological uniqueness over\nbroad spatial extents using species distribution modelling. Oikos 2022, e09063\n(2022).\nIt is clear from Dansereau et al.’s results that the strength and direction of RUR depends\non the richness of the species pool, the breadth of the richness gradient across site, the\nspatial extent of the sampling, as well as the proportion of rare species in the pool\n(Dansereau et al., 2022). Furthermore, their results clearly show that the canonical\nRUR is not linear (rather strongly concave), and therefore Pearson’s correlation is\nperhaps not a proper approximation of its strength and direction.\nResponse: Thanks very much for your constructive comments. In this revised\nversion, we carefully reviewed the study by Dansereau et al. (2022) and\nincorporated it into the initial dataset of literature review.\nFor the first comment on the incomplete analysis in our study when comparing\nwith Dansereau et al (2022), we re-constructed our synthesis to explore the\nrole of four different hypotheses in driving RUR following your suggestions\nand others (i.e. species pool, dispersal limitation, environmental filtering, and\nspatial grain size; Lines 52-100). In this case, “the regional species pool size, the\nbreadth of species richness, and the proportion of rare species” you mentioned were\nused to represent species pool in this new version. Specifically, the regional species\npool size was represented by Gamma diversity (i.e. total number of species within\neach study). For breadth of species richness, we calculated coefficient of variation\nof species richness among different communities (CV_richness) for each\nobservation, with larger CV_richness indicating wider breadth of species richness\nwithin each observation. The proportion of rare species was defined through dual\ncriteria corresponding to different data types (i.e. presence-absence (PA) and\nabundance (AB); See Section: Quantifications of different hypotheses in\nSupplementary Text 5 for details).\nAs for the spatial extent of the sampling you have pointed out, we found that spatial\nextent data were often unavailable in most studies. Moreover, it was challenging to\ndefine the overall spatial extent of the study area for certain study groups (e.g., fish\nand freshwater macroinvertebrates), which limited our ability to assess its influence\non RUR (Lines 343-350). Therefore, we alternatively considered the role of\nspatial grain size as another important component of spatial scale (See Barton\net al., 2013). It was defined as the sampling area of each study unit within each\nstudy if reported (See Section: Quantifications of different hypotheses in\nSupplementary Text 5 for details). We also highlighted that future studies would\nbenefit from more standardized and explicit definitions of sampling scales,\nenabling systematic analyses of scale-dependent RUR patterns using more\nextensive and integrative datasets (Lines 353-355).\nOur results revealed that species pool and spatial grain size were important\ndrivers of RUR (See Figures 3-4 in the main text). For species pool, the effect size\nof RUR consistently declined with increasing gamma diversity (Figure 3a-b in the\nmain text; Supplementary Table S8). By contrast, the effect size of RUR increased\nwith CV_richness and the proportion of rare species (especially in PA data; Figure\n3 in the main text; Supplementary Table S8). For spatial grain size, it negatively\naffected effect size of RUR although effects were generally weak in both single\nfactor (Figure 4a-b in the main text, Supplementary Table S8) and full model\nanalyses (Figure 5a-b in the main text; Supplementary Table S16b).\nFigure 3 Variation of the effect sizes of RUR (i.e. Fisher’s Z statistics) with gamma diversity\n(a-b), CV_richness (c-d), proportion of rare species (e-f), and dispersal limitation (g-h)\nbased on presence absence (PA) and abundance (AB) data types. The regression lines represent\nthe response to the above factors with the 95% confidence intervals shaded from hierarchical\nmixed-effects meta-regressions, and the values of the estimated coefficients (β) and 95%\nconfidence intervals are shown in each panel. The solid and dashed lines represent significant\n(P < 0.05) and non-significant (P > 0.05) relationships. The points have radii relative to the\ninverse variance of each observation (i.e., observations with larger points have lower variance\nin effect size). Gamma diversity was log-transformed before modelling. CV_richness: the\nbreadth of species richness. See Table S8 for details of model statistics. k: Number of effect\nsizes.\nFigure 4 Variation of the effect sizes of RUR (i.e. Fisher’s Z statistics) with spatial grain size\n(a-b), Clim PC1 (c-d) and Clim PC2 (e-f) based on presence absence (PA) and abundance (AB)\ndata types. The regression lines represent the response to the above factors with the 95%\nconfidence intervals shaded from hierarchical mixed-effects meta-regressions, and the values\nof the estimated coefficients (β) and 95% confidence intervals are shown in each panel. The\nsolid and dashed lines represent significant (P < 0.05) and non-significant (P > 0.05)\nrelationships. The points have radii relative to the inverse variance of each observation (i.e.,\nobservations with larger points have lower variance in effect size). Spatial grain size was log-\ntransformed before modelling. Clim PC1 and Clim PC2: The first two principal components\nof bioclimatic variables from the WorldClim database. See Table S8 for details of model\nstatistics. k: Number of effect sizes.\nFigure 5 (a-b) Explanatory power of predictor variables on effect sizes of RUR from multiple\nhierarchical mixed-effects meta-regressions, and (c-d) the relative importance of multiple\nfactors in determining the effect size of RUR across four major taxa groups (i.e. Terrestrial\nplants, freshwater macroinvertebrates, birds and terrestrial insects) based on presence absence\n(PA) and abundance (AB) data types. The effect of each factor is considered to be significant\nif the 95% confidence intervals did not overlap zero. Each factor is standardized before\nmodelling and the relative importance of each factor is expressed as the percentage of variance\nexplained. Gamma diversity and spatial grain size were log-transformed before modelling.\nCV_richness: the breadth of species richness. Clim PC1 and Clim PC2: The first two principal\ncomponents of bioclimatic variables from the WorldClim database. See Table S16a, S16b for\ndetails of model statistics. R2 : the amount of variance explained by multiple fixed factors. k:\nm\nNumber of effect sizes.\nFor the second comment regarding different formats (i.e. linear or nonlinear)\nof RUR, we identified linear (first-order) and quadratic (second-order)\nmodels as two predominant forms of RUR characterization through a\nsystematic review of previous researches (See Appendix 1) and compared these\ntwo formulations across all studies (See section: Effect size of RUR in Methods;\nLines 438-457). Our analysis revealed that 32.1% of records (616 out of 1918)\nexhibited significantly better fit with quadratic model than linear model (P < 0.05;\nSee Appendix 2), while 67.9% (1302 records) showed no statistical difference\nbetween linear and quadratic formulations (i.e. P > 0.05). For these cases, the\ndifferences in AIC or AICc between linear and quadratic models were also less than\n2 units, further supporting their statistical equivalence (Burnham & Anderson 2004;\nAppendix 2). These results suggested that the linear model provided a more\nparsimonious representation of RUR patterns for these 1302 observations (Lines\n457-464). Given the challenges in deriving standardized effect sizes for non-linear\n(quadratic) relationships of RUR, we restricted subsequent analyses to the 1302\nrecords best represented by linear models to ensure methodological consistency\nand comparability across studies (Lines 468-472).\nReferences:\nBarton, P. S. et al. The spatial scaling of beta diversity. Glob. Ecol. Biogeogr. 22,\n639–647 (2013).\nBurnham, K. P. & Anderson, D. R. Multimodel inference: Understanding AIC and\nBIC in model selection. Sociol. Methods Res. 33, 261–304 (2004).\nSpecific comments\nIt is not clear if the meta regression aims at testing the “strength” OR “direction” of\nthe RUR as both terms are used in the introduction (line 39 and line 71 of the\nintroduction) as well as in the results. There may be different processes affecting the\nstrength or the direction. These differences are not well explained (but see Dansereau\net al., 2022).\nJudging from Figure 4, there would be no association between Rate of species\nreplacement and RUR if only strength was tested (absolute Fisher Z). But then, Fisher\nZ is a metric mostly related to strength (Pearson’s r) rather than direction. This is\na problem.\nResponse: Thank you for these thoughtful and constructive comments. In the\nrevised manuscript, we focused on characterizing the overall distribution of RUR\nand identifying their ecological drivers. Specifically, we aimed to show the patterns\nthat significant negative RURs are common while positive RURs are relatively rare\n(i.e. we mainly focused on the direction in this study). To avoid ambiguity, we no\nlonger frame our analyses in terms of testing the “strength” of RUR.\nIn addition, we agree with your concerns related to the use of species replacement\nin explaining RUR in Figure 4 of former manuscript. Therefore, we have removed\nthis part of the analysis in the manuscript. We alternatively incorporated four\ndifferent hypotheses: species pool (i.e. regional species pool size; breadth of\nspecies richness and proportion of rare species), dispersal limitation, environmental\nfiltering, and spatial grain size in the current manuscript (Lines 52-100). Then, we\nmainly explored the effects of these hypotheses in determining RUR and used\nFisher's Z statistics as the effect size of RUR rather than to directly infer strength\nper se.\nWe find that variations of RUR effect sizes are associated with several factors, such\nas species pool size (i.e. gamma diversity), breadth of species richness and the\nproportion of rare species. The effect sizes of RUR tend to decrease with increasing\ngamma diversity, indicating a stronger tendency toward negative associations\n(Figure 3a-b in the main text). In contrast, the effect sizes of RUR tend to increase\nwith higher proportion of rare species, suggesting a potential shift toward more\npositive associations (Figure 3e-f in the main text). These clarifications are now\nreflected particularly in the discussions (Lines 244-247; Lines 257-259).\nIt is also not clear why should we expect changes in the strength or direction of the\nRUR across taxa and ecosystems. Are factors such as the size of the species pool, the\nbreadth of the richness gradient, the spatial extent, or the proportion of rare species,\nexpected to vary across taxa and ecosystems, thus explaining some of the differences?\nProbably. These are not detailed.\nResponse: Thank you for pointing this out. We did not pursue comparisons of\nRUR across ecosystems and taxonomic groups, as RUR is, as you rightly\npointed out, more strongly influenced by factors such as species pool size, the\nbreadth of species richness, the proportion of rare species, and spatial grain\nsize. Therefore, we have revised the manuscript to focus on how these factors\ninfluence the overall distribution and variation of RUR effect sizes, both in general\nand across different taxonomic groups (Lines 520-564), and have largely removed\ncomparisons across taxa and ecosystems.\nWe found significant negative RUR for both PA and AB data (Figure 2 in the main\ntext). Taxon-specific analyses demonstrated significantly negative RURs\nacross different taxonomic groups in both PA and AB datasets, such as birds,\nterrestrial plants, fishes and terrestrial insects (Supplementary Figure S3; Table\nS7). We also found that the drivers of RUR differed by data type: PA-based RURs\nwere primarily shaped by attributes of species pool such as gamma diversity and\nproportion of rare species, while AB-based RURs were more strongly influenced\nby dispersal limitation. Environmental filtering and spatial grain size exhibited\nrelatively weak effects and low contributions (Figures 3-5 in the main text).\nMore fundamentally, the authors partitioned the LCBD metric (beta diversity) into its\ncomponents of species replacement (LCBD_repl) and richness difference (LCBD_rich).\nThey then tested if the strength and direction of the RUR is affected by LCBD_repl.\nFirst, it is mentioned in the MS that the effect of LCBD_rich was also evaluated (line\n437), but I could not find these results. Second, since both the slope and Pearson’s\ncorrelation of the RUR depend on LCBD values, and that LCBD_repl also depends on\nLCBD values (LCBD_tot = LCBD_repl + LCBD_rich), there is a high risk of spurious\nassociation here. A full simulation should be conducted to avoid any artifactual findings.\nResponse: For the first comment on the LCBD_rich result, we partitioned\necological uniqueness (LCBD metric) into the species replacement component (the\nlocal contributions of species replacement to beta diversity; LCBD ) and richness\nrepl\ndifference (the local contributions of richness difference to beta diversity; LCBD )\nrich\n(Legendre 2014), and the proportion explained by richness difference is equal to\none minus the proportion explained by species replacement (i.e., the sum of the\nproportions explained by species replacement and richness difference equals one).\nTherefore, we only provided the proportion explained by species replacement in\nthe previous manuscript. We no longer considered this factor in the revised\nmanuscript (See response below).\nFor the second comment, we agree with the point that there is a high risk of spurious\nassociation between LCBD and its components, LCBD and LCBD .\nrepl rich\nTherefore, we have removed this part of the analysis in the revised manuscript.\nInstead, we now focused on the effects and relative importance of four different\nhypotheses: species pool, dispersal limitation, environmental filtering, and spatial\ngrain size on the effect size of RUR (See Figures 3-4 in the main text). Because\nthe general hypotheses in this study have been changed in this version\nfollowing your comments and the comments of other reviewers (Lines 52-100),\nwe believe that the full simulation you mentioned is not necessary here.\nNevertheless, we sincerely appreciate your careful suggestions.\nReferences:\nLegendre, P. Interpreting the replacement and richness difference components of\nbeta diversity. Glob. Ecol. Biogeogr. 23, 1324-1334 (2014).\nThe above concerns are major, and therefore I cannot recommend the paper for\npublications.\nResponse: We sincerely thank you for your constructive comments and valuable\nsuggestions. In response, we have carefully and thoroughly revised the manuscript.\nWe believe that the changes made in this version adequately address both the major\nand minor concerns raised.\nReviewer #2 (Remarks to the Author):\nMAJOR CONCERNS\nIn its current state, this is not a high-quality meta-analysis and does something I’ve\nnever seen before which is only use snowballing methods to find studies for the synthesis\n(e.g., screening and including only studies that cite a key paper, here Legendre & Ce\nCaceres 2013). Typically, this approach is meant to supplement an exhaustive keyword\nsearch of multiple databases -- here they skipped the exhaustive search and went\nstraight to screening and including only studies citing Legendre & Ce Caceres (2013).\nThey also supplement sources but do not justify why these sources were chosen or\nwhether they are representative of what’s available. How many studies were missed\nwith these unusual decisions? There is also no information on who screened the studies,\nor whether there were any practices to validate screening decisions. I suggest the\nauthors have a look at PRISMA guidelines to develop practices that generate high\nquality syntheses.\nResponse: Thank you very much for this helpful comment. Following your\nsuggestions, we have improved the synthesis by incorporating three distinct\napproaches to allow for a broad and representative collection of studies or\ndatasets (See Supplementary Text 1-4 for details).\n(1) Peer-reviewed literature studies: First, we conducted a systematic search of\npeer-reviewed journal articles in the ISI Web of Science using predefined keyword\ncombinations to compile a comprehensive set of studies focused on ecological\nuniqueness globally. The keyword structure and methodological framework, which\nfollow the PRISMA guidelines, are now described in detail in Supplementary Text\n1-2 and illustrated in Supplementary Figure S1.\n(2) Public data repositories: Second, we compiled additional raw datasets across\ndifferent taxonomic groups from public data repositories, as both ecological\nuniqueness and species richness can be calculated from species composition\ndatasets (Legendre & De Cáceres, 2013). The search terms and methods used for\ndataset collection are now detailed in Supplementary Text 1, Text 3 and\nSupplementary Figure S2.\n(3) Compiled datasets: Third, we compiled a series of additional raw datasets\nfrom our own data. These datasets were mainly acquired through: (1) direct\ncorrespondence with principal investigators maintaining unpublished datasets, or\n(2) using the existing data sources including: GrassPlot (Biurrun et al., 2021;\nDengler et al., 2018), sPlot (Bruelheide et al., 2019) and CESTES (Jeliazkov et al.,\n2020) by ourselves. The detailed procedures are now described in Supplementary\nText 1 and Supplementary Text 4.\nIn summary, our final dataset incorporated a total of 1918 data records from 451\nliterature studies or datasets (Figure 1 in the main text). The final dataset included\nthree types of data sources: Peer-reviewed literatures studies (100 studies, 288 data\nrecords), public data repositories (243 datasets, 1415 data records) and compiled\ndatasets (108 datasets, 215 data records; See Figure 1 in the main text and\nSupplementary Table S1). However, for consistency in description, we uniformly\nreferred to each dataset or literature study as a \"study\" in the subsequent text (Lines\n402-407).\nFigure S1 PRISMA diagram showing the process of locating Peer-reviewed literature studies\nincluded in this synthesis. N: number of studies.\nFigure S2 PRISMA diagram showing the process of locating datasets included in this synthesis\nfrom Dryad repository. N: number of datasets.\nReferences:\nBiurrun, I. et al. Benchmarking plant diversity of Palaearctic grasslands and other\nopen habitats. J. Veg. Sci. 32, e13050 (2021).\nBruelheide, H. et al. sPlot – A new tool for global vegetation analyses. J. Veg. Sci.\n30, 161–186 (2019).\nDengler, J. et al. GrassPlot – A database of multi-scale plant diversity in Palaearctic\ngrasslands. Phytocoenologia 48, 331–347 (2018).\nJeliazkov, A. et al. A global database for metacommunity ecology, integrating\nspecies, traits, environment and space. Scientific Data 7, 6 (2020).\nLegendre, P. & De Cáceres, M. Beta diversity as the variance of community data:\ndissimilarity coefficients and partitioning. Ecol. Lett. 16, 951-963 (2013).\nThey also apply a custom weighting that will have unknown consequences for the meta-\nanalysis and publication bias tests because the variance of Fisher’s Z is meant to\nstrictly model sampling error among effect sizes (and the variability among these\nweights later, in turn, impact the random-effects estimated by the meta-analysis model),\nbut here they’ve confounded this weight and its purpose with aspects of spatial extent\nwithout justification. The meta-analysis model should not be modified like this, and if\nyou believe spatial extent impacts effect sizes then treat it as a continuous predictor in\nthe model, and then explore pooled effects as marginal means corrected by spatial\nextent.\nAlso please report your full meta-analysis model; it’s unclear what random effects were\nincluded, what predictors are combined to achieve synthesis, or whether these were\nimplemented correctly with the metafor R package. Results are also incomplete as they\nemphasize P values but do not report the statistic of these P values.\nResponse: Thank you for your suggestions. We agree that including additional\nvariables in the calculation of corresponding weights is inappropriate. Therefore,\nin this revised version, we retained Fisher’s Z transformation as the effect size\nof RUR but removed other factors from the definition of weights to ensure\nconsistency and clarity in our analysis.\nWe weighted each effect size based on the sample size (i.e., the number of sites for\nwhich of ecological uniqueness and species richness were calculated) using the\nformula:\n1\n𝑣 = (2)\n𝑖 𝑁 −3\n𝑠\nWhere 𝑁 is the sample size (i.e., number of sampling units) of each estimate\n𝑠\nwithin each study. This indicates that studies had larger sample sizes would be\ngiven higher weights (Lines 495-501).\nIn accordance with your suggestion, we included spatial grain size as an\nexplanatory variable in our analysis to account for the variation in effect sizes\nof RUR. Spatial grain size was defined as the sampling area of each study unit\nwithin each study if reported (Barton et al., 2013; See Section: Quantifications of\ndifferent hypotheses in Supplementary Text 5).\nWe agree that the previous version may not provide detailed descriptions of the\nstatistical analysis. To address this issue, we have included additional details in the\nrevised version. For the inclusion of predictors, we explored the effect and\nrelative importance of four distinct hypotheses on the effect size of RUR. They are:\n(1) species pool, (2) dispersal limitation, (3) environmental filtering and (4) spatial\ngrain size (Lines 52-100). We examined the potential influence of seven predictors\nbelong to four hypotheses including: (1) gamma diversity, (2) CV_richness, (3)\nproportion of rare species, (4) dispersal limitation, (5) Clim PC1, (6) Clim PC2 and\n(7) spatial grain size on the effect sizes of RUR across different taxa groups (See\nSection: Quantifications of different hypotheses in Supplementary Text 5 for\ndetails). We first conducted single-factor analyses for each predictor. We then\nexamined relative importance of these hypotheses by conducting a full-model\nanalyses incorporating different predictors (See Statistical analysis for details,\nLines 520-564).\nFor the inclusion of random effects, we implemented an inverse-variance-\nweighted hierarchical mixed-effects meta-regression model incorporating a multi-\nlevel random effects structure when assessing the overall mean effect size across\ndata types (PA and AB data). The model specified two independent random\neffects to account for potential sources of heterogeneity: (1) nested random\neffects of repeated individual RUR estimates within each study and (2) an\nindependent random effect for taxonomic group to address between-taxa\nheterogeneity (Lines 504-509). We also conducted subgroup analyses to examine\nRURs across different taxonomic groups. we retained only the nested random\neffect (i.e., estimates nested within studies), as taxonomic group was the grouping\nvariable being analyzed (Lines 510-514).\nIt should be noted that we conducted single-factor analyses separately for each taxa\ngroup. In these analyses, the nested random effect of estimate within study was\nretained whenever applicable. However, for taxa groups where most studies\ninclude only one estimate, we specified random effects at the study level instead\n(Lines 535-540).\nRegarding the reporting of results, we agree that presenting only P-value is\ninsufficient. Therefore, in the revised manuscript, we have included the\nregression coefficients along with their 95% confidence intervals (CIs),\nvariance estimates for random effects, marginal R2 values (explained variance\nby fixed effects), overall heterogeneity (I²) and the I² values corresponding to\neach level of random effects (See Statistical analysis for details). These results\nare presented in Figures 3–5 in the main text and Supplementary Tables S5–S17.\nThe data and R code for analysis are available in Figshare at\nhttps://figshare.com/s/a974caae4f643878afc7\nReferences:\nBarton, P. S. et al. The spatial scaling of beta diversity. Glob. Ecol. Biogeogr. 22,\n639–647 (2013).\nOn Line 131, you also cannot do a meta-analysis of absolute values of effect sizes\nbecause this moves the null from zero to an unknown null value -- making significant\ntests meaningless.\nResponse: We agree that analyzing the absolute values of Fisher’s Z statistics was\nnot appropriate in the previous version of the manuscript. In the revised manuscript,\nwe have removed the corresponding results and method descriptions. To assess the\neffect of environmental filtering, we alternatively extracted 19 bioclimatic\nvariables (Bio1-Bio19) from the WorldClim database (Hijmans et al., 2005)\nwith 30 arc seconds spatial resolution using each study's mean geographic\ncoordinates. Principal Component Analysis (PCA) was applied to mitigate\nmulticollinearity among these variables. The first two principal components (Clim\nPC1 and Clim PC2) were retained as synthetic climatic predictors for explaining\nthe effect size of RUR in subsequent analyses. (See Section: Quantifications of\ndifferent hypotheses in Supplementary Text 5 for details).\nReferences:\nHijmans, R. J., Cameron, S. E., Parra, J. L., Jones, P. G. & Jarvis, A. Very high\nresolution interpolated climate surfaces for global land areas. Int. J. Climatol. 25,\n1965–1978 (2005).\nA lot of this work (some of it is quite amazing and remarkable) is simply waste because\nodd decisions were made that make it impossible to assess the robustness of what was\nachieved.\nResponse: In this revised manuscript, we have carefully updated the data sources\nand checked the statistical methods for the synthesis as well as the corresponding\ndescriptions in the Methods section. Specifically, we (1) updated the PRISMA\nflowchart and data sources (Supplementary Text 1-4; Supplementary Figures 1-2),\n(2) revised the definition of variance (Lines 495-501), and (3) provided more\ndetailed model parameter settings (especially for the specification of random\neffects; Lines 504-514) as well as model statistics and results (Supplementary\nTables S5-S17). We believe that these issues have now been properly addressed\nand that the revised manuscript provides an accurate and appropriate description of\nthe statistical methods used (Lines 504-582).\nMINOR NOTES\n1) Lines 37-38: do you mean “have yet to be synthesized”\nResponse: We have clarified this point in the text. It now reads: “Despite perceived\nimportance, the patterns of this relationship and underlying mechanisms have yet\nto be synthesized at the global scale.” (Lines 36-37).\n2) k instead of N should be used to describe number of effect sizes (N is reserved for\nnumber of data points used to estimate Z)\nResponse: Thank you for pointing this out. We have made this correction and used\nk to represent the number of effect sizes in analyses (Figures 3-5 in the main text,\nSupplementary Tables S6-S17).\n3) report the degrees of freedom for all Q tests\nResponse: We have reported the degrees of freedom corresponding to the Q\nstatistic in the supplementary table (See Supplementary Tables S5-S17 for details).\nReviewer #2 (Remarks on code availability):\nI could not review the code as it was bundled as a .rar file which my institution does\nnot permit opening or extracting files from.\nResponse: Thank you for pointing this out. In this revised manuscript, the R code\nand the associated dataset supporting the findings of this study are available in\nFigshare at https://figshare.com/s/a974caae4f643878afc7 Meanwhile, we have\nalso clarified this issue in the Code Availability statement of the revised manuscript\n(Lines 588-590).\nReviewer #3 (Remarks to the Author):\nIn this article, the authors carried out a global meta-analysis to investigate the\nrelationships between species richness and ecological uniqueness (RUR), defined as\nthe local contribution to beta diversity (LCBD), and to study the role of species\nreplacement and environmental covariates in determining the direction and strength of\nRUR across taxa and ecosystems.\nThe methodology of the meta-analysis is well explained and summarized in Figure 2.\nThe resulting dataset is very interesting for studying patterns of community composition\nacross taxa and ecosystems.\nResponse: Thank you for this positive assessment of our work.\nThe authors refer to the local beta diversity contribution (LCBD) as defined by\nLegendre and De Cáceres (2015).\nIt can be calculated using the beta.div function in the R adespatial package, which I\nwill do to illustrate my issues below.\nResponse: Thank you for this suggestion. We applied two alternative approaches\nto calculate ecological uniqueness and found these approaches are highly\ncorrelated (Lines 420-435). Therefore, we prioritized the method based on the\nbeta.div function in the current synthesis, as it is the most commonly used in studies\nof ecological uniqueness (Chen et al., 2024; da Silva et al., 2018; Vilmi et al., 2017;\nYao et al., 2021).\nReferences:\nChen, Y. et al. Multiple processes jointly determine ecological uniqueness across\nforest plant life‐forms in Northeast China. J. Biogeogr. 51, 1133-1147 (2024).\nda Silva, P. G., Medina Hernandez, M. I. & Heino, J. Disentangling the correlates\nof species and site contributions to beta diversity in dung beetle assemblages.\nDivers. Distrib. 24, 1674-1686 (2018).\nVilmi, A., Karjalainen, S. M. & Heino, J. Ecological uniqueness of stream and lake\ndiatom communities shows different macroecological patterns. Divers. Distrib. 23,\n1042-1053 (2017).\nYao, J., Huang, J., Ding, Y., Xu, Y., Xu, H. & Zang R. Ecological uniqueness of\nspecies assemblages and their determinants in forest communities. Divers. Distrib.\n27, 454-462 (2021).\nI'm not sure that the questions and hypotheses (as illustrated in Figure 1) that are\naddressed here are sufficiently relevant, and that the way of analyzing and interpreting\nRUR patterns can support the authors' conclusions.\nResponse: We have conducted a thorough revision on the manuscript based on\nyour suggestions and comments from other reviewers. We carefully checked our\nhypotheses and analytical methods. We have removed the original hypotheses and\nthe corresponding data analysis from the previous version of the manuscript. In the\ncurrent version, we focus primarily on examining the effects and relative\nimportance of four distinct hypotheses: (1) species pool; (2) dispersal limitation;\n(3) environmental filtering and (4) spatial grain on RUR (Lines 52-100).\nAccordingly, we have also revised the analytical methods to align with the updated\nhypotheses (See Statistical analysis in the main text for details). After this revision,\nwe believe the data analysis and hypotheses about RUR are clearer and more\nrelevant.\nMore specifically, my main concern is that the authors have not taken into account the\ninfluence of dispersion limitation on RUR.\nThis influence is only mentioned late in the Discussion (L227): \"Moreover, habitat\nisolation can also have non-negligible effect on the RURs since it limits the capacity of\nspecies to disperse between local communities, which are highly correlated with species\nreplacement\".\nIt is known that dispersal limitation leads to a decrease in beta diversity, e.g., as a\nfunction of distance (Condit et al. 2002), which in turn should influence LCBD variation\nin assemblages.\nHere's a specific example:\nConsider a neutral case such as the Hubbell's (2001) mainland-island metacommunity\nmodel. The following R code simulates many neutral communities connected to a\ncommon mainly, with varying migration rates.\n#devtools::install_github(\"frmunoz/ecolottery/pkg\", build_vignettes = TRUE)\nrequire(ecolottery)\n# Building a neutral mainland assemblage, as in Hubbell (2001)\nJ <- 500;\npool <- coalesc(100*J, theta = theta)$com\n# Simulate a number of neutral communities of same size J=500, for varying migration\nrate from the mainland, between 0.002 and 1, with 1/500 increment\nm <- 1/500\ncomm1 <- list()\nfor(i in 1:500)\ncomm1[[i]] <- coalesc(J, m*i, pool=pool)\n# Let us calculate the richness of simulated communities\nS1 <- unlist(lapply(comm1, functon(x) nlevels(as.factor(x$com$sp))))\nplot(S1,(1:500)/500, ylab=\"m\")\n# We note here that richness increases with decreasing dispersal limitation, as expected\n# See Reviewer #3 Review Attachment #1\n# Making a community x species table\ntab1 <- matrix(0,nrow=500,ncol=Spool)\ncolnames(tab1) <- levels(as.factor(pool$sp))\nfor(i in 1:500) {\nsad.loc1 <- table(comm1[[i]]$com$sp)\ntab1[i,names(sad.loc1)] <- sad.loc1\n}\n# Calculating the Local Contributions to Beta Diversity of the simulated communities\ncontrib1 <- adespatial::beta.div(tab1)$LCBD\nplot(contrib1, S1,log=\"xy\")\n# See Reviewer #3 Review Attachment #2\ncor.test(contrib1, S1, method=\"spearman\")\n# Spearman's rank correlation rho\n#\n#data: contrib1 and S1\n#S = 39535164, p-value < 2.2e-16\n#alternative hypothesis: true rho is not equal to 0\n#sample estimates:\n# rho\n#-0.8976955\nConsequently, there may be a strongly negative RUR (shown in Reviewer #3 Review\nAttachment #2) for purely neutral communities that differ only in immigration intensity.\nTo provide a convincing analysis of the drivers of RUR in taxa and ecosystems, the\nstudy should thus acknowledge the influence of dispersal limitation and\ndetermine/discuss the extent to which dispersal limitation and environmental variation\nmay influence RUR.\nReferences\nCondit, Richard, Nigel Pitman, Egbert G. Leigh Jr., Jérôme Chave, John Terborgh,\nRobin B. Foster, Percy Núñez, et al. « Beta-Diversity in Tropical Forest Trees ». Science\n295 (janvier 2002): 666-69.\nHubbell, S. P. The Unified Neutral Theory of Biodiversity and Biogeography. 1 vol.\nPrinceton and Oxford: Princeton University Press, 2001.\nLegendre, P. & De Cáceres, M. Beta diversity as the variance of community data:\ndissimilarity coefficients and partitioning. Ecol. Lett. 16, 951-963 (2013).\nResponse: Thank you for this insightful suggestion. Following the suggestions\nfrom you and other reviewers, we have incorporated four distinct hypotheses into\nthis synthesis: (1) species pool, (2) dispersal limitation (as you specifically\nsuggested), (3) environmental filtering, and (4) spatial grain size (Lines 52-100).\nWe explored these hypotheses and their corresponding factors to assess their\npotential impacts and relative importance on the effect sizes of RUR (See Section:\nQuantifications of different hypotheses in Supplementary Text 5 and Statistical\nanalysis for details).\nTo assess the potential influence of dispersal limitation on RUR, we employed the\nmethod of distance-based Moran’s Eigenvector Maps (dbMEM) to model spatial\nstructure across study sites (See Section: Quantifications of different hypotheses in\nSupplementary Text 5 for details). Our results showed that dispersal limitation\nexerted significantly negative effects on RUR (Figure 3g-h in the main text;\nSupplementary Table S8) especially for freshwater macroinvertebrates and\nterrestrial insects (Supplementary Table S12). The same trend persisted after\naccounting for other variables in the full multivariate models (Figure 5a-b in the\nmain text; Supplementary Table S16b), highlighting the robustness of this negative\neffect. For relative contributions, we found that dispersal limitation (22.7%) was\nthe primary determinant of RUR based on AB data while it was also an important\ndriver of PA-based RUR (Figure 5c-d in the main text; Supplementary Table S16b).\nFigure 3 Variation of the effect sizes of RUR (i.e. Fisher’s Z statistics) with gamma diversity\n(a-b), CV_richness (c-d), proportion of rare species (e-f), and dispersal limitation (g-h) based\non presence absence (PA) and abundance (AB) data types. The regression lines represent the\nresponse to the above factors with the 95% confidence intervals shaded from hierarchical\nmixed-effects meta-regressions, and the values of the estimated coefficients (β) and 95%\nconfidence intervals are shown in each panel. The solid and dashed lines represent significant\n(P < 0.05) and non-significant (P > 0.05) relationships. The points have radii relative to the\ninverse variance of each observation (i.e., observations with larger points have lower variance\nin effect size). Gamma diversity was log-transformed before modelling. CV_richness: the\nbreadth of species richness. See Table S8 for details of model statistics. k: Number of effect\nsizes.\nImportantly, dispersal limitation had greater explanatory power in AB models\n(Figure 5c-d in the main text; Supplementary Table S16b), likely because dispersal\nlimitation governs both range sizes and abundance gradients (Condit et al., 2000).\nFor instance, poor dispersers (e.g., wind-dispersed pioneer plant) may form dense,\nlocalized populations, with sharp abundance gradients along sites (Condit et al.,\n2000; Wang et al., 2010). while superior dispersers (e.g., large seeds for animal-\ndispersed trees) may achieve broad but shallow abundance distributions, thereby\nhomogenizing community abundance patterns (Wang et al., 2010). AB metrics may\ncapture these patterns more effectively than PA data (Lines 292-299).\nFigure 5 (a-b) Explanatory power of predictor variables on effect sizes of RUR from multiple\nhierarchical mixed-effects meta-regressions, and (c-d) the relative importance of multiple\nfactors in determining the effect size of RUR across four major taxa groups (i.e. Terrestrial\nplants, freshwater macroinvertebrates, birds and terrestrial insects) based on presence absence\n(PA) and abundance (AB) data types. The effect of each factor is considered to be significant\nif the 95% confidence intervals did not overlap zero. Each factor is standardized before\nmodelling and the relative importance of each factor is expressed as the percentage of variance\nexplained. Gamma diversity and spatial grain size were log-transformed before modelling.\nCV_richness: the breadth of species richness. Clim PC1 and Clim PC2: The first two principal\ncomponents of bioclimatic variables from the WorldClim database. See Table S16a, S16b for\ndetails of model statistics. R2 : the amount of variance explained by multiple fixed factors. k:\nm\nNumber of effect sizes.\nIn addition, we have added text to the Discussion to highlight how the effect of\ndispersal limitation on RUR may also depend on several dispersal-related traits (e.g.\nphysiological or morphological traits) of different taxa groups (Alzate & Onstein,\n2022) (Lines 356-362)\nReferences:\nAlzate, A. & Onstein, R. E. Understanding the relationship between dispersal and\nrange size. Ecol. Lett. 25, 2303–2323 (2022).\nCondit, R. et al. Spatial patterns in the distribution of common and rare tropical\ntree species: a test from large plots in six different forests. Science 288, 1414–1418\n(2000).\nWang, X., Ye, J., Li, B., Zhang, J., Lin, F. & Hao, Z. Spatial distributions of species\nin an old-growth temperate forest, northeastern China. Can. J. For. Res. 40, 1011–\n1019 (2010)."}, {"round": 2, "reviewer_comments": {"reviewer_2": "This new version significantly improves clarity on what was achieved, and because the approach applies more established methodology (e.g., broad search and inclusion of the primary literature, conventional weighting, PRISMA methodology for reporting), it is also easier to assess the robustness of findings. That said, there are still ways to make things more repeatable and clear.\nMore work needs to be done on how studies were found and selected for synthesis. The goal is to be clear so that the work presented is repeatable; but also that the correlations analyzed is representative of what’s been done. In particular that null and positive correlations were extracted neutrally without bias for negative correlations (publication bias tests help with this, but should be performed among the three databases to assess individual retrieval bias). It is great that the synthesis now includes correlations from primary literature; however, there are gaps on how these new data were retrieved. Given the many authors, how many people were involved in screening studies for inclusion? Was work divided among members, or was there a dual screening design used to improve repeatability of screening decisions? Surely some of the new studies did not explicitly report raw correlations, but presented results as regression coefficients, or slope t-tests, or maybe data were only available in figures. What methods were used to retrieve these more difficult, yet highly relevant correlations? There must have also been some overlap between the different approaches to find correlations—were studies identified in the web of science search that also were found in the figshare database? How were duplicate studies handled? All three PRISMA plots should be merged as a single plot; since the end-product is a combined synthesis of all these sources. Further, given that correlations can originate from three sources (e.g., primary literary, custom databases, etc.), a sensitivity analysis is necessary to assess whether they represent a homogeneous group and that the data can be combined under a single analysis.\nThe first section of the results is essentially a type of vote-counting where they enumerate the number of significant and non- significant studies (e.g., lines 120 to 122, or parts of Fig. 2). Don’t do this, this undermines the whole purpose of doing a meta-analysis—which carefully does not use significance of each study but rather the magnitude, direction, and sampling error to make inferences. It’s ok to report the composition of the dataset (e.g., number of studies per taxa group, or among categories of responses relative to richness), but not differences in significance. The rank order of effects to emphasize representation of negative effects is ok, but discussing their prevalence without the lens of meta-analysis is dangerous as many studies may be underpowered and visually may not be representative of the underlying (aggregate) effect. Having correlations as Fisher-Z transformed effects is better for clarity since results are analyzed with the metric (e.g., main inferences are derived from these) and consistence helps with understanding what was done in the paper. The results also need to define the acronyms (PA, AB) and more descriptive words should be used for CV_richness, Clim_PC1, and Clim_PC2. Since some details about the methodology are embedded outside the main manuscript, it’s important that they are defined so that the reader knows right-way on how to interpret what was reported. These include the above acronyms, but also measure of % explained (line 176 but explained using different definitions in line 568). Many of the results reported are described as “affects”, for example this “affected” that. If a moderator is being used to explain variation among effects, maybe “predict” is a better inferential description, since correlations are just associations and not causes, and similarly, moderators are predictors and not causes of variability. Finally, lines 517-518 state that significance was tested by 95% overlap, yet there are only p-values reported throughout the results; what is the test statistics these p-values, and their values with degrees of freedom need to be reported along with p-values.\nMinor comments:\n1) What is the scale of the taxonomic group? Are they Linnean rankings, or homologous groups?\n2) Line 194, also 403: “literature studies” should be “primary studies”\n3) Line 409: Define what a “unit” is.\n4) Line 439: What do you mean by “previous researches”?", "reviewer_4": "This manuscript aims to describe the relationship between, and drivers of, species richness and ecological uniqueness.\nGenerally, the manuscript reads very well and the comments of reviewer 1 and 3 are well addressed. However, when the quadratic models are tested (as suggested by reviewer 1) and found more appropriate compared to the linear models in 30% of the cases, these records are not considered as the results are harder to interpret. I think that the better fit of the quadratic relationship in such a high percentage of the records show that the relationship can not always be interpreted as linear and that these cases should be considered and not disregarded as the results are harder to interpret. In that way you are only using the results which are convenient, while this is not a good representation of what the data actually shows. Additionally, different drivers of the RURs are tested, following four different hypotheses. In the methods, I recommend a more thorough description of why and how well the variables chosen reflect these particular mechanisms underlying the proposed hypotheses. Generally, I strongly recommend to first run a fully randomised null model, to check if the obtained relationships, and therefore your key results, are obtained by chance or that ecological mechanisms drive this relationship. Lastly, from a traits perspective, beta diversity can be different between sites, but some sites are more ecologically unique then others as the niche and traits of the species differ more compared to the other sites. This is a point which should be considered in the introduction and discussion."}, "author_response": "Response to reviewers\nReviewer #2 (Remarks to the Author):\nThis new version significantly improves clarity on what was achieved, and because the\napproach applies more established methodology (e.g., broad search and inclusion of\nthe primary literature, conventional weighting, PRISMA methodology for reporting), it\nis also easier to assess the robustness of findings. That said, there are still ways to make\nthings more repeatable and clear.\nResponse: Thank you for this positive assessment of our work.\nMore work needs to be done on how studies were found and selected for synthesis. The\ngoal is to be clear so that the work presented is repeatable; but also that the correlations\nanalyzed is representative of what’s been done.\nResponse: Thank you for this valuable suggestion. In this revision, we have revised\nthe Methods section to provide a more detailed description of our study selection\nprocess. This includes specifying the number of authors involved in screening and\nhow any discrepancies were resolved to ensure consistency (See Supplementary\nText 1 for details). Furthermore, we have now listed all excluded primary studies\nand datasets under each of our pre-defined criteria, enhancing the transparency and\nreproducibility of our workflow (See Supplementary Text 5 for details).\nAdditionally, to address the representativeness of the synthesized data, we\nconducted publication bias tests on the effect sizes from different data sources (i.e.\nprimary studies, data repositories and compiled datasets). The results confirmed\nthat no significant publication bias was detected across the various data\nsources (See Supplementary Figures S5-S8), so the correlations analyzed here\nare representative.\nIn particular that null and positive correlations were extracted neutrally without bias\nfor negative correlations (publication bias tests help with this, but should be performed\namong the three databases to assess individual retrieval bias).\nResponse: Thank you for the suggestion. We have conducted publication bias tests\non the effect sizes from different data sources (i.e. primary studies, data repositories\nand compiled datasets). The results indicate that there is no significant\npublication bias in the effect sizes obtained from any of the data sources (See\nSupplementary Figures S5-S8).\nIt is great that the synthesis now includes correlations from primary literature; however,\nthere are gaps on how these new data were retrieved. Given the many authors, how\nmany people were involved in screening studies for inclusion? Was work divided among\nmembers, or was there a dual screening design used to improve repeatability of\nscreening decisions?\nResponse: Thank you for this important question regarding the screening process.\nTo ensure the reliability and reproducibility of our study screening and selection,\nwe indeed adopted a dual-independent screening approach. Specifically, two\nauthors independently screened all the identified studies and datasets against the\npre-defined exclusion criteria. This process was conducted simultaneously to\nminimize bias. These clarifications have been added to Supplementary Text 1:\n“Following the PRISMA guidelines, two authors independently conducted the\nliterature and dataset screening based on the above exclusion criteria. The two\nauthors examined 3 primary studies and 4 datasets that were only selected by one\nof the authors. Following a joint re-assessment of the studies and datasets with\ninitial screening discrepancies between the two authors, a full consensus was\nreached to exclude all of them. We have listed all unselected literature and datasets\nwith the reasons for exclusion (See Supplementary Text 5).”\nWe believe such a rigorous, dual-independent screening method significantly\nstrengthens the validity of our inclusion process and the overall findings of our\nsynthesis.\nSurely some of the new studies did not explicitly report raw correlations, but presented\nresults as regression coefficients, or slope t-tests, or maybe data were only available in\nfigures. What methods were used to retrieve these more difficult, yet highly relevant\ncorrelations?\nResponse: Thank you for this comment. As you pointed out, some studies\npresented results as regression coefficients, t-statistics, or in figures. While it is\nmethodologically possible to convert such statistics into Pearson's correlation\ncoefficients (Nakagawa & Cuthill, 2007), we made a conscious decision to\nexclude studies that provided only these summary statistics without the\nunderlying raw data or values for species richness and ecological uniqueness.\nThis decision was driven by a critical methodological consideration raised by\nReviewer #1, who noted that the RUR may not always be linear (e.g., it could be\nquadratic). Because 1) to empirically test the form of the relationship for each study,\nit was essential for us to have access to the data distributions; 2) relying solely on\nconverted summary statistics would have presupposed a linear association and\nprevented us from verifying this key assumption; in this case, we prioritized data\nintegrity and the ability to assess the relationship's form rather than simply\nmaximize sample size. This specific inclusion criterion is explicitly stated in our\nscreening protocol (See Supplementary Text 1):\n“(3) We excluded studies that did not report the distributions of both species\nrichness and ecological uniqueness values, or studies that only report related\nstatistics (e.g. regression coefficients, slope t-tests or model estimates) without\nproviding richness or uniqueness distributions. Although it is methodologically\nfeasible to convert such statistics into correlation coefficients (Nakagawa & Cuthill,\n2007), this exclusion was necessary because previous empirical evidence\ndemonstrates that the RUR may exhibit either linear or nonlinear (e.g. curve-linear)\nrelationships (Dansereau et al., 2022; Luukkonen et al., 2024; Hill et al., 2021).\nWithout detailed distributions of species richness and ecological uniqueness,\nwe could not reliably characterize the association form of RUR. To empirically test the\nform of the relationship for each study, it was essential for us to have access to the\ndata distributions. Moreover, relying solely on converted summary statistics would\nhave presupposed a linear association and prevented us from verifying this key\nassumption. Such incomplete reporting would fundamentally limit our capacity to\ndistinguish between alternative ecological models while introducing substantial\nuncertainty in the synthesis (Reason 3, 186 studies; Fig. S1)”.\nReferences:\nDansereau, G., Legendre, P. & Poisot, T. Evaluating ecological uniqueness over\nbroad spatial extents using species distribution modelling. 2022, Oikos, e09063\n(2022).\nHill, M. J. et al. Local contributions to beta diversity in urban pond networks:\nImplications for biodiversity conservation and management. Divers. Distrib. 27,\n887-900 (2021).\nLuukkonen, S., Heino, J., Hjort, J., Helm, A. & Alahuhta, J. A macroecological\nanalysis of ecological uniqueness of freshwater macrophyte assemblages across\nEurope and North America. J. Ecol. 2024, 1–15 (2024).\nNakagawa, S. & Cuthill, I. C. Effect size, confidence interval and statistical\nsignificance: a practical guide for biologists. Biol. Rev. Camb. Philos. Soc. 82, 591–\n605 (2007).\nThere must have also been some overlap between the different approaches to find\ncorrelations—were studies identified in the web of science search that also were found\nin the figshare database? How were duplicate studies handled? All three PRISMA plots\nshould be merged as a single plot; since the end-product is a combined synthesis of all\nthese sources.\nResponse: Thank you for this important comment. Indeed, there was an overlap of\nstudies identified through the literature search in Web of Science and dataset search\nfrom the Dryad repositories. All duplicate records were meticulously identified and\nremoved to ensure each unique dataset was represented only once in our synthesis.\nThis is now stated in exclusion criteria of public data repositories:\n“(5) We also excluded studies sharing identical datasets (i.e., the same dataset\nused in primary studies or multiple data publications containing the same datasets)\n(Reason 5, 11 datasets; Figure S1).” (See Supplementary Text 1 for details)\nThe number of excluded duplicates is clearly documented in the merged\nPRISMA flow diagram (See Supplementary Figure S1), which now integrates the\nsearch and selection processes for data sources including Web of Science and\nDryad repositories, into a single, cohesive figure as recommended.\nFigure S1 PRISMA diagram showing the process of locating Peer-reviewed primary studies\nfrom ISI Web of Science and datasets from Dryad repository included in this synthesis. N:\nnumber of studies.\nFurther, given that correlations can originate from three sources (e.g., primary literary,\ncustom databases, etc.), a sensitivity analysis is necessary to assess whether they\nrepresent a homogeneous group and that the data can be combined under a single\nanalysis.\nResponse: Thank you for this suggestion, we have assessed whether the overall\neffect sizes of RUR differed significantly across different data sources (i.e. primary\nstudies, data repositories and compiled datasets). The results demonstrate that\nthe overall effect sizes of RUR were consistently and significantly negative\nacross all three data sources, and no significant differences were found among\nthem either for presence absence or abundance data types (See Supplementary\nFigure S3 and Table S7a-b). This indicates that the heterogeneity originating from\ndifferent data sources is low, thereby justifying their combination into a single\ndataset for our overall analysis. We have now highlighted this point in the Results\nsection (Lines 127-131).\nTable S7b Differences in effect sizes of RUR among data sources. Non-significant\nP-values (P > 0.05) suggest no differences between the data sources in their RUR\neffect sizes.\nP-value for presence P-value for abundance\nPair-wise comparison\nabsence data data\nPrimary studies v.s. Compiled datasets 0.70 0.98\nPrimary studies v.s. Data repositories 0.072 0.48\nCompiled datasets v.s. Data repositories 0.19 0.40\nThe first section of the results is essentially a type of vote-counting where they\nenumerate the number of significant and non-significant studies (e.g., lines 120 to 122,\nor parts of Fig. 2). Don’t do this, this undermines the whole purpose of doing a\nmeta-analysis—which carefully does not use significance of each study but rather the\nmagnitude, direction, and sampling error to make inferences. It’s ok to report the\ncomposition of the dataset (e.g., number of studies per taxa group, or among categories\nof responses relative to richness), but not differences in significance.\nResponse: Thank you for raising this critical methodological point. We agree that\nvote-counting based on significance undermines the purpose of a meta-analysis. In\nthis version, we have thoroughly revised the first section of the Results to eliminate\nthis approach and drop the proportion of different RUR types (positive or negative).\nInstead, we mainly reported the overall effect size, confidence interval and\nsignificance in this revision (Lines 123-127).\nFurthermore, regarding the subgroup analyses across different taxa, we only\nreported the effect sizes and their significance for several groups as you suggested\n(Lines 132-139; See also Table S9 for details). However, our primary intent is not\nto directly compare whether the effect size of one group is significantly higher or\nlower than that of another. Rather, the key message we aim to convey is that the\noverall negative RUR pattern was found to be widespread and prevalent\nacross diverse taxonomic groups.\nThe rank order of effects to emphasize representation of negative effects is ok, but\ndiscussing their prevalence without the lens of meta-analysis is dangerous as many\nstudies may be underpowered and visually may not be representative of the underlying\n(aggregate) effect. Having correlations as Fisher-Z transformed effects is better for\nclarity since results are analyzed with the metric (e.g., main inferences are derived from\nthese) and consistence helps with understanding what was done in the paper.\nResponse: Thank you for this comment. We agree that presenting results through\nthe lens of meta-analysis is crucial for accurate interpretation. Accordingly, we\nhave revised our figures to align with this approach. Specifically, we have replaced\nthe distribution plots of raw correlation coefficients in Figure 2b with distributions\nbased on Fisher's Z transformed effect sizes. Furthermore, Figure 2a has also been\nupdated to present the overall effect sizes and confidence intervals for both\npresence absence and abundance data.\nFor the supplementary figures and tables showing the distribution of RUR\nacross different data sources and taxonomic groups (Supplementary Figure S3,\nFigure S4, Table S9), we have also converted the data to Fisher's Z transformed\neffect sizes for consistency and clarity across all visualizations and presentations.\nFigure 2 (a) The overall effect sizes of RUR (Fisher’s Z statistics) with 95% confidence intervals\nacross different studies and taxa groups based on presence absence and abundance data types. (b)\nThe global variation of RUR across different studies and taxa groups based on presence absence and\nabundance data types. Number in the left panel indicate number of observations. Details of model\noutput can be found in Table S5.\nThe results also need to define the acronyms (PA, AB) and more descriptive words\nshould be used for CV_richness, Clim_PC1, and Clim_PC2. Since some details about\nthe methodology are embedded outside the main manuscript, it’s important that they\nare defined so that the reader knows right-way on how to interpret what was reported.\nThese include the above acronyms, but also measure of % explained (line 176 but\nexplained using different definitions in line 568).\nResponse: Thank you for this suggestion. We have revised the Results section\naccordingly by replacing the acronyms with more descriptive terms. Specifically,\nwe have replaced \"PA\" with \"Presence absence\" and \"AB\" with \"Abundance\" in\nall figures of the main text (See Figures 1-5 in the main text). Also, we have\nexplicitly stated their definitions at the beginning of the Results section (Line 125).\nMoreover, \"CV_richness\" has been changed to \"Breadth of richness\" (Lines 150-\n154), and \"Clim_PC1\" and \"Clim_PC2\" have been replaced with \"PC1 of climatic\nvariables\" and \"PC2 of climatic variables,\" respectively (Lines 174-186). We have\nalso added explanatory captions to the figures to clarify the meaning of PC1 and\nPC2 (See Figures 3-5 in the main text). These changes will help readers better\nunderstand the meaning of these factors. Additionally, in response to your point on\nthe percentage of variance explained, we have clarified in the Results that this value\nrepresents the marginal R squared (as defined in Lines 625-626), which denotes the\nvariance explained by the fixed effects in the mixed-effects models (Lines 188-189;\nLine 194).\nMany of the results reported are described as “affects”, for example this “affected”\nthat. If a moderator is being used to explain variation among effects, maybe “predict”\nis a better inferential description, since correlations are just associations and not\ncauses, and similarly, moderators are predictors and not causes of variability.\nResponse: Following your suggestion, we have revised the wording throughout\nthe Results section, replacing causal terms such as \"affects\" with more appropriate\ninferential words like \"predicts\" or \"association\". (Lines 142-186)\nFinally, lines 517-518 state that significance was tested by 95% overlap, yet there are\nonly p-values reported throughout the results; what is the test statistics these p-values,\nand their values with degrees of freedom need to be reported along with p-values.\nResponse: Thank you for pointing this out. In addition to the P-values, we have\nnow reported the Q statistic and its degrees of freedom for model regression\nM\ncoefficients in the Results section (Lines 146-186). Furthermore, we have included\nthis statistical information directly in Figures 3-5 in the main text. Other statistics,\nsuch as the 95% confidence intervals and the total heterogeneity of effect sizes (i.e.,\nI²), can be found in the Supplementary Tables S5-S18.\nMinor comments:\n1) What is the scale of the taxonomic group? Are they Linnean rankings, or homologous\ngroups?\nResponse: Thank you for raising this question. Our classification system is based\non Linnean rankings. To ensure mutual exclusivity between groups, we have\nsystematically categorized each taxon according to specific taxonomic ranks (such\nas Phylum, Class, or Order). The table below (Table R1) provides detailed\ninformation on the specific finer groups included within each taxa group (where\napplicable). If insisted, we can of course give this table also as a supplement in the\nmanuscript.\nTable R1 Details of taxonomic groups in the synthesis\nTaxa group Finer groups\nAlgae Diatoms\nBenthic diatoms\nStudies only reported “Algae”\nLaminariales (Kelp)\nAmphibians Frogs\nAnurans\nTadpoles\nStudies only reported “Amphibians”\nAquatic plants Macrophyte\nMarine Plants\nStudies only reported “Aquatic plants”\nWetland plants\nArthropods Spiders\nCrabs\nScorpions\nDecapods\nPeracarida\nStudies only reported “Arthropods”\nBats No finer groups in the studies\nBirds No finer groups in the studies\nBivalves Mussels\nStudies only reported “Bivalves”\nFishes No finer groups in the studies\nFreshwater macroinvertebrates Aquatic insects\nStudies only reported “Freshwater macroinvertebrates”\nBenthic macroinvertebrates in freshwater ecosystem\nGastropods No finer groups in the studies\nInvertebrates Intertidal invertebrates\nMarine Invertebrates\nInsular invertebrates\nEarthworms\nLichens No finer groups in the studies\nMammals No finer groups in the studies\nNematodes No finer groups in the studies\nPhytoplankton No finer groups in the studies\nReptilia Snake\nLizard\nSeaweed No finer groups in the studies\nTerrestrial insects Bees\nWasps\nFireflies\nBeetles\nBombus\nAnts\nMoths\nStaphylinidae\nStudies only reported “insects” in terrestrial ecosystem\nOrthoptera\nGrasshopper\nHoverflies\nFlies\nButterflies\nTerrestrial plants Terrestrial vascular plants\nStudies only reported “plants” in terrestrial ecosystem\nMosses\nBryophytes\nZooplankton No finer groups in the studies\n2) Line 194, also 403: “literature studies” should be “primary studies”\nResponse: We have revised “literature studies” into “primary studies” throughout\nthe manuscript.\n3) Line 409: Define what a “unit” is.\nResponse: Thank you for pointing this out. In this context, a \"unit\" refers to the\nnumber of sampling points or observation sites in each case study. We have\nclarified this definition in the revised manuscript. (Line 437)\n4) Line 439: What do you mean by “previous researches”?\nResponse: Thank you for pointing this out. The term \"previous researches\" refers\nspecifically to all existing studies that have reported on the relationship\nbetween species richness and ecological uniqueness. We synthesized these\nstudies to better determine the specific form of this relationship (e.g., whether it is\nlinear or quadratic). We have now clarified this point in the revised manuscript.\n(Lines 467-468)\nReviewer #4 (Remarks to the Author):\nThis manuscript aims to describe the relationship between, and drivers of, species\nrichness and ecological uniqueness. Generally, the manuscript reads very well and the\ncomments of reviewer 1 and 3 are well addressed.\nResponse: Thank you for this positive assessment of our work.\nHowever, when the quadratic models are tested (as suggested by reviewer 1) and found\nmore appropriate compared to the linear models in 30% of the cases, these records are\nnot considered as the results are harder to interpret. I think that the better fit of the\nquadratic relationship in such a high percentage of the records show that the\nrelationship can not always be interpreted as linear and that these cases should be\nconsidered and not disregarded as the results are harder to interpret. In that way you\nare only using the results which are convenient, while this is not a good representation\nof what the data actually shows.\nResponse: Thank you for the suggestion. As you critically pointed out, these\nnonlinear relationships represent distinct patterns of variation between ecological\nuniqueness and species richness. They are equally important and have been\nrepeatedly observed in previous studies (See Appendix 1). In this revision, we have\nidentified different types for the quadratic relationships via performing a qualitative\nanalysis of the nonlinear forms of RUR (i.e., those cases in which quadratic\nrelationships exhibited better performance than linear), following the procedure\nbelow:\nFirst, we fitted quadratic models to describe the relationships between\necological uniqueness and species richness in these cases. Then, we calculated the\naxis of symmetry for each quadratic model and quantified the number of data points\ndistributed on each side of the axis. Based on these patterns, the relationships were\ncategorized as follows:\n(1) For models with R² greater than 0.1, if all data points were distributed entirely\non the left or right side of the axis of symmetry, the relationship is classified as\neither a \"negative association\" or \"positive association\" (i.e., ecological uniqueness\ndecreases or increases with increasing species richness; Fig. S2a-d).\n(2) For models with R² greater than 0.1, if more than 80% of the data points were\nlocated on one side of the axis of symmetry, the relationship is classified as a\n\"general negative association\" or \" general positive association\", indicating a\npredominant but not strictly monotonic trend (Fig. S2e-h). We also tested\nalternative thresholds (e.g., 70% and 90%) and found the results to be similar, and\n80% was accordingly selected as the distribution threshold.\n(3) When fewer than 80% of the data points were located on one side of the axis of\nsymmetry, the relationship was classified as \"no general association\", indicating\nthat both species-rich and species-poor sites could exhibit high ecological\nuniqueness.\n(4) If the model R² is less than 0.1, regardless of the distribution of points relative\nto the axis of symmetry, the relationship is also classified as \"no general\nassociation\". (Lines 501-521)\nFigure S2 Schematic illustration of classification based on the quadratic richness-uniqueness\nrelationship (RUR). (a, b) Negative association: all data points were distributed on the same\nside of the axis of symmetry, and ecological uniqueness decreases with increasing species\nrichness. (c, d) Positive association: all data points were distributed on the same side of the axis\nof symmetry, and ecological uniqueness increases with increasing species richness. (e, f)\nGeneral negative association: over 80% of data points were distributed on the same side of the\naxis of symmetry, and ecological uniqueness decreases with increasing species richness in the\nsame side of the axis of symmetry. (g, h) General positive association: over 80% of data points\nwere distributed on the same side of the axis of symmetry, and ecological uniqueness increases\nwith increasing species richness. See Appendix 3 for details.\nWe found that 22.8% (140 out of 615 records) of the quadratic relationships\nshowed \"negative association\", which indicates that ecological uniqueness\ndecreased as species richness increased (Appendix 3). 20% (122 out of 615 records)\nof the observations showed \"general negative association\", which means that we\ntend to observe higher ecological uniqueness in species-poor sites although the\nrelationship is not strictly monotonous. In contrast, fewer than 1% (3 out of 615)\nshowed “positive association” and only 3% showed “general positive association”,\nsuggesting that ecological uniqueness was generally higher in species-rich sites.\nThese results highlighted that: although RUR can be linear or quadratic depend on\ncontext, those quadratic patterns also reflect a predominant tendency for high\necological uniqueness to coincide with low species richness, reinforcing the widely\nobserved negative association between ecological uniqueness and species richness\n(Lines 220-228).\nThe results of this qualitative analysis are presented in Appendix 3 and\nSupporting data, and their implications are discussed in the revised discussion\n(Lines 228-232). However, we have to acknowledge that our analysis integrates\ndatasets spanning diverse biodiversity and environmental gradients, making\nsystematic classification challenging. Nevertheless, these cases collectively\nunderscore the point that species richness and ecological uniqueness represent\ncomplementary aspects of biodiversity.\nAdditionally, different drivers of the RURs are tested, following four different\nhypotheses. In the methods, I recommend a more thorough description of why and how\nwell the variables chosen reflect these particular mechanisms underlying the proposed\nhypotheses.\nResponse: Thank you for the suggestion. We have added an explanation in the\nSupplementary Text regarding why (and how well) these specific factors were\nselected as representative variables for each hypothesis (See Supplementary Text\n6: Quantifications of different hypotheses) as follows:\nSpecies pool attributes: These three metrics (Gamma diversity, breadth of\nrichness and proportion of rare species) provide a multi-dimensional quantification\nof the species pool hypothesis. Gamma diversity directly measures the size of the\nregional species pool, representing the total source of potential colonist species\navailable for community assembly (Lessard et al., 2012). Breadth of richness\nspecifically reflects the variability of richness within given region - a higher value\nindicates greater richness variation across sampling units, thereby illustrating how\nthe species pool maintains different community configurations across varying local\nhabitats (Dansereau et al., 2022). The proportion of rare species reveals the pool's\ninternal structure, with a higher proportion suggesting lower number of shared\nspecies or higher species turnover among communities (Yao et al., 2021).\nDispersal limitation: Distance-based Moran’s Eigenvector Maps (dbMEM)\nwere employed to infer dispersal limitation because they effectively detect spatial\nautocorrelation patterns characteristic of constrained species movement. Dispersal\nlimitation produces a distance-decay relationship in community similarity, where\nnearby sites tend to share more species (positive spatial autocorrelation). dbMEM\nvariables capture these underlying spatial structures at different scales, from broad-\nscale gradients to fine-scale patches. When these variables explain a significant\nportion of the variation in community composition, the most parsimonious\ninterpretation is that dispersal limitation plays a dominant role in shaping\ncommunities.\nEnvironmental filtering: We used 19 bioclimatic variables to represent\nenvironmental filtering, mainly because of their direct quantification of the critical\nclimatic conditions that determine species' survival and distribution. The species-\nenergy hypothesis posits that climate strongly influences primary productivity, or\nthe total energy available for partitioning within most ecological communities,\nthereby affecting species richness or species composition via its impact on the\nnumber and size of organisms that can be supported (Hutchinson 1959; Currie et\nal., 2004; Chu et al., 2016; Storch et al., 2018). These variables such as the\nminimum temperature of the coldest month (Bio6) and precipitation of the driest\nquarter (Bio17) capture the critical extremes of temperature and water availability\nthat often limit species persistence.\nSpatial grain size: Spatial grain size not only directly determines the number\nof species observed and the spatial turnover in species composition, but also\npotentially alters the influence of other ecological processes such as biotic\ninteractions and environmental heterogeneity. Larger grains tend to integrate\ngreater habitat heterogeneity and buffer fine-scale stochasticity, while smaller\ngrains better capture local biotic and environmental filtering effects.\nReferences:\nChu, C. et al. Direct effects dominate responses to climate perturbations in\ngrassland plant communities. Nat. Commun. 7, 11766 (2016).\nCurrie, D. J. et al. Predictions and tests of climate-based hypotheses of broad-scale\nvariation in taxonomic richness. Ecol. Lett. 7, 1121–1134 (2004).\nDansereau, G., Legendre, P. & Poisot, T. Evaluating ecological uniqueness over\nbroad spatial extents using species distribution modelling. Oikos, e09063 (2022).\nHutchinson, G. E. Homage to Santa Rosalia or why are there so many kinds of\nanimals? Am. Nat. 93, 145–159 (1959).\nLessard, J-P., Belmaker, J., Myers, J. A., Chase, J. M. & Rahbek, C. Inferring local\necological processes amid species pool influences. Trends Ecol. Evol. 27, 600-607\n(2012).\nStorch, D., Bohdalková, E. & Okie, J. The more-individuals hypothesis revisited:\nthe role of community abundance in species richness regulation and the\nproductivity–diversity relationship. Ecol. Lett. 21, 920–937 (2018).\nYao, J., Huang, J., Ding, Y., Xu, Y., Xu, H. & Zang R. Ecological uniqueness of\nspecies assemblages and their determinants in forest communities. Divers. Distrib.\n27, 454-462 (2021).\nGenerally, I strongly recommend to first run a fully randomised null model, to check if\nthe obtained relationships, and therefore your key results, are obtained by chance or\nthat ecological mechanisms drive this relationship.\nResponse: Thank you for this insightful suggestion. We conducted a null model\nanalysis to determine whether the observed RUR was statistically distinct from\npatterns expected by chance, thereby indicating the influence of underlying\necological processes. For presence absence data, the null model randomized\nspecies co-occurrence while preserving species incidence frequencies, site richness,\nand total occurrences (Gotelli, 2000). For abundance data, it maintained both site-\nlevel and species-level total abundances, randomizing only the assignment of\nindividuals to sites (Kraft et al., 2011; Xing & He, 2019). The details of null model\nanalysis can be found in Supplementary Text 6: Null model analysis.\nWe generated 999 randomized communities to construct the null distribution\nof RUR. Statistical significance was evaluated using a two-tailed test, where the\nobserved RUR was compared against the 2.5th and 97.5th percentiles of the null\ndistribution. A P-value < 0.05 signified that the observed pattern was\nsignificantly different from random expectation, suggesting the influence of\nnon-random ecological processes such as dispersal limitation or\nenvironmental filtering. We have provided the results (P-value, 2.5th and 97.5th\npercentiles of null distributions) of this analysis in the Supporting data.\nIn total, the null model analysis revealed that a total of 645 records (290\nof presence absence data and 355 of abundance data) showed significant\ndeviations from the null expectation (See Supporting data for details). Therefore,\nour subsequent analyses examining the predictions of different hypotheses on\nRUR and their relative importance were conducted exclusively on this subset\nof records (k = 645).\nReferences:\nGotelli, N. J. Null model analysis of species co-occurrence patterns. Ecology 81,\n2606–2621 (2000).\nKraft, N. J. B. et al. Disentangling the drivers of β diversity along latitudinal and\nelevational gradients. Science 333, 1755–1758 (2011).\nXing, D. L. & He, F. L. Environmental filtering explains a U-shape latitudinal\npattern in regional β-deviation for eastern North American trees. Ecol. Lett. 22,\n284–291 (2019).\nLastly, from a traits perspective, beta diversity can be different between sites, but some\nsites are more ecologically unique than others as the niche and traits of the species\ndiffer more compared to the other sites. This is a point which should be considered in\nthe introduction and discussion.\nResponse: Thank you for pointing out it. Following your suggestion, we have\nhighlighted the role of traits and niche differences in the introduction and\ndiscussions:\nIn the Introduction: “Moreover, sites exhibiting high ecological uniqueness\nmay also contribute disproportionately to regional biodiversity by enhancing trait\ndivergence and niche differentiation of specific species assemblages.” (Lines 31-\n34)\nIn the discussions: “Finally, while our findings highlight the conservation\nrelevance of RUR, taxonomic measures alone may not fully capture the\nmultidimensional nature of biodiversity. Variation in ecological uniqueness across\nsites can also be further interpreted through the lens of functional traits and niche\ndifferentiation (Violle et al., 2007; Violle et al., 2017; Heino et al., 2022;\nAscanio et al., 2024). For instance, sites with higher ecological uniqueness likely harbor\nspecies assemblages that are distinct not only in taxonomic identity but also in their\nfunctional characteristics and fundamental ecological niches (Heino et al., 2022;\nCasanelles-Abella et al., 2023). These communities may display greater divergence\nin key niche-related traits and occupy underrepresented or extreme regions of the\nmultivariate trait space. Critically, the observed uniqueness may reflect substantial\ndifferences in the ecological niches of the constituent species, which could vary\ngreatly in resource use strategies, environmental tolerances, or other dimensions\n(Chesson, 2000; Ascanio et al., 2024). This perspective underscores that the\necological uniqueness is intrinsically linked to both functional rarity and niche\ncomplementarity, highlighting the need to integrate trait-based and niche-based\napproaches with biodiversity partitioning to better understand the mechanisms\nunderlying spatial heterogeneity in ecological communities.” (Lines 393-407)\nReferences:\nAscanio, A., Bracken, J. T., Stevens, M. H. H. & Jezkova, T. New theoretical and\nanalytical framework for quantifying and classifying ecological niche\ndifferentiation. Ecol. Monogr. 94, e1622 (2024).\nCasanelles-Abella, J., Fontana, S., Meier, E., Moretti, M. & Fournier, B. Spatial\nmismatch between wild bee diversity hotspots and protected areas. Conserv. Biol.\n37, e14082 (2023).\nChesson, P. Mechanisms of Maintenance of Species Diversity. Annu. Rev. Ecol.\nEvol. Syst. 31, 343-366 (2000).\nHeino, J. et al. Assessing the conservation priority of freshwater lake sites based\non taxonomic, functional and environmental uniqueness. Divers. Distrib. 28, 1966–\n1978 (2022).\nViolle, C. et al. Functional rarity: The ecology of outliers. Trends Ecol. Evol. 32,\n356–367 (2017).\nViolle, C. et al. Let the concept of trait be functional! Oikos 116, 882–892 (2007).\nAdditional comment sent to the editor:\nElaborate on the relationship in Figure 3C and if they can test how this relationship is,\nor is not, driven by the outliers at the right side of the graph.\nResponse: Thank you for this careful comment, we have rigorously tested the\nrelationship in Figure 3c. We identified and removed outliers falling outside the\nrange of two standard deviations from the mean of breadth of richness and\nsubsequently re-ran the meta-regression. The results indicate that a positive trend\npersists although the relationship is no longer statistically significant after\nexcluding these outliers (See Figure R1 below). This indicates that the originally\nsignificant relationship was largely driven by a few outlier records. However, this\nfinding also suggests that if future studies encompass a broader range of\ncommunities with higher breadth of richness, the relationship could become more\nrobust and potentially reveal a clearer transition from negative to positive\nassociations. We have incorporated this interpretation and its implications into the\nrevised Discussions (Lines 278-286).\nFigure R1 The association between effect sizes of RUR and breadth of richness (a) before and\n(b) after outlier removal. The regression lines represent the response to the above factors with\nthe 95% confidence intervals shaded from hierarchical mixed-effects meta-regressions. The\nvalues of the estimated coefficients (β), 95% confidence intervals, Q statistics and degree of\nM\nfreedom (df) are shown in each panel. The solid and dashed lines represent significant (P <\n0.05) and non-significant (P > 0.05) relationships. The points have radii relative to the inverse\nvariance of each observation (i.e., observations with larger points have lower variance in effect\nsize). k: Number of effect sizes."}, {"round": 3, "reviewer_comments": {"reviewer_2": "In terms of what was reported and presented with this meta-analysis, I believe this new version has the clarity and detail that is more in line with current expectations in research synthesis. My comments below are very very minor.\nMAIN TEXT\nlines 832, 838, 851: \"Fisher’s Z statistics\" should be \"Fisher’s Z-transformed correlation coefficients\"\nlines 864-876: definition of effect sizes (see above) is missing\nSUPPLEMENTARY TEXT 1:\nline 11: where was this search done? E.g., using who's institutional access to IS Web of science?\nline 32: Which PRISMA statement did you use? cite source\nline 3863, 3872, 3879, 3890, 3898, 3906, 3978, 3990, 4005, 4018, 4030, 4033, 4037, 4054, 4070, 4086, 4102, 4118, 4134,4150,4166,4182: \"Fisher’s Z statistics\" should be \"Fisher’s Z-transformed correlation coefficients\"", "reviewer_4": "I would like to thank the authors for the thourough answers and additional analyses performed. The comments are well addressed. The only question which remains for me is about the outliers in Figure 3C. Are the datapints at the right side of the figure indeed checked and found to be valid or do these datapoints represent less robust observations? In the latter case I would disregard these datapoints and change the text accordingly, otherwise all good."}, "author_response": "Response to reviewers\nReviewer #2 (Remarks to the Author)\nIn terms of what was reported and presented with this meta-analysis, I believe this new version has the clarity and detail that is more in line with current expectations in research synthesis. My comments below are very very minor.\nResponse: Thank you for this positive assessment of our work.\nMAIN TEXT\nlines 832, 838, 851: \"Fisher’s Z statistics\" should be \"Fisher’s Z-transformed correlation coefficients\"\nResponse: As you suggested, we have revised the corresponding figure captions by replacing \"Fisher’s Z statistics\" with \"Fisher’s Z-transformed correlation coefficients\". (See captions of Figure 2-4 in the main text)\nlines 864-876: definition of effect sizes (see above) is missing\nResponse: We have included the definition of effect sizes in the corresponding figure captions. (See captions of Figure 5 in the main text)\nSUPPLEMENTARY TEXT 1:\nline 11: where was this search done? E.g., using who's institutional access to ISI Web of science?\nResponse: We have clarified the search methodology by adding the following detail in the manuscript: “The search was conducted using the Web of Science Core Collection database, accessed via the institutional subscription of the University of Chinese Academy of Sciences.” (See Supplementary Note 1 for details)\nline 32: Which PRISMA statement did you use? cite source\nResponse: We have specified the version of the PRISMA statement used in our study and added the corresponding citation in the supplementary information (See Supplementary Note 1 for details). We conducted this systematic review and meta-analysis in accordance with the Preferred Reporting Items for Systematic Reviews and Meta-Analyses (PRISMA) 2020 statement (Page et al., 2021).\nReferences:\nPage, M. J. et al. The PRISMA 2020 statement: an updated guideline for reporting systematic reviews. BMJ. 372, n71 (2021).\nline 3863, 3872, 3879, 3890, 3898, 3906, 3978, 3990, 4005, 4018, 4030, 4033, 4037, 4054, 4070, 4086, 4102, 4118, 4134,4150,4166,4182: \"Fisher’s Z statistics\" should be \"Fisher’s Z-transformed correlation coefficients\"\nResponse: We have replaced \"Fisher’s Z statistics\" with \"Fisher’s Z-transformed correlation coefficients\" as you recommended. (See captions of Supplementary Figure 3-9, Supplementary Table 3-18)\nReviewer #4 (Remarks to the Author)\nI would like to thank the authors for the thorough answers and additional analyses performed. The comments are well addressed.\nResponse: Thank you for this positive assessment of our work.\nThe only question which remains for me is about the outliers in Figure 3C. Are the data points at the right side of the figure indeed checked and found to be valid or do these data points represent less robust observations? In the latter case I would disregard these data points and change the text accordingly, otherwise all good.\nResponse: Thank you for the thoughtful suggestion. We have carefully examined the two data points on the right side of figure 3C, particularly regarding the calculation of the coefficient of variation (CV) of species richness. After verification, we confirm that these CV values were correctly calculated. Therefore, we have retained these data points in the analysis and figure. To further assess the robustness of the relationship, we have also presented the results after removing these outliers in the supplementary information (See Supplementary Figure 9). The overall pattern remains consistent, supporting the stability of our conclusion.\nFurthermore, we have expanded the Discussion to emphasize the potential influence of these extreme values on the results. We note that as future studies encompass a broader range of communities with greater breadth of richness, the estimated relationship may become more robust and potentially reveal a clearer transition from negative to positive associations. (Lines 301-307).\nSupplementary Figure 9 The association between effect sizes of RUR (Fisher's Z-transformed correlation coefficients) and breadth of richness (A) before and (B) after outlier removal. The regression lines represent the response to the above factors, with their 95% confidence intervals shown as shaded areas from hierarchical mixed-effects meta-regressions. The values of the estimated coefficients (β), 95% confidence intervals, Q statistics (one-sided chi-square test) M and degree of freedom (df) are shown in each panel. The solid and dashed lines represent significant (p < 0.05) and non-significant (p > 0.05) relationships. The points have radii relative to the inverse variance of each observation (i.e., observations with larger points have lower variance in effect size). k: Number of effect sizes."}]

Decision letter