Impact of climate change on the future quality of surface waters: case study of the Ardak River, northeast of Iran

Climate change, which is the changes in weather patterns due to increases in gas concentrations that absorb heat in the atmosphere, is changing and disrupting ecosystems and daily life. One of the consequences of this phenomenon is the catastrophic effects on water resources, which have become important in recent years due to water scarcity in many parts of the world. The effect of climate change on temperature and precipitation may influence runoff volume (Luo et al. 2013) and change the transport and dilution of contaminants (Barrow et al. 1996). Increases in temperature have a direct impact on the rate of chemical reactions, which in many cases are accompanied by reduced water quality and changes to water ecology. Increasing flow rates may affect sedimentation rates, changing the shape and structure of rivers and potentially adversely affecting drinking water resources. Climate change scenarios suggest more extreme drought and devastating floods in the summer, causing uncontrolled water flow from urban areas towards receiving waters and river estuaries (Rehana & Mujumdar 2011; Whitehead et al. 2009). Low flow rates also increase the chance of toxic algal bloom by reducing water velocity and increasing water detention time in rivers and lakes, as well as reducing the amount of dissolved oxygen. As a result, waters from upstream areas, which contain increased soluble organic matter and color, may affect the quality of water resources, including rivers. These factors affect the quality and quantity of water and may harm ecosystems and water bodies (Bagherian Marzouni et al. 2014). Despite significant attention to the impacts of climate change on water supplies, little information is available about the consequences of water quality changes caused by climate change, and current predictions are often too low. A common approach to estimate climate impacts on hydrologic behavior is the implementation of general circulation models (GCM) via downscaling methods (Rehana & Mujumdar 2011).

In this study focusing on the Ardak River (northeast of Iran), GCM and RCP scenarios are used to obtain temperature and precipitation parameters for the future; river flows are then simulated using a rainfall-runoff model, allowing assessment of the effects of climate change on various water quality parameters from the QUAL2KW one-dimensional model (Pelletier & Chapra 2008a). One of the main weaknesses of previous similar studies is the lack of accurate and reliable statistics, especially for river water quality tests (Fereidoon & Khorasani 2013). To overcoming this problem, qualitative water samples and tests were conducted to validate the QUAL2KW model, in addition to using statistics. The main objective of this study is simulation of future surface water quality, in addition to river water quality monitoring in the main river branches, simulation of future runoff, and study of the effects of climate change on river quality based on changes in rainfall and runoff temperature. On the other hand, in this study changes in land use, population and sewage discharges are considered to be constant. This made it possible to monitor only climate change effects on water quality of the river.

The Ardak Chenaran Basin is located approximately 45 km northwest of Mashhad (northeast of Iran) and has an area of 479.24 km². It is in the Kashafrud watershed and is within Mashhad and Chenaran in the central part of the region. The Ardak River consists of two main branches, the Mianmargh branch, in the northwest, and the Abghad, in the northeast, as shown in Figure 1. These branches are nearly connected at the southeastern part of the basin, and they enter the plain after about 3 km. The Ardak dam was established in recent years to supply water for drinking and agriculture. The main branch of the Ardak River (the Mianmargh) has a number of other tributaries, including Amrudak, Andishsh, and Gash, whose water finally flows into the reservoir of the Ardak dam after joining the Mianmargh River. The Abghad branch also ends at the dam and has less water due to its shorter length, high river-bed infiltration rates, and water diversion from traditional dikes to surrounding gardens. The Mianmargh and Abghad, which are the main drainage routes for surface water in this region, were the focus of this study.

In the Ardak study area, more than 90 percent of the water resources are allocated to the agricultural and horticultural sectors. However, drinking water and sanitation utilize just over 9 percent of the whole volume of available water supplies. For drinking demand, water should be clear, colorless, odorless and not contaminated with germs and pathogens. Furthermore, organic and mineral substances in water should not be harmful to the human body. Water for agricultural demand should be quite healthy with low volumes of elements and suspended materials. Given the population growth rate in this region (0.6 percent), there is expected to be a 15 percent increase in the current consumption of water at the end of the study period (2039).

As shown in Figure 2, the main uses of land are rangelands, agriculture, horticulture and woodlands. Currently, the project area has no centered and developed tourism, industrial and mining units.

Additionally, the major types of soils in the region are limestone and shale.

Information from the Abghad-Darzab and Mianmargh hydrometric stations was used to simulate the existing flow in the Abghad and Mianmargh branches (Figure 1). To simulate climate change, temperature and precipitation information from the Ardak-Saroog dike and Ghadirabad meteorological stations (Khorasan Razavi Regional Water Authority) were used (Figure 1). For calibrating the Qual2kw model, data from sampling stations indicated in Figure 1 were used. Information from the regional water company of Khorasan Razavi (Jooyab Consultant Company, November 2014) was used for calibration and for verification of the model in April 2016, with qualitative sampling and flow measurements done by a flow vane instrument. Samples were transferred to the water laboratory of the Khorasan Razavi Regional Water Authority under standard conditions (Standard Methods 2012 ), and chemical oxygen demand (COD), dissolved oxygen (DO), NO₃, and pH were analyzed for verification of the qualitative model. The Sanders method was used to design sampling stations, whose features are explained in Table 1.

One common and practical method of water quality modeling is using reciprocal exchanges between materials and quantities (Pelletier & Chapra 2008a). Qual2kw version 5.1, which was used in this study, uses this method and is a modern version of the famous stream water quality program QUAL2E, which was presented by Tao, Chapra and Pelletier in 2006 (Turner et al. 2009). This program simulates abatement and transfer of common pollution (non-toxic) in a river. Assuming steady state, non-uniform flow, Qual2kw simulates rivers in one dimension, including the effects of contaminant charge and sources (Pelletier & Chapra 2008b). Qual2kw is also capable of simulating daily changes and time steps of less than one hour (Pelletier & Chapra 2008a). The statistical evaluator root mean square error (RMSE) is used for measuring the accuracy of model predictions (Pelletier & Chapra 2008a). This model uses three forms of balance equations: hydraulic balance equation for flow rate, heat balance equation, and mass balance equation for qualitative material concentrations.

where:

• Q_in,i is the total inflow into the reach from the point and nonpoint sources [m³/d]

• Q_ab,i is the total outflow from the reach due to point and nonpoint abstractions [m³/d]

• Q_i–1 is the inflow from the upstream reach i − 1 [m³/d]

• Q_i is the outflow from reach i into the downstream reach i + 1 [m³/d].

Qual2 k simulates diffuse inputs or outputs linearly. It considers the start and end of these resources as shown in Figure 4, and weights each element according to input load. Therefore, the length of diffuse inflowing or outflowing sources must be known (Pelletier & Chapra 2008a).

The second edition of IHACRES was used, which was improved by Croke et al. (2005). This model is lumped, which means that it considers the whole catchment area as a single entity and uses only one station and its observed data. The model takes action to simulate present streamflow in the river using observed data (temperature, precipitation, and runoff) from the station. This action consists of nonlinear (loss) and linear (hydrograph) modules, which are used for making an accurate calculation of loss and transferring effective rainfall to stream discharge, respectively, as shown in Figure 6 (Croke et al. 2005).

The mentioned model requires five to seven variables for calibration. In the calibration process, values of (the catchment drying time constant) and f (the temperature modulation factor), which comes from the nonlinear module, are manually selected by the user. In this procedure, the best model is chosen by the coefficient of determination R² (Croke et al. 2005).

The atmospheric circulation models, in fact, solve the governing equations of the atmosphere based on Newton and thermodynamics laws. In other words, they are 3D models that have been widely used in recent years to simulate future climate change (Kockum et al. 2016). There are many GCMs that represent the effects of various factors, such as reflection and absorption of atmospheric water vapor, greenhouse gas concentrations, clouds, annual and daily solar heat, ocean temperatures, and ice boundaries (Wilby et al. 2009). In this study, we have used the GCMs listed in Table 2 (Yin et al. 2013) for future projection of temperature and precipitation. These parameters have direct impacts on the flow of rivers and thus affect the quality of water.

Each GCM requires initial condition data to run. In some of these models, which predict the climate patterns over long-term periods, Intergovernmental Panel on Climate Change (IPCC) emission scenarios act as boundary conditions (Wilby et al. 2009).

These scenarios estimate future greenhouse gas emissions into the atmosphere. Such estimates are based on possible projection of economic and population growth and technological development, as well as physical processes within the climate system (Trzaska & Emilie 2014). The IPCC, which is responsible for recognizing all aspects of the climate change phenomenon, has published various versions of the emission scenarios. The latest version, called RCP scenarios, was released with the Fifth Assessment Report (AR5) in 2014. It has four pathways called RCP2.6, RCP4.5, RCP6 and RCP8.5, which are named after a possible range of radiative forcing values in the year 2100. Among them, RCP2.6 and RCP8.5 with the lowest and highest emission levels (Hosseini et al. 2017) are more beneficial to reflect the changes (Wilby et al. 2009). Therefore, the study of consequences on water quality is made more tangible and possible by using them. Thus, those two scenarios were employed in this research.

Here, P_rcp and T_rcp are long-term monthly average of rainfall and temperature under RCP2.6 and RCP8.5 scenarios, P_historical and T_historical are long-term monthly average of rainfall and temperature obtained from historical data, and ΔP and ΔT are change factors of precipitation and temperature. Entering Δ values into the LARS-WG model enables it to downscale GCM outputs and generate weather data for future periods.

Data collected by the Khorasan Regional Water Company in November 2014 was used to calibrate the model. Sampling was conducted on 20 April 2016 to verify the model. Samples were analyzed in the accredited laboratory of the Khorasan Razavi Regional Water Company after transferring them under standard conditions. The parameters DO, pH, COD, and NO₃ were measured for verification. In the calibration step, reaction rates were determined based on the determination coefficient RMSE. During verification, the difference between modeled information and observed data was determined (Nuzhat & Singh 2016; Kalburgi et al. 2015).

Figure 7 shows changes in various water quality parameters along the Mianmargh branch. The low DO values observed at some points (Figures 7(a) and 7(b)) indicate a sharp decline in dissolved oxygen in some parts of the river and were mainly due to the discharge of wastewater and solid waste from villages on the river's edge. Analysis of COD (Figures 7(g) and 7(h)) confirms the location of the contamination leading to low dissolved oxygen in some parts of the river. Nitrate (Figures 7(e) and 7(f)) decreases from upstream to downstream. The high amount of nitrate in the upper part of the river was caused by rural wastewater pollution. Changes in nitrate in other parts of the river follow a slope due to diffuse sources of pollution such as agriculture and horticulture, which are located along the river. Although pH fluctuates, it is generally slightly alkaline in the Mianmargh branch (Figure 7(c) and 7(d)). This is due to dissolution of carbonates and bicarbonates from the soil and river bed (Masamba & Mazvimavi 2008). At 20 km upstream, as shown in Figures 1 and 2, the river passes through the village of Boghmech, the largest and most populated village in the region. In this village, household and agricultural wastewater is discharged directly into the river, affecting pH and acidity of the river more than at other measured points of the river. However, it should be noted that fluctuations in pH of about 0.1 occur along the river, and this is negligible. Also the difference in pH between a point 20 km upstream and other places is only 0.3, which is not a big difference.

Figure 8 shows changes of various water quality parameters in the Abghad branch, from upstream (10 km) to the point where it flows into the Ardak Dam Lake (0 km). As shown in Figures 8(a) and 8(b), DO was high because few sources of pollution such as villages exist along this branch. At some points, the concentration of dissolved oxygen is higher than saturated, which may be due to benthic plants. COD (Figures 8(g) and 8(h)) increases from upstream to downstream, which may be due to pollutants, including organic matter, from the village near the middle of this river branch. Organic matter decomposition is completely anaerobic, and as a result, there is a decrease in DO in areas where organic matter concentration is high. Nitrate (Figures 8(e) and (f)) is highest in the upstream part of the river because the main sources of nitrate production are agricultural and horticultural activities in the upstream regions of this branch of the river. Figures 8(c) and 8(d) show pH in the Abghad branch during calibration and verification periods. Here, like in the Mianmargh branch, a wide rural area is located 6 km upstream, which directly discharges untreated agricultural wastewater into the river and causes the quality indices to differ at this point more than at other points along the river. However, this difference is not great.

After selecting several periods for calibration and reviewing results, the remaining courses were selected for model verification in the Mianmargh and Abghad branches. During calibration, some values such as the catchment drying time constant and the temperature modulation factor f in the nonlinear section of the model are manually selected by the user. The best model is selected based on the Nash coefficient value, with a value closer to one indicating a better fit (Croke et al. 2005).

To simulate monthly runoff of the Mianmargh branch, IHACRES was initially calibrated for the period of 1991–2013 and then verified for 1979–1991. For 1991–2013 the model had a Nash coefficient of 0.67 (Figure 9(a)), indicating that the calibrated model simulates the flow well. Furthermore, Figure 9(b) shows runoff time series used in the model validation process, resulting in a Nash coefficient of 0.68, showing that the model simulates flows in the basin, including low and maximum flows, well.

To simulate monthly runoff in the Abghad Branch, the IHACRES model was calibrated for 1977–2002 and then verified for 2002–2013. Figure 10(a), which shows the observed and simulated flow data of the Abghad Branch, indicates that the model was successfully calibrated (Nash = 0.68). The value of Nash = 0.66 in Figure 10(b) indicates that the model simulated water flow well for verification of the Abghad branch.

Water flow in the Ardak River in both the Mianmargh and Abghad branches was inspected and simulated for the future assuming climate change (2020–2039). A 26 year (1989–2014) statistical period of atmospheric data and the LARS-WG climate change model were used to simulate climatic conditions in the coming decades. Therefore, daily climate variables, including rainfall and minimum and maximum temperature from the synoptic stations of Ghadir­abad and Band Sarooj (Figure 1), were obtained from the regional water company of Khorasan Razavi province and used for calibration and verification of Lars-WG. In the next step, with the help of the IPCC website, and using GCMs of the RCP2.6 and RCP8.5 scenarios (Fifth Report), future temperature and precipitation information were calculated. Then, using data obtained from the GCMs and a downscaling method, future temperature and precipitation change values were calculated for RCP2.6 and RCP8.5 scenarios. By implementing the LARS-WG model and applying changes to observational data, temperature and precipitation values for 2020–2039 were generated.

Rainfall and temperature simulated by the climate model fit well with observational data, as shown in Figure 11.

As seen in Figures 11(a) and 11(b), both branches have the highest rainfall in the spring and the least rainfall in the summer. Besides that, both branches have the highest temperature in July and minimum temperature in January. Modeled future climate parameters (2020–2039) are compared with observational data in Figure 12.

By entering monthly temperature and rainfall time series which are scaled down from GCM into the IHACRES model, monthly runoff time series of each hydrometric station were produced for 2020–2039 (Figure 13).

The biggest increases in runoff in the Abghad and Mianmargh are in May (RCP8.5) and June (RCP8.5), respectively (Figures 14(a) and 14(b)). The highest losses of runoff are in January (RCP2.6) for both branches of the river.

To identify qualitative changes in surface water, average runoff and long-term seasonal temperature of both branches were entered into Qual2kW, which was calibrated and verified using existing pollutant sources. Although several factors influence the quality of surface water resources, the effects of runoff and temperature are indisputable. In the following, output results of the QUAL2 K model are given for the qualitative parameters dissolved oxygen, nitrate, pH, and COD for each branch of the Ardak River.

As shown in Figure 15, DO, NO₃, COD, and pH (Figures 15(a), 15(c), 15(e) and 15(g)) in the Mianmargh branch will not change significantly in the future in the RCP2.6 and RCP8.5 scenarios. This small change poses no concerns for water quality in the future. However, the Abghad branch (Figures 15(b), 15(d), 15(f) and 15(h)) shows noticeable changes in DO, NO₃, COD, and pH. In this branch, increasing amounts of NO₃ and COD, and decreasing DO and pH values are seen in both scenarios.

In this paper, the effects of climate change on quality of surface water resources in the Ardak River were studied. To do this, we first calibrated the Qual2kw model using available information from the Khorasan Razavi Regional Water Company. For validation, in April 2016 water quality parameters and flow were measured. Insignificant difference in the RMSE results (Kannel et al. 2007) for nitrate, dissolved oxygen, pH, and COD in calibration and verification phases showed that the model is able to simulate water quality parameters along both main branches of the Ardak River. Using flow data, the current flow in both river branches was modeled by the IHACRES rainfall-runoff model. Nash coefficients for the calibration and verification phases showed that the model simulated the current river flow well. Daily temperature and precipitation data were input to the LARS-WG model, and this model was calibrated and verified. Then, using information from the IPCC website and its fifth report, information about various GCMs was collected for the RCP2.6 and RCP8.5 scenarios. Using this data, ΔP and ΔT for each model were calculated and extracted data were scaled down for the study area using change factors and the LARS-WG model. Using GCM sets which were down-scaled and applied to daily temperature and precipitation values, the LARS-WG model simulated future temperature and precipitation values. Rainfall variations in the two branches of the river are similar to each other due to their close proximity. The highest rainfall based on RCP8.5 will occur in the spring with maximum seasonal average of 42.94 mm in Abghad branch and 32.56 mm in Mianmargh branch. Besides that, the lowest rainfall values are related to RCP2.6, with average seasonal value equal to 4.24 mm in the Abghad branch and 3.84 mm in Mianmargh branch that occurs in the summer. Temperature changes in the river branches increase in almost all months of the year. While the highest mean seasonal temperature reaches up to 24.29 °C in Abghad branch and 24.03 °C in Mianmargh branch based on RCP8.5 in summer, the lowest temperature is 3.41 °C in Abghad branch and 3.17 °C in Mianmargh branch based on RCP2.6 in winter. Future stream flows were simulated given future temperatures and precipitation values. Simulated seasonal mean runoff maximum was 114.616 litres per second in Abghad branch and 1,306 litres per second in Mianmargh branch based on RCP8.5 which will occur in summer. However, the lowest runoff based on RCP 2.6 was 48.68 and 571.1 litres per second in the branches of Abghad and Mianmargh, respectively, in winter. Using future seasonal temperature and flow data allowed water quality parameters to be simulated. No significant change is predicted for water quality in the Mianmargh branch. However, due to the lower flow (Wu & Xia 2014) in the Abghad branch, water quality will change significantly in both scenarios, including amounts of nitrate, COD, and DO. This is likely due to increased temperature which activates bacteria that are responsible for nitrification. These issues make the Abghad branch highly vulnerable to changes and decline in water quality in the coming years. In this part of the river, the highest percentage variations in DO and COD are −12.19 and −35.4 in RCP8.5 and RCP2.6, respectively, while NO₃ and pH have a maximum variation of 31.25 and 0.29 percent in RCP8.5 and RCP2.6, respectively, compared with their annual average.

In general, the results showed that in both greenhouse gas emission scenarios (optimistic (RCP2.6) or pessimistic (RCP8.5)) meteorological, hydrological and water quality variables have changed and the situation would be worse for the region. Of course, it should be noted that for the time series considered for this research to 2039, the difference between the optimistic and pessimistic scenarios in terms of the amount of greenhouse gas emissions is not significant (about 100 ppm of difference), and in the case of a longer time series (for instance, by the year 2100), a large discrepancy in greenhouse gas emissions between these two scenarios (about 800 ppm) can make a clear contrast in the changes of meteorological, hydrological and water quality variables of this region.

On the other hand, it should be understood that the results of this research are presented regardless of various uncertainties, such as the uncertainty of GCMs, the uncertainty of downscaling methods, the uncertainty of rainfall-runoff simulation models and the uncertainty of the river water quality simulator models. Therefore, to provide more logical results, these uncertainties should also be taken into account in the analysis.

Finally, in this research, the lumped conceptual model of IHACRES was implemented to simulate streamflow in the main branches of Ardak River. This model covers the entire basin as an integrated entity; therefore, changes in soil, vegetation and land use of the watershed are considered the same. Hence in regions where such variations are considerable, more complex models such as SWAT (Getachew et al. 2017) perform better in the simulation process. Moreover, many qualitative models such as Qual2kw are mechanistic models (Liangliang & Daoliang 2014) that are not capable of controlling the extraordinarily complicated rules governing the nature of the aquatic environment that lead to less accurate results. Thus, using different algorithms and methods such as artificial neural networks (ANN), which can combine different model outputs, could be suggested for future studies in order to reduce the errors and better simulate water quality.