ABSTRACT
This study has used Coupled Model Intercomparison Project Phase 5 (CMIP5) and Coupled Model Intercomparison Project Phase 6 (CMIP6). Hence, the runoff simulation was done in near-future period (2030–2050) scenarios by applying climate change conditions for HadGEM2-ES model under three representative concentration pathways RCP2.6, 4.5 and 8.5 scenarios and for HadGEM3-GC31-LL model under SSP1-2.6, SSP2-4.5 and SSP5-8.5 scenarios. Examining the climatic precipitation variables and minimum and maximum temperature under climate change conditions showed a temperature increase of 1.51–2.91 °C for all models and scenarios and a precipitation decrease of 0.05–11.15% for most of them, and the SWAT model simulation showed a runoff decrease in all four stations under SSP scenarios and a runoff increase in three stations under RCP scenarios. Since the climate data for SSP scenarios have become available only recently, results of this study predict that the overall future flow will vary in the −5 to 28% range, resulting in a 5–35% decrease in the flow and, hence, a decrease in the inflow to the dam reservoir. Based on the results, there is a possibility of a 5–30% reduction of the runoff entering the dam reservoir.
HIGHLIGHTS
Three RCP scenarios from the CMIP5 report model and three SSP scenarios from the CMIP6 report model were used in predicting the future climate of Karaj Dam.
It is expected that the weather variables will decrease and increase in Karaj Dam.
Considering the uncertainties, SWAT provided a relatively good to good performance in simulating the runoff of the Karaj dam basin.
Climate change will lead to a decrease in the future runoff at the entrance of the Karaj Dam reservoir, and as a result, will decrease the volume of the dam.
INTRODUCTION
Climate change is a key water-resource management challenge because it has created problems in areas that suffer from water shortage; one such area is the Upper Neuse Basin in Carolina, America, which serves as an important water source for a large and growing region (Maghami et al. 2019). The global average greenhouse gas concentration has increased rapidly since the 1900s, leading to changes in the characteristics of climate variables and the occurrence of more extreme events (Almazroui et al. 2020). This phenomenon is catalyzed by both natural and anthropogenic forces (Karl & Trenberth 2003; NASA Earth Observatory 2023). Considering the increased greenhouse gas concentration and the earth's warming trend, climate change is expected to affect water resources and hydrology (Kim et al. 2021). Iran has not been exempted from this process. The increase in temperature and decrease in rainfall under the influence of climate change has led to frequent extreme weather conditions such as drought, especially in the Karaj Dam basin (Zohrabi et al. 2014; Naseri et al. 2020). Therefore, it is necessary to evaluate the consequences of climate change on water resources (Meraj et al. 2022). The basin hydrology is affected by the increased temperature and rainfall distribution variations, resulting in changes in the water amount and availability (Aryal et al. 2019). Severe climatic events and the related variations can, in general, have high effects on society and the ecosystem and cause great economic losses every year. In this regard, the Intergovernmental Panel on Climate Change (IPCC) has warned, in its 5th report (AR5), about such increased natural hazards as drought and flood that may occur due to climate change (Houshmand Kouchi et al. 2019). Therefore, highly accurate predictions of, especially, the future changes in the hydrological behavior are necessary for water resource policies and planning to provide accurate information for adaptation and relief (Wang et al. 2020). So far, different climate models have been developed as effective tools for past and future climate simulations. CMIP5 includes more than 40 models that use a new set of RCP emission scenarios to produce valuable climate information for policymakers and the scientific community (Chen & Sun 2015), and CMIP6 continues the CMIP's previous evolution pattern and adaptation features through newly organized global climate modeling scenarios designed to understand various climate mechanisms (Eyring et al. 2016). In general, CMIP6 models are more accurate, involve improved dynamical processes and use joint SSP/RCP socio-economic emission scenarios to simulate future climate changes (O'Neill et al. 2016). Identifying the hydrological conditions, predicting future trends and planning the existing water resources requires a correct understanding of the climatic changes on a spatial-temporal scale (Houshmand Kouchi et al. 2019). Although general circulation models (GCMs) are the most up-to-date tools for predicting climate changes and their effects at the global level (Sachindra et al. 2018), their low spatial resolution (250–600 km) makes them unreliable for research works (Taylor et al. 2012); however, the ‘downscaling’ method can improve the resolution of GCMs (Wilby & Wigley 1997; Maraun et al. 2010). In a study investigating the impacts of climate change on hydrological and meteorological parameters in the Naroud watershed (part of the Hyrcanian forests of the Caspian Sea region) in northern Iran, the researchers utilized outputs from seven global climate models (GCMs) from CMIP6 (including Tmin, Tmax and precipitation) under two scenarios, SSP2-4.5 and SSP5-8.5. The results indicated that precipitation and maximum and minimum temperatures are projected to increase in three future time periods, leading to an increase in runoff under the influence of these climate factors (Lotfirad et al. 2022). In order to better understand climate change, there is a need to examine the Standardized Precipitation Index (SPI) and Standard Precipitation Evapotranspiration Index (SPEI) to calculate short-term, medium-term and long-term meteorological droughts of the base and future periods. Therefore, in research to investigate the impact of climate change on the water level and shrinkage of Lake Urmia, they selected the top 10 GCMs from among the 23 GCMs CMIP5. Based on the results, the temperature will increase in all seasons of the future period, and according to the RCP4.5 scenario, the amount of precipitation will decrease in spring and autumn and increase in summer and winter. In addition, the RCP8.5 scenario reduces precipitation in spring, autumn and winter, while increasing it in summer (Radmanesha et al. 2022). Applying management policies for sustainable development and optimum use of surface water resources requires a better understanding of the environment and the climatic conditions of the region because preserving and storing such resources is vital in arid and semi-arid climates. Today, the increased greenhouse gases, reduced rainfalls, overuse of water and environmental issues have all doubled the food and water security (Ashofteh et al. 2015), and the existing methods use different exploitation and optimization scenarios, based on the dam reservoir runoff inflow, to apply management policies, optimize consumption and meet the water requirements needed for drinking, industry and agriculture (Khanjari Sadati et al. 2014). In many regions that have problems supplying water for drinking and economic and environmental development, human health is at risk, and productivity and maintaining a clean environment and healthy ecosystem are quite far-reaching. According to Eyring et al. (2016), climate change in regions with uncertain future behavior will result in fast population growth, urbanization and social/economic globalization. Today, with the advancement of programming and computer science to understand the hydrological behavior of the basin, many machine learning algorithms have been widely used. For example, conventional models include artificial neural networks, least squares support vector machines (LSSVMs), K-nearest neighbor (KNN), M5 model trees, random forests, multiple adaptive regression splines and multivariate nonlinear regression (Anaraki et al. 2023). Hydrological models are essential tools used to understand basin-scale natural processes and many of them have been developed for water resource studies and related modeling. They are increasingly used in flow quantity/quality analyses, flood forecasting, reservoir/water distribution systems, underground water-resource development, surface/underground water preservation, water utilization, climate/land studies, environment and a wide range of water management activities (Wurbs 1998; Singh & Gosain 2013). SWAT, which is a conceptual, semi-distributive, soil and water assessment model developed by a group in the United States Department of Agriculture (USDA) (Arnold et al. 1998; Arnold & Fohrer 2005), has been used, in recent years, for water resource planning and management under climate change conditions (Houshmand Kouchi et al. 2019; Kim et al. 2021; Rabezanahary Tanteliniaina et al. 2021). In a watershed study in Tabriz, Iran, the Salman Dam inflow variations were investigated by three GCM models – HadGem2-ES, IPSLCM5A-LR and NorESM1-M – and the SWAT hydrological model under the RCP2.6 and RCP8.5 climate scenarios. According to the results, there was a future increase in the frequency of heavy rains and a decrease in that of light rains compared to the past, as well as an increase in the basin's min/max temperatures, leading to an increase in runoff (Houshmand Kouchi et al. 2019).
In research, the runoff and climate changes of Task-e Bakhtaran watershed were examined and analyzed with two reports using SWAT, and CMIP5 and CMIP6 climate models under RCP2.6 and RCP8.5 scenarios (in GFDL-ESM2M and IPSL_CMA5_LR) and SSP1-2.8.5 and SSP1-2.8 SSPen scenarios (in GFDL-ESM4 and IPSL_CMA6_LR). Results of R2, NSH and RMSE indices revealed that the SWAT model calibration and validation ranges were, respectively, 0.70–0.99, 0.51–0.98 and 0.9–14.4 m3/s, indicating high adjustment accuracy, and validation of the model results, rainfall climatic variables and min/max temperature under climate change conditions indicated a temperature increase of 1.51–2.91 °C (for all models and scenarios), and a rainfall reduction of 0.15–0.5% (for most of them). SWAT simulation of the models and scenarios under climate change conditions showed a runoff reduction in all four stations, under SSP scenarios, and a runoff increase in three stations, under RCP scenarios (Fallah Kalaki et al. 2021). Since SWAT is suitable for simulating the climate change effects, it was used in Lebanon to simulate the Al-Kalb River flow. The model performance was studied for monthly flow in two stations under a 9-year calibration period (2003–2011) and a 4-year validation period (2012–2015) using different statistical indices, and the results indicated good accuracy in fitting the observed and simulated flows. The near-future (2021–2040) climate change predictions for RCP2.6, 4.5 and 8.5 showed that the River's average annual runoff would drop by 28–29%, and the end-of-century (2081–2100) predictions revealed that the drops would be, respectively, 23, 28, and 45% for the mentioned RCPs (Saade et al. 2021). In a systematic study based on the latest climate simulations, Wang et al. (2022) used SWAT and 14 GCMs under two common socio-economic pathways (SSP2-4.5 and SSP5-8.5) to investigate the climate change effects on runoff and floods in Fujiang River Basin (a major Yangtze River tributary in China). Results showed a warmer and wetter climate and, hence, an overall increase in average annual and monthly runoffs, as well as in the high and low monthly flows (respectively, Q05 and Q95) in 2021–2060 and 2061–2100 periods, with a long-term period more important than the near future, especially for SSP5-8.5. However, the changes predicted in the monthly runoff showed a high increase in the GCM, with the greatest value occurring mainly in the early rainy season.
Although various studies have been done to detect and attribute the changes of climate variables, until now, the detection and attribute of runoff as a hydrological variable that has not been affected by the two climate models CMIP5 and its comprehensive and more complete versions CMIP6 has not been done. Therefore, in this research, for the first time, it has been tried to detect and quantify the runoff simulated with the SWAT semi-distributive model with the aim of investigating the effects of climate change on the runoff upstream of the Karaj Dam basin. Also, in line with the management measures, it was done by examining the number of changes in the runoff entering the dam reservoir under the influence of climate change and taking into account the water needs of users downstream of the dam in order to reduce water consumption to adapt the conditions of the region to the future climate changes of the basin. Also, efforts were made to analyze climate changes and inflow and take measures to reduce its effects on the environment, which was necessary to identify these effects and their relationships in the basin and downstream of the dam.
MATERIALS AND METHODS
Study area
Inset (a) shows the distribution of weather and hydrological stations, main rivers and the reservoir across the Karaj Dam, overlaid on DEM. Inset (b) exhibits the land use/cover map of the Karaj Dam.
Inset (a) shows the distribution of weather and hydrological stations, main rivers and the reservoir across the Karaj Dam, overlaid on DEM. Inset (b) exhibits the land use/cover map of the Karaj Dam.
Controlling the amount of inflow to the dam reservoir and checking the volume of water stored in the Karaj Dam reservoir, especially in times of water shortage, is of great help in managing water crisis and water shortage (Hasan & Pradhanang 2017). Therefore, in this study, it was necessary to investigate the flow behavior under the influence of climate change. On the other hand, the simulation of mountain watersheds is associated with modeling challenges, and SWAT software was used in this research. One of the effective environmental factors in the hydrological simulation of the basin is the spatial changes of climatic parameters, including rainfall and temperature changes in altitudes, which make the modeling process difficult in the studied area, which can be solved by using SWAT software.
Data collection
This research has used: (1) a hydrological, semi-distributive SWAT model to simulate and investigate the flow behavior and rainfall effects upstream of the Karaj Dam, (2) a DEM (digital elevation model) map, with 12.5 m accuracy, for sub-basin dividing and extracting the waterway network, (3) Landsat 2018 satellite-image land-use map to prepare the basin's hydrological-response-group map from the soil (mainly loam and sand) map and (4) meteorological station data to estimate the runoff. The land use defined according to the SWAT model format is equal to the urban area which is about 1.5%, garden and agricultural land which is about 3%, water about 3%, bushland about 81.5%, barren land about 4%, pasture about 5% and rocky land about 2%. It has been introduced as general land use. Regarding the slope, the mountainous study area was classified into 0–15°, 15–30°, 30–60° and above 60° groups, and elevation classes were defined for each sub-basin so that the precipitation and snow factors were distinguishable. As this research was not aimed to evaluate varied land use effects, only one Landsat 8 image was selected because of its good quality and, finally, the data of the related meteorological and water stations were used to extract the runoff. It is worth mentioning that the SWAT model was implemented by the Arc GIS 10.3 Software, the total number of sub-basins was 45 and that of HRU was 383.
SWAT hydrological simulation
In this research, in order to investigate the flow behavior of the Karaj Dam watershed, a simulation process has been carried out for the period from 1990 to 2018, in which the period from 1990 to 1993 was considered as the model preparation period in order to perform the calibration and extraction process. The effective parameters on the flow from 1994 to 2012 were selected, and finally, they were considered for the purpose of validating the period from 2014 to 2018.
Sensitivity analysis of the SWAT model
In order for the proper performance of the model under the influence of specific watershed processes, for the accurate development of the hydrological simulation process, especially in areas with complex physiological characteristics, taking into account the parameters affecting the distribution of temperature and precipitation under the influence of elevation changes, along with snow accumulation and snowmelt in complex mountain watersheds are processes that can significantly affect model performance (Hussainzada & Lee 2021). Hence, this research has used the SWAT-CUP software along with the SUFI-2 algorithm to obtain the effective parameters and their effective ranges. Table 1 lists some effective flow parameters to check the SWAT results and do their sensitivity analyses. To calibrate the model, this algorithm generates, based on the Latin square sampling method, a large number of random sets of the desired parameters within the range defined for each of them. These sets are then placed in the model separately to find the related value of the objective function based on the desired variables; they are quite important in model outputs and the optimal selection of the calibration parameters. Among various methods used for sensitivity analyses, uncertainty analyses and calibration of the SWAT model, although GLUE, ParaSol, SUFI-2, MCMC and BIS are quite common, different researchers believe that SUFI-2 is the best (Yang et al. 2008). This study has used the NSH, R2 and RMSE indices (Equations (2)–(4)) to evaluate the SWAT model results. It should be noted that the period from 1994 to 2012 was chosen for calibration and from 2013 to 2018 for validation.
Parameters affecting the SWAT runoff in the SWAT-CUP model
Range of changes (provided by the SWAT-CUP model) . | Name of the parameters . | Parameter . | Row . | |
---|---|---|---|---|
−0/3 | 0/2 | Curve number | R__CN2.mgt | 1 |
0 | 1 | Coefficient α of groundwater | V__ALPHA_BF.gw | 2 |
230 | 500 | The delay time of water transfer from the last profile of the soil layer to the underground water level (days) | V__GW_DELAY.gw | 3 |
180 | 549 | Minimum required depth of stagnation level in shallow aquifers for flow to occur | V__GWQMN.gw | 4 |
0/29 | 0/88 | Coefficient of determination of infiltration into deep underground water or capillary rise from shallow water table | V__GW_REVAP.gw | 5 |
0 | 180 | The delay time of subsurface flows to the river in each HRU | V__LAT_TTIME.hru | 6 |
0/018 | 0/025 | Manning coefficient for the main river | V__CH_N2.rte | 7 |
34 | 39 | Effective hydraulic conductivity of the main river bed (mm/h) | V__CH_K2.rte | 8 |
0 | 0/06 | Coefficient α of base water for canal coastal storage | V__ALPHA_BNK.rte | 9 |
0 | −56/5 | Factor affecting precipitation | V__PLAPS.sub | 10 |
0 | 5/97 | Factor affecting temperature | V__TLAPS.sub | 11 |
0/05 | 0/12 | Slope in each hydrological group | R__HRU_SLP.hru | 12 |
0/02 | 0/07 | Manning coefficient | R__OV_N.hru | 13 |
−2/54 | −1/3 | Average air temperature for turning rain into snow | V__SFTMP.bsn | 14 |
1/5 | 1/51 | V__SMTMP.bsn | 15 | |
0/64 | 0/66 | Plant maintenance compensation factor | V__ESCO.hru | 16 |
Range of changes (provided by the SWAT-CUP model) . | Name of the parameters . | Parameter . | Row . | |
---|---|---|---|---|
−0/3 | 0/2 | Curve number | R__CN2.mgt | 1 |
0 | 1 | Coefficient α of groundwater | V__ALPHA_BF.gw | 2 |
230 | 500 | The delay time of water transfer from the last profile of the soil layer to the underground water level (days) | V__GW_DELAY.gw | 3 |
180 | 549 | Minimum required depth of stagnation level in shallow aquifers for flow to occur | V__GWQMN.gw | 4 |
0/29 | 0/88 | Coefficient of determination of infiltration into deep underground water or capillary rise from shallow water table | V__GW_REVAP.gw | 5 |
0 | 180 | The delay time of subsurface flows to the river in each HRU | V__LAT_TTIME.hru | 6 |
0/018 | 0/025 | Manning coefficient for the main river | V__CH_N2.rte | 7 |
34 | 39 | Effective hydraulic conductivity of the main river bed (mm/h) | V__CH_K2.rte | 8 |
0 | 0/06 | Coefficient α of base water for canal coastal storage | V__ALPHA_BNK.rte | 9 |
0 | −56/5 | Factor affecting precipitation | V__PLAPS.sub | 10 |
0 | 5/97 | Factor affecting temperature | V__TLAPS.sub | 11 |
0/05 | 0/12 | Slope in each hydrological group | R__HRU_SLP.hru | 12 |
0/02 | 0/07 | Manning coefficient | R__OV_N.hru | 13 |
−2/54 | −1/3 | Average air temperature for turning rain into snow | V__SFTMP.bsn | 14 |
1/5 | 1/51 | V__SMTMP.bsn | 15 | |
0/64 | 0/66 | Plant maintenance compensation factor | V__ESCO.hru | 16 |
As shown, the vegetation cover, underground water-feed coefficient and the delay time of water transfer from the last soil layer profile to the underground water level are also among the 16 flow-affecting parameters of this research.
Results: Evaluation indices
Climate change models and scenarios
In this article, two climate models CMIP5 and CMIP6 were used to investigate the effect of climate change on runoff. Due to the fact that the use of different models of climate reporting reduces the uncertainty in the estimation of climate changes, it has been used. At the beginning of the research, among the HadGEM2-ES, CNRM-CM5, EC-Earth3, MPI-ESM, FIO-ESM and GFDL-CM3 models according to Table 2, observational data were evaluated with basic data based on the statistical period from 1990 to 2014 based on R2 and RMSE criteria until the most appropriate model, HadGEM2-ES, was selected based on R2 and RMSE criteria, which were accepted as 0.97 and 0.2, respectively. The future data were generated for the Karaj Dam basin using RCP2.6, RCP4.5, RCP8.5 and HadGEM3-GC31-LL models under SSP1-2.6, SSP2-4.5 and SSP5-8.5 scenarios. The data were found from Climate Change Services at https://cds.climate.copernicus.eu/, the global-scale files of which are available in NetCDF format, but for local-scale data (e.g., this study), the data can be extracted first in the ACSII format and then converted to the input format of the SWAT model (done similarly in this research). This project aimed to quantify uncertainties in climate change scenarios and models that examine the climate change effects on different sectors such as water. As the climate data found from GCM models are usually biased compared to observational ones, they should be pre-modified. Another limitation of GCM models is their large-scale outputs, which are converted to smaller scales by scale-reducing methods.
Summary of GCMs used in CMIP6 and CMIP5 experiments
![]() | ![]() | Horizontal resolution (Longitude × Latitude) . | CMIP6 . | ![]() | ![]() | Horizontal resolution (Longitude × Latitude) . | CMIP5 . | Row . |
---|---|---|---|---|---|---|---|---|
0.2 | 0.97 | 1.875° × 1.25° | HadGEM3-GC31 | 0.2 | 0.92 | 1.9° × 2.5° | HadGEM2-ES | 1 |
0.83° × 0.56° | ||||||||
0.35° × 0.23° | ||||||||
0.2 | 0.91 | 1.406° × 1.406° | CNRM-CM6-1 | 0.3 | 0.89 | 1.4° × 1.4° | CNRM-CM5 | 2 |
0.5° × 0.5° | ||||||||
0.31 | 0.87 | 0.703° × 0.703° | EC-Earth3P | 0.2 | 0.85 | 1.1° × 1.1° | EC-Earth3 | 3 |
0.32 | 0.86 | 0.563° × 0.563° | MRI-AGCM3-2 | 0.3 | 0.82 | 1.9° × 1.9° | MPI-ESM | 4 |
0.188° × 0.188° | ||||||||
0.33 | 0.83 | 1.25° × 1° | FGOALS-f3 | 0.31 | 0.8 | 2.8 × 2.8° | FIO-ESM | 5 |
0.34 | 0.8 | 0.625° × 0.5° | GFDL-CM4C192 | 0.35 | 0.78 | 2.5° × 2.0° | GFDL-CM3 | 6 |
![]() | ![]() | Horizontal resolution (Longitude × Latitude) . | CMIP6 . | ![]() | ![]() | Horizontal resolution (Longitude × Latitude) . | CMIP5 . | Row . |
---|---|---|---|---|---|---|---|---|
0.2 | 0.97 | 1.875° × 1.25° | HadGEM3-GC31 | 0.2 | 0.92 | 1.9° × 2.5° | HadGEM2-ES | 1 |
0.83° × 0.56° | ||||||||
0.35° × 0.23° | ||||||||
0.2 | 0.91 | 1.406° × 1.406° | CNRM-CM6-1 | 0.3 | 0.89 | 1.4° × 1.4° | CNRM-CM5 | 2 |
0.5° × 0.5° | ||||||||
0.31 | 0.87 | 0.703° × 0.703° | EC-Earth3P | 0.2 | 0.85 | 1.1° × 1.1° | EC-Earth3 | 3 |
0.32 | 0.86 | 0.563° × 0.563° | MRI-AGCM3-2 | 0.3 | 0.82 | 1.9° × 1.9° | MPI-ESM | 4 |
0.188° × 0.188° | ||||||||
0.33 | 0.83 | 1.25° × 1° | FGOALS-f3 | 0.31 | 0.8 | 2.8 × 2.8° | FIO-ESM | 5 |
0.34 | 0.8 | 0.625° × 0.5° | GFDL-CM4C192 | 0.35 | 0.78 | 2.5° × 2.0° | GFDL-CM3 | 6 |
Investigating climate models in terms of evaluating observational data with basic data to choose the best GCM model
Hydrological and climatic drought indices
![](https://iwa.silverchair-cdn.com/iwa/content_public/journal/jwcc/15/7/10.2166_wcc.2024.721/1/m_jwc-d-23-00721if03.gif?Expires=1739896111&Signature=jPleYFex2H5ExKjnicyyl5fDHuPJDyPX9UNFci9KxX0f2vi9wL1lWzjXfWjg0O~3PAPF4wvzxXoa858DikstyshV3fTm~7jSxwIdg6lGOOYej7v3NjEaPA7BK~LnmQVL62B8PCkXmiIA5SO8cpqa4LqS2vfM100AODCo18CRsyA~xAcjTfCiOPdM~vfHbRZCz6MBBYdM3mhCv4ILscb2w3d3nxV2d9Ksym5HijMZKn1FQhHBTN028eplv516E701M14APeOvxDd8xxOx2jzcegMp15h43NuVvyx6BqrALTC4OF3t3uKmmrgiuRs~LqL-HGq3-AZZqWbX~F2TfVah2g__&Key-Pair-Id=APKAIE5G5CRDK6RD3PGA)
RESULTS AND DISCUSSION
Calibrating SWAT model
To check the results of the SWAT model and analyze the sensitivity of the effective flow parameters, use was made of the data of four hydrometric stations at the basin level along with the SUFI-2 calibration algorithm in the SWAT-CUP software. In the field of flow calibration of the studied area upstream of Karaj Dam, taking into account the various topographic characteristics of the region due to its mountainous nature and also the high areas connected to urban and pasture areas, it creates several challenges in the direction of calibration. In order to solve this problem, to consider the effects of topographical conditions, the importance of accurate representation of the parameter can be emphasized. Hence, when models are not accurately calibrated to describe these processes, predictions can be inaccurate, leading to potential misunderstandings of watershed behavior and the effectiveness of interventions. This step is very important, considering that topographic gradients can cause many changes in precipitation, snow, freezing and melting events. Therefore, this research applies the calibration of snow parameters (SMTMP, SFTMP) separately from other parameters, following the recommendations of Abbaspour et al. (2017). Parameters that played great roles in the formation of surface runoff and base flow in the model were selected (Table 3). Considering the model setting with the existing user maps and the climatic and hydrological data accessibility, 1991–2019 was considered for the model simulation period, data of the first 3 years were used for model preparation, and the model was calibrated for the first 19-year period (1994–2012). To calibrate the Table 3 parameters separately for different HRUs, the first parameters upstream of sub-basins of the hydrometric stations of each basin were calibrated. In the SUFI-2 program, the sensitivity of the parameters was determined using two criteria: T-stat and P-value. In the case of P-value, the closer the value of the parameter is to zero, the greater the sensitivity of the parameter, and also in the case of T-stat, its absolute value is greater than one. The sensitivity of the parameter is higher (Abbaspour et al. 2017).
Parameters affecting the amount of Karaj watershed runoff in the SWAT-CUP model
P-value . | T-stat . | Optimal value of the parameter . | Range of changes (provided by the SWAT-CUP model) . | Name of the parameters . | Parameter . | Row . | |
---|---|---|---|---|---|---|---|
0 | 18 | −0/02 | −0/3 | 0/2 | Curve number | R__CN2.mgt | 1 |
0 | −4 | 0/50 | 0 | 1 | Coefficient α of groundwater | V__ALPHA_BF.gw | 2 |
0.04 | 2 | 384 | 230 | 500 | The delay time of water transfer from the last profile of the soil layer to the underground water level (days) | V__GW_DELAY.gw | 3 |
0.08 | 1.2 | 250 | 180 | 549 | Minimum required depth of stagnation level in shallow aquifers for flow to occur | V__GWQMN.gw | 4 |
0.03 | 2 | 0/45 | 0/29 | 0/88 | Coefficient of determination of infiltration into deep underground water or capillary rise from shallow water table | V__GW_REVAP.gw | 5 |
0.1 | 1.2 | 45 | 0 | 180 | The delay time of subsurface flows to the river in each HRU | V__LAT_TTIME.hru | 6 |
0.5 | 3 | 0/024 | 0/018 | 0/025 | Manning coefficient for the main river | V__CH_N2.rte | 7 |
0 | 12 | 35 | 34 | 39 | Effective hydraulic conductivity of the main river bed (mm/h) | V__CH_K2.rte | 8 |
0.05 | 3 | 0/004 | 0 | 0/06 | Coefficient α of base water for canal coastal storage | V__ALPHA_BNK.rte | 9 |
0 | 4.5 | 56/6 | 0 | −56/5 | Factor affecting precipitation | V__PLAPS.sub | 10 |
0 | 2.8 | 5/99 | 0 | 5/97 | Factor affecting temperature | V__TLAPS.sub | 11 |
0.2 | 1 | 0/11 | 0/05 | 0/12 | Slope in each hydrological group | R__HRU_SLP.hru | 12 |
0.05 | 2 | 0/04 | 0/02 | 0/07 | Manning coefficient | R__OV_N.hru | 13 |
0.02 | 1.2 | −1/68 | −2/54 | − 1/3 | Average air temperature for turning rain into snow | V__SFTMP.bsn | 14 |
0 | 2.5 | 1/51 | 1/5 | 1/51 | V__SMTMP.bsn | 15 | |
0.2 | − 1 | 0/65 | 0/64 | 0/66 | Plant maintenance compensation factor | V__ESCO.hru | 16 |
P-value . | T-stat . | Optimal value of the parameter . | Range of changes (provided by the SWAT-CUP model) . | Name of the parameters . | Parameter . | Row . | |
---|---|---|---|---|---|---|---|
0 | 18 | −0/02 | −0/3 | 0/2 | Curve number | R__CN2.mgt | 1 |
0 | −4 | 0/50 | 0 | 1 | Coefficient α of groundwater | V__ALPHA_BF.gw | 2 |
0.04 | 2 | 384 | 230 | 500 | The delay time of water transfer from the last profile of the soil layer to the underground water level (days) | V__GW_DELAY.gw | 3 |
0.08 | 1.2 | 250 | 180 | 549 | Minimum required depth of stagnation level in shallow aquifers for flow to occur | V__GWQMN.gw | 4 |
0.03 | 2 | 0/45 | 0/29 | 0/88 | Coefficient of determination of infiltration into deep underground water or capillary rise from shallow water table | V__GW_REVAP.gw | 5 |
0.1 | 1.2 | 45 | 0 | 180 | The delay time of subsurface flows to the river in each HRU | V__LAT_TTIME.hru | 6 |
0.5 | 3 | 0/024 | 0/018 | 0/025 | Manning coefficient for the main river | V__CH_N2.rte | 7 |
0 | 12 | 35 | 34 | 39 | Effective hydraulic conductivity of the main river bed (mm/h) | V__CH_K2.rte | 8 |
0.05 | 3 | 0/004 | 0 | 0/06 | Coefficient α of base water for canal coastal storage | V__ALPHA_BNK.rte | 9 |
0 | 4.5 | 56/6 | 0 | −56/5 | Factor affecting precipitation | V__PLAPS.sub | 10 |
0 | 2.8 | 5/99 | 0 | 5/97 | Factor affecting temperature | V__TLAPS.sub | 11 |
0.2 | 1 | 0/11 | 0/05 | 0/12 | Slope in each hydrological group | R__HRU_SLP.hru | 12 |
0.05 | 2 | 0/04 | 0/02 | 0/07 | Manning coefficient | R__OV_N.hru | 13 |
0.02 | 1.2 | −1/68 | −2/54 | − 1/3 | Average air temperature for turning rain into snow | V__SFTMP.bsn | 14 |
0 | 2.5 | 1/51 | 1/5 | 1/51 | V__SMTMP.bsn | 15 | |
0.2 | − 1 | 0/65 | 0/64 | 0/66 | Plant maintenance compensation factor | V__ESCO.hru | 16 |
SWAT model simulation results in the flow calibration interval in selected stations of the Karaj Dam watershed according to Table 4, based on the criteria R2, NSE, P-factor and R-factor. According to the proposed criteria, the value of range extracted in the evaluation of observational data and simulation based on Table 5 according to Congalton (1999), Houshmand Kouchi et al. (2017). The accuracy of simulation in the calibration period from 2012 to 2013 at the Sierra-Karaj hydrometric station according to the R2, NSE range value is 0.72 and 0.73, Good, and in the Sleep Bridge hydrometric station, according to the R2, NSE range value for the station equal to 0.58 and 0.64 was in the Satisfactory condition and for Gachsar it was equal to 0.54 and 0.59 in the Satisfactory condition. Also, in Sierra-Kelvan hydrometric station, according to the R2, NSE range of 0.45 and 0.52, it is in the Satisfactory condition.
Statistical criteria to evaluate the model accuracy in runoff simulation of hydrometric stations
Verification . | Calibration . | Elevation . | Station name . | Row . | ||||||
---|---|---|---|---|---|---|---|---|---|---|
R-factor . | P-factor . | R2 . | NSE . | R-factor . | P-factor . | R2 . | NSE . | |||
0/65 | 0 | 0/54 | 0/43 | 0/46 | 0/73 | 0/59 | 0/54 | 2263 | Gachsar | 1 |
0/13 | 0 | 0/48 | 0/83 | 0/41 | 0/66 | 0/52 | 0/54 | 1833 | Sierra-Kelvan | 2 |
0/15 | 0 | 0/68 | 0/6 | 0/24 | 0/72 | 0/73 | 0/72 | 1831 | Sierra-Karaj | 3 |
0/04 | 0 | 0/61 | 0/54 | 0/35 | 0/57 | 0/72 | 0/68 | 1790 | Sleep bridge | 4 |
Verification . | Calibration . | Elevation . | Station name . | Row . | ||||||
---|---|---|---|---|---|---|---|---|---|---|
R-factor . | P-factor . | R2 . | NSE . | R-factor . | P-factor . | R2 . | NSE . | |||
0/65 | 0 | 0/54 | 0/43 | 0/46 | 0/73 | 0/59 | 0/54 | 2263 | Gachsar | 1 |
0/13 | 0 | 0/48 | 0/83 | 0/41 | 0/66 | 0/52 | 0/54 | 1833 | Sierra-Kelvan | 2 |
0/15 | 0 | 0/68 | 0/6 | 0/24 | 0/72 | 0/73 | 0/72 | 1831 | Sierra-Karaj | 3 |
0/04 | 0 | 0/61 | 0/54 | 0/35 | 0/57 | 0/72 | 0/68 | 1790 | Sleep bridge | 4 |
General performance ratings for model simulation (Congalton 1999; Houshmand Kouchi et al. 2017)
Performance rating . | R2 . | NSE . |
---|---|---|
Very good | 0.75 < R2 ≤ 1 | 0.75 < NSE ≤ 1 |
Good | 0.65 < R2 ≤ 0.75 | 0.65 < NSE ≤ 0.75 |
Satisfactory | 0. 5 < R2 ≤ 0.65 | 0. 5 < NSE ≤ 0.65 |
Unsatisfactory | R2 ≤ 0.5 | NSE ≤ 0.5 |
Performance rating . | R2 . | NSE . |
---|---|---|
Very good | 0.75 < R2 ≤ 1 | 0.75 < NSE ≤ 1 |
Good | 0.65 < R2 ≤ 0.75 | 0.65 < NSE ≤ 0.75 |
Satisfactory | 0. 5 < R2 ≤ 0.65 | 0. 5 < NSE ≤ 0.65 |
Unsatisfactory | R2 ≤ 0.5 | NSE ≤ 0.5 |
The accuracy of simulation in the validation period from 2013 to 2018 at the Sierra-Karaj hydrometric station according to the R2, NSE range value is 0.58 and 0.5, Satisfactory, and in the Sleep Bridge hydrometric station, according to the R2, NSE range value for the station equal to 0.51 and 0.44 was in the Satisfactory condition and for Gachsar it was equal to 0.54 and 0.34 in the Satisfactory condition. Also, in the Sierra-Kelvan hydrometric station, according to the R2, NSE range of 0.43 and 0.38, it is in Unsatisfactory status. The main reason for the low estimation of the model in the validity period was due to the unauthorized withdrawals along the route of the waterway network, the information of which was not available to be included in the model to increase the accuracy of the simulation.
Uncertainty and validating SWAT model
Another criterion for evaluating parameters affecting the flow rate of the Karaj Dam watershed in this research is to consider the P-factor and R-factor to quantify the existing uncertainty, parameters and input data in the conceptual model. The value of R-factor smaller than one generally shows a good calibration result of Ran and also the value of P-factor is at maximum equal to one indicating the confidence of 95% placement of all data (Faramarzi et al. 2009). In general, the results showed that P-factor and R-factor criteria obtained in the investigation of station runoff are within the acceptable range.
As observed in the hydrograph at the Gachsar station, the model has not been able to simulate the maximum flow well within a certain range, and according to the evaluation results, the modeling accuracy was in an acceptable condition, indicating the model's weakness in estimating the flow at this station. Another factor affecting the modeling accuracy is the precipitation and temperature at the station near the Gachsar sub-basin, which, due to missing data in the specified intervals, means the model has not been able to estimate the flow properly. In the hydrograph of the Sierra-Kelvan station, we also see this effect. But the simulation process in the sleep bridge and Sierra-Karaj stations is in a better condition, so it is considered good based on the evaluation index, and the model was able to estimate the maximum discharge, considering the fact that the precipitation and temperature statistics in the stations related to the basins have been more complete and there have been fewer missing numbers. In general, it can be concluded that the higher the number of climate stations in the region, the more appropriate statistical quality and the smaller the number of missing statistics, the more they have a positive effect on simulating the flow process in the model.
As shown in the table, the observational and simulation results of the NSE, r2, P-factor and R-factor evaluation indices reveal that the model has good simulation accuracy.
Climate change model
The rainfall and temperature estimation is possible by different atmospheric GCMs, but the process involves high uncertainties because the difference in location makes the physical process parameterization a difficult task (Reichler & Kim 2008). To reduce the estimation error in climate change studies in a region, use can be made of the results of several different models (Dessai et al. 2005). In this research, according to the assessment made in Section 2–4 and the selection of the HadGEM model, CMIP6 and CMIP5 precipitation and temperature changes were carried out under RCP and SSP emission scenarios. Three release scenarios of the CMIP5 model are divided into optimistic RCP2.6, moderate RCP4.5 and pessimistic RCP8.5. In CMIP6, the three emission scenarios SSP1-2.6, SSP2-4.5 and SSP5-8.5 are all upgraded from RCP2.6, RCP4.5 and RCP8.5, which is an important improvement over the CMIP6 scenarios. The SSP1-2.6 scenario in CMIP6 is optimistic, indicating the combined effects of low social vulnerability, low downsizing pressure and low radiative forcing. The SSP2-4.5 scenario is moderate and shows a combination of moderate social vulnerability and moderate radiative forcing. The SSP5-8.5 scenario is pessimistic, showing a combination of high social vulnerability and high radiative forcing, and the only path to achieving an anthropogenic radiative forcing level of 8.5 Wm2 by 2100.
Average annual minimum temperature changes of GMC models (2050–2030) and past (1991–2019) under RCP.
Average annual minimum temperature changes of GMC models (2050–2030) and past (1991–2019) under RCP.
Average annual maximum temperature changes of GMC models (2050–2030) and past (1991–2019) under RCP.
Average annual maximum temperature changes of GMC models (2050–2030) and past (1991–2019) under RCP.
Average annual minimum temperature changes of GMC models (2050–2030) and past (1991–2019) under SSP.
Average annual minimum temperature changes of GMC models (2050–2030) and past (1991–2019) under SSP.
Average annual maximum temperature changes of GMC models (2050–2030) and past (1991–2019) under SSP.
Average annual maximum temperature changes of GMC models (2050–2030) and past (1991–2019) under SSP.
Average annual precipitation changes of GMC models (2050–2030) and past (1991–2019) under RCP.
Average annual precipitation changes of GMC models (2050–2030) and past (1991–2019) under RCP.
Average annual precipitation changes of GMC models (2050–2030) and past (1991–2019) under SSP.
Average annual precipitation changes of GMC models (2050–2030) and past (1991–2019) under SSP.
Results found by the GCM model show an increase in the min and max temperature in 2030–2050 compared with those in 1991–2019; the increase has been more in the max and min temperatures for SSP and RCP scenarios, respectively, and pessimistic scenarios (RCP8.5, SSP5-8.5) have increased these temperatures more than optimistic scenarios (RCP2.6, SSP1-2.6). This research has used ‘lars-wg6’ for the GCM micro-scaling and measured the daily rainfalls at a total 11 meteorological and regional water stations during 1991–2019 (min and max daily temperatures belong only to four synoptic stations of the ‘Iran Meteorological Organization’), and to extract the future (2030–2050) climate data use has been made of ‘lars-wg6’ under the HadGEM (global environmental model).
Temperature changes can cause an increase in surface evaporation and transpiration in the region (Choi et al. 2016; Mei et al. 2021). According to the IPCC report, it is likely that the average evaporation in the world has an indirect relationship with temperature and causes an increase in evaporation of 1–3% for every 1 °C increase in temperature (IPCC 2021). The multi-year average maximum and minimum air temperatures in the base period (1990–2018) are −13 and 10 °C, respectively. According to the extracted data, the maximum and minimum temperature for each of the release scenarios of the climate model of the fifth and sixth reports in the scope of the study of changes in the average air temperature in the period of 2030–2050 compared with the base period for all three release scenarios of SSP and RCP as shown in Figure 3. Up to six are provided. In the future, the multi-year average temperature in the river basin upstream of the Karaj Dam will be accompanied by a rapid warming trend. This may be due to the intensification of global warming caused by the continuous increase in global greenhouse gas emissions, which is reflected in the increase in the rate of warming with an equivalent increase in carbon dioxide emissions. Figures 3–6 show that all three SSP and RCP release scenarios estimate that the maximum and minimum temperatures will increase to different degrees in the future.
Taking into account Figures 5 and 6, which show the average percentage of changes in maximum air temperature in the RCP4.5, RCP2.6 and RCP8.5 scenarios, the values are 2.38, 2.11 and 2.73 respectively. These values indicate an increasing trend in the maximum temperature of the study area. Additionally, the average percentage of changes in minimum air temperature are 2.58, 2.59 and 2.91, indicating an increase. It can be concluded that based on the output of all three RCP scenarios, the climatic conditions in the upper Karaj River basin will continue to show warm and humid trends in the future.
Considering Figures 5 and 6, which represent the average percentage changes in maximum air temperature under the SSP1-2.6, SSP2-4.5 and SSP5-8.5 scenarios, the values are 2.4, 2.63 and 1.94, respectively. These values indicate an increasing trend in the maximum temperature of the study area. Additionally, the average percentage changes in minimum air temperature are 1.62, 2.03 and 2.19, respectively, showing an increase. It can be concluded that based on the output of all three SSP scenarios, the hydrothermal conditions in the upstream basin of Karaj River will continue to demonstrate warm and humid trends in the future.
One of the signs of climate change in the study area is the variability of rainfall events, which has a direct effect on the process of extracting runoff and the hydrological cycle. In order to investigate the change in precipitation events according to the base period from 1990 to 2018, the average annual rainfall in the base period is 600 mm. Rainfall extraction was done from 2030 to 2050. As shown in Figures 7 and 8, the average percentage of rainfall changes in the three release scenarios SSP1-2.6, SSP2-4.5 and SSP5-8.5 shows a decrease in rainfall compared with the base period. In the SSP release scenarios, the highest amount of rainfall reduction has been predicted so that the average of 10–20% rainfall reduction in the future is shown. For the RCP emission scenarios, the amount of precipitation events from 2030 to 2050 has decreased by 5–10% compared with the base period. On average, the average annual rainfall in the river basin upstream of the Karaj Dam changes slightly under each SSP scenario and shows a more decreasing trend than the RCP scenarios.
Climate change impact on runoff
According to the results obtained from the investigation of the average flow rate in the hydrometric stations of the study area, especially the Sierra-Karaj hydrometric station, which is located near the entrance to the Karaj Dam reservoir, it can be concluded that the changes in the average flow rate under the release scenarios of each of the CMIP5 and CMIP6 models indicate that there is a possibility of a 5–30% decrease in the runoff entering the dam reservoir, which will reduce the water volume storage behind the dam. Table 6 shows the percentage of changes in each of the stations.
The average annual changes in flow under different scenarios
HadGEM3-GC31-LL . | HadGEM2-ES . | Station name . | Row . | ||||
---|---|---|---|---|---|---|---|
SSP8.5 . | SSP4.5 . | SSP2.6 . | RCP8.5 . | RCP4.5 . | RCP2.6 . | ||
−26.84 | −24.42 | −24.76 | −35.93 | −29.54 | −28.39 | Gachsar | 1 |
−23 | −20.62 | −20.71 | −33.28 | −19.32 | −22.56 | Sierra-Kelvan | 2 |
−8.88 | −5.97 | −6.16 | −18.92 | −9.74 | −7.9 | Sierra-Karaj | 3 |
−2.18 | −1.1 | −0.96 | −1.1 | −1.33 | −0.1 | Sleep bridge | 4 |
HadGEM3-GC31-LL . | HadGEM2-ES . | Station name . | Row . | ||||
---|---|---|---|---|---|---|---|
SSP8.5 . | SSP4.5 . | SSP2.6 . | RCP8.5 . | RCP4.5 . | RCP2.6 . | ||
−26.84 | −24.42 | −24.76 | −35.93 | −29.54 | −28.39 | Gachsar | 1 |
−23 | −20.62 | −20.71 | −33.28 | −19.32 | −22.56 | Sierra-Kelvan | 2 |
−8.88 | −5.97 | −6.16 | −18.92 | −9.74 | −7.9 | Sierra-Karaj | 3 |
−2.18 | −1.1 | −0.96 | −1.1 | −1.33 | −0.1 | Sleep bridge | 4 |
The negative sign in the above table indicates a decrease in the flow rate in percentage terms.
Drought effects
The SPI and SRI values, found from long-term rainfall data for a certain time period, follow a normal distribution. Different values of these indices are shown in Tables 7 and 8 for the present and future time period simulations based on effective rainfalls in the study area. As shown, on an annual scale, 5 years of the 29-year simulation period are normal, 19 years are semi-arid and the rest are dry; however, after the 20-year statistical period, 6 years are normal, 13 years are semi-arid and the rest are dry.
SPI evaluation of CMIP5 and CMIP6 climate models
Future period 2030–2050 . | Again 1991–2019 . | SPI index . | Drought classes . |
---|---|---|---|
2031–2037–2038–2040–2042–2047 | 1993–1994–2005–2006–2011 | > 1 | Normal and above |
2030–2032–2033–2034–2035–2039–2041–2043–2044–2045–2046–2049–2050 | 1991–1992–1995–1996–1997–1998–2000–2001–2002–2003–2004–2008–2009–2010–2013–2014–2015–2016–2018 | −1 to 1 | Close to normal |
2036 | 1999–2007–2012 | −1.5 to −1 | Moderate drought |
2048 | 2017 | −2 to −1 | Severe drought |
– | – | ≤ −2 | Very severe drought |
Future period 2030–2050 . | Again 1991–2019 . | SPI index . | Drought classes . |
---|---|---|---|
2031–2037–2038–2040–2042–2047 | 1993–1994–2005–2006–2011 | > 1 | Normal and above |
2030–2032–2033–2034–2035–2039–2041–2043–2044–2045–2046–2049–2050 | 1991–1992–1995–1996–1997–1998–2000–2001–2002–2003–2004–2008–2009–2010–2013–2014–2015–2016–2018 | −1 to 1 | Close to normal |
2036 | 1999–2007–2012 | −1.5 to −1 | Moderate drought |
2048 | 2017 | −2 to −1 | Severe drought |
– | – | ≤ −2 | Very severe drought |
Evaluation of the hydrological SRI
Future period 2030–2050 . | Again 1991–2019 . | SPI index . | Drought classes . |
---|---|---|---|
2031–2037–2038–2040–2042–2047 | 1993–1994–2005–2006–2011 | > 1 | Normal and above |
2030–2032–2033–2034–2035–2039–2041–2043–2044–2045–2046–2049–2050 | 1991–1992–1995–1996–1997–1998–2000–2001–2002–2003–2004–2008–2009–2010–2013–2014–2015–2016–2018 | −1 to 1 | Close to normal |
2036 | 1999–2007–2012 | −1.5 to −1 | Moderate drought |
2048 | 2017 | −2 to −1 | Severe drought |
– | – | ≤ −2 | Very severe drought |
Future period 2030–2050 . | Again 1991–2019 . | SPI index . | Drought classes . |
---|---|---|---|
2031–2037–2038–2040–2042–2047 | 1993–1994–2005–2006–2011 | > 1 | Normal and above |
2030–2032–2033–2034–2035–2039–2041–2043–2044–2045–2046–2049–2050 | 1991–1992–1995–1996–1997–1998–2000–2001–2002–2003–2004–2008–2009–2010–2013–2014–2015–2016–2018 | −1 to 1 | Close to normal |
2036 | 1999–2007–2012 | −1.5 to −1 | Moderate drought |
2048 | 2017 | −2 to −1 | Severe drought |
– | – | ≤ −2 | Very severe drought |
As shown in both tables, the actual and future rainfall statistics of the region reveal that the drought index is almost normal for most of the time period.
CONCLUSION
The results of this study can be summarized as follows:
- Results of the sensitivity analyses, calibration and validation of the SWAT model showed that it is quite sensitive to such controlling parameters of the basin runoff as the curve number (CN). The model calibration and validation results were also evaluated by R2, NSH and RMSE indices, which showed that the model simulation of the water components had acceptable accuracy.
- Results of the evaluation of the climatic components (min and max daily temperature) with GCMs under different scenarios indicated an increase in the future temperature of the basin; RCP2.6, RCP4.5 and RCP8.5 scenarios reported an increase in the minimum temperature, but SSP1-2.6, SSP2-4.5 and SSP5-8.5 scenarios showed an increase in the maximum temperature, which increased the evapotranspiration and thus reduced the future runoff.
- Contrary to GCMs' temperature results, their descaling rainfall results were contradictory; the rainfall prediction uncertainty range varied from a 1.19% increase to a 7.66% decrease, in RCP scenarios, and from 8.31 to 11.15% decrease, in SSP scenarios.
- The SWAT model's runoff simulation with GCM's climate data (rainfall and min and max temperature) had different results under CMIP5 and CMIP6 scenarios. The future runoff under climate change conditions was reduced in scenarios 2.6, 4.5 and 8.5 for both models and four hydrometric stations due, maybe, to the increased max temperature compared to the min temperature and, hence, the increased evaporation.
- It is suggested that the process of developing the management and planning of the upstream basin of the dam reservoir should be reviewed for each of its programs and topics based on a proper understanding of the relevant problems.
DATA AVAILABILITY STATEMENT
All relevant data are included in the paper or its Supplementary Information.
CONFLICT OF INTEREST
The authors declare there is no conflict.