Patterns of nitrogen and phosphorus pools in terrestrial ecosystems in China

Recent increases in atmospheric carbon dioxide (CO₂) and temperature relieve their limitations on terrestrial ecosystem productivity, while nutrient availability constrains the increasing plant photosynthesis more intensively. Nitrogen (N) and phosphorus (P) are critical for plant physiological activities and consequently regulate ecosystem productivity. Here, for the first time, we mapped N and P densities and concentrations of leaves, woody stems, roots, litter, and soil in forest, shrubland, and grassland ecosystems across China based on an intensive investigation at 4868 sites, covering species composition, biomass, and nutrient concentrations of different tissues of living plants, litter, and soil. Forest, shrubland, and grassland ecosystems in China stored 6803.6 Tg N, with 6635.2 Tg N (97.5 %) fixed in soil (to a depth of 1 m) and 27.7 (0.4 %), 57.8 (0.8 %), 71.2 (1 %), and 11.7 Tg N (0.2 %) in leaves, stems, roots, and litter, respectively. The forest, shrubland, and grassland ecosystems in China stored 2806.0 Tg P, with 2786.1 Tg P (99.3 %) fixed in soil (to a depth of 1 m) and 2.7 (0.1 %), 9.4 (0.3 %), 6.7 (0.2 %), and 1.0 Tg P (< 0.1 %) in leaves, stems, roots, and litter, respectively. Our estimation showed that N pools were low in northern China, except in the Changbai Mountains, Mount Tianshan, and Mount Alta, while relatively higher values existed in the eastern Qinghai–Tibetan Plateau and Yunnan. P densities in vegetation were higher towards the southern and north-eastern part of China, while soil P density was higher towards the northern and western part of China. The estimated N and P density and concentration datasets, “Patterns of nitrogen and phosphorus pools in terrestrial ecosystems in China” (https://doi.org/10.5061/dryad.6hdr7sqzx), are available from the Dryad digital repository (Zhang et al., 2021). These patterns of N and P densities could potentially improve existing earth system models and large-scale research on ecosystem nutrients.

1 Introduction

Nitrogen (N) and phosphorus (P) play fundamental roles in plant physiological activities and functioning, such as photosynthesis, resource utilization, and reproductive behaviours (Fernández-Martínez et al., 2019; Lovelock et al., 2004; Raaimakers et al., 1995), ultimately regulating plant growth and carbon (C) sequestration efficiency (Terrer et al., 2019; Sun et al., 2017). Under the background of global warming, the limiting factors for the plant growth, such as carbon dioxide (CO₂) and temperature, are becoming less restrictive for terrestrial ecosystem productivity (Norby et al., 2009; Fatichi et al., 2019), while nutrient availability tends to constrain the increasing plant photosynthesis more intensively (Cleveland et al., 2013; Du et al., 2020). As the key nutrients for plant growth, N and P independently or jointly limit biomass production (Elser et al., 2007; Finzi et al., 2007; Hou et al., 2020). N influences CO₂ assimilation in various ways (Vitousek and Howarth, 1991; Campany et al., 2017). For example, N is a critical element in chlorophyll (Field, 1983), and plant metabolic rates are also regulated by N content (Elser et al., 2010). P is crucial in RNA and DNA construction, and its content is associated with water uptake and transport (Carvajal et al., 1996; Cheeseman and Lovelock, 2004) as well as energy transfer and exchange (Achat et al., 2009). P shortage could lower photosynthetic C-assimilation rates (Lovelock et al., 2006).

In spite of the key importance of N and P for plants, knowledge on the patterns of their storage in terrestrial ecosystems is limited. With additional CO₂ entering the atmosphere, more N could be allocated to plant growth and soil organic matter (SOM) accumulation, which may lead to less available mineral N for plant uptake (Luo et al., 2004). Direct and indirect evidence shows that N limits productivity in temperate and boreal areas (Bonan, 1990; Miller, 1981; Vitousek, 1982). P originates from bedrock weathering and litter decomposition in terrestrial ecosystems, and it experiences long-term biogeochemical processes before it is available to plants (Föllmi, 1996), which consequently makes P a more predominant limiting factor to ecosystem productivity (Reed et al., 2015). Additionally, P decomposition rates are constrained by limited soil labile P storage, especially in tropical forests, where soil P limitation is extreme (Fisher et al., 2012).

Ecosystem models based on Amazon forest free air CO₂ enrichment (FACE) experiments consistently showed that biomass C positively responded to simulated elevated CO₂, but the models incorporating N and P availability showed lower plant growth than those that did not (Wieder et al., 2015). Moreover, a recent study suggested that the inclusion of N and P availability into the earth system models (ESMs) remarkably improved the estimation accuracy of C cycles over previous models (Fleischer et al., 2019). Hence, understanding and predicting the patterns and mechanisms of global C dynamics require good characterization of N and P conditions.

N and P pools in ecosystems consist of several components that cast different influences on ecosystem C storages and fluxes. For example, N and P in plants directly affect C sequestration (Thomas et al., 2010), but their activities differ among organs (Elser et al., 2003; Parks et al., 2000); the soil pools are the source of plant nutrition, and the litter pools act as a transit link that returns nutrients from plants to soil (McGrath et al., 2000). Thus, an accurate estimation of ecosystem N and P pools involves calculating specific nutrient densities in all these components.

Terrestrial ecosystems in China play a considerable part in the continental and global C cycles. Satellite data verified that China contributed to one-fourth of the global net increase in leaf area from 2000 to 2017 (Chen et al., 2019). The total C pool in terrestrial ecosystems in China is 79.2 Pg C, and this number is still growing because of the nationwide ecological restoration and construction, which accounted for 56 % of the total C sequestration in the restoration area in China from 2001 to 2010 (Lu et al., 2018). N and/or P limitations are ubiquitous in natural ecosystems in China (Augusto et al., 2017; Du et al., 2020; Elser et al., 2007; LeBauer and Treseder, 2008; Hou et al., 2020). Understanding the distribution and allocation of N and P in ecosystems is of great significance for a precise projection of the C cycle in China. Although there are a few studies on the spatial patterns of soil nutrient storages in China (Shangguan et al., 2013; Xu et al., 2020; Yang et al., 2007; Zhang et al., 2005), a thorough study on the distribution of N and P pools of all the ecosystems is still lacking as vegetation (living or dead biomass) composes the most active part of the nutrient stocks.

To fill this knowledge gap, here we identified N and P density patterns in China based on an intensive field investigation, covering all components of the entire ecosystem, including different plant organs, litter, and soil. The present study aims to provide high-resolution maps of nutrient densities in different ecosystem components and to answer the following questions.

• How much N and P is stored in different components, i.e. leaf, stem, root, litter, and soil, of terrestrial ecosystems in China?

• How are N and P pools in different components spatially distributed in China?

2 Material and methods

2.1 Field sampling and nutrient density calculation

Forest, shrublands, and grasslands constitute major vegetation type groups in China. Focusing primarily on these three groups, a nationwide, methodologically consistent field investigation was conducted in June and September 2011–2015.

In total, 4868 sites, including 3022 forest, 1123 shrubland, and 723 grassland sites, were investigated (Fig. S1a in the Supplement). At each site, one 20×50 m² plot was set for forests, three replicated 5×5 m² plots were set for shrublands, and ten 1× 1 m² plots were established for grasslands. Species composition and abundance were investigated in plots. Height (for trees, shrubs, and herbs), diameter at breast height (DBH; at height 130 cm) (for trees), basal diameter (for shrubs), and crown width (for shrubs and herbs) were measured for all plant individuals in the plots (X. Tang et al., 2018).

Leaves, stems (woody stems), and roots (without distinguishing coarse and fine roots) were sampled for the top five dominant tree and shrub species, and above- and belowground parts were sampled for dominant herb species. Soil was sampled to a depth of 1 m or to bedrock at depths of 0–10, 10–20, 20–30, 30–50, and 50–100 cm, with at least five replications per site to measure nutrient concentrations and bulk density after removing roots and gravels. Litter was sampled in at least three 1×1 m² quadrats per site (for detailed survey protocol, see X. Tang et al., 2018).

All samples were transported to the laboratory, dried, and measured. N concentrations of all samples were measured by a C–N analyser (PE-2400 II; Perkin-Elmer, Boston, USA), while P concentrations were measured using the molybdate–ascorbic acid method after H₂SO₄–H₂O₂ digestion (Jones Jr., 2001). For the three organs, the community-level N or P density was the cumulative sum of the products of the corresponding biomass density (i.e. biomass per area, Mg ha⁻¹) and community-level concentrations for each co-occurring species, as shown by Eq. (1).

$$\begin{array}{}\text{(1)}& \mathrm{N}\left(\mathrm{P}\right)={\sum }_{i=\mathrm{0}}^n{B}_i\times {\mathit{\theta }}_i\end{array}$$

N (P) represents the community-level N or P density (Mg ha⁻¹); n is the total number of plant species at one site; B_i is the biomass density of a specific organ of the ith plant species in that site, where the plant organ biomass was estimated by allometric equations or harvesting; θ_i represents the N or P concentration (g kg⁻¹) of the same organ of the ith plant species at that site. Allometric equation methods were adapted to trees and some shrubs (tree-like shrubs and xeric shrubs) for biomass estimation, while the biomass of grass-like shrubs and herbs was obtained by direct harvesting. Litter N or P density was litter biomass density (by harvesting) multiplied by litter N or P concentration of each sampling site. The soil N or P density was calculated to a depth of 1 m. Soil N or P concentration and bulk density were measured at different depths (0–10, 10–20, 20–30, 30–50, and 50–100 cm) to determine the community-level soil N or P density using Eq. (2):

$$\begin{array}{}\text{(2)}& \mathrm{SND}\left(\mathrm{SPD}\right)={\sum }_{i=\mathrm{0}}^n(\mathrm{1}-{\mathit{\delta }}_i)\times {\mathit{\rho }}_i\times C_i\times T_i/\mathrm{10},\end{array}$$

where SND(SPD) is the total N or P density of the soil within the top 1 m (Mg ha⁻¹), n is the total number of soil layers (ranging from one to five) at one site, δ_i is the volume percentage of gravel with a diameter > 2 mm, ρ_i is the bulk density (g cm⁻³), C_i is the soil N or P concentration (g kg⁻¹), and T_i is the depth (cm) of the ith layer. For detailed calculations of species biomass and community-level concentrations at each site, please refer to previous studies (X. Tang et al. 2018; Z. Tang et al. 2018).

2.2 Climatic and vegetation data

The daily meteorological observation data from 2400 meteorological stations across China were averaged over the 2011–2015 period to generate a spatial interpolation dataset of mean annual temperature (MAT) and precipitation (MAP) using a smooth spline function (McVicar et al., 2007) with a spatial resolution of 1 km. MAT and MAP of each site were extracted from this dataset.

Elevation was extracted from GTOPO30 with a spatial resolution of 30 arc seconds (https://www.usgs.gov/centers/eros/science/, last access date: 29 October 2021). The mean enhanced vegetation index (EVI) from June to September during the 2011–2015 period was calculated based on MOD13A3 data with a resolution of 1 km (https://lpdaac.usgs.gov/products/mod13a3v006/ last access date: 29 October 2021).

The ranges of these variables of our field sites (EVI: 0.03–0.7; elevation: −137–5797 m; MAP: 19.8–2316.3 mm; MAT: −5.2–26.0 ^∘C) could generally cover the ranges of corresponding variables in the focused vegetation types across China (99 % ranges of EVI: 0.03–0.6; of elevation: 24–5628 m; of MAP: 50.6–2956.5 mm; of MAT:−6.6–22.8 ^∘C).

Based on the level II vegetation classification of ChinaCover (Land Cover Atlas of the People's Republic of China Editorial Board, 2017), we classified the vegetation type groups into the following 13 vegetation types: five forest types, i.e. evergreen broadleaf forests, deciduous broadleaf forests, evergreen needle-leaf forests, deciduous needle-leaf forests, and broadleaf and needle-leaf mixed forests; four shrubland types, i.e. evergreen broadleaf shrublands, deciduous broadleaf shrublands, evergreen needle-leaf shrublands, and sparse shrublands; and four grassland types, i.e. meadows, steppes, tussocks, and sparse grasslands.

2.3 Prediction the nationwide nutrient pools and distribution patterns

We used random forests to predict the nutrient densities and concentrations across China. The predictors included MAT, MAP, longitude, latitude, elevation, EVI, and vegetation types (as dummy variables). We established one random forest model for N or P in each component (three plant organs, litter, and five soil layers), respectively. In each model, six variables were randomly sampled at each split, and 500 trees were grown. Larger values of these parameters did not increase validation R² obviously. Model prediction was repeated 100 times to obtain the average results. When modelling the nutrient densities in woody stems, we excluded the four grassland types. All densities were log-transformed based on e, and explanatory variables were transformed using the following equation to ensure they were in the same range before modelling.

$$\begin{array}{}\text{(3)}& x_i{}^{\prime }={\displaystyle \frac{x_i-\min \left(x\right)}{max\left(x\right)-\min \left(x\right)}},\end{array}$$

where x_i means the ith value of the environmental variables x, and max(x) and min(x) represent the maximum and minimum values of x, respectively. We estimated the relative importance of predictors using the increase in node purity for the splitting variable, which was measured by the reduction in residual sum of squares. The same procedures were repeated for the prediction of N and P concentrations in different components across China. The spatial pattern of the N:P ratio was calculated from the predicted N and P density datasets of the corresponding component.

The vegetation N or P density was the sum of all plant organs, the soil N or P density was the sum of all soil layers, and the ecosystem N or P density was the sum of all components. The soil depth data across China were obtained from Shangguan et al. (2017). The N and P pools in 13 vegetation types were estimated, respectively. The N and P pools were calculated from the predicted nationwide densities. The predicted N and P densities were in 1 km spatial resolution, so the nutrient stock is the density multiplied by the grid area (1 km²) for each grid. The nutrient pools of a given vegetation type equal the sum of stocks of the grids belonging to that type.

2.4 Model validation and uncertainty

To evaluate the model performance, we calculated the linear relationship between the observed validation data (10 % of the dataset by random sampling) and predicted data that were estimated based on training data (90 % of the dataset by random sampling) 100 times with the models for every component. We then calculated means of validation R², slopes, and intercepts of the 100 relationships. We also calculated the standard deviations (SDs) of the 100 predictions of each component in each map grid to show the uncertainty in the models.

All statistical analyses were performed using R 3.6.1 (R Core Team, 2019); random forests were built using the randomForest package (Liaw and Wiener, 2002).

3 Results

3.1 Allocation of nutrients among ecosystem components

The mean N and P densities varied among forest, shrubland, and grassland sites and among different tissues (Figs. 1 and 2) according to the measured data. On average, leaves and woody stems in forests stored more N than those in shrublands (0.1±0.1 (mean ± SD) Mg N ha⁻¹ vs. $\mathrm{4.2}\pm \mathrm{10}\times$ 10⁻² Mg N ha⁻¹ for leaves and 0.3 ± 0.6 Mg N ha⁻¹ vs. $\mathrm{5.1}\pm \mathrm{20}\times {\mathrm{10}}^{-\mathrm{2}}$ Mg N ha⁻¹ for woody stems). Similarly, P densities were higher in leaves and woody stems for forests than for shrublands ($\mathrm{1.3}\pm \mathrm{1.5}\times$ 10⁻² Mg P ha⁻¹ vs. $\mathrm{3.1}\pm \mathrm{6.5}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg P ha⁻¹ for leaves and $\mathrm{5.6}\pm \mathrm{11}\times {\mathrm{10}}^{-\mathrm{2}}$ Mg P ha⁻¹ vs. $\mathrm{4.7}\pm \mathrm{19}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg P ha⁻¹ for woody stems). The root N and P densities for forests (0.1±0.2 Mg N ha⁻¹ and $\mathrm{2.1}\pm \mathrm{3.9}\times {\mathrm{10}}^{-\mathrm{2}}$ Mg P ha⁻¹) and grasslands (0.2±0.2 Mg N ha⁻¹ and $\mathrm{1.5}\pm \mathrm{1.6}\times {\mathrm{10}}^{-\mathrm{2}}$ Mg P ha⁻¹) were remarkably higher than for shrublands ($\mathrm{6.6}\pm \mathrm{11}\times {\mathrm{10}}^{-\mathrm{2}}$ Mg N ha⁻¹ and $\mathrm{5.6}\pm \mathrm{8.8}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg P ha⁻¹).

Figure 1 Frequency distributions of N densities in soil, roots, leaves, litter, and woody stems in forests (a–e), shrublands (f–j), and grasslands (k–n) in China.

Figure 2 Frequency distributions of P densities in soil, roots, leaves, litter, and woody stems in forests (a–e), shrublands (f–j), and grasslands (k–n) in China.

The mean litter N densities for forest, shrubland, and grassland sites were $\mathrm{6.1}\pm \mathrm{7.6}\times {\mathrm{10}}^{-\mathrm{2}}$ Mg N ha⁻¹, $\mathrm{3.8}\pm \mathrm{4.6}\times$ 10⁻² Mg N ha⁻¹, and $\mathrm{5.5}\pm \mathrm{9.3}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg N ha⁻¹, respectively. The mean litter P densities for forest, shrubland, and grassland sites were $\mathrm{5.3}\pm \mathrm{9.3}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg P ha⁻¹, $\mathrm{2.5}\pm \mathrm{2.3}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg P ha⁻¹, and $\mathrm{4.1}\pm \mathrm{7.1}\times {\mathrm{10}}^{-\mathrm{4}}$ Mg P ha⁻¹, respectively.

The mean soil N densities for forest, shrubland, and grassland sites were 12.1±10.8 Mg N ha⁻¹, 8.8±7.4 Mg N ha⁻¹, and 9.9±8.9 Mg N ha⁻¹, respectively. The mean soil P densities were 4.9±6.5 Mg P ha⁻¹ for forest sites, 3.9±3.7 Mg P ha⁻¹ for shrubland sites, and 4.4±2.8 Mg P ha⁻¹ for grassland sites.

Belowground vegetation N and P densities were higher than aboveground in grasslands and sparse shrublands. By contrast, this condition was reversed in forests and the other three shrubland types (Fig. 3). Among various forest types, deciduous broadleaf forests and deciduous needle-leaf forests held the highest aboveground N and P densities, respectively. Evergreen needle-leaf forests held the lowest vegetation N density, and evergreen broadleaf forests had the lowest P density. For grassland types, meadows held higher N and P densities in belowground biomass than the other three grassland types, whereas these four grasslands types had relatively approximate nutrient densities in aboveground biomass. Shrublands possessed the lowest vegetation N and P densities among three vegetation groups. Sparse shrublands had the lowest vegetation nutrient densities and soil N density but the highest soil P density among the four shrubland types.

Figure 3 N and P density allocations among leaf, stem, and root (a, b) and between vegetation and soil (c, d) in 13 vegetation types. See Table 1 for abbreviations. The error bar represents standard error. Notice that the y axes above and below zero are disproportionate.

3.2 Mapping of N and P densities in China's terrestrial ecosystems

All models of the N and P densities of different components performed well, with the validation R² ranging from 0.55 to 0.78 for plant organs and litter (Fig. 4) and from 0.47 to 0.62 for soil layers (Fig. 5). As for the concentration models, the validation R² varied from 0.45 to 0.63 for plant organs and litter (Fig. S2) and from 0.53 to 0.70 for soil layers (Fig. S3). Prediction results of 100 repetitions were quite stable, as shown by SDs of the predictions close to zero in all components (Figs. S4 and S5).

Figure 4 Fitting performance of random forest models for nutrient densities of leaves (a, b), woody stems (c, d), roots (e, f), and litter (g, h) of terrestrial ecosystems in China based on 100 replications with the 10 % validation data. Solid lines represent all the fitting lines, and the displayed parameters stand for the average conditions. The dashed line denotes the 1:1 line.

Figure 5 Fitting performance of random forest models for nutrient densities of 0–10 cm (a, b), 10–20 cm (c, d), 20–30 cm (e, f), 30–50 cm (g, h), and 50–100 cm (i, j) soil layers of terrestrial ecosystems in China based on 100 replications with the 10 % validation data. Solid lines represent all the fitting lines, and the displayed parameters stand for the average conditions. The dashed line denotes the 1:1 line.

Leaf N density was high in southern and eastern China but low in northern and western China. It was especially high in the Changbai Mountains, southern Tibet, and the south-eastern coastal areas (Fig. 6a; see Fig. S1b for the topographic map of China), while it was low in northern Xinjiang and northern Inner Mongolia. The woody stem and litter N densities showed similar patterns to that of the leaves (Fig. 6c and g), whereas root N density was high in Mount Tianshan, Mount Alta, the Qinghai–Tibetan Plateau, the north-eastern mountainous area, and eastern Inner Mongolia (Fig. 6e). The vegetation N density was relatively higher in eastern China, the eastern Qinghai–Tibetan Plateau, Mount Tianshan, and Mount Alta (Fig. 7a). The soil and ecosystem N densities were low in northern China except the Changbai Mountains, Mount Tianshan, and Mount Alta but high in the eastern Qinghai–Tibetan Plateau and Yunnan Province (Fig. 7c and e).

Figure 6 Predicted spatial patterns of N and P densities with a resolution of 1 km in leaves (a, b), woody stems (c, d), roots (e, f), and litter (g, h) of terrestrial ecosystems in China.

Figure 7 Predicted spatial patterns of N and P densities with a resolution of 1 km in vegetation (a, b: the sum of leaves, stems, and roots), soil (c, d: the sum of five layers), and ecosystems (e, f: the sum of vegetation, litter, and soil) of terrestrial ecosystems in China.

The P densities in leaves, woody stems, roots, litter, and the whole vegetation showed similar patterns to the N densities in the corresponding components, respectively (Figs. 6b, d, f, and h and Fig. 7b). However, soil and ecosystem P densities were high in western and northern China but low in eastern and southern China (Fig. 7d and f).

The N and P concentrations in plant organs and litter were generally higher in northern and western mountain regions, but larger values of the former often occur in the north-western part of China, while those of the latter often occur in the north-eastern part of China (Fig. S6a–h). The spatial patterns of soil nutrient concentrations at different depths were consistent with those of soil nutrient densities (Fig. S6i–r).

The N:P ratio of plant organs and litter showed similar distribution patterns, with higher values occurring in south-eastern and north-western China and the Qinghai–Tibetan Plateau (Fig. S7a–d). Soil N:P ratio was higher in north-eastern and southern China but lower in north-western China (Fig. S7e).

3.3 N and P pools in China's terrestrial ecosystems

In total, the terrestrial ecosystems in China stored 6803.6 Tg N, with 2634.9, 873.0, and 3295.8 Tg N stored in the forests, shrublands, and grasslands, respectively (Table 1). Vegetation, litter, and soil stored 156.7 (2.3 %), 11.7 (0.2 %), and 6635.2 Tg N (97.5 %), respectively (Table 1).

Table 1 N and P stocks of vegetation, litter, soil, and total ecosystem in forests, shrublands, and grasslands in China.

EBF, evergreen broadleaf forest; DBF, deciduous broadleaf forest; ENF, evergreen needle-leaf forest; DNF, deciduous needle-leaf forest; MF, broadleaf and needle-leaf mixed forest; EBS, evergreen broadleaf shrubland; DBS, deciduous broadleaf shrubland; ENS, evergreen needle-leaf shrubland; SS, sparse shrubland; ME, meadow; ST, steppe; TU, tussock; and SG, sparse grassland.

China's terrestrial ecosystems stored 2806.0 Tg P, with 981.1, 381.8, and 1443.0 Tg P stored in the forest, shrublands, and grasslands, respectively. Vegetation, litter, and soil accounted for 18.8 (0.7 %), 1.0 (<0.1 %), and 2786.1 Tg P (99.3 %), respectively (Table 1).

Meanwhile, N and P stocks among plant organs showed different allocation patterns (Table 2). Compared with the other two vegetation type groups, forests allocated the majority of N and P to the stem pool (55.5 Tg N and 9.2 Tg P), followed by the root pool (23.4 Tg N and 3.3 Tg P) and leaf pool (21.0 Tg N and 2.1 Tg P). However, the root pools in shrublands and grasslands held the most N and P (3.8 Tg N and 0.3 Tg P for shrublands, 44.1 Tg N and 3.1 Tg P for grasslands) (Table 2).

Table 2 N and P stocks of plant organs (leaf, stem, and root) in forests, shrublands, and grasslands in China.

See Table 1 for abbreviations.

Among the four grassland types, steppes had the largest N and P stocks (1370.1 Tg N and 574.6 Tg P), taking the ecosystem as a whole. Deciduous broadleaf shrublands had the largest N and P stocks considering the whole ecosystem (577.6 Tg N and 234.2 Tg P) as well as vegetation (5.5 Tg N and 0.5 Tg P) compared with the other three shrubland types. The largest ecosystem N and P stocks across all five forest types appeared in evergreen needle-leaf forests (984.0 Tg N) and deciduous broadleaf forest (353.8 Tg P) (Table 1).

4 Data availability

The datasets of N and P densities and concentration of different ecosystem components, “Patterns of nitrogen and phosphorus pools in terrestrial ecosystems in China”, are available from the Dryad digital repository along with the geographic coordinates of field sites and layer files of environmental factors for prediction (https://doi.org/10.5061/dryad.6hdr7sqzx) (Zhang et al., 2021).

5 Discussion

5.1 Performance of density models

The accuracy of the density models varied among different components. Models for soil showed relatively poorer accuracy than models for plant organs and litter (Figs. 4 and 5), partly because the soil N and P were largely influenced by geological conditions, soil age, and parent material (Buol and Eswaran, 1999; Doetterl et al., 2015), which were not included in our analysis because of the limited data availability. This can be evidenced by the decreasing validation R² of the models for soil N and P concentrations as well as N densities with soil depths (Figs. 5 and S3). The models preformed best for the stem N and P because woody stems occupied the most biomass in the forest and shrublands (stem biomass and vegetation biomass accounted for 68 % and 48 % of vegetation biomass for forest and shrublands, respectively). Climate variables could affect vegetation growth and biomass accumulation, and the variation in stem biomass could be the most direct reflection (Kirilenko and Sedjo, 2007; Jozsa and Powell, 1987; Poudel et al., 2011).

It is also noteworthy that the values of the validation R² of the density models were higher than those of the concentration models for plant organs and litter (Figs. 4 and S2), which was opposite for soil layers (Figs. 5 and S3). They might reflect that biomass was more constrained by the selected factors in this study than nutrient concentrations in vegetation, while bulk density was less affected than nutrient concentrations in soil.

5.2 Nutrient pools in terrestrial ecosystems in China

Previous research has estimated N and P stocks in soil across China. For example, Shangguan et al. (2013) estimated that the storage of soil total N and P in the upper 1 m of soil in China were 6.6 and 4.5 Pg. Yang et al. (2007) estimated China's average density of soil N at a depth of 1 m, which was 0.84 kg m⁻², and the soil N stock was 7.4 Pg. Zhang et al. (2005) investigated the soil total P pool at a depth of 50 cm in China and concluded that the soil stock was 3.5 Pg, with the total P density of soil 8.3×10² g m⁻³. Our estimation of the soil N pool in China (6.6 Pg) agreed with Shangguan et al. (2013), but the estimated soil P pool (2.8 Pg) was lower than the results of the aforementioned studies. The mean soil N:P ratio in our study (2.5 for the predicted outcome and 2.1 for the sampling sites) was lower than the result of Tian et al. (2010), 5.2, while the spatial patterns in both studies are similar. Other than the research focusing on soil, Xu et al. (2020) estimated China's N storage by calculating the mean N densities of vegetation and soil from different ecoregions and reported that there was 10.43 Pg N in China's ecosystems, 10.14 Pg N in the top 1 m of soil and 0.29 Pg N in vegetation, both higher than our results (6.6 Pg N in soil and 0.16 Pg N in vegetation).

5.3 Potential driving factors of the N and P densities in various components

The distribution and allocation of N and P pools in ecosystems were largely determined by vegetation types and climate. The difference in the spatial patterns of nutrient pools could reflect the spatial variation in local vegetation. For example, it is obvious that the regions covered by forests tend to have higher aboveground nutrient densities than those covered by other types, while the regions covered by sparse shrublands tend to have the lowest nutrient densities (Fig. 3). Despite its decisive influences on vegetation types, climate also greatly impacts the nutrient utilization strategies of vegetation (Kirilenko and Sedjo, 2007; Poudel et al., 2011). For example, in south-eastern China, with higher precipitation and temperature, forests tend to allot more nutrients to organs related to growth, such as leaves that perform photosynthesis and stems that are related to resource transport and light competition (Zhang et al., 2018). These influences were reflected in our models (Figs. S8–S11). In the models of densities for plant organs and litter, vegetation types and climate variables showed higher relative importance. Heat and water are usually limited in the plateau and desert regions in western China, where shrublands and grasslands are dominant vegetation type groups. More nutrients are allocated to root systems by dominant plants in such stressful habitats to acquire resources from soil (Eziz et al., 2017; Kramer-Walter and Laughlin, 2017). Spatial variables, longitude and latitude, also held high importance, especially in the models for soil nutrients. On the one hand, it may result from their tight links with climate conditions. On the other hand, it may imply the influence of spatial correlation on nutrient pools. The effects of elevation and spatial variables were obvious from the prediction maps. There were relatively larger values of soil nutrient densities in the plateau and mountainous area in western China, possibly because of the lower rates of decomposition, mineralization, and nutrient input as well as less leaching loss in high-altitude regions (Bonito et al., 2003; Vincent et al., 2014). However, the distribution patterns of soil nutrient densities in eastern China were generally consistent with the soil substrate age hypothesis that the younger and less leached soil in temperate regions tends to be more N-limited but less P-limited than the elder and more leached soil in tropical and subtropical regions (Reich and Oleksyn, 2004; Vitousek et al., 2010; Walker and Syers, 1976). Additionally, such patterns reflect that the factors not investigated in this study, such as soil age and parent material, could contribute to the patterns of nutrient pools, which should be considered in future research as potential drivers (Augusto et al., 2017; Porder and Chadwick, 2009).

5.4 Potential applications of the data

An atmospheric CO₂ enrichment trend was undoubtable, but how this procedure will develop is still unclear (Fatichi et al., 2019). A number of previous studies proved that global carbon cycle models would produce remarkable bias if the coupled nutrient cycle is overlooked (Fleischer et al., 2019; Hungate et al., 2003; Thornton et al., 2007). However, high-resolution and accurate ecosystem nutrient datasets were unattainable and hard to be modelled without an enormous field investigation basis. This study relied on nationwide field survey data, providing comprehensive N and P density datasets of different ecosystem components. Based on the present dataset, enhancement could be made in various ecosystem research aspects.

First and foremost, the dataset could facilitate the improvement in the prediction of large-scale terrestrial C budget, thereby allowing us to better understand patterns and mechanisms of the C cycle as well as the future trend of climate change (Le Quéré et al., 2018). Numerous projections of future C sequestration overestimated the amount of C fixed by vegetation due to the neglect of nutrient limitation (Houghton et al., 2001; Cramer et al., 2001). Global C cycling models coupled with the nutrient cycle may make more accurate predictions of carbon dynamics. Moreover, our dataset illustrated N and P densities of major ecosystem components and vegetation types at a high spatial resolution for the first time, which could help identify C and nutrient allocation patterns from the tissue level to the community level, especially for vegetation organs which still lack large-scale nutrient datasets.

In addition, large-scale N and P pool spatial patterns could provide the data references for vegetation research using remote sensing (Jetz et al., 2016). Vegetation nutrient densities were important traits but hard to be extracted and detected remotely. With the development of hyperspectral remote sensing technology and the theory of spectral diversity, foliar nutrient traits can be successfully predicted (Skidmore et al., 2010; Wang et al., 2019). However, previous studies still focused on finer-scale patterns and were constrained by the lack of large-scale field datasets for uncertainty assessment (Singh et al., 2015). Our nationwide nutrient dataset offers an opportunity to enlarge the generality of remote sensing models and algorithms at large scales.