Distribution and Attribution of Gross Primary Productivity Increase Over the Mongolian Plateau, 2001-2018

Gross primary productivity (GPP) over the Mongolian Plateau (MP) is a vital component of the global terrestrial carbon cycle. Using the latest MODIS GPP estimates at the best achievable spatial resolution along with several ancillary datasets, we investigated GPP variations in the MP region during 2001–2018 and attributed these changes to land-surface temperature (<inline-formula> <tex-math notation="LaTeX">$T_{\mathrm {s}}$ </tex-math></inline-formula>), total precipitation (<inline-formula> <tex-math notation="LaTeX">$P_{\mathrm {t}}$ </tex-math></inline-formula>), land-cover change (LCC), and atmospheric carbon dioxide (CO<sub>2</sub>) concentrations. The 18-year-averaged annual cumulative GPP in the MP region was 357.02 ± 24.76 gC <inline-formula> <tex-math notation="LaTeX">$\text{m}^{-2}$ </tex-math></inline-formula> yr<sup>−1</sup> during the study period, ranging from 60.51 ± 6.10 gC m<sub>−2</sub> yr<sup>−1</sup> in deserts to 596.41 ± 35.49 gC m<sub>−2</sub> yr<sup>−1</sup> in forests. A linear regression analysis indicated a significant overall increase in GPP, at a rate of 3.91 gC m<sub>−2</sub> yr<sup>−1</sup> (<inline-formula> <tex-math notation="LaTeX">$p < 0.01$ </tex-math></inline-formula>). In comparison, GPP increased at a rate of 0.79 gC m<sub>−2</sub> yr<sup>−1</sup> in deserts (<inline-formula> <tex-math notation="LaTeX">$p < 0.01$ </tex-math></inline-formula>), 4.79 gC m<sub>−2</sub> yr<sup>−1</sup> in forests (<inline-formula> <tex-math notation="LaTeX">$p < 0.01$ </tex-math></inline-formula>), and 5.76 gC m<sub>−2</sub> yr<sup>−1</sup> in grasslands (<inline-formula> <tex-math notation="LaTeX">$p < 0.01$ </tex-math></inline-formula>). Our detailed attribution analysis indicates that GPP is positively sensitive to surface air temperature (0.15 gC °C<sub>−1</sub>) and total precipitation (0.25 gC mm<sup>−1</sup>) but negatively sensitive to atmospheric CO<sub>2</sub> concentrations (−0.20 gC mol<sup>−1</sup>) and LCC (−0.93 gC class<sup>−1</sup>). Furthermore, we reported large differences in the spatial patterns and magnitudes among individual variables in the GPP attribution analysis, with LCC proving to be the dominant factor followed by CO<sub>2</sub> fertilization effects; climatic factors had comparatively little influence on GPP variations during the study period. Although MODIS GPP does not take CO<sub>2</sub> fertilization effect into account, the close relationship between MODIS GPP and atmospheric CO<sub>2</sub> concentrations still pose referencing value in attributing the GPP increase in this period. Overall, the findings of this study contribute to our understanding of the responses of sensitive ecosystems to the competing effects of climate change and human disturbance at regional scales.


I. INTRODUCTION
Terrestrial gross primary productivity (GPP), defined as the carbon uptake by terrestrial ecosystems through plant photosynthesis, is the largest component of the global carbon cycle and a key indicator of land ecosystem dynamics [1]- [3]. GPP is essential for several ecosystem applications, such as crop growth and yield prediction [4], vegetation disturbance monitoring [5], and solar-induced fluorescence retrieval [6] at regional and global scales. As an important part of the The associate editor coordinating the review of this manuscript and approving it for publication was John Xun Yang .
East Asian ecosystem, the Mongolian Plateau (MP) not only represents an important ecological barrier in China but also plays an important role in the global carbon cycle [7]. Therefore, accurate information on GPP in this region is vital for both local ecosystem monitoring and global carbon budgeting [8], [9] against the backdrop of a changing climate.
The increasing attention toward patterns and attribution of GPP change in past few years has been driven by the global need for accurate peak carbon and carbon neutral calculations [1]- [3], [8], [9]. Statistical interpolation [10], [11], satellite-based estimation [12], [13], and carbon cycle models [14]- [17] are three widely employed approaches in GPP VOLUME 10, 2022 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ estimation. Compared with statistical interpolation and model simulations, estimation using satellite data has distinct advantages for GPP monitoring, including large spatial coverage, fine spatial resolution, high temporal resolution, and economic and practical in application. Several published GPP datasets are currently available [12], [13], [18], e.g., the Moderate Resolution Imaging Spectroradiometer (MODIS) [12] and Global Land Surface Satellite (GLASS) [13] datasets. Although GPP estimates from carbon models and satellite observations show similar patterns at a global scale, there are large differences in the temporal trends, seasonality, and interannual variability of GPP calculated from individual regional and local datasets because of varying inputs and forcing mechanisms [1], [19]. In addition, due to the lack of observational GPP data over the MP, accurate estimation in this region remains difficult, and high-quality GPP datasets are crucially needed. Among the published satellite-retrieved estimates, the good performance of MODIS GPP, based on a light-use efficiency algorithm, has been highlighted. The validation of GPP calculated from ground sites and MODIS-based estimates shows good spatiotemporal correlation [10], [11], [20]. Moreover, the newly released gap-filled MODIS GPP estimates (MOD17A2HGF) Version 6 [21] provide spatially complete GPP data at a 500-m spatial resolution, offering new opportunities for regional-scale studies.
In addition to studies of GPP variability, the attribution of GPP anomalies has recently gained significant attention. Several putative drivers have been attributed to GPP changes including surface air temperature (T s ) [22], precipitation (P t ) [23], [24], atmospheric carbon dioxide (CO 2 ) concentrations [25], [26], aerosol optical depth (AOD) [27], plant phenology and physiology [28], and land-cover change (LCC) [25]. Generally, increases in atmospheric CO 2 concentrations dominate GPP changes via CO 2 fertilization effects; the T s and P t produce similar but lower-magnitude positive GPP effects; solar radiation (R s ) is associated with an overall decrease in GPP; and LCC is positively correlated with GPP [25]. However, the mechanisms underlying regional-scale GPP attribution remain debated. The putative driving factors of large-scale variability in GPP are determined by location and biome type. For example, AOD is positively correlated with GPP in forests but negatively correlated with GPP in grasslands over China [27]. Furthermore, although the sensitivity of GPP to LCC is positive at a global scale, a negative effect is observed in the rainforests of South America and Eurasia [25]. Although attribution analyses have been carried out at different scales, the relative contributions of different drivers of GPP change remain highly uncertain [25], [29], especially at regional scales.
Affected by the Siberia-Mongolia High in winter, the East Asian Monsoon in summer, and the westerly circulation, the MP is extremely fragile and ecologically sensitive to climate variations [30], [31]. Under the competing effects of climate change and human activity, the MP has experienced drastic ecosystem shifts over the past several decades including grassland degradation [32], forest decline [33], cooling effects from re-vegetation [34], and the rapid loss of lakes [35]. Nevertheless, evidence of GPP variability in this region remains limited. To address this gap, here we investigate the pattern of GPP changes in the MP region and its drivers using the latest MODIS GPP simulations from 2001 to 2018. We specifically focus on the relationships between GPP and four major putative drivers (T s , P t , CO 2 concentration, and LCC) and their underlying forcing mechanisms. Results in this study provide valuable evidence of GPP change and its underlying mechanisms in this globally important region, having implications for global carbon budgets as well as understanding ecosystem sensitivities and responses to climate change and human disturbance.

A. STUDY AREA
Due to its arid and semi-arid climate, habitats of the MP transition from deserts to grasslands and forests, from southwest to northeast [31]. As rapid climate and land-cover changes in this region may have a marked impact on GPP, we attributed GPP variability in grasslands, forests, and deserts according to the global terrestrial ecoregions classification developed by the Nature Conservancy (TNC) [36], as shown in Figure 1. According to the TNC classification, the MP is covered by nine ecoregions, which were re-classified into forests, grasslands, and deserts to facilitate statistical analysis. The ''forests'' classification included temperate conifer forests, temperate broadleaf and mixed forests, boreal forests/taiga, and tundra; the ''grasslands'' category included montane grasslands and shrublands, flooded grasslands and savannas, as well as temperate grasslands, savannas, and shrublands; and the ''deserts'' category included deserts and xeric shrublands. Forests are mostly distributed in mountainous areas, including the Altay Mountains in the west, the Khangai Mountains in the west, and the Khentii Mountains in the north, which account for 32.58% of the total study area. Grasslands are distributed in the middle of the MP and account for 38.84% of the total area. Finally, deserts are located in the southwest of the MP and account for 28.04% of the total area. Although the deserts in the north of the region are more arid, they receive some precipitation in summer and support sparse vegetation that contributes to GPP.

B. DATASETS
We used the latest satellite-based MODIS GPP dataset and datasets of major putative drivers covering the period from 2001 to 2018.

1) GPP DATASETS
To capture the spatial pattern and assess interannual variability and long-term trends in GPP, eight-day composite gap-filled MOD17A2HGF data [21] with a 500-m spatial resolution were downloaded from the National Aeronautics and Space Administration (NASA) Land Processes Distributed Active Archive Center (https://lpdaac.usgs.gov/products/ mod17a2hgfv006/).
MODIS GPP was the first satellite-based modeled GPP dataset for monitoring global vegetation productivity with complete spatial coverage and high spatial resolution.
To ensure data quality, poor-quality inputs from the eight-day Fraction of Photosynthetically Active Radiation and Leaf Area Index have been removed from the MOD17A2HGF dataset based on the quality control flag for each grid cell.

2) LAND-COVER DATASET
To attribute GPP changes to LCC, we employed the annual International Geosphere-Biosphere Programme (IGBP) classification from the Terra and Aqua combined MODIS yearly land-cover dataset (MCD12Q1) with a 500-m spatial resolution. This dataset covers 17 IGBP categories encompassing 11 natural vegetation categories, three land-use and land-mosaic categories, and three vegetation-free land categories [37]. For the purpose of our study, the land-cover data were analyzed in a 0.10 • grid, in which all 500-m pixels encompassed in the 0.10 • cells were used to calculate the proportions of the dominant land-cover types.

3) ATMOSPHERIC CO 2 CONCENTRATION DATASET
To attribute GPP changes to atmospheric CO 2 concentrations, we used the Carbon Tracker CT2019B product, providing global monthly continuous spatial surface-atmosphere flux of atmospheric CO 2 concentrations. CT2019B is distributed by the Global Monitoring Laboratory (GML) of the National Oceanic and Atmospheric Administration (NOAA), and is based on observations from the NOAA Earth System Research Laboratories (ESRL) greenhouse gas observational network and collaborating institutions (spatial resolution = 1.0 • ) [38]. In this dataset, terrestrial biosphere, wildfire, fossil fuel emissions, atmospheric transport, and other factors are assimilated to estimate atmospheric CO 2 mole fractions.

4) REANALYSIS TEMPERATURE AND PRECIPITATION DATASET
Previous studies have reported that climatic factors are closely related to GPP changes at a regional scale [22]- [24]. To attribute GPP changes to local meteorological variables, we employed monthly averaged 2-m surface air temperature (T s ) and total precipitation (P t ) datasets with a resolution of 0.10 • (approximately 9 km) from the European Center for Medium Range Weather Forecasts (ECMWF) Reanalysis v5-Land (ERA5-Land) [39]. These data are widely considered suitable for ground-level modeling [40]- [42]. Compared with other ERA5 reanalysis data, ERA5-Land offers more accurate classifications of land parameters and land status.

5) SOLAR RADIATION DATASET
Incoming solar radiation (R s ) is the energy source for photosynthesis and causes changes in T s , relative humidity, and evaporation, and thus, indirectly affects plant productivity [43]. Therefore, R s is the one of the most important environmental factors affecting terrestrial ecosystem productivity and carbon budgets. Here, to explore whether changes in R s are associated with GPP anomalies over the MP, we employed the Clouds and the Earth's Radiant Energy System (CERES) Energy Balanced and Filled (EBAF) dataset [44] (spatial resolution = 1.0 • ).

6) DATA PREPARATION
All of the used datasets are summarized in Table 1. Prior to analysis, all the gridded datasets were resampled to gridcells with a spatial resolution of 0.10 • in geographic projection. Using ''gdalwarp'' software (https://gdal.org/programs/gdalwarp.html), datasets with a spatial resolution coarser than 0.10 • were regridded using the ''cubic-spline'' function; datasets with spatial resolution finer than 0.10 • were regridded using the ''average'' function.

C. RELATIVE CONTRIBUTION CALCULATION
The influences of four putative drivers on GPP (T s , P t , atmospheric CO 2 concentration, and LCC) were analyzed in detail using multiple linear regression analysis [45], [46]. To avoid multicollinearity among the four drivers, we employed ridge regression. First, each driver was normalized to facilitate cross-comparison between variables with different units and scales. For variable X i , the z-score (X iz ) was calculated using Equation (1) [47]: where X iz is the normalized variable X i , µ x is the mean value of variable X i , and δ x is the standard deviation of variable X i . Second, ridge regression analysis was carried out to calculate the sensitivity of GPP to the four putative drivers using Equation (2): where GPP z is the z-score-normalized GPP; X iz is the normalized putative driver X i ; β i is the standard ridge regression coefficient for putative driver X i ; α is the residual error, representing the contribution of unknown factors to GPP, such as fire, pests, wind, and disease; and n is the number of putative drivers. Third, the relative contributions of the putative drivers to GPP variability were obtained using the ridge regression coefficients and the z-score series of each driver, as follows: where η ci is the contribution of putative driver X i to the z-score-normalized GPP variation, X iz is the z-scorenormalized X i , and β i is the linear slope of X iz . Finally, the relative contribution (η rc ) of each individual putative driver to GPP was confirmed using Equation (4): where η rci is the relative contribution of X i to GPP change, and n is the number of putative drivers. Ridge regression analysis is an effective approach for solving collinearity problems between independent variables [48] and has been successfully applied in GPP trend analysis of the 'Three North' region of China (northeastern, northern, and northwestern regions) [46] and Yellow River Basin [49].

III. RESULTS
In the subsequent subsections, changes in GPP over the MP are first described for the period 2001-2018. Then, concurrent changes in the four putative drivers are described. Finally, the results of attribution analysis are presented.

A. CHANGES IN GPP OVER THE MP
The spatial pattern of the 18-year-averaged (2001-2018) annual cumulative GPP over the MP is displayed in Figure 2, as calculated from the MOD17A2HGF dataset.
With the transition from desert to forest, GPP generally increases from the southwest to northeast parts of the MP (Figure 2a  (in 2018). Over the entire study period, the linear regression analysis indicated a significant increase in annual cumulative GPP at a rate of 3.91 gC m −2 yr −1 (p < 0.01); in deserts, forests, and grasslands, the rate of increase was 0.79 (p < 0.01), 4.79 (p < 0.01), and 5.76 gC m −2 yr −1 (p < 0.01).
The observed GPP increase was strongest in the southeastern part of the MP, coinciding with the grasslands in this region [50]. Moreover, forests, grasslands, and deserts accounted for 39.94%, 53.72%, and 6.44% of the overall change in GPP during the study period.

B. CHANGES IN METEOROLOGICAL VARIABLES OVER THE MP
Changes in T s and P t are shown in Figure 3 for the 2001-2018 study period, as calculated from the latest ERA5-land atmospheric reanalysis dataset.
The spatial pattern of annual mean T s exhibited a decreasing trend from south to north over the MP (Figure 3a). The changes in T s exhibited a conspicuous warming trend in grasslands and deserts areas (Figure 3b) compared to a cooling trend in the northwest part of the MP, which is consistent with the reported Eurasian cooling caused by a strengthening Siberian High [51].
Compared with T s , P t ( Figure 3c) and associated overall changes during the study period (Figure 3d) exhibited larger spatial differences within the study area. The 18-year averaged annual-mean P t showed an increasing trend from the deserts in the southwest to the forests in the northeast. P t showed an increasing trend in the eastern MP but decreased at the southern and northern margins of the MP. Notably, changes in P t were opposite to the Palmer Drought Severity Index (PDSI) anomalies in this region during the 2000s [52].
The detailed changes in annual-mean T s and P t are summarized in Table 2.  In most regions, the changes in T s and P t were insignificant, except for the desert regions. In deserts, P t decreased by −3.21 mm between 2001 and 2018, which was positively correlated with the increase in GPP in this region (r = 0.52, p < 0.05). In addition, although grassland P t showed an insignificant increasing trend between 2001 and 2018, this was positively correlated with GPP changes (r = 0.45, p < 0.10).

C. CHANGES IN ATMOSPHERIC CO 2 CONCENTRATION AND LAND COVER
The spatial pattern of 18-year-averaged atmospheric CO 2 concentrations and the dominant land-cover types over the MP are shown in Figure 4. Influenced by the intensity of human activity in the MP region, atmospheric CO 2 concentrations are high in Inner Mongolia relative to lower values in the southwestern deserts ( Figure 4a). The dominant land-cover types exhibit a remarkable spatial transition from the deserts in the southwest to the grasslands in the central area and the forests in the northeast (Figure 4b).
The linear regression analyses revealed a significant increase in atmospheric CO 2 concentrations in the MP region between 2001 and 2018 (1.99 mol m −2 , p < 0.05). Notable regional differences were revealed, however, with significant decreases in the northeast forest region and increases in the grassland regions, especially at the southern margin. These changes are similar to the decrease in CO 2 concentrations previously reported in this region for the period between 2009 and 2018 [53]. The dominant landcover types exhibited a decreasing trend during the study period (−0.11 class −1 , p < 0.05), indicating an overall transition from desert to grassland and from grassland to forest. As shown in Figure 4d, the most remarkable changes have occurred near the transition regions between deserts and grasslands, which corroborates the findings reported by Jiang et al. [34].
The observed changes in annual-mean atmospheric CO 2 concentrations and the dominant land-cover types during the study period are detailed in Table 3.   Table 3, atmospheric CO 2 concentrations over grassland and desert areas increased at rates of 4.83 (p < 0.05) and 0.68 mol m −2 (p < 0.05) between 2001 and 2018. Moreover, atmospheric CO 2 was highly correlated with annual GPP augmentation in the forest (r = −0.54, p < 0.05), grassland (r = 0.72, p < 0.05), and desert (r = 0.68, p < 0.05) areas.

D. ATTRIBUTION ANALYSIS
The changes in annual cumulative GPP over the MP during the period 2001-2018 and their sensitivities to individual putative drivers, calculated using Equation (3), are shown in Figure 5.
Annual cumulative GPP shows a conspicuous rising trend in the central and southeastern regions of the MP (Figure 2b), and this is positively correlated with increases in T s (Figure 5a) and atmospheric CO 2 concentrations (Figure 5c) in these regions. As shown in Figure 5d, the evident LCC decrease is also associated with the increase in GPP at the eastern margin of the MP.
Detailed ridge regression coefficients reveal that GPP change was positively correlated with T s (0.1495 gC • C −1 ) and P t (0.2424 gC mm −1 ) but negatively correlated with atmospheric CO 2 concentration (−0.1993 g C mol −1 ) and LCC (−0.9348 gC class −1 ) during the 2001-2018 study period. However, there were large differences among the desert, grassland, and forest land-cover categories.  For example, annual cumulative GPP was negatively correlated to atmospheric CO 2 concentration (−0.3958 gC mol −1 ) in forests but significantly and positively correlated with grasslands (0.3247 gC mol −1 ) and deserts (0.1101 gC mol −1 ) over the same period.
The relative contributions of each putative driver to GPP anomalies as calculated from Equation (4) are shown in Figure 6; the contributions of T s , P t , atmospheric CO 2 concentration, and LCC over the entire MP were 9.16%, 14.87%, 12.64%, and 63.32% for the period 2001-2018, respectively.

A. INTERANNUAL VARIATION IN GPP
According to Sun et al. [25] and Alexandrov [54], an increase in both the amplitude of the GPP seasonal cycle and in growing season length can lead to elevated GPP.
To comprehensively analyze the reasons for observed GPP increases in the MP region, we mapped eight-day seasonal variations in cumulative GPP for the entire region based on four sub-periods (2001-2005, 2006-2010, 2011-2015, and 2016-2018), as shown in Figure 7.
Notably, the multi-year averaged annual GPP for the earlier (2001-2010) and latter (2011-2018) sub-periods showed a remarkable change in amplitude (height of peaks in Figure 7a) whereas changes in growing season length were negligible (width of peaks in Figure 7a). We also found similar patterns for the amplitudes of the GPP seasonal cycle in forest (Figure 7b), grassland (Figure 7d), and desert ( Figure 7e). Based on this, we suggest that changes in the length of the growing season can be excluded from the attribution of GPP changes during our study period.

B. CHANGES IN SOLAR RADIATION
Previous studies have shown that R s changes significantly in response to AOD anomalies, which may lead to GPP changes [27] and [55]. To explore whether R s changes Contributed to the observed increase in GPP in the MP region, we mapped the 18-year annual mean R s and its associated overall change in Figure 8.
Between 2001 and 2018, the 18-year averaged annual mean R s in the MP region was 293.31 (± 0.12) W m −2 , with a clear latitudinal gradation from south to north (Figure 8a). However, overall changes during the study period were minimal (Figure 8b). Detailed statistical analysis showed that the regions with R s changes between −0.05 and 0.05 W m −2 and between −0.10 and −0.05 W m −2 account for 89% and 11% of the study area, respectively. This indicated that R s had little influence on annual GPP variability during the study period, and thus, could be excluded from the attribution analysis.

C. UNCERTAINTY ANALYSIS
Global terrestrial GPP has increased by 31 ± 5% since 1900, driven by the effects of a changing climate and increasing atmospheric CO 2 concentrations on several ecosystem processes including vegetation productivity, harvesting, deforestation, and secondary forest regrowth [29]. However, the patterns of GPP changes in terrestrial ecosystems demonstrate high spatial variability because of the coupled interactions between vegetation characteristics and environmental variables, e.g., rising atmospheric CO 2 concentrations, changing land cover, AOD, and climatic variability [2], [25], [27].
Satellite-based GPP estimates have been considered highly reliable because of the advantages they offer, with respect to consistent spatial and temporal information on vegetation dynamics. Compared with other satellite-based and model-simulated GPP estimates, the MOD17A2HGF product has proven a feasible and useful means of evaluating GPP in the MP region. We observed a remarkable increase in GPP in this region between 2001 and 2018, and this has partly coincided with enhanced net CO 2 uptake by terrestrial and marine ecosystems, as reported for the period 1959-2010 [56]. Nevertheless, our results differ from recent reports suggesting that there has been no proportional augmentation in terrestrial gross carbon sequestration from enhanced greening of the Earth's surface [57]. Moreover, in contrast to the findings of Zhou et al. [28], our analysis reveals a remarkable increase in the amplitude of the seasonal GPP cycle rather than an increase in the length of the growing season. Therefore, we excluded phenological shifts from our attribution analysis.
Previous work has indicated that atmospheric CO 2 concentration is the dominant factor driving global GPP increases between 1982 and 2015 [25]. However, our attribution analysis revealed that LCC dominated the GPP changes in the MP region over the last few decades, which contrasts with VOLUME 10, 2022 the findings of regional-and global-scale GPP attribution studies by Sun et al. [25] and Zhang et al. [27]. Since MOD17A2HGF does not take CO 2 fertilization effect into account, the MOD17A2HGF-based analysis on CO 2 fertilization effect is not valid. However, the close relationship between MOD17A2HGF and atmospheric CO 2 concentrations still pose referencing value in attributing the GPP increase in this period. Furthermore, as shown in Figure 4d, land-cover change has dominated the transition between desert and grassland, and between grassland and forest, in the MP region, with GPP showing a significant and sensitivity to LCC in these regions. In addition, although we found that P t has not been the dominant factor affecting GPP change overall, it has played a key role in the grassland area of the MP region. This is corroborated by findings of Yuan et al. [24], suggesting that P t plays an important role in GPP dynamics in the grasslands of northern China. Moreover, different from the general global rising trend, atmospheric CO 2 concentrations have decreased in the forested areas of MP, and this has been negatively correlated with the increase in GPP across the region. This is not unexpected, as CO 2 fertilization only partially accounts for changes in GPP, as reported by Keenan et al. [58].

V. SUMMARY AND CONCLUSION
The estimation and attribution of GPP changes are of great significance for understanding regional terrestrial carbon cycling, the response of vegetation to climate change and human activity, and for assessing ecosystem health. Using satellite-based and model ensemble GPP estimates and several ancillary datasets, we mapped the pattern of annual GPP over the MP and explored its associations with both climate change and human activity indicators. Our results help fill the current gap in GPP studies in the MP region and can inform future regional-scale terrestrial ecosystem studies.
Using the latest MODIS GPP dataset, we found a general increase in GPP in the MP region from southwest to northeast, with remarkable spatial variability among forest, grassland, and desert areas between 2000 and 2018. Overall, our linear regression analysis indicates that annual GPP has shown a significant increase during this period, at a rate of 3.91 gC m −2 yr −1 (p < 0.01). Based on multiple datasets, we attribute these trends to climatic factors, LCC, and atmospheric CO 2 concentrations, with LCC accounting for 63.32% of the observed increase followed by atmospheric CO 2 concentration (12.64%). In comparison, we suggest that climatic factors have had a limited influence on the GPP increase in the MP region during this period.
Distinguishing the effects of climate, LCC, and atmospheric CO 2 concentrations not only provides valuable insight into the impacts of future climate change on productivity, but also reveals the ecological consequences of ongoing LCC and rising CO 2 concentrations. As such, in addition to being relevant to the MP region, our results can inform localscale GPP studies in inland arid and semi-arid regions, as well as provide additional motivation for mitigating the negative impacts of ongoing human activity and future climate change.