induced biodiversity loss - Wiley Online Library

6 downloads 91 Views 1MB Size Report
Mar 9, 2015 - 5. Synthesis. As nitrogen (N) enrichment typically occurs at scales much larger than individual plots, understanding .... 2004b), which will reduce immigration between different communities and ...... E.B. et al. (2011) Effects of ...
Journal of Ecology 2015, 103, 750–760

doi: 10.1111/1365-2745.12395

Testing the scaling effects and mechanisms of N-induced biodiversity loss: evidence from a decade-long grassland experiment Zhichun Lan1,2, G. Darrel Jenerette3, Shuxia Zhan1,2, Wenhuai Li1,2, Shuxia Zheng1 and Yongfei Bai1* 1

State Key Laboratory of Vegetation and Environmental Change, Institute of Botany, The Chinese Academy of Sciences, Beijing 100093, China; 2University of Chinese Academy of Sciences, Beijing 100049, China; and 3 Department of Botany and Plant Sciences, University of California, Riverside, CA 92521-0124, USA

Summary 1. Although extensive studies demonstrate that nitrogen (N) enrichment frequently reduces plant diversity within small quadrats (0.5–4 m2), only a few studies have evaluated N effects on biodiversity across different spatial scales. 2. We conducted the first experimental test of the scale dependence of N effects on species richness from a 10-year N treatment (1.75–28 g N m2 year1) in a typical steppe. We used species–area relationship (SAR) to analyse the scale dependence of species loss with power model S = cAz (S is species number, A is area, c is intercept, and z is slope). 3. Absolute species loss decreased at sampling area > 8 m2. Proportional species loss (compared to control) decreased and critical threshold (Ncrit) for biodiversity losses increased with sampling areas. 4. These scale dependences were quantified as increasing slope (z-value) of SAR with N addition. Through SAR decomposition, we found that this overall positive effect was in response to positive effects of changes to the species abundance distribution over negative effects of overall species richness losses. 5. Synthesis. As nitrogen (N) enrichment typically occurs at scales much larger than individual plots, understanding how N enrichment affects the scaling patterns of biodiversity is necessary for biodiversity conservation and ecosystem management in response to anthropogenic N deposition. Key-words: critical threshold, determinants of plant community diversity and structure, N enrichment, number of individuals, SAR decomposition, scaling effect, species abundance distribution, species aggregation, species–area relationship (SAR)

Introduction While terrestrial primary production is widely limited by nitrogen (N) availability (Vitousek & Howarth 1991), the increasing rate of nitrogen pollution is changing biodiversity and ecosystem functioning world-wide (Vitousek et al. 1997; Phoenix et al. 2006; Galloway et al. 2008; Bobbink et al. 2010; Pardo et al. 2011). Decreases in biodiversity are commonly observed following N input (Suding et al. 2005; Clark & Tilman 2008; Hautier, Niklaus & Hector 2009; Clark et al. 2013), and N deposition is recognized as one of the major causes of biodiversity decline in terrestrial ecosystems (Sala et al. 2000; Stevens et al. 2004a). However, these findings are limited to quadrat scale responses and are not necessarily representative of the regional *Correspondence author. E-mail: [email protected]

N deposition. Only a few studies have evaluated scaling relations of N-induced species loss (Crawley et al. 2005; Hillebrand et al. 2007; Chalcraft et al. 2008; Leps 2014). These findings suggested that the patterns of N-induced biodiversity loss may be changed at larger scales, and few studies have tested mechanisms underlying these patterns (Chalcraft et al. 2008; Leps 2014). The scale dependence of biodiversity change also suggests that a scale-dependent critical threshold (Ncrit) for N enrichment may exist, below which N might have minimal effects on biodiversity (Clark & Tilman 2008; Bai et al. 2010; Lan & Bai 2012). It is still unclear how the sensitivity and magnitude of biodiversity loss after N addition vary with spatial scales. A mechanistic understanding of scale-dependent N-induced species loss is needed to extend findings from plot level and develop more efficient management strategies for biodiversity conservation (Chalcraft et al. 2008; Cleland & Harpole 2010).

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society

Scale dependence of N-induced species loss 751 Species–area relationship (SAR) is one of the most frequently used approaches to analyse scale dependence of biodiversity changes (Keil et al. 2011; Sandel & Corbin 2012; Powell, Chase & Knight 2013). The scale dependence of biodiversity is commonly modelled as a power function (Preston 1960; Tjørve 2003), S = cAz, where S is species number, A is area, and c and z are constants (Arrhenius 1921). The slope (z-value) and intercept (c-value) of the SAR in log–log space are frequently used to quantify how species richness varies across multiple spatial scales (Tjørve 2003; Powell, Chase & Knight 2013). For example, the commonly observed reduction of species richness at the plot level with increased N may be associated with decreases in the intercept (c-value) and/or increases in the slope (z-value) of the SAR (Crawley et al. 2005; Chalcraft et al. 2008). If N increases the slope (z-value) of SAR, the rate of species loss would be smaller at broader scales, whereas if N decreases the slope of SAR, then species reduction will be larger at broader scales (Fig. 1). However, the majority of research into understanding variation in SARs has been observational and experimental approaches, though more challenging, are needed to directly test potential mechanisms (Jenerette & Shen 2012). Recent SAR decomposition approaches provide an analytical framework to quantify changes in the c- and z-values as caused by a combination of changing species pool, number of individuals, species abundance distribution and/or species aggregation (spatial distributions) (He & Legendre 2002; Tjørve et al. 2008; Scheiner et al. 2011; Powell, Chase & Knight 2013). N addition may affect overall SAR and consequently regulate the scale dependence of biodiversity changes through its distinct influences on each of the four individual factors above. First, N addition may lead to decreased species pool through loss of rare species (H1: reduced richness)

(a)

(b)

(Crawley et al. 2005; Chalcraft et al. 2008). Through the reduced richness hypothesis, we expect N additions to decrease both the c- and z- value of the SAR, and N-induced species loss would increase with sampling area. Secondly, N addition increases biomass production, causes assemblagelevel thinning and thereby decreases the number of individuals because of increased plant size and intensified light competition (H2: reduced individuals) (Goldberg & Miller 1990; Oksanen 1996). Through reductions in the number of individuals, we expect N addition would decrease the c-value but increase the z-value of SAR (Oksanen 1997), and thus, the effects of N enrichment on species richness would diminish with increasing sampling area. Thirdly, N addition increases the abundance of a few N-favoured species at the expense of less favoured species (H3: altered abundance distributions) (Suding et al. 2005; Bai et al. 2010), which will decrease the evenness of community (Gilliam 2006; Hillebrand et al. 2007). Through decreased evenness of species abundance distribution, we expect N additions would decrease the c-value but increase the z-value, and thereby the Ninduced reduction in species richness would become increasingly moderated with sampling area. Finally, N addition intensifies recruitment limitation (Stevens et al. 2004b), which will reduce immigration between different communities and increase species aggregations (H4: species aggregation). Through increases in species aggregation, we expect N addition would decrease the c-value but increase the z-value of SAR; hence, the effects of N enrichment on species richness would decrease with increasing sampling area (Table S1 in Supporting Information). Our objective in this study was to experimentally quantify the scale dependence of N-induced species loss and test the underpinning mechanisms. We established a 10-year field

(c)

Fig. 1. Three scenarios of how N addition affects the slope of species–area relationship (SAR) and the corresponding scale effects on proportional species loss (PSL). In each scenario, the discrepancy between the two SARs represent N-induced species loss, that is DS = log (Scontrol) – log (SN) = log (Scontrol/SN), or exp (DS) = Scontrol/SN. Proportional species loss (PSL) was defined as proportion of species loss in N treatments compared to control, that is PSL = 1  SN/Scontrol = 1 – 1/[exp (DS)]. Therefore, PSL was positively correlated with DS. (a) N addition has no effects on the slope of SAR, and DS and PSL are independent of spatial scale; (b) N addition increases the slope of SAR, and DS and PSL decrease with increasing sampling area; (c) N addition decreases the slope of SAR, and DS and PSL increases with sampling area. © 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

752 Z. Lan et al. experiment with six levels of N addition rate in a semi-arid Inner Mongolia grassland, China. The Inner Mongolia grassland both historically and currently received little N deposition (Galloway et al. 2008) and has not yet experienced a rapid loss of plant species induced by chronic N input. This low N deposition background allows testing N effects on biodiversity. Specifically, we address two questions: first, how do the N-induced biodiversity loss (i.e. plant species richness, absolute and proportional species loss and critical threshold) change with spatial scales? Secondly, what are the relative contributions of alternative mechanisms in terms of species pool, number of individuals, species abundance distributions and species aggregation, mediating the scale dependence of N-induced species loss?

Materials and methods STUDY SITE

The study was conducted at the Inner Mongolia Grassland Ecosystem Research Station (IMGERS; 116°420 E, 43°380 N and 1250 m a.s.l.), which is located in the Xilin River Basin, Inner Mongolia Autonomous Region of China (Bai et al. 2010). The mean annual temperature in the study area is 0.3 °C, with mean monthly temperature ranging from 21.6 °C in January to 19.0 °C in July. The mean annual precipitation is 346.1 mm, with about 80 per cent occurring in the growing season from May to September. The field experiment was carried out on a natural steppe community, which was fenced to prevent grazing since 1979. At the beginning of the experiment, there were 48 vascular plant species and 80% of them were rare species (relative abundance < 1%) (Lan & Bai 2012). The dominant and common species in plant community include Leymus chinensis (perennial rhizomatous grass, PR), Stipa grandis (perennial bunchgrass, PB), Agropyron cristatum (PB), Achnatherum sibiricum (PB), Carex korshinskyi (sedge), Kochia prostrate (perennial forbs, PF), Allium tenuissimum (PF), Chenopodium glaucum (annual forbs, AF), Axyris amaranthoides (AF) and Salsola collina (AF) (Bai et al. 2004).

N ADDITION EXPERIMENT

To explore the effects of N addition on plant species diversity and ecosystem functioning, an N addition experiment was initiated in 2000 and is still ongoing (Bai et al. 2010; Lan & Bai 2012). The experiment has six treatments and nine replicates for each treatment, with a total number of 54 5 m 9 5 m plots. Experimental plots were arranged in a randomized block design, with 1-m-width walkways between plots. N treatments included six levels of N addition rate, 0, 1.75, 5.25, 10.5, 17.5 and 28.0 g N m2 year1 as commercial granular NH4NO3 (Bai et al. 2010). The low N input rates (e.g. 0, 1.75 and 5.25 g N m2 year1) allow us to approximate N deposition and estimate critical thresholds for species loss (Bobbink et al. 2010). The high N input rates (i.e. 10.50, 17.50, 28.00 g N m2 year1) allow us to approximate management techniques (e.g. N fertilization) and assess the effects of N enrichment on plant diversity (Schellberg et al. 1999; Crawley et al. 2005; Socher et al. 2012). For all treatments, to ensure that N was the only limiting nutrient, based on the soil census data, we also added the same amount of P (10 g P2O5 m2 year1), S (0.2 mg m2 year1) and trace elements (Zn, 190 mg m2 year1, Mn, 160 mg m2 year1 and B, 31 mg m2 year1). Each year,

nutrients were uniformly applied to each treatment plot with manual broadcasting in the mid-growing season (July 1–5) coinciding with high temperature and precipitation, as in other N addition studies (Tilman 1987).

SPECIES–AREA RELATIONSHIP

To test the scale dependence of N effects on species richness, we collected species abundance data from fifty-four 5 m 9 5 m plots for constructing species–area relationship (SAR). Each 25 m2 plot was divided into twenty-five spatially explicit 1 m 9 1 m quadrats. In each quadrat, we determined the abundance (number of individuals) for each species in August 2010. For non-clonal species with no physical connections among different stems (e.g. C. glaucum), an individual was defined as a separately rooted above-ground shoot. For caespitose species with densely packed ramets, such as bunchgrasses (i.e. S. grandis), an individual was defined as a basal. For rhizomatous species with widely spaced ramets (e.g. L. chinensis), an individual was defined as a ramet (Perez-Harguindeguy et al. 2013). Based on our data sets, the empirical SAR in each plot was constructed with type IIA design (Scheiner 2003), that is using a set of spatially explicit quadrats arrayed in a continuous grid (Fig. S1 in Supporting Information). The set of quadrat sizes in each 5 m 9 5 m plot includes 1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 20 and 25 m2, respectively. For a given quadrat (sampling area), we aggregated contiguous 1 m2 quadrats and counted the number of species. For a given quadrat size within each plot, the species richness was the average of number of species across all quadrats with the same size. For example, the species richness of 1 m2 quadrats within each plot was the average of number of species across twenty-five 1 m 9 1 m quadrats. In the same way, we determined the species richness of 2, 3, 4 m2 and other quadrat sizes within each plot. SARs were fitted and model parameters were estimated using MMSAR R-package (Guilhaumon, Mouillot & Gimenez 2010). Specifically, we fitted all the eight candidate models of SARs in the MMSAR R-package (Fig. S2). These models include power, exponential, negative exponential, Monod, rational function, logistic, Lomolino and cumulative Weibull functions of SAR (Guilhaumon, Mouillot & Gimenez 2010). The preliminary analysis indicated that, among the eight models, power, logistic and Monod functions were three best fit models with the highest R2 value but lowest AIC, AICc and BIC values. However, based on the experimental design and sampling scheme of this study, the shape of the species–area curve is expected to be convex upward, which can be best fitted by the power model (Tjørve 2003). It can also be used to determine the predicted Ncrit and extrapolate the results to larger scales. Therefore, we used the power model S = cAz (S, species; A, sampling area) to assess the scale dependence of N effects on species richness.

ESTIMATING EFFECTS OF N ADDITION ON SPECIES LOSS

Absolute species loss was calculated as the difference of species richness between control and each N treatment (i.e. Richnesscontrol  Richnesstreat). To allow direct comparison across sampling areas, proportional species loss was defined as loss proportion compared to control and calculated as one minus species richness response ratio (i.e. 1  Richnesstreat/Richnesscontrol). For each N treatment, the percentage of overall species richness was calculated as the ratio of species richness at a given sampling area to the overall species richness of the plot (25 m2), with the values ranging from 0% to 100%.

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

Scale dependence of N-induced species loss 753 Table 1. Parameters of 25 m2 plots used for each simulation to partition mechanisms contributing to changes in the slope and intercept of species–area relationship (SAR). Parameters were calculated from control and N addition plot. Simulations 1 and 2 represent observed SAR of control and N addition plot, respectively. Simulations 3, 4 and 5 change parameters to test the contribution of each mechanism. The contribution of the four potential mechanisms to changes in the z- and c-value of SAR was calculated as follows: overall species richness = S2–S3; total number of individuals = S3–S4; species abundance distribution = S4–S5; species aggregation = S5–S1; and total effect = S2–S1 (S1–S5 represent the z-value of simulated SAR from simulations 1–5, respectively) Control community parameters

Simulation 1 2 3 4 5

N addition community parameters

Overall species richness

Number of individuals

Species abundance distribution

Species aggregation









✓ ✓ ✓

✓ ✓



Overall species richness



Number of individuals

✓ ✓

Species abundance distribution

✓ ✓ ✓

Species aggregation

✓ ✓ ✓ ✓

For a given simulation, parameters with check marks (‘✓’) were included in the simulation, whereas parameters without check marks were not included in the simulation.

Based on the relationships between proportional species loss and N addition rate, we extrapolated our fitted curve and calculated the predicted critical threshold of N rate (Ncrit), below which there are no reductions in species diversity (Clark & Tilman 2008). Therefore, for a given sampling area, when N addition rates < Ncrit, species loss would not occur as the extant community is capable to sequester the added N; when N rates > Ncrit, species richness would be decreased in response to N addition (Bai et al. 2010).

SAR DECOMPOSITION

To partition the relative contribution of the four potential mechanisms responsible for N-induced changes in the slope (z-value) and intercept (c-value) of SAR, we created a set of null models based on sampled data using methods of Collins & Simberloff (2009) and Powell, Chase & Knight (2013). For each 25 m2 plot of each treatment, SARs were determined and simulated by four factors: species aggregation overall, species abundance distribution, total number of individuals and species richness (Table 1). We used rarefaction bias to quantify the effect of species aggregation (i.e. inter- and intraspecies aggregation) on changes in the z- and c-values of SAR, and a greater bias indicated higher spatial aggregation (Collins & Simberloff 2009; Powell, Chase & Knight 2013). For each N treatment, the null model was constructed using spatially explicit species abundance data collected in nine 25 m2 plot. For each 25 m2 plot, we randomized individuals within and among the twentyfive 1 m2 quadrats, determined species richness at each sampling area (quadrat size) and constructed the rarefaction curve based on the set of rarefaction species richness and sampling area (Powell, Chase & Knight 2013). The mean rarefaction curve for each plot was the average of 1000 simulated rarefaction curves, while the mean rarefaction curve for each treatment was the average of 9000 simulated rarefaction curves. These processes were performed using the software EstimateS. Rarefaction bias (RB) was calculated as the difference in area between the mean rarefaction curve (non-spatially explicit) and empirical SAR (spatially explicit) divided by the area under the empirical SAR, that is RB = (Areararefaction – Areaempirical)/Areaempirical (Collins & Simberloff 2009; Powell, Chase & Knight 2013). For a given sampling area, per cent rarefaction bias (PRB) was calculated as the difference of species richness between mean rarefaction and empirical

SAR divided by the species richness of empirical SAR, that is PRB = (SRrarefaction – SRempirical)/SRempirical (Collins & Simberloff 2009). To examine the effect of species abundance distribution on shifts in the z- and c-values of SAR, we fitted the ranked species abundance distributions for each 25 m2 plot by using a log-series model: Y = (1/log(1c)) * c X/X, where Y is species’ relative abundance, X is rank of species’ relative abundance, and c is coefficient (Powell, Chase & Knight 2013). In the null models, species abundance was calculated by multiplying total number of individuals collected in each 25 m2 plot and the relative abundance of each species obtained from the log-series distribution of the same plot. We randomly filled twenty-five 1 m2 quadrats based on the number of individuals in each plot and reassigned all individuals randomly in space for each 25 m2 plot. We determined the average species richness at each sampling area (i.e. 1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 20 and 25 m2). Then, the average species richness was divided by (1+ per cent rarefaction bias) to get the simulated species richness at each sampling area (Collins & Simberloff 2009). By the set of simulated species richness and sampling areas, we constructed the simulated SAR using power model (S = cAz). We replicated each simulation 1000 times and then calculated the mean simulated SAR. All simulations were conducted with R 2.14 (R Development Core Team). As our experiment was a random block design with nine blocks, we compared simulated SARs between N treatment plots and control plot (zero N) from the same block for examining the N effects on the slopes and intercepts of SARs in log–slog space. For a given pair of N treatment and control plots, we generated five simulations using the parameters as shown in Table 1. Simulations 1 and 2 represent observed SARs. Simulations 3, 4 and 5 changed the parameters to test the contribution of each mechanism. The contribution of each mechanism was calculated by subtracting the z- or c-values of the simulations (Powell, Chase & Knight 2013). Specifically, the contribution of the four potential mechanisms to changes in the z-value of SAR was calculated as follows: overall species richness = S2–S3; total number of individuals = S3–S4; species abundance distribution = S4–S5; species aggregation = S5–S1; and total effect = S2–S1 (S1–S5 represent the z-value of simulated SAR from simulations 1–5, respectively). We repeated these processes to partition the mechanisms contributing to the changes in the c-value of SAR.

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

754 Z. Lan et al. STATISTICAL ANALYSIS

Statistical analyses were performed using SAS v. 8.0 (SAS Institute, Cary, NC, USA). Two-way ANOVAs were performed on the effects of N treatment, sampling area and their interactions on species richness and proportional species loss. Regression models were used to determine the relationships between N addition rate and biodiversity responses (adjustments for nonlinearity were made as needed).

Table 2. Results of two-way anovas on the effects of N treatment (N), sampling area (A) and their interactions on species richness, absolute species loss and proportional species loss Response

Term

d.f.

F

P value

Species richness

N A N9A N A N9A N A N9A

5 13 65 4 13 52 4 13 52

201.92 89.24 0.54 201.95 3.73 0.54 214.86 1.25 0.07

< 0.0001 < 0.0001 0.9986 < 0.0001 < 0.0001 0.9986 < 0.0001 0.2375 1.0000

Absolute species loss

Results N EFFECTS ON SPECIES RICHNESS ACROSS

Proportional species loss

DIFFERENT SAMPLING AREAS

Two-way ANOVAs on the effects of N treatment (N), sampling area (A) and their interactions on species richness (SR) and absolute species loss (ASL) showed that the effects of N addition and sampling area on SR and ASL were highly significant (P < 0.0001), with their interactions being non-significant (Table 2). Effects of N treatment on proportional species loss were also highly significant (P < 0.0001), while the effects of sampling area were not significant (P = 0.2375; Table 2). After 10 years of treatment, N addition significantly (P < 0.05)

decreased species richness at each sampling area for most treatments except the treatment of 1.75 g N m2 year1 (Fig. 2a). There was no significant (P > 0.10) difference in species richness between control and treatment of 1.75 g N m2 year1 across all sampling areas in our study (Fig. 2a; Table S2). At the treatment of 5.25 g N m2 year1, species richness was not significantly different (P > 0.10) from that of control when sampling area ≥ 20 m2 (Fig. 2a; Table S2).

(a)

(c)

(b)

(d)

Fig. 2. Effects of N addition on (a) species richness, (b) percentage of overall species richness, (c) absolute species loss and (d) proportional species loss across different sampling areas. For a given sampling area, percentage of overall species richness was the ratio of species richness to overall species richness of 25 m2 plot. Absolute species loss was calculated as the difference of species richness between control and each N treatment (i.e. Richnesscontrol  Richnesstreat). Proportional species loss was calculated as the ratio of species loss in each N treatment to the mean species richness in control (i.e. 1  Richnesstreat/Richnesscontrol). The values for each treatment were the average of nine replicates (error bars denote SEM). Error bands in Fig 1c and d were standard errors at each sampling area in control plots. Arrows indicated the least sampling area that species richness of N addition plots was not significantly different (P > 0.10) from control. © 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

Scale dependence of N-induced species loss 755 (a)

(c)

(b)

(d)

Fig. 3. Proportional species loss (PSL) across N treatments at sampling area (a) 5 m2 and (b) 25 m2. (c) Slope of the relationship between PSL and N addition rate and (d) predicted critical threshold across different sampling areas. Solid lines are plot averages for each treatment fitted to a logarithmic curve excluding the control (error bands correspond to 95% confidence interval). Arrow indicates the predicted critical threshold (Ncrit) of N rate, below which there is no species loss at each sampling area. The x-axis is the logarithm of N addition rate for Fig. 3a and b.

Absolute species loss (i.e. Richnesscontrol  Richnesstreat) responded nonlinearly to sampling area (P < 0.05; Fig. 2c). At sampling area < 8 m2, the absolute species loss increased with sampling area; but when sampling area > 8 m2, species loss was relatively constant or even decreased (Fig. 2c). For example, at the sampling area of 1 m2, mean species richness of control and 5.25 g N m2 year1 treatment was 8.36 and 7.04, respectively, and absolute species loss was 1.32. At the sampling area of 8 m2, mean species richness of these two treatments was 16.41 and 13.36, respectively, and absolute species loss increased to 3.05. At the sampling area of 25 m2, mean species richness of these two treatments was 20.00 and 17.89, respectively, and absolute species loss decreased to 2.11 (Fig. 2c). Proportional species loss (PSL, i.e. 1  Richnesstreat/Richnesscontrol) significantly decreased with increasing sampling area for all N treatments (mean R2 = 0.82, P < 0.0001; Fig. 2d). Relationships between PSL and N rate were well fitted by loglinear models at all sampling rates, that is PSL = Intercept + Slope 9 Log (N) (R2 > 0.90, P < 0.01; Fig. 3a and b). Slope of PSL–N relationship significantly decreased with increasing sampling area (R2 = 0.54, P = 0.0017; Fig. 3c). We extrapolated the fitted curve of PSL–N relationship and calculated the predicted critical threshold of N rate (Ncrit) in which no species loss occurred at given sampling area. Ncrit significantly increased with sampling area (R2 = 0.48; P = 0.0035; Fig. 3d). The predicted

Ncrit was 1.14 (95% CI 0.52–2.45) and 1.74 (95% CI 0.52– 5.44) g N m2 year1 at sampling areas of 5 and 25 m2, respectively (Fig. 3d). SPECIES–AREA RELATIONSHIPS (SARS)

N addition did not change the shape of SARs, and the SARs were best fitted by power-law model (mean R2 = 0.97; Fig. 2a). Compared to control, SARs in N addition treatments saturated more slowly (Fig. 2b). For example, to encounter 90% of overall species richness in each plot, the minimum sampling area was about 13, 16 and 19 m2 for treatment of control, 5.25 and 28.00 g N m2 year1, respectively (Fig. 2b). N addition significantly increased the slope (z-value) of SAR, which saturated at 10.50 g N m2 year1 (R2 = 0.92, P = 0.0236; Fig. 4a). However, N addition decreased the intercept (c-value) of SARs with c = 9.08 in control and c = 6.15 in 10.50 g N m2 year1 (Fig. 4b). Relationships between c-value and N rate were best fitted by a loglinear model (R2 = 0.97, P = 0.0002; Fig. 4b). PLANT COMMUNITY COMPOSITION AND STRUCTURE

Compared to control plots, total number of individuals were slightly decreased at low N addition rates (1.75 and

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

756 Z. Lan et al. PARTITIONING MECHANISMS CONTRIBUTING TO

(a)

CHANGES IN SARS

Results from the simulations showed that, after N addition, changes in species abundance distribution linearly increased the slope (z-value) of SAR (R2 = 0.89, P = 0.0106), but this positive effect was moderated by loss of overall species richness (R2 = 0.99, P = 0.0002; Fig. 6a). However, the response of the intercept (c-value) of SAR to N were mainly due to changes in species abundance distribution (R2 = 0.82, P = 0.0222), and loss of overall species richness had no significant effects on c-value (R2 = 0.01, P = 0.7067; Fig. 6b). Compared to changes in species abundance distribution, changes in total number of individuals and species aggregation had much less effects on either z- or c-value (Fig. 6).

(b)

Discussion

Fig. 4. Effects of N addition on (a) slope (z-value) and (b) intercept (c-value) of the species–area relationship (SAR). For each treatment, SAR was fitted by the power-law model (S = cAz). The z- and c-values for each treatment were the average of nine replicates (error bars denote SEM).

5.25 g N m2 year1), but increased at high N addition rates (10.50, 17.50 and 28.00 g N m2 year1) by 36% (Fig. 5b). At the low N addition rates, decreases of individual numbers were mainly due to loss of small-sized species, Carex korshinskyi (perennial sedge). This species had the largest number of individuals among all species in control, but it was competitively excluded after N addition, leading to the decreased total individuals in the community. At the high N addition rates, the increase of individual numbers was due to the increase of a few N-favoured annual species (i.e. C. glaucum and A. amaranthoides). N addition decreased the abundance of rare species and altered species abundance distribution (Fig. 5c). Shannon evenness firstly increased at lower N addition rates and then declined at higher N addition rates (R2 = 0.79, P = 0.0442; Fig. 5c). To evaluate the effects of species aggregation on SAR, we calculated the rarefaction bias, with greater rarefaction bias indicating greater species aggregation. After N addition, rarefaction bias was significantly increased from 6.48% in control to 10.33% at the three highest N addition rates (R2 = 0.91, P = 0.0021; Fig. 5d), indicating that plant species in N-fertilized communities were more aggregated.

Our study provides a critical test of the scale dependence of N effects on species richness with a gradient of N addition rates. Although N addition decreased species richness across all scales, the proportional loss of species richness decreased with increasing sampling area. Meanwhile, the predicted critical threshold (Ncrit) of species loss increased with sampling area and ranged from 1.1 to 1.7 g N m2 year1, which is greater than a previous study based on 0.5 m2 quadrat (0.85 g N m2 year1) (Lan & Bai 2012). Although the predicted Ncrit is still within the range of N crit (0.5–2.0 g N m2 year1) reported from North American and European grasslands (Clark & Tilman 2008; Bobbink et al. 2010), our study clearly demonstrates that both the sensitivity (Ncrit) and magnitude of species loss in response to N addition are scale dependent. The regional species pool is less sensitive to N addition, and the Ncrit could be greater than this range at much larger scales. Our results are consistent with two other studies that N-induced species loss depends on scales (Crawley et al. 2005; Chalcraft et al. 2008). Our findings from long-term experimental studies in the Eurasian steppe, combined with those of related studies in European and North American terrestrial herbaceous plant communities (Crawley et al. 2005; Chalcraft et al. 2008), clearly suggest that N-induced species loss was generally mitigated at larger spatial scales. In addition, previous studies showed that N-induced species loss increased with experiment duration (Clark & Tilman 2008; Lan & Bai 2012). A recent study also found that the N effects on slope of SAR increased with experiment duration (Leps 2014). Therefore, the effects of spatial scale on biodiversity loss may also depend on the duration of N enrichment. Our results demonstrated that the scale dependence of Ninduced biodiversity loss was characterized by changes in the slope and intercept of species–area relationship (SAR). Through SAR decomposition, we found that contributions to changes in SAR were highly uneven among four potential mechanisms, including overall species richness, number of individuals, species abundance distribution and species

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

Scale dependence of N-induced species loss 757 (a)

(c)

(b)

(d)

Fig. 5. Effects of N addition on (a) overall species richness, (b) number of individuals (or total abundance of all species), (c) species abundances distribution and (d) rarefaction bias for each plot (25 m2). Evenness of species abundance distribution was calculated as Shannon evenness. A greater bias indicated that the community was more aggregated. For each treatment, the values were the average of nine replicates (error bars denote SEM).

aggregation. Changes in overall species richness and species abundance distribution significantly affected SAR after N addition, supporting the findings that these two factors were key to SAR (He & Legendre 2002; Tjørve et al. 2008; Chase & Knight 2013). The observed increases in the slope of SAR were mainly caused by the positive effect of changes in species abundance distribution over the negative effect of decline in overall species richness. Contrary to other studies (He & Legendre 2002; Tjørve et al. 2008; Chase & Knight 2013; Powell, Chase & Knight 2013), changes in total number of individuals had no significant effect on SAR in response to N addition. It seems likely that the total number of individuals in our study were so high (> 100 individual m2) that increasing individual numbers after N addition had a few effects on SAR. At low N addition rates (≤ 5.25 g N m2 year1), the total number of individuals decreased because of the decline of C. korshinskyi, a small and short sedge species. However, at the high N addition rates (≥ 10.5 g N m2 year1), the total number of individuals was increased. The increase in the total number of individuals at the high N addition rates was mainly caused by two reasons. First, in the highest N treatments, bunchgrasses (i.e. S. grandis, A. cristatum) were replaced by a large number of tall annual species (e.g. A. amaranthoides, C. glaucum); thus, the total individual number increases.

Secondly, without N enrichment, plant above-ground cover was discontinuous in semi-arid grasslands because of strong below-ground competition (e.g. N and water) (Burke et al. 1998; Bai et al. 2008). The high levels of N inputs removed N limitation for plants, reduced below-ground competition and thereby increased the total number of individuals, which decreased community openness and increased above-ground vegetation cover (Bai et al. 2010; Lan & Bai 2012). Therefore, changes in species composition and decreased belowground competition were attributable to the increased total number of individuals at the three highest N addition rates. Moreover, the increases in species aggregation in response to N addition also had little effect on the slope SAR, which is consistent to previous findings (Tjørve et al. 2008). Effect size of species aggregation on SAR has been suggested to increase with spatial heterogeneity and the extent of study (Tjørve et al. 2008; Scheiner et al. 2011). The experimental area in our study is comparatively homogeneous, and the extent of study is relatively small (25 m2) for each SAR. These factors could contribute to the limited effects of changes in species aggregation on SAR in our study. The Ninduced decreases in the intercept (c-value) of SAR were mainly due to negative effects of changes in species abundance distribution and species aggregation. The decline in

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

758 Z. Lan et al.

(a)

(b)

Fig. 6. Contribution of changes in overall species richness, number of individuals, species abundance distribution and species aggregation to the changes in (a) slope (z-value) and (b) intercept (c-value) of the species–area relationship (SAR) after N addition. For each treatment, the values were the average of nine replicates (error bars denote SEM).

overall species richness, in contrast, had no significant effects on the intercept of SAR. The result was inconsistent with the predictions of species pool hypothesis that species richness at small scales (associated with intercept of SAR) was determined by species pool size (Rosenzweig 1995; P€artel et al. 1996; Zobel 1997). Following N addition, the decline of overall species richness was mainly caused by the loss of rare species, which had limited contribution to community species richness (or c-value) at small scales. Given that N addition intensified interspecific competition, the results support the hypothesis that strong species interactions would lead to curvilinear local–regional relationships, that is regional species richness changes had little effect on the local assemblage (Cornell & Lawton 1992). Our results suggest that without identifying the scaling effects, we may substantially overestimate biodiversity losses associated with regional N deposition. For example, at addition rate of 5.25 g N m2 year1, the proportional species loss was 16–20% at the sampling area of 1–4 m2; however, it was only 11% at the sampling area of 25 m2. Increase in the slope (z-value) of SAR was in response to the positive effect of

changes in species abundance distribution over the negative effect of losses in overall species richness. The effects of N on overall species richness and species abundance distribution differ among ecosystems (Hillebrand et al. 2007; Chalcraft et al. 2008), and thus, slopes of SARs show context-dependent responses to N addition (Crawley et al. 2005; Leps 2014). Meanwhile, as sampling increases to include more area and spatial heterogeneity, changes in SAR will be driven more by species aggregation at larger scales (Keil et al. 2011). Future efforts should include more ecosystem types and larger ranges of scales to identify how site-specific characteristics (i.e. moisture, temperature and spatial heterogeneity) influence N effects on slope of SAR, which is the key to understand the scaledependent effects of N on biodiversity. Will regional species pool be maintained as z-value increases after N addition? Our findings showed that, at the treatments of 5.25 g N m2 year1, species richness decreased at small scales, but it was not significantly different from control when sampling area ≥20 m2. Above this rate of N addition, species richness declined at larger scales. Chalcraft et al. also found that, among 11 experiments which suffered N-induced plant diversity loss at small scales, only three of them maintained diversity at large scales after N addition (Chalcraft et al. 2008). In our study, the proportional species loss (i.e. 1  Richnesstreat/Richnesscontrol) and sampling area were significantly negatively correlated. By extrapolating the fitted curve of proportional species loss and sampling area, we predict that the smallest sampling area with no species loss compared to control is 41, 54, 85, 125 and 165 m2 at the N addition rate of 1.75, 5.25, 10.50, 17.50 and 28.00 g N m2 year1, respectively. Large-scale observational studies along N deposition gradients should be conducted as a first test of these predictions. Our study has important implications for the biodiversity conservation and ecosystem management of the vast semi-arid grasslands in the Eurasian Steppe region. First, our results suggest that both the critical threshold and magnitude of N-induced biodiversity loss are scale dependant. These findings highlighted the urgency of reassessing the previously estimated values of critical loads, which have been widely used as the backbone of biodiversity conservation and environmental legislation on transboundary pollution in Europe and United States (Posch, Vries & Hettelingh 1995; Sala et al. 2000; Pardo et al. 2011). Secondly, our findings demonstrate that the scale dependence of N-induced biodiversity loss is mainly triggered by reductions in overall species richness and changes in species abundance distribution. It has also been predicted that N deposition will continue to increase in this century (Galloway et al. 2008). Therefore, future biodiversity conservation should put high priorities on the grassland management practices that could efficiently mitigate the negative consequences of Ninduced decreases in species richness and shifts in species dominance, such as mowing or grazing (Collins et al. 1998).

Acknowledgements We thank Dunmei Lin, Xihunghua Sui, Kristin Powell, Xiujuan Qiao, Zhiyao Tang, Lorraine Clarke and Liyin Liang for help with data analysis, and Xue

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

Scale dependence of N-induced species loss 759 Bai, Xiaoliang Wang and many others at the Inner Mongolia Grassland Ecosystem Research Station, Chinese Academy of Sciences for their help with fieldwork. This project was supported by the National Natural Science Foundation of China (31030013, 31320103916) and Strategic Priority Research Program of the Chinese Academy of Sciences (XDA05050400).

Data accessibility Species richness and species abundance data are available from the Dryad Digital Repository (Lan et al. 2015).

References Arrhenius, O. (1921) Species and area. Journal of Ecology, 9, 95–99. Bai, Y.F., Han, X.G., Wu, J.G., Chen, Z.Z. & Li, L.H. (2004) Ecosystem stability and compensatory effects in the Inner Mongolia grassland. Nature, 431, 181–184. Bai, Y.F., Wu, J.G., Xing, Q., Pan, Q.M., Huang, J.H., Yang, D.L. & Han, X.G. (2008) Primary production and rain use efficiency across a precipitation gradient on the Mongolia plateau. Ecology, 89, 2140–2153. Bai, Y.F., Wu, J.G., Clark, C.M., Naeem, S., Pan, Q.M., Huang, J.H., Zhang, L.X. & Han, X.G. (2010) Tradeoffs and thresholds in the effects of nitrogen addition on biodiversity and ecosystem functioning: evidence from inner Mongolia Grasslands. Global Change Biology, 16, 358–372. Bobbink, R., Hicks, K., Galloway, J., Spranger, T., Alkemade, R., Ashmore, M. et al. (2010) Global assessment of nitrogen deposition effects on terrestrial plant diversity: a synthesis. Ecological Applications, 20, 30–59. Burke, I.C., Lauenroth, W.K., Vinton, M.A., Hook, P.B., Kelly, R.H., Epstein, H.E., Aguiar, M.R., Robles, M.D., Aguilera, M.O., Murphy, K.L. & Gill, R.A. (1998) Plant-soil interactions in temperate grasslands. Biogeochemistry, 42, 121–143. Chalcraft, D.R., Cox, S.B., Clark, C., Cleland, E.E., Suding, K.N., Weiher, E. & Pennington, D. (2008) Scale-dependent responses of plant biodiversity to nitrogen enrichment. Ecology, 89, 2165–2171. Chase, J.M. & Knight, T.M. (2013) Scale-dependent effect sizes of ecological drivers on biodiversity: why standardised sampling is not enough. Ecology Letters, 16, 17–26. Clark, C.M. & Tilman, D. (2008) Loss of plant species after chronic low-level nitrogen deposition to prairie grasslands. Nature, 451, 712–715. Clark, C.M., Bai, Y., Bowman, W.D., Cowles, J.M., Fenn, M.E., Gilliam, F.S., Phoenix, G.K., Siddique, I., Stevens, C.J., Sverdrup, H.U. & Throop, H.L. (2013) Nitrogen deposition and terrestrial biodiversity. Encyclopedia of Biodiversity, 2nd edn. (ed. S.A. Levin), pp. 519–536. Elsevier Publishing, Waltham, MA. Cleland, E.E. & Harpole, W.S. (2010) Nitrogen enrichment and plant communities. Year in Ecology and Conservation Biology 2010, Vol. 1195 (eds R.S. Ostfeld & W.H. Schlesinger), pp. 46–61. Blackwell Publishing, Oxford. Collins, M. & Simberloff, D. (2009) Rarefaction and nonrandom spatial dispersion patterns. Environmental and Ecological Statistics, 16, 89–103. Collins, S.L., Knapp, A.K., Briggs, J.M., Blair, J.M. & Steinauer, E.M. (1998) Modulation of diversity by grazing and mowing in native tallgrass prairie. Science, 280, 745–747. Cornell, H.V. & Lawton, J.H. (1992) Species interactions, local and regional processes, and limits to the richness of ecological communities: a theoretical perspective. Journal of Animal Ecology, 61, 1–12. Crawley, M.J., Johnston, A.E., Silvertown, J., Dodd, M., de Mazancourt, C., Heard, M.S., Henman, D.F. & Edwards, G.R. (2005) Determinants of species richness in the Park Grass Experiment. American Naturalist, 165, 179–192. Galloway, J.N., Townsend, A.R., Erisman, J.W., Bekunda, M., Cai, Z.C., Freney, J.R., Martinelli, L.A., Seitzinger, S.P. & Sutton, M.A. (2008) Transformation of the nitrogen cycle: recent trends, questions, and potential solutions. Science, 320, 889–892. Gilliam, F.S. (2006) Response of the herbaceous layer of forest ecosystems to excess nitrogen deposition. Journal of Ecology, 94, 1176–1191. Goldberg, D.E. & Miller, T.E. (1990) Effects of different resource additions on species-diversity in an annual plant community. Ecology, 71, 213–225. Guilhaumon, F., Mouillot, D. & Gimenez, O. (2010) mmSAR: an R-package for multimodel species-area relationship inference. Ecography, 33, 420–424. Hautier, Y., Niklaus, P.A. & Hector, A. (2009) Competition for light causes plant biodiversity loss after eutrophication. Science, 324, 636–638. He, F. & Legendre, P. (2002) Species diversity patterns derived from speciesarea models. Ecology, 83, 1185–1198.

Hillebrand, H., Gruner, D.S., Borer, E.T., Bracken, M.E., Cleland, E.E., Elser, J.J., Harpole, W.S., Ngai, J.T., Seabloom, E.W. & Shurin, J.B. (2007) Consumer versus resource control of producer diversity depends on ecosystem type and producer community structure. Proceedings of the National Academy of Sciences of the United States of America, 104, 10904–10909. Jenerette, G.D. & Shen, W. (2012) Experimental landscape ecology. Landscape Ecology, 27, 1237–1248. Keil, P., Biesmeijer, J.C., Barendregt, A., Reemer, M. & Kunin, W.E. (2011) Biodiversity change is scale-dependent: an example from Dutch and UK hoverflies (Diptera, Syrphidae). Ecography, 34, 392–401. Lan, Z. & Bai, Y. (2012) Testing mechanisms of N-enrichment-induced species loss in a semiarid Inner Mongolia grassland: critical thresholds and implications for long-term ecosystem responses. Philosophical Transactions of the Royal Society B: Biological Sciences, 367, 3125–3134. Lan, Z.C., Jennerette, G.D., Zhan, S.X., Li, W.H., Zheng, S.X. & Bai, Y.F. (2015) Data from: Testing the scaling effects and mechanisms of N-induced biodiversity loss: evidence from a decade-long grassland experiment. Dryad Digital Repository, doi: 10.5061/dryad.h7893. Leps, J. (2014) Scale- and time-dependent effects of fertilization, mowing and dominant removal on a grassland community during a 15-year experiment. Journal of Applied Ecology, 51, 978–987. Oksanen, J. (1996) Is the humped relationship between species richness and biomass an artefact due to plot size? Journal of Ecology, 84, 293–295. Oksanen, J. (1997) The no-interaction model does not mean that interactions should not be studied. Journal of Ecology, 85, 101–102. Pardo, L.H., Fenn, M.E., Goodale, C.L., Geiser, L.H., Driscoll, C.T., Allen, E.B. et al. (2011) Effects of nitrogen deposition and empirical nitrogen critical loads for ecoregions of the United States. Ecological Applications, 21, 3049–3082. P€artel, M., Zobel, M., Zobel, K. & van der Maarel, E. (1996) The species pool and its relation to species richness: evidence from Estonian plant communities. Oikos, 75, 111–117. Perez-Harguindeguy, N., Dıaz, S., Garnier, E., Lavorel, S., Poorter, H., Jaureguiberry, P. et al. (2013) New handbook for standardised measurement of plant functional traits worldwide. Australian Journal of Botany, 61, 167–234. Phoenix, G.K., Hicks, W.K., Cinderby, S., Kuylenstierna, J.C.I., Stock, W.D., et al. (2006) Atmospheric nitrogen deposition in world biodiversity hotspots: the need for a greater global perspective in assessing N deposition impacts. Global Change Biology, 12, 470–476. Posch, M., Vries, W.d. & Hettelingh, J.P. (1995) Critical loads of sulfur and nitrogen. Calculation and Mapping of Critical Thresholds in Europe: Status Report 1995 (eds M. Posch, P.A.M. de Smet, J.P. Hettelingh & R.J. Downing), pp. 31–41. RIVM, Bilthoven. RIVM Rep. 259101004. Powell, K.I., Chase, J.M. & Knight, T.M. (2013) Invasive plants have scaledependent effects on diversity by altering species-area relationships. Science, 339, 316–318. Preston, F.W. (1960) Time and space and the variation of species. Ecology, 41, 611–627. Rosenzweig, M.L. (1995) Species Diversity in Space and Time. Cambridge University Press, Cambridge, UK. Sala, O.E., Chapin, F.S., Armesto, J.J., Berlow, E., Bloomfield, J., Dirzo, R. et al. (2000) Global biodiversity scenarios for the year 2100. Science, 287, 1770–1774. Sandel, B. & Corbin, J.D. (2012) Scale-dependent responses of species richness to experimental manipulation of productivity and disturbance in Californian coastal grasslands. Journal of Vegetation Science, 23, 906–918. Scheiner, S.M. (2003) Six types of species-area curves. Global Ecology and Biogeography, 12, 441–447. Scheiner, S.M., Chiarucci, A., Fox, G.A., Helmus, M.R., McGlinn, D.J. & Willig, M.R. (2011) The underpinnings of the relationship of species richness with space and time. Ecological Monographs, 81, 195–213. Schellberg, J., Moseler, B.M., Kuhbauch, W. & Rademacher, I.F. (1999) Longterm effects of fertilizer on soil nutrient concentration, yield, forage quality and floristic composition of a hay meadow in the Eifel mountains, Germany. Grass and Forage Science, 54, 195–207. Socher, S.A., Prati, D., Boch, S., M€uller, J., Klaus, V.H., H€olzel, N. & Fischer, M. (2012) Direct and productivity-mediated indirect effects of fertilization, mowing and grazing on grassland species richness. Journal of Ecology, 100, 1391–1399. Stevens, C.J., Dise, N.B., Mountford, J.O. & Gowing, D.J. (2004a) Impact of nitrogen deposition on the species richness of grasslands. Science, 303, 1876–1879. Stevens, M.H., Bunker, D.E., Schnitzer, S.A. & Carson, W.P. (2004b) Establishment limitation reduces species recruitment and species richness as soil resources rise. Journal of Ecology, 92, 339–347.

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760

760 Z. Lan et al. Suding, K.N., Collins, S.L., Gough, L., Clark, C., Cleland, E.E., Gross, K.L., Milchunas, D.G. & Pennings, S. (2005) Functional- and abundance-based mechanisms explain diversity loss due to N fertilization. Proceedings of the National Academy of Sciences of the United States of America, 102, 4387– 4392. Tilman, D. (1987) Secondary succession and the pattern of plant dominance along experimental nitrogen gradients. Ecological Monographs, 57, 189–214. Tjørve, E. (2003) Shapes and functions of species-area curves: a review of possible models. Journal of Biogeography, 30, 827–835. Tjørve, E., Kunin, W.E., Polce, C. & Tjørve, K.M.C. (2008) Species-area relationship: separating the effects of species abundance and spatial distribution. Journal of Ecology, 96, 1141–1151. Vitousek, P.M. & Howarth, R.W. (1991) Nitrogen limitation on land and in the sea – how can it occur. Biogeochemistry, 13, 87–115. Vitousek, P., Aber, J., Howarth, R., Likens, G., Matson, P., Schindler, D., Schlesinger, W. & Tilman, D. (1997) Human alteration of the global nitrogen cycle: sources and consequences. Ecological Applications, 7, 737–750. Zobel, M. (1997) The relative of species pools in determining plant species richness: an alternative explanation of species coexistence? Trends in Ecology & Evolution, 12, 266–269.

Supporting Information Additional Supporting Information may be found in the online version of this article: Figure S1. Comparisons of species area relationships (SARs) built from two different sampling schemes (Scheiner 2003). Figure S2. Statistics of eight candidate models for the species area relationships (SARs). Table S1. Hypotheses regarding the effects of N addition on the intercept (c- value) and slope (z-value) of species area relationship. Table S2. Duncan’s multiple range tests reported from one-way ANOVA among different N addition rates at each sampling area.

Received 12 October 2014; accepted 9 March 2015 Handling Editor: Will Cornwell

© 2015 The Authors. Journal of Ecology © 2015 British Ecological Society, Journal of Ecology, 103, 750–760