Introduction

Elemental cycles are strongly interrelated through coupled biotic and abiotic processes (Gao et al. 2013; Soussana and Lemaire 2014; Ochoa-Hueso et al. 2019). In the ecological literature, the coupling or decoupling of chemical elements, many of which are essential nutrients for living organisms, often refers to variations in nutrient proportions/ratios (Sterner and Elser 2003; Ochoa-Hueso et al. 2019). However, the coupling of elemental cycles can also be defined based on the co-variation of the different chemical elements found in ecosystems (Ochoa-Hueso et al. 2021a), which has been shown to depend on the atomic properties of elements, such as atomic mass and ionization energy, as well as on biotic (e.g., plant cover and fungi:bacteria ratios) and abiotic (e.g., soil pH, soil texture, temperature, precipitation) drivers (Kabata-Pendias and Pendias 2000; Schlesinger et al. 2011; Ochoa-Hueso et al. 2021c). Within this new paradigm, elemental coupling can be defined based on the level of spatial association of multiple soil elements, which provides an approximation of how biogeochemical cycles are interrelated through biotic and abiotic processes (Ochoa-Hueso et al. 2021a). Greater elemental coupling emerges when processes, both biotic and abiotic, controlling the release and immobilization of multiple elements in the soil solution are coordinated (i.e. they covary, either positive or negatively), while reductions in the degree of coupling occur when these processes are out of sync.

Ecosystem processes such as plant growth or microbial activity alter the degree of coupling among elemental cycles (Delgado-Baquerizo et al. 2013; Soussana and Lemaire 2014). For example, environmental factors triggering plant growth may temporarily lead to a decoupling of elemental cycles in soils because plants release and capture different nutrients in different proportions depending on their needs (Rumpel et al. 2015). However, under normal circumstances, this break down of the spatial co-occurrence of soil elements is temporary, as elemental cycles will be recoupled again through litter deposition, facilitating their dynamic spatial reorganization in soils (Ochoa-Hueso et al. 2019). In contrast, environmental factors that impair plant growth, such as reduced water availability, should lead to more coupled associations among multiple elements because there would be no active uptake of nutrients, thus maintaining the orderly distribution of elements in the form of minerals, organic matter, and salts (Ochoa-Hueso et al. 2019). Further, global environmental change drivers that alter resource availability, such as nitrogen (N) deposition, are expected to influence the degree of association among elemental cycles by distinctly affecting the bioavailability as well as biological demand and immobilization of N and other nutrients in various biotic and abiotic pools (Tian et al. 2019). Given that plants need at least 17 chemical elements in stoichiometrically conserved proportions to carry out their metabolic functions adequately (Sterner and Elser 2003; Kaspari and Powers 2016), understanding the simultaneous response of the chemical elements that are indispensable for plants is of critical importance in the current context of global environmental change, including that of increased N deposition (Finzi et al. 2011; Schlesinger et al. 2011).

Nitrogen is an essential nutrient for plants as it is a key constituent of proteins, nucleic acids, and other biomolecules such as photosynthetic pigments (Chapin et al. 2002). However, increased N deposition has now doubled the amount of reactive N entering into terrestrial ecosystems, with dramatic consequences in terms of biodiversity loss and community reconfigurations (Fowler et al. 2013). Increased N deposition can also trigger changes in plant production and alter the physicochemical properties of soils, including soil acidification, loss of buffering capacity, solubilization of toxic metals like aluminum (Al), and leaching of essential cations like potassium (K) and magnesium (Mg) (Gruber and Galloway 2008; Fowler et al. 2013; Schlesinger and Bernhardt 2013). Moreover, changes in soil nutrient bioavailability via fertilization can directly impact the ability of plants and microbes to maintain their chemical composition within the stoichiometric proportions that are optimal for their growth (Bowman et al. 2008; Tian et al. 2019). These soil nutrient alterations can then result in a cascade of effects, potentially altering the spatial co-variation of multiple soil elements (Finzi et al. 2011). However, despite the considerable amount of information available regarding the impacts of increased N deposition on plant productivity and soil nutrient pools and bioavailability, we still do not know whether and how increased N deposition alters the coupling of soil elemental cycles. Moreover, we know even less about the role that the atomic properties of atoms such as those previously mentioned (i.e., mass, valence orbital group, or electronegativity) may play in modulating the impacts of increased N deposition on soil element coupling. Yet, investigating this will be a critical step towards developing a more predictive framework for understanding the impacts of N deposition on the coupled biogeochemistry of our planet.

Here, we evaluated the effect of four years of experimental N deposition (0, 10, 20, and 50 kg N·ha−1·yr−1 over the background, 6–7 kg N·ha−1·yr−1) on the coupling of soil elements (sensu Ochoa-Hueso et al. 2021a) in a semiarid Mediterranean shrubland. For this, we measured the bioavailability of primary and secondary soil macronutrients (N, K, sodium [Na], calcium [Ca], and Mg), micronutrients (iron [Fe] manganese [Mn], copper [Cu], zinc [Zn]), and non-nutrient elements (Al). In addition, we investigated the role of the atomic properties of elements in driving the coupling of elements as well as the potential of these properties to modulate the response of different elemental cycles to increased N deposition (Ochoa-Hueso et al. 2021c). We calculated the degree of soil elemental coupling based on the mean of Spearman rank correlation coefficients, in absolute value, of all possible pairwise interactions among elemental cycles (Ochoa-Hueso et al. 2021c). Based on previous investigations, we predicted (H1) a decrease in soil elemental coupling under high N loads (Ochoa-Hueso and Manrique 2011), with (H2) such reductions being more evident during the growing season because of a non-coordinated release and uptake of soil essential nutrients by organisms (Ochoa-Hueso et al. 2019). In addition, given the potential of the atomic properties of chemical elements to determine the reactivity and mobility of elements (Kabata-Pendias and Pendias 2000; Schlesinger et al. 2011), we also predicted that (H3a) these properties will be tightly related to the coupling of soil elements (Ochoa-Hueso et al. 2021a), and that (H3b) the effects of N on soil elemental coupling will be modulated by these properties.

Material and methods

Study site and experimental design

This study was conducted in a semiarid Mediterranean shrubland at the Nature Reserve “El Regajal-Mar de Ontígola” in Aranjuez, central Spain (40°00′ N, 3°36′ W; 580 m above the sea level). The site has a semiarid Mediterranean climate, with a mean annual temperature and precipitation of ~ 15 °C and ~ 425 mm, respectively. Rainfall occurs mostly between October and May, and summers are dry. Soils are quaternary limestones, well-buffered, and rich in bases. Soil texture ranges from loamy to sandy-loam and nitrate (NO3−) is the dominant form of inorganic N in soil (Ochoa-Hueso et al. 2013). Soil pH is alkaline (~ 8), but river siliceous stone deposits from the Tajo valley contribute to slightly reduce the soil pH in some parts of the area (Ochoa-Hueso et al. 2013).

The experiment was established in 2007 and consisted of six independent blocks on open areas dominated by shrubs, mainly rosemary (Rosmarinus officinalis) and kermes oaks (Quercus coccifera). At each block, four plots of 2.5 × 2.5 m were set up. Each plot received fertilization loads of 0, 10, 20, and 50 kg N ha−1 yr−1 in the form of ammonium nitrate (NH4NO3). During the study period, two liters of fertilizer solution were added to plots once a month from September to June, when most of the rain occurs, resulting in those plots receiving N doses of 0, 0.019, 0.037 and 0.093 M, respectively (Ochoa-Hueso and Manrique 2011). In order to simulate the increase of N bioavailability after first rains, a triple dose of N was applied at the end of September (Ochoa-Hueso et al. 2013; Ochoa-Hueso 2016). These experimental treatments were designed to mimic the future N deposition scenarios for the Mediterranean Basin by the year 2050 (Phoenix et al. 2006; Fowler et al. 2013) and are within the range of N deposition measurements reported in other Mediterranean areas (Fenn et al. 2003; Ochoa-Hueso et al. 2013, 2014a). Plot selection was based on their representativeness of the study area and included shrubs (23% cover on average), well-developed biocrusts (14% cover on average), and a diverse community of annual plants during the spring growing season (Ochoa-Hueso and Stevens 2015).

Soil sampling and chemical analysis of soil elements

Soil sampling was conducted seasonally in order to capture the temporal variability of the N effects. Soil sampling was carried out on a seasonal basis from autumn 2008 to summer 2009, and also in spring 2010, spring 2011, and autumn 2011. During each sampling event, eight soil samples from the mineral horizon were collected randomly from each plot (0–4 cm depth and 2 cm wide; Bowker et al. 2011) and then combined into a composite sample. Once in the laboratory, samples were dried at room temperature and sieved through a 2-mm mesh prior to measuring soil elemental bioavailability.

Alkali metals (Na, K, Mg and Ca) were extracted by ion exchange procedures, as detailed in Tan (2005). Briefly, 2.5 g of soil was extracted in 25 ml 1 M ammonium acetate (CH3COONH4) adjusted to pH 7. The soil solution was orbitally shaken for 30 min and then filtered using quantitative ashless filter paper. Samples were analyzed by inductively coupled plasma optical emission spectrometry (ICP–OES; Perkin Elmer 4300 DV). Transition metals (Fe, Mn, Zn, Cu and Al) were extracted using the Lakanen solution (0.5 M CH3COONH4, 0.5 M CH3COOH and 0.02 M EDTA) at pH 4.5 (Ochoa-Hueso et al. 2014b). We shook 2.5 g of soil in 25 ml for 60 min and the extractant was then filtered and analyzed by ICP–OES (Perkin Elmer 4300 DV). Soil NO3−–N and NH4+–N were extracted by orbitally shaking 10 g of soil in 50 ml of 0.01 M CaCl2 for 30 min at 160 rpm and then filtered (Ochoa-Hueso et al. 2011, 2013, 2014b). The filtered extracts were colorimetrically assessed in a Kjeltec-auto 1030 analyzer (Tecator, Sweden) for autumn and winter 2008 soil samples. Samples from spring and summer 2009 were analyzed by inductively coupled plasma–optical emission spectroscopy (ICP–OES; Perkin Elmer 4300 DV). Samples from spring 2010–autumn 2011 were analyzed in an auto analyzer (AutoAnalyzer 3, High Resolution Digital Colorimeter, SEAL). Soil mineral N was calculated as the sum of NO3−–N and NH4+–N.

Statistical analyses

For all numerical calculations and statistical analyses, we used R version 4.3.0 (R Development Core Team 2023). Soil elemental coupling (Eq. 1) for each N treatment level and sampling event was calculated as the mean of the absolute values of pairwise Spearman correlations (n = 6) among all elements (Ochoa-Hueso et al. 2021a).

$$\mathrm{Elemental\;coupling }= \frac{{\sum }_{i=1}^{n}|\rho i|}{n}$$
(1)

where n = number of unique pairwise Spearman correlations, ρi = Spearman correlation coefficient between two elements.

We also calculated the coupling for each soil element separately by following the same procedure, but only accounting for correlations involving a given element. In order to obtain an integrated measure of nutrient bioavailability at the study site (Ochoa-Hueso et al. 2019), we calculated the average of z-scores (i.e., standardized and centered coefficients) of all elements, as routinely undertaken in multifunctionality studies (Maestre et al. 2012).

The effects of N deposition on overall and single nutrient availability were tested via mixed-effects models using the lme function of the nlme package in R (Pinheiro et al. 2017). In this analysis, N treatment (four levels), sampling season, and their interactions were included as fixed-effects terms, while plot and block were included as random-effects term to account for the block design and repeated sampling over seasons (Zuur et al. 2009). The effects of N deposition on soil elemental coupling, either for all elements simultaneously or the coupling of single elements, were also tested with linear mixed-effects, where N treatment and sampling season were included as fixed-effects terms, and pairwise interactions nested within N treatments were used as random effects (Zuur et al. 2009). All analyses were checked for residual normality and variance homoscedasticity, and log-transformations were applied when necessary to meet assumptions. When interactions were significant at P < 0.05, post-hoc group differences were examined with Tukey’s Honest Significant Difference test for multiple comparisons using the lsmeans function within the lsmeans package in R (Lenth 2016).

We also used a null modeling approach to explore the extent to which coupling departed from purely random associations among elemental cycles (Ochoa-Hueso et al. 2019). For this, we ran 999 null model randomizations within each N treatment level and season and, for each randomization, we calculated coupling as previously explained. Based on this approach, elemental cycles can be found in three coupling states (Ochoa-Hueso et al. 2021a): (i) coupled, when measured coupling observations are above the 97.5% quantile of random observations (i.e., P < 0.05 for a two-tailed test); (ii) decoupled, when observations fall within the 2.5–97.5% quantile envelope; and (iii) anticoupled, when measured coupling observations are below the 2.5% quantile of random observations.

To evaluate the role of atomic properties of elements (mass, electronegativity, density, and valence orbitals [i.e., whether they belong to s, p, d, or f blocks]) on soil elemental coupling, we carried out linear mixed effects models as described above. In this analysis, N treatment, atomic properties, and their interaction were included as fixed terms, while random terms included the interaction between soil elements and N treatment.

Results

Nutrient bioavailability

Soil nutrient bioavailability varied seasonally and over the years, with very different peaks in bioavailability depending on the nutrient studied (Fig. 1; Table S1). With the exception of N, which increased with simulated N deposition treatments, none of the nutrients, nor the overall nutrient bioavailability index, responded to N addition over the study period (Fig. 1; Table S1).

Fig. 1
figure 1

Simulated N deposition effects on nutrient bioavailability under the four simulated N deposition treatments (x axis): 0 (red), 10 (green), 20 (blue), 50 (purple) kg N·ha−1·yr−1

Full size image

Soil elemental coupling

Soil elemental coupling increased in response to N addition, but this effect was only marginally significant (F3,176 = 2.248, P = 0.085; Fig. 2). Under control conditions, elemental cycles were found in a decoupled state (i.e., coupling did not deviate from null model expectations), while the three N addition treatments were more coupled than expected by chance. We used undirected correlation networks to visualize associations among elements for each N level (Fig. 3). We found a greater number of significant associations (P < 0.05) among elemental cycles with the addition of 20 kg N ha−1 yr−1, and all these associations were positive. Moreover, positive associations between nutrients were common in response to N, while negative associations only appeared under control conditions. Despite this general trend, coupling responses to N also varied across time (F18,1056 = 3.523, P =  < 0.001; Fig. S1), although coupling was maintained across seasons (F6,1056 = 1.272, P = 0.268). Moreover, when analyzed as single elements, the coupling of three elements, Mg, K and to a lesser extent Zn, was significantly affected by N (Table 1). We identified three groups of elements depending on the shape of their response to N: (i) a first group of essential macronutrients like Mg, K, and N that became increasingly coupled in response to N; (ii) a second group of metallic cations like Fe, Mn, Zn, Cu and Al that showed a critical load type of response to simulated N deposition, increasing with the addition of 10 and 20 kg N ha−1 yr−1, and then decreasing with 50 kg N ha−1 yr−1; and (iii) a third group of elements like Ca and Na that did not show a clear response pattern (Fig. 2). Similar to the overall coupling pattern, the coupling of some elemental cycles, namely Mg, Ca, Na, Cu, and Al, also responded to N in a season-dependent manner (Table 1 and Figure S1).

Fig. 2
figure 2

Effects of simulated N deposition treatments (x axis, 0, 10, 20, 50 kg N·ha−1·yr−1) on global coupling (a) and on elemental coupling (b–k) over time. Open circles indicate P < 0.05 and solid circles indicate P > 0.05, as compared to null models based on 1000 random permutations of our own data set (grey shaded area). Error bars are the standard error of the mean. Colors show the orbital group to which the elements belong (blue = d-block; green = p-block; blue = s-block)

Full size image
Fig. 3
figure 3

Network graphs showing the associations among elements for each N level. Dark grey lines show positive associations, and light grey ones show negative associations. The width of connecting lines is proportional to the strength of the relationship (Spearman rho, in absolute values). Colors show the orbital group to which the elements belong (blue = d-block; green = p-block; blue = s-block)

Full size image
Table 1 Simulated N deposition effects on the coupling of single elements
Full size table

Predicting soil element coupling from atomic properties

Regarding individual elements, Ca was the least coupled element, while Fe was the most coupled one (Fig. 4a), with variations in coupling between these two elements representing a shift in the membership of elements to blocks of the periodic table (Fig. 4b). In this sense, elements belonging to the s-block (i.e., those with valence electrons in the most internal s orbital) were the least coupled, while elements from the d-block (i.e., those with valence electrons in the more external d orbital) were the most coupled ones (Fig. 4b). Elements from the p-block showed intermediate coupling values. Soil elemental coupling was positively associated with atomic density, mass, and electronegativity (Figs. 4c–e), but negatively associated with ionic radius (Fig. 4f). These relationships were particularly evident when excluding N and Ca from the analyses. In contrast, we did not find any significant interaction between N addition and the atomic properties of elements for soil elemental coupling (Table 2), which means that, in contrast to our predictions, the effects of N on soil elemental coupling were not modulated by the atomic properties of those elements.

Fig. 4
figure 4

Elemental coupling depending on atomic properties. a Ranked coupling. Error bars are the standard error of the mean. b Coupling of elements based on their orbital group. Relationships between biogeochemicaelemental coupling and atom weight (c), electronegativity (d), density (e) and ionic radii (f). Nitrogen (N) and calcium (Ca) are represented in yellow triangles. The solid line is the fitted regression line for all elements excluding N and Ca and the shaded area represents the 95% confidence interval. Associated statistics are in Table 2

Full size image
Table 2 Summary of regression analyses testing the effects of atomic properties and N treatment on soil biogeochemical coupling, with and without N and Ca
Full size table

Discussion

Previous studies have documented that elemental cycles are strongly interrelated in natural ecosystems (Gao et al. 2013; Soussana and Lemaire 2014; Ochoa-Hueso et al. 2019), which means that their pools and fluxes tend to be more spatially and temporally associated than what is expected by chance (Schlesinger et al. 1996; Ochoa-Hueso et al. 2021c). However, given the expected increases in N deposition worldwide (Galloway et al. 2003; Josep et al. 2012), the main global elemental cycles may have become out of sync (Vourlitis et al. 2009; Ochoa-Hueso et al. 2021c). Here, we used the concept of coupling developed by Ochoa-Hueso et al. (2021a) to assess the level of spatial association among multiple soil elements under both manipulated (i.e., increased N) or unmanipulated experimental conditions in a semiarid shrubland from central Spain.

Contrary to our expectations, we found that soil elemental cycles became more coupled in response to four years of continuous N addition at 10, 20, 50 kg N ha−1 yr−1 as compared to control plots. Despite this general trend, N deposition impacts varied widely among elements and N doses. Interestingly, the coupling of several elements either increased in response to N (K and Mg) or reached a threshold at intermediate levels of N addition (Zn). Moreover, the effects of N addition on coupling were detected despite the lack of effects of N addition on the bioavailability of most soil elements evaluated, with the exception of N, which increased due to obvious reasons. This supports the usefulness of evaluating changes in associations among ecosystem components, in this case among elemental cycles, to detect subtle effects of global change drivers that may later on cascade to more easily observable, and perhaps more difficult to revert, ecosystem reconfigurations (Fridley et al. 2011; Ochoa-Hueso 2016; Ochoa-Hueso et al. 2021a).

The observed increases in soil elemental coupling in N fertilized plots and the greater number of positive and significant associations among elemental cycles in response to N suggest a more coordinated response between the processes that control the release, immobilization, and removal of multiple nutrients from the soil solution in response to N deposition. Noteworthy, simulated N deposition resulted in a concomitant reduction of cover of the dominant shrubs (Cabal et al. 2017; Benvenutto-Vargas and Ochoa-Hueso 2020) along with lower production and diversity of the annual plant community and soil microbial communities (Ochoa-Hueso et al. 2013). These effects were mostly attributed to the toxic effects of excessive NH4+ and nutritional imbalances (Ochoa-Hueso et al. 2013; Cabal et al. 2017; Benvenutto-Vargas and Ochoa-Hueso 2020). As such, a reduction of the potential selective uptake of nutrients by the vegetation might have contributed to reduce the spatial coherence of nutrient cycles under N fertilization. Thus, the increase of elemental coupling in N fertilized plots may be linked to a lower control of the biotic components on soil nutrient cycling driven by the presence of a low productivity and more homogeneous plant community (Lambers et al. 2008; Schlesinger et al. 2011). Supporting this statement, Ochoa-Hueso et al. (2019) showed that N-driven disruptions of biotic interactions (i.e., those between microbes, vegetation, and edaphic fauna) were greater than those involving biotic components and soil properties (i.e., biotic-abiotic interactions) at this study site. The simultaneous disruption of biotic interactions and greater elemental coupling also suggest that, once biotic controls are relaxed, soil nutrient cycles and their coupling might become more controlled by abiotic processes (e.g., adsorption-desorption onto organo-mineral surfaces) (Delgado-Baquerizo et al. 2013). Recent research shows that the physical and chemical properties of soil affect the degree of biogeochemical coupling, which is positively correlated with the amount of particulate organic matter and pH levels across various ecosystems and soil conditions. However, during our experiment, there were only minor changes in soil pH and organic C in response to N addition (Ochoa-Hueso et al. 2013), suggesting a minor role of N-driven changes in these soil properties in driving increased soil elemental coupling in response to N additions. Alternatively, the observed increases in soil elemental coupling under N deposition are less likely to be explained by N deposition relieving natural N limitation at our study site given the previously reported deleterious effects caused by N deposition on vegetation and microbial communities.

In order to understand the temporal dynamics of soil elemental coupling, we analyzed how the coupling of soil elements varied across seasons. We found that the effects of treatments on soil element coupling evolved over time, varying seasonally and across years, which is consistent with the high seasonality in the bioavailability of chemical elements at the experimental site (Ochoa-Hueso et al. 2013). The coupling of single elements also varied widely across time and in response to N, but without a common response pattern. This lack of a common seasonal pattern can be attributed to the fact that not all elements have the same mobility and reactivity in the soil, are not required by plants in the same proportions across seasons, and do not behave equally in response to increased N deposition (Schlesinger et al. 2011). For instance, N deposition can increase the loss of mobile cations such as K+ and Mg2+ (Ochoa-Hueso et al. 2013; Sardans and Peñuelas 2015), thus reducing their availability in soil, particularly during rainy periods (such as autumn 2011 in our case). Moreover, autumn is usually a nutrient-demanding period for plant growth at our study site (Ochoa-Hueso and Manrique 2011; Ochoa-Hueso et al. 2013), further contributing to plant-driven changes in the bioavailability of secondary macronutrients such as K and Mg (Bowman et al. 2008), thus triggering changes in the relationships among elements. Considering that plant nutrition allows a selective depletion of nutrients from soils, resulting in the spatial decoupling of soil elements (Rumpel and Chabbi 2019), reductions in plant biomass at our study site may explain why specific soil elements tend to become more coupled than others in response to increased N deposition (Ochoa-Hueso et al. 2021b). Calcium, which was highly abundant in our soils due to their calcareous nature and thus unlikely to be limiting for plant growth, was the least coupled element, perhaps indicating that this element was not under biological control. Instead, Ca may be redistributed across the soil due to chemical weathering of limestone, thus resulting in the loss of spatial coherence with the rest of elements. The potential leaching of this ion due to the slight acidification of the soil in this study area following N deposition may be another reason behind this pattern (Stevens et al. 2009; Ochoa-Hueso et al. 2013), implying that the decoupling of Ca may not be under biological control.

Given the differences in the capacity of soil elements to interact with one another naturally, we predicted that the spatial co-occurrence of different soil elements could be influenced by their atomic properties (i.e., mass, electronegativity, valence orbitals, etc.) and also that the effects of N addition on the coupling of individual elements may depend on these properties. In our study, we found partial support for our first prediction that atomic properties are important to determine soil elemental coupling, but not for the second one, indicating that the way in which the coupling of chemical elements responds to N is independent of the atomic properties of those elements. Except for ionic radius, mass, electronegativity, density, and separation of valence orbitals from the atomic nucleus (as represented by s-, p-, and d-blocks), were all significantly and positively correlated with soil elemental coupling, particularly when two elements (N and Ca) were removed from the analyses. This suggests that the ability of elements to react with another, and the potential energy involved in those relationships (valence orbitals, electronegativity), together with the potential of elements to move across space (density, mass), may play an underlying and often unappreciated role in the way in which elements are spatially distributed in topsoils (Kabata-Pendias and Pendias 2000).

Recently, Ochoa‐Hueso et al. (2023) suggested that deviations of specific elements from predictable coupling-atomic property relationships may serve as indicators of important reorganizational processes involving those elements. In our study, the fact that the coupling of Ca and N tended to deviate from the predictable relationships with all atomic properties evaluated may be related to two independent aspects. First, the coupling of N deviated from the predictable relationship with atomic properties because its bioavailability, but not that of other elements, was altered by our experimental additions. In contrast, we attributed the deviation of Ca to the biogenic and evaporitic origin of the bedrock at our study site, altering the typically co-occurring distribution of Ca with other elements such as Mg, K, Fe, and Na, which, together with silicon (not evaluated in this study), are typically found in silicates.

Our results are highly comparable to Ochoa‐Hueso et al. (2023), who found a positive relationship between the coupling of the bioavailability of 13 elements and their atomic mass across 15,000 globally distributed observations. Also similar to our results, this previous study found that one element, in that case P, was more decoupled than what the coupling-atomic mass predicted. The decoupling of P in relation to the rest of elements was interpreted in the context of the generally low bioavailability of P in relation to demand across terrestrial ecosystems (Ochoa‐Hueso et al. 2023). Our results were, however, opposite to previous results found by Ochoa-Hueso et al. (2021c), who found a negative relationship between soil element coupling and atomic mass across 16 globally distributed chronosequences. Two main reasons could account for this difference between our results and Ochoa‐Hueso et al. (2023) on one hand, and Ochoa-Hueso et al. (2021c), on the other hand. (i) First of all, Ochoa-Hueso et al. (2021c) considered total element concentrations, and not extractable forms. Total elements are an integrated representation of soil mineral composition and organic matter content and may represent the result of processes that occur at longer timescales (years to decades to centuries) (Benton 2009; Sommer et al. 2014; Ochoa-Hueso et al. 2021c) In contrast, soil element coupling based on bioavailability, like we did in this study, is a more dynamic process that takes place at much shorter temporal scales (days to weeks to seasons) (Schlesinger and Bernhardt 2013), and is highly affected by various biological processes leading to changes in the release and capture of soil elements (Lambers et al. 2008; Tian et al. 2019). (ii) Secondly, in Ochoa-Hueso et al. (2021c) soil samples were collected across a much wider spatial scale, while data from this study is based on a local experiment. Given that the coupling of elements that occur across spatial scales may be driven by different processes (Delgado-Baquerizo et al. 2013; Soussana and Lemaire 2014), it is not surprising that results are not fully equivalent, or even contrasting. What is highly relevant though is the fact that in both cases atomic mass consistently contributed to explain soil element coupling, implying an underlying fingerprint of atomic properties in the way in which elements are related to one another, even in the presence of element-redistributing living organisms (Kabata-Pendias and Pendias 2000; Ochoa-Hueso et al. 2021c).

Conclusions

In this study, we found that four years of N deposition at realistic levels resulted in an increase in soil elemental coupling in a semi-arid Mediterranean shrubland from central Spain. In general terms, we found three patterns of response to treatments when we analyzed the elements separately: (i) the coupling of N, K and Mg increased as compared to the control plots; (ii) the coupling of metal cations (Fe, Mn, Zn, Cu and Al) showed a critical load type response to simulated N deposition; and (iii) the coupling of Na and Ca showed no clear pattern. While the temporal dynamics of soil elemental coupling varied widely between elements in response to treatments, atomic properties and the electronic configuration of atoms (as represented by s-, p-, and d-blocks) seemed to control the degree of coupling of individual elements, but not their response to simulated N deposition. Indeed, differences in soil elemental coupling in response to N deposition were more likely controlled by biological processes that determine the release, demand, and immobilization of soil elements. Importantly, our study opens up new and exciting questions linked to the importance of considering the way in which soil elements are linked to one another in topsoils (i.e., their coupling) in order to capture subtle responses of terrestrial ecosystems to global change. Considered together, our study also strongly suggests that by integrating information on the biotic, abiotic, and anthropogenic controls on soil elemental coupling together with the role of atomic properties of elements, we may be able to develop a more unified way to understand the elemental reorganization that is currently taking place in the face of widespread global environmental change.