Introduction

Grasslands cover approximately 40% of the terrestrial ice-free surface of the Earth and provide diverse ecosystem services including climate regulation, soil carbon storage, plant diversity maintenance and support for pollinators while contributing to human nutrition (Lamarque et al. 2011). In particular, legumes are one of the key plant functional groups in grasslands for their capacity to increase the nitrogen (N) availability by symbiotic N2 fixation, which, in turn, enhances the grassland net primary productivity, mitigates environmental pollution, and increases forage quality and productivity, critical for livestock production and stockless organic cropping systems (Lüscher et al. 2014; Suter et al. 2015). However, anthropogenic N inputs (in the form of fertilizers, manure, and atmospheric deposition) are changing the supply of N relative to phosphorus (P) in grasslands (Peñuelas et al. 2013) which can affect the symbiotic N2 fixation by legumes (Høgh-Jensen et al. 2002; Carlsson and Huss-Danell 2003; Stevens et al. 2004). The N and P availability can affect symbiotic N2 fixation by changing the legume biomass production in grasslands and the contribution of N derived from symbiotic N2 fixation to the total N content of legumes (percentage of legume N derived from atmosphere, %Ndfa) (Høgh-Jensen et al. 2002; Nyfeler et al. 2011; Peoples et al. 2012). Therefore, understanding the effects of N and P inputs on symbiotic N2 fixation by legumes is crucial to maintain grassland herbage productivity and additional ecosystem services such as plant biodiversity.

Nitrogen inputs can affect the legume biomass production and the %Ndfa and thus, symbiotic N2 fixation by legumes in grasslands (West et al. 2005; Nyfeler et al. 2011; Oberson et al. 2013). The %Ndfa often declines with availability of both ammonium and nitrate in soil (Leidi and Rodríguez-Navarro 2000; Peoples et al. 2012), because symbiotic N2 fixation is energetically expensive, and legumes take up reactive N, if available in the soil. In addition, grasses and non-leguminous forbs often displace legumes at high N availability due to their higher competitiveness for light (Soussana and Tallec 2010; Tognetti et al. 2021). Thus, anthropogenic N inputs to grasslands, can lead to decreased N2 fixation by legumes affecting the grassland functioning.

Symbiotic N2 fixation by legumes in grasslands can also be influenced by P inputs (Høgh-Jensen et al. 2002; Edwards et al. 2006). Low soil P levels reduce the activity of N2-fixing legume-associated bacteria (Edwards et al. 2006) due to the high ATP requirements of N2 fixation (Valentine et al. 2017). At the plant level, long-term P deprivation decreases nodular P concentration and reduces the energy status of the nodules and their capacity to fix N2, shifting the source of legume N nutrition from atmospheric N2 fixation towards uptake of reactive soil N (Valentine et al. 2017). Low P availability could reduce the rate of N2 fixation per unit legume biomass, removing the advantage that legumes might have over non-N2 fixing plants in N-poor conditions and reducing legume biomass production (Edwards et al. 2006). Thus, it can be expected that P addition increases legume biomass production and symbiotic N2 fixation. Furthermore, it has been observed that simultaneous additions of N and P (NP) can offset the negative effect of N addition on N2 fixation resulting in higher symbiotic N2 fixation compared to only N addition in tropical leguminous trees (Zheng et al. 2016). However, a recent study revealed that P addition enhances legume abundance in grasslands but does not mitigate the negative N effect when both elements are added simultaneously (Tognetti et al. 2021). Nevertheless, the extent to which the NP addition affects the N2 fixation rates in grasslands remains to be studied.

Most studies on the effect of N and P addition on N2 fixation in grasslands differ in experimental design, element addition rate, type of fertilizer used, and sampling procedure which leads to biases and uncertainties hampering our understanding of the main drivers of symbiotic N2 fixation at the global scale (Zheng et al. 2019). Therefore, a standardized and globally replicated experiment is needed to gain insight into N2 fixation. Recent work using a standardized and globally replicated experiment has shown that legume biomass production declines with N addition (Tognetti et al. 2021), but the extent to which N and P additions affect symbiotic N2 fixation rates was not evaluated.

Here, we use a standardized evaluation of symbiotic N2 fixation in a globally coordinated grassland experiment, replicated at sites spanning a wide range of climatic and edaphic characteristics, to shed light on the response of symbiotic N2 fixation in grasslands to N and P inputs. To evaluate the influence of single and combined N and P additions under various environmental conditions, we studied symbiotic N2 fixation in 17 natural and semi-natural grassland sites with natural abundance of legumes (i.e. legumes were not deliberately introduced for the study) located on four continents. We determined symbiotic N2 fixation based on the natural abundance of 15N in shoot plant biomass. This approach relies on the distinct isotopic N signature of atmospheric N2 and reactive soil N, which affects the plant N isotopic signature, depending on the source from which plants take up N (Amarger et al. 1979; Hoegberg 1997). We hypothesized that (i) N addition leads to a reduction in symbiotic N2 fixation since it decreases legume biomass and N2 fixation per unit legume biomass (i.e., the content of N derived from the atmosphere per unit of legume biomass), (ii) P addition enhances symbiotic N2 fixation by increasing legume biomass and N2 fixation per unit legume biomass, and (iii) the combined application of N and P increases N2 fixation by offsetting the N-induced P deficiency caused by N application.

Material and methods

Study sites

The 17 study sites (Table 1, Figure S1) are experimental grasslands part of the Nutrient Network Global Research Cooperative (NutNet, https://nutnet.org) (Borer et al. 2014, 2017) and were selected according to the criterion that legumes were recorded in a minimum of six out of the 12 experimental plots (Table S1). The study sites are natural and semi-natural grasslands representing the regional flora with natural abundance of legumes (i.e. legumes were not deliberately introduced for the study). The selected sites are distributed across four continents (Table 1, Figure S1) covering a wide range of climatic conditions: Mean annual temperature (MAT) ranged between -3.1 and 27.3 °C and the mean annual precipitation (MAP) from 243 to 1222 mm. Sites were located between 87 and 2320 m above sea level (Table 1).

Table 1 Continent, country, region/state, site name, latitude and longitude, elevation, mean annual precipitation (MAP), mean annual temperature (MAT), potential evapotranspiration (PET), estimated atmospheric nitrogen deposition (Ndep), year of study establishment and seasons of element addition (i.e., seasons from the establishment to the sample collection) of the 17 NutNet sites included in the study

An identical experiment design is replicated at each site with four treatments: control (Ctrl; no element addition), N addition (N; 100 kg N ha−1 yr−1 as slow-release urea with δ15N close to 0‰ (Choi et al. 2017)), P addition (P; 100 kg P ha−1 yr−1 as triple superphosphate (Ca(H2PO4)2.H2O)), and combined N and P addition (NP) (100 kg ha−1 yr−1 of both N and P). All treatments are replicated three times (n = 3) at each site, and the experiments are organized in a randomized block design with 25 m2 plots (5 × 5 m). All sites follow the same experimental design and sampling protocol. Further details about the experimental design and sampling can be found in Borer et al. (2014)

The climatic data were derived from Hijmans et al. (2005) based on the location of each site (Table 1). Nitrogen deposition (as kg N ha−1 yr−1) was estimated based on the location of each site (longitude and latitude) using the model output of Ackerman et al. (2019) for the year 2016 (Table 1). The soil properties of the control treatment at the time of establishment of the experiment are summarized in the Table 2. The methods used for soil analysis have been described in Seabloom et al. (2021).

Table 2 Selected soil properties (0–10 cm) at the time of experiment establishment of the 17 Nutrient Network (NutNet) sites evaluated in the study

Plant sampling

Aboveground biomass was sampled annually at the time of peak biomass using a standardized protocol (Borer et al. 2014). For the present study, we used the aboveground biomass of a single sampling between the years 2015 and 2020 depending on the site (as detailed in Table 1). Two 10 × 100 cm strips (covering area of 0.2 m2) of vegetation were clipped directly above the soil surface in a subplot of 1 × 1 m within each plot. The clipped plant biomass was sorted into the three functional groups: grasses, non-leguminous forbs and legumes, and oven-dried at 60 °C to a constant mass prior to weighing. Hereafter, we will refer to aboveground plant biomass as plant biomass. Representative subsamples of the biomass of the three plant functional groups from all plots were sent to the University of Bayreuth (Germany) for further analyses. The two most abundant grass, forb and legume species based on the cover estimates in the control plots for each year of sampling at each site are shown in Table 3. The cover estimate of each vascular plant was estimated visually in a 1 m2 subplot within each treatment plot when plant biomass was highest (peak biomass). In addition, the sites were classified as grassland with perennial or annual legumes according to the type (perennial or annual) of the two most abundant legume species (Table 3).

Table 3 Names of the two most abundant grass, forb and legume species and the type of dominating legumes (annual or perennial) at each site according to the cover estimates in the three control plots of each site for the year of sampling. In brackets, mean values of cover estimates (%) of each species in the control plots is given

Plant C, N and P concentration and stable isotope determination (δ15N)

In total, 553 dried plant samples were processed (201 grasses, 185 forbs and 167 legumes, Table S1). Plant biomass of each functional group was cut with scissors, homogenized and ground in a ball mill. The total C and N concentration and the isotopic composition were analyzed using continuous-flow isotope ratio mass spectrometry (NA 1108 elemental Analyzer, CE Instruments, Milano, Italy) coupled via ConFlo III open-split interface (Finnigan MAT, Bremen, Germany) to a delta S isotope ratio mass spectrometer (Finnigan MAT, Bremen, Germany) at the University of Bayreuth. The isotopic composition of N was expressed in δ notation, which represents the ‰ of variation compared to the international standard for natural 15N abundance measurements (atmospheric N isotope ratio). In addition, the total P concentration of plant biomass was determined by ICP-OES (Vista-Pro radial, Varian, Aschaffenburg, Germany) after nitric acid digestion.

Calculations

The percentage of legume N in aboveground biomass derived from atmospheric N2 fixation (%Ndfa), via N2 fixation by legumes, was calculated following the approach described by Amarger et al. (1979) and Hoegberg (1997) using Eq. 1:

$$\mathrm{Ndfa}\, (\%)=({\updelta }^{15}{\mathrm{N}}_{\mathrm{reference}}-{\updelta }^{15}{\mathrm{N}}_{\mathrm{legume}}\, )/({\updelta }^{15}{\mathrm{N}}_{\mathrm{reference}} -\mathrm{B})\times 100$$
(1)

where δ15Nreference is the δ15N of a non-fixing reference plant, δ15Nlegume is the δ15N of the legume aboveground biomass (called legume δ15N hereafter) in the evaluated plot, and B is the legume δ15N fully relying on atmospheric N2 fixation which accounts for any internal isotopic fractionation of the legume (Hoegberg 1997). We considered the mean of the δ15N of the two non-fixing functional groups (grasses and forbs) as the reference, which was calculated separately for each plot at each site (and is called reference δ15N hereafter).

We used the lowest legume δ15N value of all plots at each site as the B value (Eq. 1), similar to previous studies (Hansen and Vinther 2001; West et al. 2005; Roscher et al. 2011; Oberson et al. 2013). This approach relies on the assumption that the legumes in the plot with the lowest legume δ15N receive 100% of the N from symbiotic N2 fixation. We used the lowest legume δ15N from all four treatments because the legume δ15N was not significantly affected by element addition (see below in the Results section), similar to previous studies that considered different element addition treatments (Oberson et al. 2013; Tzanakakis et al. 2017).

We detected a relative 15N-depletion in reference δ15N compared to legume δ15N (i.e., lower reference δ15N than legume δ15N) at seven sites [Koffler (Canada), Hopland, Spindletop and Bunchgrass (USA), Bad Lauchstädt and Bayreuth (Germany) and Val Mustair (Switzerland)] which challenged the estimation of %Ndfa using Eq. 1 because it resulted in negative %Ndfa estimates. We observed a decrease in reference δ15N (mean of grasses and forbs) with increasing elevation of the study site (r2 = 0.301, p = 0.024, Figure S2), and no significant relationship between elevation and legume δ15N (p = 0.337). This observation is consistent with previous studies, showing that δ15N of non-fixing plants decreases with increasing elevation (Jacot et al. 2000; Craine and Lee 2003; Huber et al. 2007; Zhou et al. 2016). We adjusted the reference δ15N for the effect of elevation (elevation adj. δ15Nreference), assuming that all sites would be located at an elevation of 0 m a.s.l., using the slope of the regression line describing the relationship between the δ15N of the reference and elevation (see Figure S2), as:

$$Elevation\,adj.\;\mathrm\delta^{15}N_{reference}\,\left(\permil \right)=\mathrm\delta^{15}N_{reference}\,\left(\permil \right)-\lbrack-0.002\ast Elevation(\mathrm m)\rbrack$$
(2)

Next, Ndfa (%) was calculated based on the Elevation adj. δ15Nreference as follows:

$${\text{Ndfa }}\left( \% \right) \, = \, ({\text{Elevation adj}}.\delta^{{{15}}} {\text{N}}_{{{\text{reference}}}} - \, \delta^{{{15}}} {\text{N}}_{{{\text{legume}}}} ) \, /({\text{Elevation adj}}.\delta^{{{15}}} {\text{N}}_{{{\text{reference}}}} - {\text{ B}}){\text{ x 1}}00$$
(3)

Further details about the elevation adjustment and %Ndfa calculation are presented as Supporting Information.

Legume N uptake was calculated based on the legume N concentration at peak biomass as:

$$\mathrm{Legume}\;N\;\mathrm{uptake}\;(\text{mg N g}^{-1}\;\mathrm{yr}^{-1})=\mathrm{legume}\;N\;\mathrm{concentration}\;(\text{mg N g}^{-1})\times\mathrm{yr}^{-1}$$
(4)

N2 fixation per legume biomass was calculated for the plots with legumes as:

$${\mathrm N}_2\;\mathrm{fixation}\;\mathrm{per}\;\mathrm{unit}\;\mathrm{legume}\;\mathrm{biomass}\;(\text{mg N g}^{-1} \;\mathrm{yr}^{-1})=\mathrm{Legume}\;N\;\mathrm{uptake}\;(\text{mg N g}^{-1}\;\mathrm{yr}^{-1})\times\;\mathrm{Ndfa}\;(\%)\times0.01$$
(5)

Symbiotic N2 fixation per area grassland was calculated as follows:

$${\mathrm N}_{2\;}\mathrm{fixation}\;\mathrm{per}\;\mathrm{area}\;(\text{g N m}^{-2} \;\mathrm{yr}^{-1})=\mathrm{legume}\;\mathrm{biomass}\;(\text{g m}^{-2})\times\;{\mathrm N}_2\;\mathrm{fixation}\;\mathrm{per}\;\mathrm{legume}\;\mathrm{biomass}\;(\text{mg N g}^{-1} \;\mathrm{yr}^{-1})$$
(6)

Symbiotic N2 fixation per total grassland biomass (with grassland biomass being the sum of the biomasses of all three functional groups) was calculated for all plots as:

$${\text{N}}_{2\;}\text{fixation per total biomass }(\text{mg N g}^{-1}\text{yr}^{-1})={\text{ N}}_2\text{fixation per area }\left(\text{g N m}^{-2}\text{yr}^{-1}\right)\text{ x 1}000\,/\text{total biomass }\left(\text{g m}^{-2}\right)$$
(7)

We calculated means and standard errors of all four treatments across all 17 sites. If legumes were absent in a plot, the legume N stock (in g ha−1) or legume N uptake of this plot was assumed to be zero. If the legume biomass was zero, we assumed that symbiotic N2 fixation per unit biomass or area was also zero. When calculating the mean of the N or P concentration of the biomass of all three functional groups, we considered only plots with biomass of the respective functional group. Similarly, when calculating the mean of symbiotic N2 fixation per legume biomass, only plots with legume biomass > 0 and valid %Ndfa were considered. In contrast, when calculating symbiotic N2 fixation per unit area of grassland or total grassland biomass, all plots with legumes and valid %Ndfa as well as plots without legumes were considered.

Statistical analyses

Data were analyzed using linear mixed models with the software SPSS 27 (IBM SPSS, Inc., Chicago, USA). Because data were not normally distributed (Shapiro–Wilk-test, p > 0.05), all variables except the δ15N values (including both positive and negative values) were log-transformed. The different element addition treatments (Ctrl, N, P, NP) and site- and plot-level covariates were used as fixed factors, and block as a random factor where block was nested within site. The site-level covariates included in the linear mixed model were MAP, MAT, water availability index (MAP/potential evapotranspiration), the estimates of N deposition and the legume proportion of biomass in the control treatments at each site. Because sites were set up in different years, the number of years of element addition was considered as a site-level covariate (redundant in most of the evaluated parameters with no significant effect). Soil properties summarized in Table 2 (except total organic carbon (TOC), which was highly correlated with total nitrogen (TN)) were included in the linear mixed model as plot-level covariates. The interactions between treatment and the plot-level covariates were initially considered in the model, although after a selection based on Akaike Information Criterion only the interactions ‘treatment x TN’ and ‘treatment x soil pH’ remained in the final model as covariates. When a significant treatment effect (p < 0.05) was found, LSD post hoc test (p < 0.05) was used for comparison of means of the element addition treatments.

Additionally, we evaluated how the type of dominating legumes (perennial or annual) at each site affected legume biomass, legume N concentration, legume δ15N, %Ndfa and symbiotic N2 fixation per area of grassland. A linear mixed model was used with treatment (Ctrl, N, P and NP), type of dominating legumes (annual or perennial), their interaction and the site- and plot-level covariates as fixed factors, and block as a random factor nested within site. The model was performed as previously described.

We calculated the response to nutrient addition of legume biomass, %Ndfa and N2 fixation per unit area as:

$$\mathrm{Response}=\mathrm{Ln}\;(({\mathrm Y}_{\mathrm{treatment}}+1)/({\mathrm Y}_{\mathrm{control}}+1))$$
(8)

where Ytreatment is the value of legume biomass, %Ndfa or N2 fixation per unit area in the N, P or NP addition treatment and Ycontrol is the mean value of legume biomass, %Ndfa or N2 fixation per unit area in the control. The response was calculated separately for each site. We added 1 to the numerator and denominator to remove zeros before the logarithmic transformation.

We performed stepwise multiple regression analyses to evaluate the impact of site- and plot-level covariates on the response of legume biomass, %Ndfa and N2 fixation per unit area to nutrient addition. Stepwise multiple regressions analyses were performed using the site-scale factors (MAP, MAT, water availability, N deposition and the legume proportion of biomass in the control treatment at each site) and plot-scale soil properties (TN, carbon-to-nitrogen ratio (C:N ratio), available P, and soil pH). Collinearity was evaluated based on the variance inflation factor. The multiple regression analyses were performed separately for the three different element addition treatments (i.e., N, P and NP).

Results

Plant aboveground biomass

On average across sites, N and P addition increased total plant biomass (the sum of the biomasses of all three functional groups) by 32 and 28%, respectively, compared to the control (Fig. 1). The combined addition of N and P increased the total plant biomass by 72%, from a mean of 3040 kg ha−1 in the control to 5222 kg ha−1 in the NP treatment. The biomass of both the grasses and the forbs increased significantly with NP addition (Fig. 1) but did not respond significantly to N or P addition alone.

Fig. 1
figure 1

Aboveground biomass at peak biomass of grasses, forbs, and legumes as well as total biomass (grasses, forbs and legumes together) as affected by single and combined N and P addition. Bars show means (n = 51) with standard errors. Ctrl, control; N, nitrogen addition; P, phosphorus addition; NP, nitrogen and phosphorus addition. All plots were included in this calculation irrespective of amount of biomass of the functional groups. Different letters indicate significant differences (p < 0.05) among treatments, tested separately for each plant functional group

Legume biomass was highly variable among sites and ranged from 0.1 kg ha−1 (Potrock, Argentina) to 1082 kg ha−1 (Bad Lauchstädt, Germany) in the control treatment (Table S2). Compared to the control, biomass of legumes was reduced by 65% and 45% in the N and NP treatments, respectively, while in the P treatment was increased by 77% (although this difference was not statistically significant) (Fig. 1). Compared to the P treatment, biomass of legumes was significantly reduced by 81% and 69% by N and NP addition. The percent of legumes in the total plant biomass (the sum of the biomasses of all three functional groups) was 9.9% in the control, 4.7% in the N, 12.0% in the P, and 3.3% in the NP treatment (Fig. 1). No significant interaction between the treatments and the type of dominant legumes (annual or perennial) at each site was observed on legume biomass, although legume biomass was significantly higher in the grasslands dominated by annual than by perennial legumes (Figure S3A). Similar differences between P and N and NP treatments were observed in the response ratio to nutrient addition of legume biomass (Figure S4A).

Nitrogen and phosphorus concentrations in aboveground biomass

The addition of N significantly increased plant N concentrations (Fig. 2A), while the addition of P significantly increased plant P concentrations (Fig. 2B) in both grasses and forbs compared to the control. In contrast, legume N concentration was not affected by N addition, whereas legume P concentration was enhanced by both the P and NP treatments (Fig. 2A, B). No significant interaction between the treatments and the type of dominating legumes (annual or perennial) was observed in legume N concentration (Figure S3B). The biomass N:P ratio in the control treatment was 8.5 for grasses, 7.3 for forbs, and 13.9 for legumes (Fig. 2C). The N:P ratio of legumes was significantly reduced by the addition of P (N:P = 8.7) and NP (N:P = 7.9) compared to the control and the N treatment (Fig. 2C). Similarly, P addition decreased the N:P ratio of grasses and forbs (Fig. 2C).

Fig. 2
figure 2

Plant aboveground nitrogen (A) and phosphorus (B) concentrations and plant aboveground nitrogen to phosphorus ratio (mass: mass) (C) of grasses, forbs and legumes as affected by single and combined N and P addition. Bars show means with standard errors. White numbers at the base of each column indicate the number of replicates. Only plots with biomass of the functional group > 0 were included in the calculation of the element concentrations of each functional group. Ctrl, control; N, nitrogen addition; P, phosphorus addition; NP, nitrogen and phosphorus addition. Different letters indicate significant differences (p < 0.05) among treatments, tested separately for each plant functional group

The legume N stock was significantly higher in the control (0.71 g N m−2) and P treatment (1.34 g N m−2) than in the N (0.24 g N m−2) and NP (0.36 g N m−2) treatment (Figure S5A). The P stock of grasses and forbs was increased by P and NP addition in comparison to the control, while the P stock of legumes was only increased in the P treatment compared to the control (Figure S5B).

Plant isotopic composition

We observed a decrease in δ15N of the reference functional groups (grasses and forbs) in the control treatment with increasing elevation of the study site (r2 = 0.301, p = 0.024, Figure S2). Therefore, we adjusted the reference δ15N of all plots and treatments for elevation (see Sect. Calculations). After the recalculation of the reference δ15N, the mean of the elevation-adjusted δ15Nreference was + 0.81‰ in the control treatment (Table S3). The δ15Nreference before and after the elevation-adjustment was significantly higher in the N and NP treatments than in the control and P treatment (Table S3). The legume δ15N across all sites was unaffected by treatments (Table S3).

Symbiotic N2 fixation per unit legume biomass

The mean %Ndfa in the control treatment was 65.8% across all 17 sites (Fig. 3A). No significant difference in %Ndfa among treatments was found. However, %Ndfa was slightly higher in the P and NP treatments (69.8 and 70.9%, respectively) than in the control and N treatments (65.8 and 64.2%, respectively). Similarly, no significant difference was observed in the %Ndfa response to nutrient addition (Figure S4B). The response of %Ndfa to NP addition was positively related to the proportion of legumes in the total biomass in the control treatment, initial soil N and water availability index, while it was negatively related to N deposition (Table 4, Table S4). There were no significant linear regression models (p > 0.05) for the %Ndfa response to N or P addition. We observed a significantly higher %Ndfa in the ten sites dominated by perennial legumes (74.0%) compared to the seven sites dominated by annual legumes (58.0%) (Fig. 4). In addition, the interaction between treatment and type of dominating legumes revealed that the single N addition limited the differences in %Ndfa between annual and perennial legumes (Fig. 4).

Fig. 3
figure 3

Percentage of legume biomass-nitrogen derived from atmospheric N2 fixation (%Ndfa) (A), symbiotic N2 fixation by legumes per legume biomass (B) as well as symbiotic N2 fixation by legumes per unit area of grassland (C) and per total grassland biomass (D) as affected by single and combined N and P addition. White number at the base of each bar shows the number of replicates. Only plots with legume biomass > 0 were included in panels A and B, whereas panels C and D also include plots with no legume biomass. Bars represent means with standard errors. Ctrl, control; N, nitrogen addition; P, phosphorus addition; NP, nitrogen and phosphorus addition. Treatments with different letters represent significant differences at p < 0.05

Table 4 Regression models of the response of legume biomass, the proportion of N derived from atmosphere (%Ndfa), and symbiotic N2 fixation per unit area to element addition (nitrogen (N), phosphorus (P) and their combined application (NP)) as a function of site-scale environmental factors (MAP, MAT, water availability (Aw), N deposition (Ndep) and legume proportion of biomass in the control treatment at each site (prop)) and plot-scale soil properties (total nitrogen (TN), soil carbon to nitrogen ratio (C:N), available phosphorus (P) and soil pH). A dash (-) indicates that no significant (p < 0.05) model was found
Fig. 4
figure 4

Percentage of legume biomass-nitrogen derived from atmospheric N2 fixation (%Ndfa) as affected by single and combined N and P addition (Treatment), the type of dominating legumes (annual or perennial) (Type) and their interaction (T x T). Bars show means with standard errors. Numbers at the base of each column indicate the number of replicates. Ctrl, control; N, nitrogen addition; P, phosphorus addition; NP, nitrogen and phosphorus addition. * and *** indicate significant effect at p < 0.05 and p < 0.001, respectively; n.s. not significant. Different letters indicate significant differences (p < 0.05) between annual and perennial within each treatment

Mean N2 fixation per unit legume biomass in the control treatment across all 17 study sites was 18.1 mg N g−1 yr−1 legume biomass (Fig. 3B). No significant difference among treatments in N2 fixation per legume biomass was observed due to the lack of element addition effect on legume N concentration and %Ndfa. However, N2 fixation per unit legume biomass was slightly higher in the NP treatment (20.2 mg N g−1 legume biomass yr−1) than in the other three treatments (Fig. 3B).

Symbiotic N2 fixation per unit area grassland

Symbiotic N2 fixation per unit area in the control treatment across all 17 study sites was 3.5 kg N ha−1 yr−1 (Fig. 3C). It ranged from 0.002 kg ha−1 yr−1 (Potrock, Argentina) to 11.9 kg N ha−1 yr−1 (Bad Lauchstädt, Germany). Across all sites, N2 fixation per area was 1.39 kg N ha−1 yr−1 in the N treatment and 2.13 kg N ha−1 yr−1 in the NP treatment. Thus, N2 fixation was significantly reduced in the N and NP treatment by 60 and 39%, respectively, compared to the control. Similarly, N2 fixation per unit area was also 63% lower in the N and 43% lower in the NP treatment than in the P treatment (3.71 kg N ha−1 yr−1) (Fig. 3C). In contrast, P addition had no significant effect on N2 fixation per area compared to the control. The number of years of element addition was a redundant site-level covariate in the linear mixed model, indicating that the different ages of the sites did not influence the estimation of N2 fixation in the present study. In addition, no significant effect of the type of dominating legumes (annual or perennial) on N2 fixation per area was observed (Figure S3D). Similarly, the response to N and NP addition of the N2 fixation per unit area was significantly lower than the response to P addition (Figure S4C). Soil pH, N deposition, MAT and MAP were negatively related with the response of N2 fixation per unit area to N addition (r2 = 0.550, p = 0.001, Table 4, Table S4). In addition, the response of N2 fixation per unit area to P addition was negatively related with soil N (r2 = 0.162, p = 0.009, Table 4, Table S4), while the response of N2 fixation per unit area to NP addition was positively related with soil N, MAT and N deposition (r2 = 0.406, p = 0.001, Table 4, Table S4).

The mean N2 fixation per unit total biomass in the control treatment was 1.23 mg N g−1 yr−1 across all 17 study sites (Fig. 3D). The N2 fixation per total biomass was reduced by 60% under N and 73% under NP addition compared to the control (Fig. 3D). Similarly, N2 fixation per total biomass was lower in the N and NP treatments compared to the P treatment.

Discussion

Our results reveal that the addition of N decreased the rate of N2 fixation per grassland area compared to the control, and this effect was not reversed by additional P amendment (Fig. 3C). The reduced N2 fixation per area grassland was due to the reduction in legume biomass, and not an altered N2 fixation rate per unit legume biomass. Tognetti et al. (2021) recently showed that legume biomass was negatively affected by N addition in grasslands on several continents. Our study goes further, demonstrating that this reduction in legume biomass causes the N2 fixation rate per unit area to decrease from 3.50 kg N ha−1 yr−1 to 1.39 kg N ha−1 yr−1 due to N addition across all sites. Further, we found that N2 fixation per legume biomass was not significantly affected by N or P addition, which is an important finding since it suggests that the expression of the enzyme that fixes atmospheric N2 is not downregulated in legumes in response to high availability of reactive N in the soil as observed by Menge and Hedin (2009).

δ15N patterns in plant functional groups

We observed a significant negative correlation between site elevation and the δ15N of the reference plants (grass and forbs; Figure S2), which suggests low soil δ15N at high elevations. Similar observations were described in previous global reviews (Amundson et al. 2003) and in studies along altitudinal gradients (Vitousek et al. 1989; Jacot et al. 2000; Craine and Lee 2003; Huber et al. 2007; Zhou et al. 2016). The underlying reason is the relationship between elevation and MAT and MAP, two of the main drivers of plant δ15N as described by Craine et al. (2009) and Zhou et al. (2016). The elevation-dependent δ15N of non-fixer plants (reference plants) likely caused difficulties when applying the 15N natural abundance method to determine symbiotic N2 fixation about altitudinal gradients in previous studies (Vitousek et al. 1989; Jacot et al. 2000). The unique global design of this study allowed us to correct this elevation effect on reference plants and to determine symbiotic N2 fixation across a large number of sites on different continents. The relationship between elevation and δ15N of grasses and forbs identified here will likely be of use also in future studies. However, external inputs of 15N-depleted N, such as cattle urine, large inputs of legume-derived N or atmospheric N deposition cannot be dismissed as other factors affecting δ15N of grasses and forbs (Jacot et al. 2000; Hansen and Vinther 2001; Gehring and Vlek 2004). Further details about the elevation adjustment and 15N natural abundance method are presented as Supporting Information.

Symbiotic N2 fixation per area decreased by N addition

We found that N addition significantly reduced the rate of symbiotic N2 fixation per area of grassland due to a reduction in legume biomass without affecting the N2 fixation rate per unit legume biomass. The negative effect of N addition on symbiotic N2 fixation is consistent with previous single-site studies observing a reduction in symbiotic N2 fixation by N addition (West et al. 2005; Nyfeler et al. 2011; Burchill et al. 2014; Tzanakakis et al. 2017) and a recent study showing that also non-symbiotic N2 fixation in the soil is decreased by N addition at sites of the Nutrient Network experiment (Schleuss et al. 2021). Our results suggest that continuous anthropogenic N enrichment of grasslands can lead to a decrease in legume biomass production, which in the long-term can limit symbiotic N2 fixation. The most plausible explanation for reduced rates of N2 fixation per unit area due to N addition is that higher soil N availability allows grasses and forbs to outcompete legumes (via competition for light) and reduce legume biomass (Suding et al. 2005; Soussana and Tallec 2010; Tognetti et al. 2021). The size of the response of symbiotic N2 fixation per area to N addition was affected by several abiotic factors. Addition of N reduced N2 fixation per area more strongly at sites with higher soil pH, atmospheric N deposition, MAT, and MAP (Table 4, Table S4). The reasons for this could be that i) legumes and rhizobium strains from neutral and alkaline sites are less tolerant to soil acidification caused by urea addition (Hungria and Vargas 2000), that ii) N addition has a larger effect at sites where the antrophogenic N input through atmospheric N deposition is already large, and iii) N2 fixation is more sensitive to N addition at sites where the N2 fixation is not constrained by temperature or water availability (Houlton et al. 2008; Tognetti et al. 2021).

In contrast to our second hypothesis, we observed similar rates of symbiotic N2 fixation per area of grassland in the control and P treatment, which can be attributed to a lack of P limitation of N2 fixation in the control treatment as further indicated by the low N:P ratio of legume biomass (13.9) (Güsewell 2004). The reason for the lack of effect to P addition is likely that legumes have evolved very effective mechanisms to increase their P uptake from different soil P pools through the release of phosphatases or organic acids into the rhizosphere (Nuruzzaman et al. 2006). The negative relationship between the response of N2 fixation per area and soil N to P addition (Table 4, Table S4) might indicate that at N-limited sites, legumes could invest the added P in symbiotic N2 fixation to overcome the N limitation (McKey 1994; Houlton et al. 2008; Soussana and Tallec 2010).

The application of P in combination with N did not counterbalance the negative impact of N addition on symbiotic N2 fixation per area, in contrast to our third hypothesis. Our results indicate that the negative impact of N addition on symbiotic N2 fixation per area is not a result of N-driven P deficiency because in this case, combined addition of NP would have offset the negative impact of single N addition. Thus, the addition of P does not seem to be a suitable strategy to enhance symbiotic N2 fixation in a scenario of anthropogenic N enrichment of grasslands.

We observed no significant difference in the symbiotic N2 fixation per area between the sites with annual or perennial legumes, because the higher biomass at sites dominated by annual legumes was counterbalanced by the higher %Ndfa of the sites dominated by perennials, irrespective of treatment (Fig. 4, Figure S3). This shows that there are no substantial differences in N2 fixation on an area basis between grasslands dominated by annual and perennial legumes.

Across treatments, we found relatively low symbiotic N2 fixation per area grassland (3.4 kg N ha−1 yr−1 in the control) compared to other studies (Carlsson and Huss-Danell 2003; Nyfeler et al. 2011; Peoples et al. 2012; Oberson et al. 2013). The reason for this seems to be our focus on natural and semi-natural grasslands with natural abundance of legumes (i.e. legumes were not deliberately introduced for the study), and the inclusion of some sites with low overall biomass production. The rate of N addition used in the present study (100 kg N ha−1 yr−1) exceeds any present and even projected atmospheric N deposition levels (Ackerman et al. 2019). However, considering the common fertilization rates used in managed grasslands ranging from 20–30 up to 400 kg N ha−1 yr−1 (Oenema et al. 2012; Klaus et al. 2018), the experimental rate used in our study resembles a realistic situation for many grasslands. The reduction of symbiotic N2 fixation in grasslands by N addition increases the dependence of grassland biomass productivity on fertilization, which has several economic and environmental drawbacks (Lüscher et al. 2014).

Symbiotic N2 fixation per unit legume biomass not affected by element addition

Although relatively high rates of N and P were added, we found no significant response of N2 fixation per unit legume biomass to N, P or NP addition (Fig. 3B). This is due to the lack of element addition effect on legume N concentration and %Ndfa. Our finding is in disagreement with previous field studies in grasslands reporting a positive effect of P addition and a negative effect of N addition on %Ndfa (Høgh-Jensen et al. 2002; West et al. 2005; Burchill et al. 2014; Tzanakakis et al. 2017). The lack of response of N2 fixation per unit legume biomass to P addition could indicate that symbiotic N2 fixation was not limited by P in the majority of grasslands included in the present study, as P addition should increase N2 fixation per unit biomass under strong P limitation due to the high ATP requirements of N2 fixation (Almeida et al. 2000; Høgh-Jensen et al. 2002; Edwards et al. 2006). Similarly, the lack of differences between N and NP addition on N2 fixation per unit biomass indicates that N addition did not induce a P limitation of N2 fixation as previously described in pot experiments or tree plantations (Leidi and Rodríguez-Navarro 2000; Zheng et al. 2016). Otherwise, combined NP would have increased the N2 fixation per unit biomass compared to single N addition. The reason why N addition did not cause a P limitation of N2 fixation is likely that soil P availability is relatively high since the sites are located in the temperate zone which is dominated by relatively young soils (Figure S1).

The lack of response of N2 fixation per unit biomass to N addition contrasts with previous results (Carlsson and Huss-Danell 2003) including experiments using urea as N source (Burchill et al. 2014), as in our study. We speculate that N addition had no significant effect on N2 fixation per unit biomass because grasses and forbs were N limited (as indicated by the low N:P ratio of grasses and forbs in the control treatment), and their efficient uptake of additional N reduced the availability of added N to legumes, as described in previous studies (Nyfeler et al. 2011; Peoples et al. 2012; Oberson et al. 2013). In contrast, in other studies where legumes were deliberately introduced leading to a higher legume abundance (West et al. 2005; Burchill et al. 2014), non-fixing plants may not take up the added N which, in turn, could reduce the N2 fixation per unit biomass. Another explanation might be that most legume species are permanent, rather than facultative N2 fixers and cannot shift their N source in spite of increased soil N availability (Menge and Hedin 2009).

We observed a higher N2 fixation per unit biomass in the grasslands dominated by perennial legumes compared to the sites dominated by annual legumes. The reason for this could be that perennials can build up a symbiosis with N2 fixation that last for several years, whereas annuals have to establish a new symbiosis with N2 fixing microorganisms every year which makes this symbiosis likely less effective (Primieri et al. 2022). The differences in N2 fixation per unit biomass between annual and perennial legume sites disappeared in the N treatment (Fig. 4D). This finding suggests that N addition has a very similar effect on N2 fixation per unit biomass irrespective of whether the grassland is dominated by annual or perennial legumes.

Conclusions

Our results show that soil N enrichment by anthropogenic activities significantly reduces N2 fixation in grasslands, and these effects cannot be reversed by additional P amendment. The negative effect of N addition on symbiotic N2 fixation by legumes per area was caused exclusively by the negative effect of N addition on legume biomass, and not by an effect on the N2 fixation rate per unit biomass, indicating that the legumes did not alter the rate at which they fixed N2 per unit biomass. This reduction in symbiotic N2 fixation per area increases the dependence of grassland productivity on fertilization and can ultimately change the ecological functioning of grasslands, affecting their net primary productivity as well as their above and belowground biodiversity, forage quality and provision of additional ecosystem services. Further, the unique global design of this study allowed us to derive an equation to correct for the effect of elevation on the isotope signature of N in grasses and non-fixing forbs which will be useful in future studies.