Identifying the runoff variation in the Naryn River Basin under multiple climate and land-use change scenarios

Hydrological processes are intricate, involving precipitation, interception, infiltration and many other factors (Zhang et al. 2014). Climate and land-use changes are considered to be the main two drivers affecting hydrological processes, which bring great challenges to the accurate prediction of future runoff amount (Hung et al. 2020; Zango et al. 2021). In detail, climate change is dominated by a continuous increment of temperature and nonstationary fluctuation of precipitation, direct/indirectly modifying local evaporation and precipitation (Ahmadi et al. 2019; Yao et al. 2019). Land-use change (e.g., arable land expansion) affects surface and subsurface flows by influencing infiltration rate, evapotranspiration, albedo and canopy interception (Shrestha et al. 2018). Consequently, runoff processes would become more complicated by the individual and combined impacts of changing climate and land-use (Sharma et al. 2019; Yu et al. 2020). Therefore, assessing the impacts of climate and land-use changes on hydrological processes is crucial to determine the potential scale of these changes and enhance the precision of runoff simulation and prediction.

Over the past decades, the impact of climate change on runoff has been discussed worldwide by coupling the global climate model (GCM) and the hydrological model. Climate change has significant impacts on runoff and exhibits different characteristics in various regions (Bhatta et al. 2019), such as the Guajoyo River Basin (Blanco-Gómez et al. 2019), Baitarani River Basin (Padhiary et al. 2019), Southern Quebec (Berthot et al. 2021) and upper Godavari River Basin (Sanyal et al. 2021). These studies suggest that the GCM contributes to simulating climate change. Owing to the complexity of the GCM in terms of model structure and parameter estimation, the impact of uncertainty cannot be ignored when using a single GCM (Zhuang et al. 2018; Nie et al. 2021). Conversely, the use of multiple GCMs cannot only provide more accurate climate projections but also reduce errors between simulation and reality. A number of studies have proved that a variety of GCMs can effectively deal with uncertainty in climate change projection caused by heterogeneity of structures and parameters of the GCM as well as improve climate simulation and prediction precision (Tegegne et al. 2019; Mei et al. 2021).

Land-use change is another important factor influencing the filtration and routine processes of runoff. Degradation of ecological land cover (e.g., forest, pasture and cropland) would lead to runoff volume changes (Vicente-Serrano et al. 2019). Rational prediction of future land-use change would be useful for runoff simulation. For instance, Newman et al. (2018) employed the land transformation model to forecast future land-use change in Florida, where three future land-uses were developed by increasing the built-up area and the associated impact on flood risk was analyzed. Ni et al. (2021) used the linear regression method for predicting the land-use of the Big Sunflower River Watershed, where the map in 2025 was forced into the Soil and Water Assessment Tool (SWAT) model to explore its effect on runoff generation. Getachew et al. (2021) combined Cellular Automata–Markov (CA–Markov) and SWAT models to assess the changes of runoff and lateral flow based on the future land-use prediction in the Lake Tana Basin. In general, CA–Markov is a popular tool for land-use prediction due to its computational effectiveness in randomly presenting the current state depending only on the nearest previous period (Yang et al. 2019). It is also easily conducted as it needs limited constraints (e.g., slope and water body), which is prior to regression methods, especially in the data-scarce region.

Analysis of the combined impacts of both climate and land-use changes on runoff generation has bloomed in recent years. Shrestha et al. (2018) quantified the integrated effect of changes to climate and land-use in runoff in the Songkhram River, which showed that combined factors created a net effect of decreased future runoff (16% on average). Chen et al. (2019) explored potential changes of runoff in the Loess Plateau region using the SWAT model, which revealed that climate and land-use changes accounted for 36.48 and 63.52% of the runoff reduction, respectively. Dosdogru et al. (2020) evaluated the runoff in urbanizing watersheds in Alabama; results found that more frequent flood and drought events would occur considering the changes of climate and land-use. All the existing studies suggest that runoff varies spatially under different climate and land-use scenarios. However, limited future climate data and simulated land-use scenarios could bring uncertainties in runoff impact assessment. More in-depth studies are needed to disclose the impacts of climate and land-use changes through incorporating multiple GCMs, CA–Markov and SWAT within a general framework, especially for the Naryn River Basin (an important subregion of the Syr Darya River in Central Asia), where it is still unclear how the runoff responds to both climate and land-use changes (Liu et al. 2019). Due to its cold-arid climate, the Naryn River Basin is very sensitive to global warming and precipitation fluctuation (Jia et al. 2020). The lower reach of the Naryn River Basin is also an important agricultural development region. The expansion of arable lands would destroy the land surfaces, thus influencing the runoff generation processes. Furthermore, the Naryn River Basin is a mountainous region with scarce meteorological gauged stations (Gan et al. 2015). Multiple climate and land-use scenarios would be more flexible for providing plausible information on runoff prediction in this basin.

Therefore, this study aims to develop a multiple scenario-based ensemble prediction (MSEP) method for analyzing the impact of climate and land-use changes on runoff of the Naryn River Basin. MSEP consists of multiple downscaled GCMs under two representative concentration pathways (RCPs), CA–Markov method and a distributed hydrological model (i.e., SWAT). The novelty and contribution of the MSEP framework lie in (i) analyzing combined impacts of climate and land-use changes under varying conditions, (ii) reflecting the uncertainties derived from future diverse climate and land-use scenarios under ensemble prediction and (iii) identifying future runoff trends and supporting regional water resources management. This is also the first application of the MSEP method to the Naryn River Basin to simulate and predict its runoff variation and disclose the impacts of climate and land-use changes.

The Naryn River Basin is located in the upper reaches of the Syr Darya River in Central Asia, with an area of about 59,100 km² (Figure 1). The river is about 807 km long and has an annual flow of 13.7 billion m³. It provides 34% of the total flows of the Syr Darya River, which plays a significant role in downstream human activities (especially irrigation and hydropower) (Radchenko et al. 2017). The average temperature is −18 °C in January and 10 °C in July. The annual precipitation is spatiotemporally heterogeneous varying regionally between 280 and 450 mm. Precipitation is mainly concentrated in the mountainous areas (Kriegel et al. 2013). Annual precipitation is primarily contributed by that in spring and early summer. The maximum runoff occurs in June, which is mainly discharged by both melted snow and precipitation (Gan et al. 2015). Climate change and high demand for irrigation water have gradually changed the underlying surface conditions throughout the region (Wang et al. 2020). The relationship between water and humans in the Naryn River Basin is multifaceted and influenced by water volume and reliance degrees of livelihoods on water (Liu et al. 2019).

In this study, the MSEP method was developed through integrating multi-GCMs, CA–Markov model and SWAT model into a framework (as shown in Figure 1). Each component of this framework has its unique contribution to runoff assessment in the Naryn River Basin. In detail, eight GCMs under two RCPs are employed to extract future climatic conditions (e.g., temperature and precipitation) for reflecting the uncertainties in the climate. The CA–Markov model can generate future potential land-use change scenarios, where six land-use scenarios related to different agricultural developments are produced. The SWAT model is applied to predict the runoff under the total 96 combinations (8 × 2 × 6) of climate and land-use change scenarios during the period of 2020–2099, which can be divided into three slices, i.e., near-term (the 2030s: 2020–2039), mid-term (the 2050s: 2040–2069) and long-term (the 2080s: 2070–2099). 96 scenarios are defined as GiSj, where i ranges from 1 to 8 representing the 16 climate scenarios and j ranges from 1 to 6 representing the six land-use maps.

Two sets of data were required in this study: (i) remote sensing data and (ii) hydrometeorological data. In detail, the remote sensing data includes digital elevation model (DEM), land-use maps and soil types. The DEM data are obtained from the CGIAR Consortium for Spatial Information (CGIAR-CSI, https://cgiarcsi.community) with a resolution of 90 m×90 m. The soil data are abstracted from the Harmonized World Soil Database (HWSD) provided by the Food and Agriculture Organization of the United Nations. The 1-km resolution land-use in the years 1990, 2000 and 2010 are obtained from the Data Sharing and Service Portal (http://data.casearth.cn).

The hydrometeorological data for this study includes daily flow, historical and future precipitation, maximum temperature and minimum temperature. Due to data limitations, the daily runoff data from 1951 to 1974 are derived from the Kurgan hydrological station (http://www.bafg.de/GRDC). In this study, this period is used as the baseline period. Daily gridded precipitation from 1951 to 1974 is extracted from the dataset of Asian Precipitation-Highly Resolved Observational Data Integration Towards Evaluation (APHRODITE, http://aphrodite.st.hirosaki-u.ac.jp). APHRODITE is a long-term continental daily precipitation product. It is based on an intensive rain gauge network from thousands of stations in Asia, with a spatial resolution of 0.25° and the time period from 1951 to 2007. It has been proved that this dataset is better than other models and can be used directly (Ghimire et al. 2015). Historical daily temperature data from 1951 to 1974 are downloaded from the website of Princeton Global Meteorological Forcing Dataset (PGMFD, http://hydrology.princeton.edu). PGMFD is established by combining a series of datasets based on global observations with the National Center for Environmental Prediction-National Center for Atmospheric Research (NCEP-NCAR) reanalysis. The spatial resolution of PGMFD is 0.5° and the dataset reflects the real temperature exactly (Ahmed et al. 2019). For keeping consistent with precipitation data subsequently, the daily temperatures are bilinearly interpolated into 0.25°.

The future precipitation and temperature are abstracted from eight GCMs under RCPs of 4.5 and 8.5 (from 2020 to 2099). These GCMs are CCSM4 (defined as G1), CNRM-CM5 (G2), MPI-ESM-MR (G3), IPSL-CM5A-MR (G4), MRI-CGCM3 (G5), MIROC-ESM (G6), GFDL-ESM2M (G7) and BUN-ESM (G8) from different institutes with various initial and feedback mechanisms. These GCMs are downloaded from the NASA Earth Exchange Global Daily Downscaled Projections (NEX-GDDP) website (https://www.nccs.nasa.gov). NEX-GDDP is a dataset converted from daily maximum temperature, minimum temperature and precipitation data to a resolution of 0.25°×0.25° by means of bias-correction spatial disaggregation, and consistent with the resolution of historical data, which is available for the entire globe (Musie et al. 2020).

Figure 2 shows the projected annual mean temperature and precipitation over the Naryn River Basin under all GCMs and RCPs. Results also reveal that both the maximum temperature (T_max) and the minimum temperature (T_min) would experience increasing trends in the future. In detail, T_max would reach [10.4, 15.7] °C by the 2080s, with an increasing rate of [0.1, 0.6] °C/decade. T_min would increase to [−3.3, −1.8] °C, with a rate at [0.1, 0.5] °C/decade. The increasing T_max would lead to evaporation intensification, and the increasing T_min would lead to snow/ice melt acceleration. Precipitation would change from [329.1, 384.7] mm in the 2030s to [296.6, 443.8] mm in the 2080s under RCP4.5 and [301.6, 397.0] to [292.1, 479.4] mm under RCP8.5. By the 2080s, the mean precipitation under all GCMs would increase by 12.6% (RCP4.5) and 16.6% (RCP8.5) more than values under the baseline (i.e., 329 mm). The increasing amount may lead to precipitation events more frequently and complicatedly. Results indicate that different GCMs and RCPs would generate changed projections of temperature and precipitation. These differences are associated with the climate model structure, inherent model parameterization and emission policies (Stephens et al. 2018). Therefore, results based on multiple GCMs can improve the accuracy of climate simulation and improve the credibility of prediction result.

The results of the CA–Markov model indicate that the κ between the observed and simulated land-use maps in 2010 is 0.80, revealing a good fitness and suitability for land-use predictions in the Naryn River Basin based on the CA–Markov model. Through driving the validated CA–Markov model with the defined scenarios, land-use maps from the 2030s, the 2050s (as shown in Supplementary Material, Figure S1) and the 2080s (as shown in Figure 3) are obtained. The land areas under all scenarios are summarized in Table 3, which shows that under the natural development scenario (i.e., S1), the proportion of pasture would increase from 56.8 to 65.3%, and the forest would increase from 5.2 to 7.9% from the 2030s to the 2080s. The natural expansion of pasture and forest would consume much bare land, which is consistent with that of the past decades. In the future, agricultural development scenarios (i.e., S2, S3 and S4 related to intensified human activities with different degrees) show that the maximum agricultural area would expand to 32.1% by the 2080s. Results related to eco-environment protection scenarios (i.e., S5 and S6) suggest that the forest area would increase to 19.7% by the 2080s, and the agricultural area would decline to 2.8%. The land-use patterns under different scenarios highly related to human activities would impact the runoff generation processes, such as the canopy, evapotranspiration and filtration.

The optimal 13 parameters calibrated by the SUFI-2 algorithm in SWAT are presented in Table 4. As shown in Figure 4, the values of NSE are 0.80 and 0.63, and the values of R² reach 0.80 and 0.64 in calibration and validation periods, respectively. This means that the simulated runoff sequence matches well with the observed, revealing the applicability of SWAT in the Naryn River Basin's runoff prediction (Moriasi et al. 2007; Zhang et al. 2016). Through introducing GCMs' outputs and land-use prediction values into the validated SWAT model, 96 runoffs under different climate and land-use scenarios were generated (as shown in Figure 5). Results show that the multi-year average flows would range between 454.1 m³/s (G4S6) and 614.5 m³/s (G7S4) under RCP4.5 and range between 410.6 m³/s (G4S6) and 649.8 m³/s (G5S4) under RCP 8.5. These differences are caused by the various climate and land-use scenarios, implying that large model dependence may exist in the results from individual GCM or land-use simulations. Besides, the intra-annual runoffs show a huge undulate due to the heterogeneousness and nonstationary precipitation and temperature. For example, under S1, the maximum annual runoff would be 1,044.8 m³/s under MPI-ESM-MR in 2052, and the minimum value would be 183.1 m³/s under IPSL-CM5A-MR in 2055. Results disclose that the adoption of multiple GCMs can reflect the uncertainties derived from different climate institutes due to various initial mechanisms and physical characteristics.

Figures 6 and 7 describe the combined impacts of climate and land-use changes on the peak, average and low flows in the 2030s, the 2050s and the 2080s. Results show that the peak, average and low flows would change across all scenarios, and such changes would be enlarged along with time. For example, by the 2080s, under RCP4.5, the peak flow would fluctuate between 1,115.9 m³/s (G4S6) and 1,648.9 m³/s (G2S4) with a deviation of 174.6 m³/s. The low flow would vary between 37.2 m³/s (G4S5) and 72.7 m³/ (G6S2) with a deviation of 12.3 m³/s. Obviously, peak flow is sensitive to climate and land-use changes. This is because peak flow would occur accompanied by plentiful precipitation. Due to the steep slope of local ground surface and river bed, the changes in peak precipitation would easily lead to variations in the confluences. Besides, the maximum peak flow is generated under the agricultural development scenario. Excessive development of agriculture and reduction of forests would cause a decline in the soil's ability to retain water (Hussain et al. 2021). Consequently, more water would be gathered and routed to the nearby channels, leading to an increment in peak flow and a potential in floods. Comparably, the low peak would be generated under the combination of eco-environment protection scenario and IPSL-CM5A-MR. This is because IPSL-CM5A-MR corresponds to the lowest precipitation and the highest temperature. Besides, eco-environment protection would increase the soil water retention rate and the infiltration of runoff.

Figure 8 shows the probabilistic distributions of peak, average and low flows in the 2030s, the 2050s and the 2080s. The flow difference between RCP 4.5 and RCP 8.5 would be intensified along with time. For example, in the 2030s, the low flow is around 41.8 m³/s under RCPs; however, in the 2080s, the low flows would climb to 52.2 m³/s (RCP4.5) and 69.8 m³/s (RCP8.5). Such a trend is similar to those of temperature and precipitation, illustrating that uncertainty caused by climate change would be amplified along with time. Compared with the baseline period (1951–1974), future peak and average flows would increase. For instance, in the 2080s, under RCP4.5, the 75th percentile of future peak flows would be 20.1–53.1% higher than that in the baseline period (1,058.1 m³/s), and the 75th percentile of future average flows would be 23.7–58.6% higher than that in the baseline period (397.7 m³/s). Conversely, the 75th percentile of projected low flows would be 26.3–47.2% lower than that in the baseline period (83.3 m³/s). This is because peak flow mainly occurs in early summer discharged by both heavy precipitation and snowmelt, which would be promoted by the argument of temperature and precipitation. However, low flows are mainly generated in winter discharged by scarce precipitation, which may be underestimated. The magnitude of peak flow increase would be larger than that of low flow, leading to a raised average flow. Generally, the variation of runoff in the future is seriously differentiated, so it is necessary to pay attention to both floods and droughts.

Figure 9 describes the monthly runoffs under all scenarios. Results show that the runoff during December and May of next year (i.e., winter and spring) would gradually increase in the 21st century. For example, under RCP8.5, the average flow in December would be 104.5 m³/s in the 2030s, which would increase to 138.4 m³/s in the 2080s. This may be related to the uptrend of the future minimum temperature and corresponded intensive snowmelt-discharged runoff. Besides, the runoffs of June–August (i.e., summer) would vary between 147.7 and 3,381.0 m³/s with a large distinction under all scenarios. Results reveal that runoff in summer is more sensitive to climate and land-use changes. In summer, the runoff generation processes are complicated with both plentiful precipitation and melted snow/ice; the precipitation projections by multiple GCMs would be plagued with more complicated variation and uncertainty. Therefore, attention should be paid to such uncertainties and complexities in order to provide flexible information for water resources management.

To further explore how climate change impacts monthly runoff, the monthly runoffs generated by all climate scenarios under S1 were analyzed via comparing with baseline runoffs. As shown in Figure 10(a), results show that the runoff in August would be the most sensitive to GCMs, which would increase by 28.4–155.3% compared to the value in the baseline period (i.e., 531.5 m³/s). March would be the least sensitive to GCMs with a changing rate between −37.9 and 32.5% to the baseline. This may be because the summer temperature exceeds the snowmelt base temperature so that the summer runoff is affected by a number of factors such as precipitation, snowmelt and evapotranspiration. The early spring temperature is generally lower than the snowfall temperature, and the precipitation is easy to form snow, so that the runoff is less affected. From the perspective of seasonal scale, future climate change would lead to a lower amount of runoff in winter and higher amount of runoff in spring, summer (especially) and fall compared to the baseline (Figure 10(b)). Although precipitation decreases only in summer, runoff increases in summer. This suggests that variations in runoff are not only dependent on precipitation, and changes in vegetation and ET can also influence the runoff.

Figure 11 depicts the impacts of land-use change on the temporal distribution of monthly runoffs. By comparing monthly runoffs with that of S1 (i.e., natural development scenario), results show that runoff change is large (i.e., varying from −31.7 to 28.6 m³/s from March to October) under the five scenarios. From November to February of the following year, the monthly runoff change is small (i.e., varying from −5.7 to 1.5 m³/s). Results reveal that monthly runoffs in spring, summer and fall are more easily impacted by land-use changes. This is possibly caused by the combined effects of intensive agricultural activities (e.g., reclamation and tillage), sufficient precipitation and melted snow/ice in spring and fall. Besides, compared to the base scenario (S1), arable land expansion under S4 would lead to 29.0 m³/s of runoff increment, and afforestation under S6 would result in runoff declination (5.3 m³/s). In fact, land-use patterns are highly related to a hydrological process parameter (i.e., CN) that is positive to surface runoff. For example, when the land-use patterns switch from S1 to S4 and S6, the area-weighted average values of CN would change from 63.55 to 63.83 and 63.22, respectively. As a result, these phenomena may also be caused by the agricultural operations (e.g., plowing and irrigating) during this period that may lead to the soil surface layer crusted and the infiltration rate reduced.

From the perspective of the annual scale, intensive agricultural development would increase runoff (corresponding to increased water supply). For example, by the 2080s, the generated annual average runoff under S4 would be [1.8, 6.8] × 10⁹ m³ higher than that in the baseline period. However, the increasing water supply would be far below the soaring water consumption by agricultural irrigation. According to the Kyrgyzstan Statistical Yearbook (2014), the agricultural water consumption quota is about 12.0 × 10³ m³/ha; this means that by the 2080s, the total agricultural water consumption would rise by 19.8 × 10⁹ m³ under S4; and there would be [13.0, 18.0] × 10⁹ m³ of water resources deficit (about averagely 10.9 × 10³ m³/ ha). Conversely, S6 would lead to the highest water resources surplus, because agricultural water consumption would be 0.9 × 10⁹ m³ lower than that in the baseline period and water supply would be on average 4.3 × 10⁹ m³ higher. Results indicate that rapid agriculture development would cause intensified water shortage problems in the study basin. Results provide new insights into the linkages between agricultural development and water resources at regional and local scales for the effective management of water resources.

In this study, an MSEP method has been developed for assessing the response of future runoff to the climate and land-use changes. The MSEP method has been applied to a data-scarce region (i.e., the Naryn River Basin in Central Asia). Some major findings can be summarized as follows: (i) in the 21st century, the maximum and minimum temperatures (T_max and T_min) would increase by [0.1, 0.6] and [0.1, 0.5] °C/decade, respectively. Precipitation would increase by an average of [12.6, 16.6]% compared with the baseline; (ii) future multi-year average runoff would range between 410.6 and 649.8 m³/s under all scenarios, and climate change plays a leading role in runoff variation; (iii) by the 2080s, peak flow would on average be 36.6% higher than that in the baseline period (1,058.0 m³/s), and low flow would on average be 36.8% lower than that the baseline period (83.3 m³/s), indicating that extreme hydrological events (e.g., flooding and drought) would get frequent in the future; (iv) summer runoff (especially in August) is the most sensitive to climate and land-use changes, which would be about 155.3% higher than the baseline period; (v) expansion of arable land would increase the CN of water yield, which would lead to the increased runoff; and (vi) every additional hectare of arable land would increase the water deficit by an average of 10.9×10³ m³. The results obtained can effectually identify future runoff change trends. To keep local water and environmental security, rational planning of agricultural development and improvement of water use efficiency are crucial to narrowing the gap between water and water supply.

This study is the first attempt to apply the MSEP method to analyze the sensitivity of runoff to the combined effects of climate and land-use changes in the Naryn River Basin. However, this method also has some potential extensions. First, the study basin (i.e., the Naryn River Basin) is located in the mountainous region of Central Asia, and little observation and lack data could increase the uncertainty of the model performance and result. Therefore, in order to improve the simulation accuracy of the model, the runoff sequence can be increased by means of remote sensing inversion in future work. Secondly, the glacier hydrology module can be added to the SWAT code to improve the accuracy of glacier hydrology in alpine basins. These methods should be considered in future research to improve the accuracy of runoff prediction and reduce the total uncertainty of hydrological modeling.

This research was supported by the Strategic Priority Research Program of the Chinese Academy of Sciences (No. XDA20060302). The authors are grateful to the editors and the anonymous reviewers for their insightful comments and suggestions.

All relevant data are available from an online repository or repositories.

DEM: https://cgiarcsi.community

Land-use: http://data.casearth.cn

Daily runoff: http://www.bafg.de/GRDC

Daily precipitation: http://aphrodite.st.hirosaki-u.ac.jp

Daily temperature: http://hydrology.princeton.edu

GCM: https://www.nccs.nasa.gov