Assessment of Drought Impact on Net Primary Productivity in the Terrestrial Ecosystems of Mongolia from 2003 to 2018

Abstract

Drought has devastating impacts on agriculture and other ecosystems, and its occurrence is expected to increase in the future. However, its spatiotemporal impacts on net primary productivity (NPP) in Mongolia have remained uncertain. Hence, this paper focuses on the impact of drought on NPP in Mongolia. The drought events in Mongolia during 2003–2018 were identified using the Moderate Resolution Imaging Spectroradiometer (MODIS) normalized difference vegetation index (NDVI). The Boreal Ecosystem Productivity Simulator (BEPS)-derived NPP was computed to assess changes in NPP during the 16 years, and the impacts of drought on the NPP of Mongolian terrestrial ecosystems was quantitatively analyzed. The results showed a slightly increasing trend of the growing season NPP during 2003–2018. However, a decreasing trend of NPP was observed during the six major drought events. A total of 60.55–87.75% of land in the entire country experienced drought, leading to a 75% drop in NPP. More specifically, NPP decline was prominent in severe drought areas than in mild and moderate drought areas. Moreover, this study revealed that drought had mostly affected the sparse vegetation NPP. In contrast, forest and shrubland were the least affected vegetation types.

1. Introduction

The global surface temperature has been increasing since the pre-industrial period [1,2,3]. It is undeniable that climate change has a notable impact on natural ecosystems and indirectly affects human life and the economy across the world [4]. Rangeland ecosystems, however, are usually second- or third-class quality land in terms of production potential and are therefore vulnerable to climate change. Worldwide, 40–50 million km² of rangelands constitute about 40% of the total land area [5] and have been chiefly used in fodder production. Moreover, rangelands play an essential role in regulating the carbon cycle of forests due to their contribution to a large part of regional, continental, and global net primary productivity (NPP).

The climatic and multiple environmental attributes have been considered as an important factor driving the global terrestrial NPP in previous studies. For example, Nemani et al. [6] reported that climate change has resulted in increasing global terrestrial NPP during the period between 1982 to 1999. In contrast, climate change events like drought have resulted in a decline in global terrestrial NPP over the past decades according to Zhao and Running [7]. Thus, the influence of climate change on NPP has largely remained ambiguous owing to its complexity in nature [8,9]. Whilst the influence of climate change on the variability of NPP has been studied extensively in Mongolia [10,11], the impact of droughts either singly or in combination with other environmental factors, such as fires and land use change, on the Mongolian terrestrial NPP is inadequately explored.

Mongolia has widespread drought-sensitive land, which merits the assessment of NPP response to droughts because half (52%) of the land area comprises arid and semiarid environments [12]. An earlier study revealed that drought affects one-fourth of the country every 2 to 3 years, while half of it is affected every 4 to 5 years [13]. More recent findings have shown that the average drought-affected areas accounted for one-third of the total area every year from 2000 to 2016 [14]. Thus, the impact of droughts in Mongolia is more severe than in other regions due to which it has been identified as one of the extremely vulnerable countries to climate change [15,16].

Drought is one of the most common catastrophes induced by climate change throughout the world. Drought can be simply defined as long-term environment dryness usually resulting from a low precipitation at the specific location [17]. As both drought and aridity characterize environment dryness, some authors use the term drought and aridity/dryness interchangeably [18]. Unfortunately, its occurrence is predicted to increase in the context of increasing climatic variability in the future [19]. Drought inadvertently reduce vegetation productivity. Recently, the notable effects of drought events on NPP at both regional and global levels have been observed [7,20,21,22]. Studies have shown that both observed and simulated droughts resulted in reduced NPP worldwide [22,23]. However, some studies have shown that drought increased the NPP of the natural ecosystem [7,24]. Moreover, some studies have found that productivity remained stable under long-term drought events [25]. Global-scale research based on Moderate Resolution Imaging Spectroradiometer (MODIS) imagery revealed that NPP decreased by 50% during 2000–2014 in Mongolia [22]. Although the MODIS NPP products (MOD17) at the annual scale are excellent at the global scale, they may contain systematic errors in arid and semiarid regions such as Mongolia because drought is a major disastrous event during the growing season for this country. Therefore, daily to seasonal scale modeling efforts need to be tested considering the independent contribution of drought on NPP across Mongolia. Moreover, the mechanism and magnitude of the impact of drought on terrestrial NPP variability at a local to regional scale have largely remained uncertain. Thus, evaluating the impacts of drought on terrestrial NPP in Mongolia, particularly for the rangelands, is a critical challenge.

Historical records show a high precipitation variability during June–August in Mongolia, which makes up the majority of the annual total rainfall. Therefore, a drought may have long-lasting effects on vegetation production in the growing seasons. Droughts especially have a far-reaching effect on growing season NPP (NPP_gs) dominant vegetation. The response of NPP to drought conditions is one of the notable indicators to understand the effects of climate change on regional or global carbon cycling. This study intends to evaluate the impacts of droughts on rangelands in Mongolia. In the light of the anticipated increase in drought frequency and severity in the future (5–15%) [26], it is also imperative to evaluate the status of vegetation NPP over the next century.

Previously, in situ experiments [27], remote sensing imageries [28], and ecological modeling [29] have been employed to assess the responses of NPP to drought conditions. The station-based drought index (e.g., the standardized precipitation evapotranspiration index, SPEI [30]; the Palmer drought severity index, PDSI; [31], and the standardized precipitation index, SPI; [32]) has been used. However, such indicators have only considered temperature and precipitation to assess the impacts of drought on NPP. Moreover, these drought indices can effectively evaluate the impact of drought on NPP. However, sparse distribution of field and station data limits the assessment of the impact of drought on terrestrial NPP at the broader spatial and multi-temporal coverage. Remote-sensing indices and methods, such as NPP and NDVI, have made it possible to understand vegetation response to drought from a broader perspective and larger spatial extent [33,34,35,36,37].

We used the MODIS NDVI and a process-based ecosystem model, namely the Boreal Ecosystem Productivity Simulator ((BEPS), Liu et al. [38]) to determine the response of NPP on drought events in Mongolia during 2003–2018. The objectives of this study were to (1) quantify the spatial and seasonal (growing season) variability of drought events, (2) examine the spatial and temporal variations in terrestrial NPP, (3) assess the impact of drought on NPP change, and (4) identify uncertainties and knowledge gaps associated with NPP estimation in Mongolia.

2. Data and Methods

2.1. Study Area

Mongolia spans 41°35′–52°09′ N in latitude and 87°44′–119°56′ E in longitude, with an area of 1.565 million km². It has an average elevation of 1580 m AMSL, which gradually increases from east to west (Figure 1). The country has a variety of land features, including vast semi-desert and desert plains, grassy steppe, and mountains. Most of the places (52%) experience arid climatic condition, which can be further categorized into hyper-arid (5.13%), arid (6.74%), and semiarid (40.13%) environments [12]. Mongolia receives relatively low precipitation, although there is high variation in its spatial distribution. As per the station data, the annual rainfall ranges from 50 to 150 mm in desert and semi-desert areas, while precipitation in other regions ranges from 250 to 400 mm. The region receives most of its rainfall (about 85%) in the growing season (May–September) [39]. The annual mean temperature in Mongolia is around 0 °C, with a north–south gradient ranging from −8 to 6 °C. However, it experiences highly fluctuating temperatures across different seasons in the year. The maximum temperature in winter (November–April) remains between −30 and −40 °C in many regions and as low as −48 °C in the mountains. In summer (June–August), temperatures rise sharply and reach as high as 40 °C in some places such as the South Gobi. The major land use types include alpine tundra (3.0%), mountain taiga (4.1%), mountain steppe and forest steppe (25.1%), grass steppe (26.1%), desert steppe (27.2%), and desert (14.5%) [40]. Generally, forest steppe has the highest productivity, followed by the steppe, desert steppe, and desert zones. There are 3191 taxa of vascular plants belonging to 684 genera of 108 families documented in the flora of Mongolia [41].

2.2. Data Processing Method

The overview of the methodology for this study is outlined in Figure 2. We identified several drought events that occurred in Mongolia during the past 16 years (2003–2018) using the NDVI anomaly as a drought indicator. Meanwhile, NPP in the growing season of 2003 to 2018 across Mongolia was estimated using the BEPS model. The daily NPP was computed using the BEPS model. Then, the growing seasonal NPP was obtained summating the daily values between days 113 and 289, which corresponded to the drought monitoring period. The relationship between drought severity and NPP was examined by calculating separately the NDVI anomaly and NPP for various vegetation types based on the ESA land cover map (2015) of Mongolia. Later, we analyzed the time series analysis of the growing season total NPP and mean NDVI anomaly in the drought-affected areas.

2.2.1. NDVI Anomaly as the Indicator of Droughts

Drought was defined as relative to the statistical average for one or several seasons or years. The relative NDVI anomaly renders one of the simplest NDVI-based approaches for monitoring drought [42]. We used Terra MODIS NDVI to monitor the drought variation in Mongolia during the vegetation growing season (May to September) from 2003 to 2018. To derive NDVI anomaly, the individually mean NDVI ($NDV{I}_{\mathrm{mean}i}$) was determined as the annual NDVI using Equation (1) following Nanzad et al. [14]:

$$NDV{I}_{\mathrm{mean}i}=\left(NDV{I}_1+NDV{I}_2+\dots \dots ,NDV{I}_n\right)/n$$

(1)

where $NDV{I}_{\mathrm{mean}i}$ is mean NDVI value of the growing season of the ith year, and $NDV{I}_1$ is the first 16-day NDVI composite to $NDV{I}_n$, which is the last 16-days NDVI composite during the growing season of the ith year.

In addition to that, the growing season long-term mean NDVI (2003–2018) was calculated using Equation (2) following Chopra et al. [43]:

$$\overline{NDVI}={\displaystyle \sum }_{i=1}^n\frac{NDV{I}_{\mathrm{mean}i}}{n}$$

(2)

where $n$ is the number of years, which is equal to 16 years.

We derived the spatial patterns of mean growing season NDVI anomaly for each grid cell using Equation (3) in the study area [44]:

$$NDV{I}_{anomalyi}=\frac{NDV{I}_{\mathrm{mean}i}-\overline{NDVI}}{\overline{NDVI}}\times 100$$

(3)

where $NDV{I}_{anomalyi}$ is the percentage of the growing season at a grid cell during the ith year, $NDV{I}_{\mathrm{mean}i}$ is the individual mean growing season, and $\overline{NDVI}$ is the mean growing season from 2003 to 2018.

Here, we defined drought conditions as a negative variation of the NDVI, while positive values indicated normal or wet conditions [45]. Based on the NDVI anomaly classification scheme, the drought severity was classified into mild drought, moderate drought, severe drought, and extreme drought. Details about the categorization of dryness/wetness grade of NDVI anomaly are described in

Table 1. Categorization of dryness/wetness of normalized difference vegetation index (NDVI) anomaly values.

Category	NDVI Anomaly Values (%)	Category	NDVI Anomaly Values (%)
Extreme drought	Below −50	Mild wet	0 to 10
Severe drought	−25 to −50	Moderate wet	10 to 25
Moderate drought	−10 to −25	Severe wet	25 to 50
Mild drought	0 to −10	Extreme wet	Above 50

.

2.2.2. BEPS for the Calculations of NPP

We used the BEPS model to evaluate the terrestrial NPP in Mongolia from 2003 to 2018. The model was initially built and validated for daily boreal ecosystems [38] and has been commonly used for calculating global and regional NPP. Later, it was employed to monitor the dynamic changes in global and regional NPP of various vegetation types, such as shrubland, farmland, and grassland ecosystems [46]. The model computes NPP based on the difference between total carbon uptake from the air through photosynthesis (GPP) and the carbon lost due to respiration ($R_a$) [38,47]. Equations (4)–(6) were used to run the NPP [48]. The model was operated on a daily step throughout the year from 2003 to 2018, pixel by pixel.

$$NPP=GPP-R_a$$

(4)

$$GPP=A_{canopy}{L}_{day}{F}_{GPP}$$

(5)

$$R_a=R_m+R_g$$

(6)

where NPP refers to net primary productivity (g C m⁻²), GPP (gross primary productivity; g C m⁻²) is the carbon flux during photosynthesis, $A_{canopy}$ is the carbon assimilation rate, $L_{day}$ is the day length, and $F_{GPP}$ is the coefficient of conversion from umol m⁻² s⁻² to g C m⁻² day⁻¹. Daily $R_a$ is the autotrophic respiration, while $R_m$ and $R_g$ are the maintenance respiration and growth respiration by all living parts including fine roots and leaves, respectively [49,50].

The assessment of drought-induced NPP anomalies was done by utilizing an anomaly index based on the growing seasonal NPP from 2003 to 2018 as defined in Equation (7) following Lai et al. [51]:

$$NP{P}_{anomalyi}=\frac{NP{P}_{\mathrm{mean}i}-\overline{NPP}}{\overline{NPP}}\times 100$$

(7)

where $NDV{I}_{anomalyi}$ is the NPP anomaly of an individual grid cell for a drought event, $NP{P}_{\mathrm{mean}i}$ is the accumulated NPP value of the ith drought event, and $\overline{NPP}$ is the long-term average growing seasonal value of NPP corresponding to the period of the ith drought. For example, if the drought occurred from May 2004 to September 2004, the $NP{P}_{\mathrm{mean}i}$ represents the accumulated NPP value for this period, while the $\overline{NPP}$ represents the mean of NPP values of the corresponding drought period (May–September across 2003–2018).

2.2.3. Data Acquisition and Preparation

The Terra MODIS NDVI product (MOD13Q1, collection v006) with 250 m ground resolution and 16-day revisit cycle was downloaded from the National Aeronautics and Space Administration Earth Observing System Dataset and Information System (https://ladsweb.modaps.eosdis.nasa.gov, accessed on 17 October 2017) for drought index. The daily climatic data, including precipitation, air temperature, net surface radiation, wind speed, snow depth, and dew point temperature, for the period of 2003–2018 were retrieved from the European Centre for Medium-Range Weather Forecast (ECMWF) products. The ECMWF provides global gridded three-hour time-varying climatology products (http://www.ecmwf.int/research/era, accessed on 10 December 2018) with a spatial resolution of 0.125° × 0.125°. The data provided in the NetCDF-4 (network Common Data Form) file format was processed and converted to a common img. format in their respective units as mentioned in

Table 2. List of major parameters of BEPS model, their source, and spatiotemporal extent. ECMWF, European Centre for Medium-Range Weather Forecasts; MODIS, Moderate Resolution Imaging Spectroradiometer; ESA, European Space Agency; SRTM, Shuttle Radar Topography Mission; DAAC, Distributed Active Archive Center.

S.N.	Complete Name	Format	Units	Source	Spatial Resolution	Temporal Resolution	Reference
1	Precipitation	img	mm	ECMWF	12.5 km	Twice per day	[52]
2	Net solar radiation	img	W·m−2	ECMWF	12.5 km	Twice per day	[52]
3	Net thermal radiation	img	W·m−2	ECMWF	12.5 km	Twice per day	[52]
4	Temperature	img	°C	ECMWF	12.5 km	Twice per day	[52]
5	Snow depth	img	mm	ECMWF	12.5 km	Twice per day	[52]
6	Net solar radiation Downward	img	W·m−2	ECMWF	12.5 km	Twice per day	[52]
7	Dew point temperature	img	°C	ECMWF	12.5 km	Twice per day	[52]
8	Wind speed	img	ms−1	ECMWF	12.5 km	Twice per day	[52]
10	Leaf area index	tiff	NA	MODIS	500 m	Twice per day	[53]
11	Vegetation type (C3)	tiff	NA	ESA	300 m	Annual	[54]
12	Land cover	tiff	NA	ESA	300 m	Annual	[54]
13	DEM	tiff	m	SRTM	30 m	One	[55]
14	Latitude	tiff	°(Degree)	SRTM	30 m	One	[55]
15	Tree cover	tiff	%	MODIS	250 m	Yearly	[56]
16	Field capacity	tiff	KPA kg−1	DAAC	8 km	One	[57]
17	Wilting point	tiff	KPA kg−1	DAAC	8 km	One	[57]

.

The vegetation index is an important input variable for the BEPS model. We used collection 6 (C6) leaf area index (LAI) data from the MODIS as input for vegetation dynamics that constrain BEPS for a more realistic estimation of NPP at a regional scale. The MODIS LAI product (MCD15A3H) with 4-day temporal and 500 m spatial resolution from NASA’s Terra and Aqua satellites was obtained through the Earth data search client website (https://0-doi-org.brum.beds.ac.uk/10.5067/MODIS/MCD15A3H.006, accessed on 1 January 2019). Data gaps in the 4-day temporal MODIS LAI were interpolated to daily intervals from the original 4 days raster using the direct linear interpolation method for 2003–2018.

For the BEPS model, land cover information was obtained from the European Space Agency (ESA) Climate Change Initiative Land Cover (CCI-LC) product through the project’s website (http://maps.elie.ucl.ac.be/CCI/viewer/download.php, accessed on 17 June 2019). In addition, an additional layer of percent tree cover at a spatial resolution of 250 m [56] and the available soil water-holding capacity (AWC), the 30 m digital elevation data (DEM), and the latitude data were also major input parameters for the BEPS model (Figure 2). Wilting point and field capacity at a resolution of 0.083° × 0.083° were downloaded from the International Geosphere–Biosphere Program, Data and Information System (IGBP-DIS; https://daac.ornl.gov/cgi-bin/dsviewer.pl?ds_id=569, accessed on 15 January 2019) to calculate the AWC. Shuttle Radar Topography Mission (SRTM) 30 m DEM was downloaded from an online portal (https://dwtkns.com/srtm30m/, accessed on 15 June 2019), whereas the latitude data was derived from the DEM data based on its y-coordinates. Finally, the measured aboveground biomass (AGB) was utilized to verify the BEPS model (Figure 1). The observations of AGB data were collected from 205 sample points through 41 stations across the country covering the year 2003–2018, provided by the Information and Research Institute of Meteorology, Hydrology, and Environment (IRIMHE) of Mongolia [10]. A monthly (every 25th day) measurement was made in the growing season during the study period. The observed AGB provided in the unit of dry matter (g/m⁻²), was converted to NPP value (g·C m⁻²year⁻¹) by applying a carbon conversion factor of 0.475 for grass and foliage components as C per gram dry mass [58].

3. Results

3.1. Validation of the NPP Calculation

Two methods were used to validate the estimated NPP. First, we compared estimated NPP with observation data, followed by corroboration with the results estimated by other models. The correlation analysis between BEPS-simulated NPP with the observed data gathered from May to September for the study period showed a moderate positive correlation (R² = 0.59), and the corresponding root mean square error (RMSE) was 13.22 g·C m⁻²year⁻¹ for the overall dataset (Figure 3).

In this study, the growing season total NPP over the Mongolian rangeland was determined to be 0.1713 Pg C year⁻¹ from 2003 to 2018, which is less than what was reported by Bao et al. (0.3 Pg C year⁻¹) [11] and Lin and Dugarsuren (0.71 Pg C year⁻¹) [59]. They both used the CASA (Carnegie–Ames–Stanford Approach, a light-use efficiency) model to estimate NPP in Mongolia for the different periods. Such variation in NPP estimates between the current study and previous studies is due to differences in the study periods, resolution of remote sensing imageries, and the fundamental models employed for the studies. Bao et al. [11] used NDVI products of 8 km spatial resolution for 30 years (1982–2011), while Lin and Dugarsuren [59] used the standard product of MOD13Q1 vegetation indices with 250 m spatial and 16-day temporal resolutions from 2000 to 2004. Meanwhile, this study used LAI of 500-m spatial resolution. Moreover, the foliage clumping effect in the BEPS reduces sunlit leaves and increases shaded leaves, which may be another reason.

3.2. Spatiotemporal Trends of NPP

The Mann–Kendall test and Sen’s slope estimator were applied to generate the time series (2003–2018) NPP and its driving climatic variables during the vegetation growing season over the Mongolian rangeland (Figure 4 and

Table 3. Statistic summary of the percentages of pixels showing different trends in growing season NPP and climatic variables across Mongolia during 2003–2018.

Variable	Nonsignificant Increase (%)	Significant Increase (%, p < 0.05)	Nonsignificant Decrease (%)	Significant Decrease (%, p < 0.05)
NPP	56.37	7.23	46.63	2.60
Precipitation	51.15	0.07	48.85	1.39
Temperature	99.88	20.60	0.12	-
Net radiation	9.13	-	90.86	59.14

) [60]. The growing season NPP across Mongolia increased over 56.37% of vegetated areas (up to 13.46 g C m⁻²), and the significant increasing trends accounted for 7.23% of vegetated areas at a 95% confidence level. Areas with decreasing trends were observed in 40.94% areas (up to −14.54 g C m⁻²), and only 2.60% of them exhibited a significant decreasing trend. The dynamics of the growing season total NPP were found parallel to that of the annual NPP [10] as the growing season total NPP accounted for 93.14% of the annual total NPP.

3.3. The Impact of Drought on NPP

3.3.1. Characteristics of Drought

Based on the NDVI anomaly, the spatial distribution of drought across Mongolia from 2003 to 2018 was identified, and the drought years and their affected area percentage are illustrated in

Table 4. Drought-affected area in Mongolia over the period 2003–2018. The percentage number (unit: %) indicates the drought-affected area.

Year	Normal	Mild	Moderate	Severe	Extreme
2003	58.95	26.23	13.24	1.41	0.16
2004	39.99	33.59	20.20	6.12	0.09
2005	31.63	32.21	29.34	6.70	0.12
2006	33.01	40.19	21.58	5.08	0.14
2007	36.58	20.96	25.04	17.01	0.41
2008	45.48	31.25	18.47	4.67	0.14
2009	28.48	35.57	27.43	8.34	0.19
2010	49.62	32.87	15.15	2.24	0.12
2011	66.86	25.70	6.83	0.50	0.11
2012	72.81	22.06	4.64	0.37	0.13
2013	73.64	20.64	5.19	0.37	0.16
2014	71.27	21.43	6.89	0.30	0.11
2015	58.81	28.42	11.74	0.92	0.12
2016	74.69	14.97	7.94	2.30	0.11
2017	27.56	36.32	31.33	4.72	0.06
2018	64.59	21.39	11.98	1.92	0.11

. On average, drought-affected areas accounted for 47.87% every year. To explore the effects of drought on NPP, drought events with larger affected areas in 2004 (60.01%), 2005 (60.37%), 2006 (66.99%), 2007 (63.42%), 2009 (71.52%), and 2017 (72.44%) were further investigated by performing several statistical analyses.

We assessed the effects of extensive droughts on vegetation productivity across the Mongolian rangeland by comparing the spatial patterns of growing season total NPP anomaly and drought severity. The most influential drought occurred between May to September 2017 and affected 1.13 million km², covering 72.43% of the country’s land areas (Figure 5). During this drought event, severe and extreme drought areas covered approximately 4.71 and 0.06% and moderate and mild drought areas covered 31.33 and 36.32% of land areas of the entire country, respectively. The drought in 2004 mainly occurred in the eastern regions of Mongolia, and the drought areas covered nearly 0.94 million km², which accounts for about 60% of the total land area of the country. Most of the country experienced mild (35.5%) and moderate (27.4%) drought during the 2009 drought events. The drought from May to September in 2005 mainly affected the central and southeast regions. In 2006, drought was predominant in the southwest and the eastern regions. The country experienced mild (20.95%), moderate (25.04%), and severe (17%) drought during the year 2007. Our result was consistent with previous findings [61,62], thereby indicating that the NDVI anomaly can also be used for long-term drought monitoring in arid and semiarid regions.

3.3.2. NPP Variation during Drought Events

We used growing season total net primary productivity (NPP_gs; Pg C m⁻²year⁻¹) value to determine the NPP response to different drought conditions. On average, the growing season total NPP_gs accounted for 93.14% of the annual total NPP_a in Mongolian rangeland from 2003 to 2018 and exhibited a near-perfect positive correlation with the fluctuations of annual total NPP (r = 0.98, p < 0.001). Thus, small changes in NPP_gs caused by droughts can significantly affect annual NPP dynamics. Figure 6 shows the spatial patterns of percent NPP anomalies across Mongolia for growing seasons in the studied period. Our result showed a good match between widespread drought-affected areas and NPP anomalies. A total of 62.33% of the national area showed a decreasing NPP trend during the drought of May to September 2004. A 30–50% reduction of NPP was detected in severe drought areas, mainly in the eastern regions, and the NPP in a large fraction of the region decreased by 50–75% in the same areas. A decrease of NPP was observed over 60.55% of the total vegetated area during the drought period of May to September 2005. The country-wide NPP reduced by 10–30%, while a 30–50% decrease in NPP areas was observed for moderate and severe drought areas in Mongolia (Figure 6b). From May to September 2006, a 10–30% decline in NPP was found in mild and moderate drought conditions, while the NPP of severe drought areas, which were scattered across the country, decreased by 30–40% (Figure 6c). About 84.85% of vegetated pixels showed a decline in NPP to a different extent in drought areas in Mongolia from May to September 2007, while vegetation not affected by drought indicated an NPP increment of 15.15% (Figure 6d). A considerable decline in NPP (by 10 to 50%) was observed in mild to severe drought areas from May to September 2009, which constitutes 76.43% of the vegetated lands. Approximately 87.75% of the entire vegetated lands showed a 10–50% decrease of NPP in mild to severe drought areas during May to September 2017, although the NPP reduced by more than 50% in some areas (Figure 6f).

3.3.3. Relationship between NPP and Drought

We ascertained the decrease in NPP to different levels during the drought events. Further, we focused on the relationship between NPP and the droughts expressed as NDVI anomalies. Figure 7 shows the spatial correlation between the growing season mean NDVI anomaly and total NPP in Mongolia during the period 2003–2018. The majority of areas in Mongolia suffered NPP loss due to drought (indicated through positive correlations), accounting for 92.95% of the vegetation. Negative correlations were mainly distributed throughout the forest and high mountain regions (covering 7.05% of the total area).

The temporal variations of NPP and droughts were also explored by Pearson’s correlation analysis between NPP and NDVI anomaly concerning accumulated and average values in the growing season for 2003–2018 (Figure 8). A slightly increasing trend in both NDVI anomaly and BEPS-derived NPP was observed across Mongolia. The calculated correlation coefficient between the growing season total NPP and mean NDVI anomaly was 0.61 (p < 0.001, n = 16). Meanwhile, the Pearson’s correlation coefficient of the NPP and the drought-affected area was −0.55 (r = −0.55, p < 0.001, n = 16), suggesting that drought was one of the factors for the decline in terrestrial NPP in Mongolia. We also quantified the contribution of drought to NPP in different vegetation types. Figure 8 illustrates the influence of droughts on vegetation productivity and accumulated NPP over different vegetation types during 2003–2018. The droughts reduced the grassland NPP in Mongolia (r = 0.64, p < 0.001) in the growing seasons of 2007, 2009, 2010, 2015, and 2017. For the forest, the correlation between NPP and NDVI anomaly was lowest (r = 0.39, p < 0.001) because trees have deeper roots and have access to groundwater. However, the droughts in the growing seasons of 2003, 2009, 2010, 2015, 2016, and 2017 reduced the forest NPP (Figure 8). The droughts reduced the growing seasonal NPP for the wetland (r = 0.70, p < 0.001) in 2008, 2009, and 2012; for sparse vegetation (r = 0.74, p < 0.001) in 2004, 2005, 2007, 2009, 2015, and 2017; and for shrubland vegetation (r = 0.52, p < 0.001) across Mongolia in 2009, 2015, and 2017 (Figure 8).

4. Discussion

4.1. Impact of Climatic Variables on NPP Trend

During the 16-year study period, the terrestrial NPP in Mongolia increased from 0.16 Pg C year⁻¹ in 2003 to 0.18 Pg C year⁻¹ in 2018, at the rate of 0.00053 Pg C year⁻¹ (p-value = 0.065). The increasing NPP trend in the region indicates that Mongolia plays a crucial role in the global carbon budget. The variability of growing season NPP was probably caused by the combined effects of climatic variables. In order to investigate the influence of the main climatic parameters, namely temperature, precipitation, and radiation, on the inter-annual variability of NPP, we calculated the trends of these parameters from 2003 to 2018 using daily reanalysis data from ECMWF ERA-Interim. The results based on the Mann–Kendall test showed that 48.85% of the country experienced a decreasing trend in precipitation (up to −2.83 mm per year), while 51.15% of the area indicated an increasing trend (up to 2.33 mm per year) (Figure A1a). In the regions, increased precipitation may inhibit plant growth by decreasing solar radiation. While analyzing the growing season temperature change between 2003 and 2018, the increased rate (up to 0.18 °C) accounted for 99.88% of the country (Figure A1b). In the case of net radiation, about 90.86% of the areas exhibited a decreasing trend. Of these, 59.14% were statistically significant, and they mainly occurred in the southern regions of the country (Table 3). The increasing precipitation and temperature trends resonated with the overall trends of average annual temperature and precipitation worldwide from 2000 to 2014 [63]. Typically, the combination of warmer temperature, higher precipitation, and elevated CO₂ concentration during the early growing season possibly increased NPP, partly by lengthening the growing season [64]. However, previous studies have indicated that drought is the main primary influential factor for NPP in Mongolia [65].

In addition, we investigated land use change based on the ESA land cover map from 2003 to 2018, which showed that barren land had declined over the study period. Earlier, Wang et al. [66] reported that barren lands were being restored to non-degraded lands. A decline in barren land is another possible reason for an increase in total NPP in Mongolia. Therefore, our study suggests that further studies should focus on the potential effects of drought in a more historical manner by considering attributes such as overgrazing, urbanization, deforestation, land use change, and wildfire, which can limit the NPP variation. However, this study provides an insight into how droughts affect the terrestrial NPP in Mongolia, a region that is critical to global climate change.

4.2. NPP Response to Drought Events

Mongolia experienced six major drought events in 2004, 2005, 2006, 2007, 2009, and 2017 during which more than 60% of the total area suffered from drought, resulting in a decline in NPP to different extents. For the temporal dimension, there was a slight increasing trend in both NDVI anomaly and BEPS-derived NPP. The growing season NPP variance was moderately (r = 0.61 p < 0.001) correlated with the variability in wet and dry conditions as expressed by the drought index NDVI anomaly at the country scale, which is consistent with the finding from previous studies [67,68].

Spatially, the growing season NDVI anomaly was positively correlated with total NPP for 92.95% of the land, and this association was more pronounced in the western and southeastern regions. This suggests that the growing season productivity dynamics are parallel with the drought conditions in these regions. This also agrees well with the findings of a recent study that used SPEI as an indicator of drought [37]. However, for forested and mountain areas (which constitute 7.05% of the country’s land area), the growing season mean NDVI anomaly was negatively correlated with the total NPP. These regions experience a low temperature and relatively high precipitation and do not suffer from drought conditions.

Droughts have significantly affected the NPP, which was also ascertained in previous studies on a regional to a global scale by Wang et al. [20], Zhao and Running [7], and Peng et al. [22]. A recent study based on the MODIS NPP product and a drought index (SPEI anomalies) found that the NPP decreased by 50% in Mongolia at a country level from 2000 to 2014 [22]. However, we found a decreasing trend in NPP during the drought events, and NPP anomalies were the highest when drought severity reached its peak (decreased by 75%). This discrepancy is possibly because of differences between the regional and global scales. As a representative of sparse vegetation, reduced NPP in the eastern regions by severe drought was approximately 10–75% during those events. For mild and moderate drought, a 10–50% decrease in NPP was observed, and they even reduced 50–75% in some parts in the growing season of 2017 (Figure 6).

More specifically, the impacts of drought on vegetation productivity were intertwined with the severity of drought events and the types of vegetation. The severe droughts had a larger impact on vegetation growth than moderate and mild droughts. This finding agrees with that of Chen et al. [69], who found that NPP decreased significantly with the increase in drought severity. Finally, we found that the correlation between the annual total NPP and NDVI anomaly varied across vegetation types, from 0.39 in the forest to 0.79 in sparse vegetation (Figure 8). This means that the NPP of sparse vegetation is more susceptible to drought conditions than other types of vegetation. This was also consistent with a similar previous study, which showed that NPP of sparse vegetation was highly sensitive to temperature, whereas forest NPP showed a weak response to precipitation [10]. Additionally, it was indicated that the NPP value was highest in the wettest years, while the lowest NPP value was witnessed in the drought years (Figure 6). Thus, drought events that typically lasted for 5 months have largely impacted NPP in Mongolia.

5. Conclusions

In this study, the BEPS-derived NPP and MODIS-derived NDVI anomaly for 2003–2018 were used to examine the response of NPP to drought events. We affirmed that NPP and NDVI anomalies are good indicators of drought and its impacts on vegetation productivity. The following provides a summary of the main findings of this study:

• The study indicated that the BEPS model could predict terrestrial NPP in Mongolia. The growing season NPP was increased in 56.37% of the vegetated lands, whereas 43.63% of areas showed a decreasing trend from 2003 to 2018. Overall, the variation of NPP exhibited a slightly increasing trend over the entire study period.
• Whilst the NPP trend was partly influenced by various climatic variables, drought events played a key role in NPP variations in the Mongolian rangeland. A 10–50% drop in NPP was observed for mild and moderately drought-affected areas of the county, while the NPP in severe drought areas, mainly in the eastern regions, declined by approximately 10–75% during the six major drought events.
• The response of the NPP to drought conditions for different vegetation types was diverse due to variation in their drought resistance ability. The forests were generally more drought resistant than shrubland, while sparse vegetation exhibited extreme sensitivity to drought followed by wetland and grassland.
• The model and data used in this study identified the major drought events in Mongolia and highlighted vegetation types that would most likely be sensitive to future drought events. Therefore, our results suggest that this methodology could be applied to assess the impacts of drought on NPP for arid and semiarid regions, which is the basis for drought preparedness strategies.