Materials and Methods
Field Experiment
The field experiment was conducted in Duolun County (42°02’N, 116°70’E,
1324 m a.s.l.), a semiarid temperate steppe in Inner Mongolia, northern
China. The multiple-year mean annual precipitation and mean annual
temperature are approximately 385 mm and 2.1°C, respectively. The
topography is characterized by low foothills at elevations of 1150-1800
m. Soil is Calcis-orthic Aridisol (the US Soil Taxonomy classification),
with 62.75 ± 0.04% sand, 20.30 ± 0.01% silt and 16.95 ± 0.01% clay
(Wu et al. 2010). The dominant plant species include Stipa
krylovii , Artimesia frigida , Potentilla acaulis ,Cleistogenes squarrosa , Allium bidentatum andAgropyron cristattum .
To assess the nitrogen (N) and phosphorus (P) impacts on
CH4 flux, a block design experiment with different
combinations of N and P additions was established in early 2013 and run
to 2016. A complete random design, with three block replicates, was
adopted to address the high spatial heterogeneity. There were four
experimental treatments (i.e., control (CK), N addition (N), P addition
(P), both N and P additions (N+P)), which were randomly assigned to four
6 m × 6 m plots in each block. Two chambers for greenhouse gas
measurement were set up in each plot. All blocks were separated with a
3-m buffer zone. The N and P were applied twice per month from May to
September during 2013-2016; we sprayed the fertilizer solution to ensure
that the application was evenly distributed in the plots. The dose of N
addition was 100 kg N ha-1 y-1 as
NH4NO3 solution; this dose was selected
because it meets the N required to sustain the local maximum vegetation
productivity (Bai et al. 2010). We monitored the atmospheric N
deposition and estimated it to be 20.4 kg N ha-1y-1 in our experimental site in 2012. The dose of P
was 100 kg P ha-1 y-1, equivalent to
previous nutrient addition experiments in grasslands (Phoenix et
al. 2003). The plots with both N and P additions received the same
amounts of N and P as in the N-only or P-only addition treatments.
Control plots were not fertilized, but rather watered with the same
amount of water as used in the fertilizer solutions; the water used to
dissolve N and P was approximately 800 ml, which is equivalent to 0.8 mm
for a 1 m2 plot, a pre-treatment field experiment
found that this small amount of water addition did not cause any
significant changes in water supply, thus unlikely to have altered the
ecosystem functions.
Over the study period, soil and plant properties were measured
(Supplementary Online Materials). Soil samples were collected to a depth
of 10 cm once per year during 2013-2016. Six soil cores were randomly
taken in each plot and mixed completely. The measured soil variables
included ammonium (NH4+) and nitrate
(NO3-), soil pH, soil organic carbon
(SOC), soil total nitrogen (STN), soil total phosphorus (STP), soil
microbial biomass carbon (MBC) and soil microbial biomass nitrogen
(MBN). It should be noted that no nitrite was detected at our field
site. The measured vegetation variables included aboveground biomass,
plant total carbon (PTC) and total nitrogen (PTN).
Field Measurements of CH4Flux
A static chamber technique was used to measure the CH4flux (Song et al. 2009; Zhang et al. 2017). Stainless
steel permanent bases (50 × 50 × 12 cm) with a 3-cm-deep groove for
water sealing were inserted into soil down to a depth of 12 cm in the
plots approximately one month before the experiment started in the first
year. The chamber base was left in the field for one month before any
flux measurement; this was to let the system evolve into an equilibrium
state, avoiding any potential disturbances to methanogens and
methanotrophs and thus the CH4 flux. An opaque
stainless-steel top chamber (50 cm height) with heat-isolation was
installed over the base, with bottom in the groove. The grooves were
filled with water for sealing the chamber. Two electric fans were
installed inside the chamber to mix the air in the headspace
continuously and thoroughly during the measurements (Fan et al.2007). Sixty-ml gas samples were collected at 10-min intervals for 30
min using sixty-ml syringes with airtight stopcocks. Simultaneously, the
air temperature and pressure in the chamber were measured, and the soil
temperature and moisture were measured at a 5 to 10-cm depth using a
thermometer (Spectrum Technologies, Inc. East - Plainfield, Illinois)
and a portable soil moisture measuring kit ML2x (ThetaKit, Delta-T
Devices, Cambridge, UK), respectively. The CH4concentrations of the gas samples were analyzed using a gas
chromatograph (Agilent 7890A, Agilent Technologies Limited Co., USA)
equipped with flame ionization detector (FID) using 60–80 mesh 13 XMS
column (2 mm inner diameter and 2 m long), with an oven temperature of
55°C. Pure nitrogen was used as carrier gas at a flow speed of 30 ml
min-1, and the CH4 flux was determined
from the linear slope of the mixing ratio changes in four samples taken
at 0, 10, 20 and 30 min after chamber closure. The CH4flux was measured weekly from May to September in 2013, 2014, 2015 and
2016, respectively.
The CH4 fluxes were linearly interpolated and
sequentially cumulated to estimate the total flux over growing seasons
(Zhang et al. 2017). The linear interpolation was adopted due to
two reasons: (1) it has been widely used in previous studies and has
been suggested to be effective for seasonal interpolation of ecological
variables (Song et al. 2009; Nikiema et al. 2011); (2) the
auxiliary data to assist annually interpolation was lacking. Throughout
this analysis, positive fluxes represent gas uptake by the grassland
ecosystems.
Meta-analysis
A meta-analysis was conducted to investigate the N and P impacts on
CH4 across global grasslands. We collected publications
by searching for “nitrogen and phosphorus”, “methane”,
“grassland”, and “upland” in Google Scholar in July 2016, and later
updated in September 2018. The search returned 6230 publications. A few
criteria were then used to screen these publications for our purpose.
The criteria applied to determine whether or not to use the data were:
(1) it must be manipulation experiments with either external N and/or P
addition; (2) the field measurements must cover at least one full
growing season, which makes the estimation of annual budget more
reliable; (3) the studies report clear information for the field site
that is useful when extracting edaphic and meteorological data from
global datasets. For sites without clear latitude and longitude, we
Googled and found their geographic coordinates with country and site
names. Finally, 35 papers were selected for data extraction. When the
data were presented in figures, we extracted mean values and standard
errors using GraphClick (http://www.arizona-software.ch/graphclick/).
For studies with measurements from different N and P addition levels in
one paper, they were extracted and treated as independent data points.
Finally, we compiled a comprehensive dataset of 382 data points for
meta-analysis (Figs. S4 and S5 ). The dataset covers the major
grassland types across the globe and all experiments were carried out
between 1980-2017. For all field experiments, the atmospheric
depositions of N and P from the global dataset were treated as
background rate and were added to the reported N and/or P addition
rates. Since a majority of the data points from Yu et al. did not
contain N or P deposition (Yu et al. 2017), we used the extracted
contemporary nutrient deposition (Mahowald et al. 2008), based on
latitudes and longitudes, as nutrient inputs. Every CH4flux rate corresponds to one site with auxiliary information including
latitude, longitude, and factors such as N deposition, P deposition,
soil temperature (ST), soil moisture (SM), soil pH, soil organic carbon
(SOC), bulk density (BD), and clay content (CL). Across the studies in
our data set, the N and P treatments fell within the ranges of
0~200 kg N ha-1 and
0~50 kg P ha-1, respectively. The
non-growing season CH4 flux is normally not available
for most field experiments, therefore, the annual rate of
CH4 uptake (kg C-CH4ha-1 yr-1) thereafter was calculated
by the ratio of the growing season to non-growing season
CH4 uptake as reported by a few studies (Li et
al. 2012; Yue et al. 2016).
The soil properties for each point data site were extracted from global
datasets according to their latitude and longitude. Specifically, the
soil pH, SOC, BD and clay were retrieved from the Re-gridded Harmonized
World Soil Database v1.2 in the Oak Ridge National Laboratory
Distributed Active Archive Center for Biogeochemical Dynamics (available
online: https://daac.ornl.gov/SOILS/guides/HWSD.html); soil
temperature and moisture for the top 10 cm were from the NCEP/DOE
AMIP-II Reanalysis (Reanalysis-2):
http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.derived.surfaceflux.html).
Global Extrapolation
The data for the grassland distribution were derived from the Global
Mosaics of the standard MODIS land cover type data product (MCD12Q1)
from http://glcf.umd.edu/data/lc/index.shtml. The spatial coverage is
-180.0° ~ 180.0° in longitude, and -64.0°
~ 84.0° in latitude. The global land cover data sets are
available as ESRI ASCII Grid format files and was re-projected to be
consistent with soil and climate datasets (Fig. S4 ). The global
simulations were carried out at a spatial resolution of 0.5° × 0.5°.
According to the definition of Food and Agriculture Organization of the
United Nations (FAO), any geographic area dominated by natural
herbaceous plants including grasslands, prairies, steppes, and savannahs
with coverage of at least 10% was designated as grassland (Lathamet al. 2014). The annual CH4 uptake was
determined by summing up all grassland grids (equation 1). For the model
simulation, we set N and P depositions to pre-industrial level to
represent ambient condition, N set to pre-industrial level and P set to
contemporary level represent P addition, P set to pre-industrial and N
set to contemporary represent N addition, both N and P were set to
contemporary level to represent N and P concurrent treatment.
Both climate and soil data were resampled and re-projected to be
consistent through NCL (NCAR Command Language, current version 6.4.0).
The data with higher resolution were uniformly transformed to the lowest
resolution at 0.5° × 0.5°. Based on the calculated CH4uptake rate in every grid and grid area of grassland, we scaled up the
results from this analysis by multiplying them for target
CH4 uptake with the corresponding grid areas currently
summarized:
\(F(CH4)=\sum_{k=0}^{n}{\left(\frac{n}{k}\right){\text{CH}_{4}}_{k}A_{k}}\)(1)
where F(CH4) is the sum of the global grassland
CH4 uptake expressed as Tg C-CH4year-1, and CH4k is the
CH4 uptake rate in kth grid as kg C-CH4ha-1 y-1, and Ak is
the area of the kth grid (Fig. S8 ).
Statistical Analysis and Uncertainty
Analysis
We used a one-way ANOVA analysis (ANOVA; R package) followed by Duncan’s
post-hoc test to examine the N and P impacts on annual
CH4 uptake rate, soil
NH4+ and
NO3-, and plant N content. A stepwise
multivariate analysis was adopted to evaluate the response of the soil
plant and microbial communities to environmental factors as well as
their correlations with each other. A number of key variables, including
N input rate, P input rate, soil temperature, soil pH, soil organic
carbon density, bulk density, and clay content, are kept in the multiple
linear equation to quantify their impact on CH4 flux
(Equation 2).
A structure equation model (SEM) was used to explore how nutrient (N, P
and N + P)-induced changes in CH4 uptake considering
meteorological, plant, edaphic factors, and microbial biomass. The
“lavaan” package
(https://cran.r-project.org/web/packages/lavaan/lavaan.pdf) for R
program (version 3.4) was used for the structural equation modeling.
Data used in the SEM and linear regression analyses were calculated as
the mean of every year during the four years (Figs. S2 and 3 ).
Prior to the SEM procedure, we conducted principal component analysis
(PCA) to remove the redundant variables for meteorological factors, soil
edaphic parameters, plant and microbial parameters. Meanwhile, due to
the strong multicollinearity among other variables, the variance
inflation factor (VIF) was used to quantitatively select the right
factor for SEM model (VIF < 5). One SEM model (Fig. S3) was
used to disentangle the direct and indirect environmental impacts on the
CH4 uptake, and the other (Fig. 3) was used to analyze
the N and P impacts on CH4 uptake. A conceptual model
was created to represent the key linkages and connecters between
different parameters. The fitted conceptual model was further modified
with “mod_ind” (within the R package “lavaan”) to yield a reliable
multivariate causal network. The SEM results were evaluated with the
comparative fit index (CFI), the normed fit index (NFI), and the
chi-square test (χ2). The chi-square value is the
traditional measure for evaluating the overall model fit and accessing
the magnitude of discrepancy between the sample and fitted covariance
matrices; NFI is an incremental fit index that assesses the model by
comparing the chi-square value of the model to the null model that
assumes that all variables are uncorrelated; CFI is a revised form of
the NFI that takes into account sample size; the model fits with NFI
> 0.95 and CFI > 0.95, indicating a good SEM
(Hooper et al. 2008).
The Latin Hypercube Sampling based Monte Carlo method was adopted to
quantify the uncertainties of the global grassland CH4uptake under different treatments (Control, N, P and N+P conditions)
estimated by the empirical models. The LHS approach has been used in our
previous study (Xu 2010); detailed procedure can be found inSupplementary Online Material .