Climate change and hydrological response of Megech catchment, Upper Blue Nile River Basin, Ethiopia

Abstract

Climate change scenarios of precipitation and temperature were divided into three time windows of 30 years each from 2011 to 2099. The new version of Soil and Water Assessment Tool (QSWAT v2.6.1) was used to simulate the hydrological response and it was first calibrated and validated using observed data as an input for assessing the hydrological responses. The results showed that the QSWAT calibration and validation reveals a good agreement (R²=0.77) during calibration and validation. Whereas nash and Sutcliffe simulation efficiency (NSE) was found to be 0.76 during calibration and 0.73 during validation. Based on changes of precipitation, maximum and minimum temperature, the monthly flow volume did not show systematic trends i.e. increases non-rainy months and decreases in the wet-months. Thus, the hydrology of Megech River is highly vulnerable to climate change which causes water stress for domestic use and environmental flow.

Keywords:

• Climate change
• hydrology
• rainfall
• downscaling
• global circulation
• models
• Blue Nile river basin

Introduction

The global average temperature has been steadily rising due to climate change. NASA reports that Earth’s average surface temperature in 2023 was the warmest on record since recordkeeping began in 1880. The planet was approximately 2.45 degrees Fahrenheit (1.36 degrees Celsius) warmer in 2023 compared to the late 19th-century preindustrial average (Mertz 2023). Climate change is a pressing global issue that significantly impacts the hydrological cycle.

It involves the movement of water through land, oceans, and the atmosphere. As climate changes, so does this intricate dance of water. Evaporation and precipitation are key components of the hydrological cycle. Climate change affects both these processes. Warming global temperatures accelerate evaporation, leading to more water vapor in the atmosphere. The Hydrological Responses to Climate Change intensifies hydrological extremes, altering not only their magnitude but also their timing and floods and droughts may occur at different times than historically observed. Understanding these potential changes is crucial for policy decisions and adaptation measures (Mertz 2023). Global average temperature would rise between 1.4 and 5.8 °C by 2100 with the doubling of the CO₂ concentration in the atmosphere (Gregory et al. 2007). Sea level rise, change in precipitation pattern (up to ± 20%), and change in other local climate conditions are expected to occur as a consequence of rising global temperature (Pachauri et al. 2014). Being one of the very sensitive sectors, climate change can cause significant impacts on water resources by resulting changes in the hydrological cycle. The IPCC finding also indicates that those developing countries, such as Ethiopia will be more vulnerable to climate change. Because of the less flexibility to adjust the economic structure and being largely dependent on the rain-fed agriculture (Adger et al. 2003). From the point of view of the design and management of water resource systems, experts are required to make accurate predictions of the impacts of climate change on the intensity, amount, and spatial and temporal variability of rainfall. Furthermore, it has been also examined how the stream flow regime (e.g. stream flow hydrographs, peak flow, etc.) at different spatial and temporal scales is affected by rainfall variability and by the expected changes in that variability as a result of climate change (Cooper et al. 2008). One of the most important impacts on society of future climatic changes will be changes in regional water availability. Such hydrologic changes will affect nearly every aspect of human well-being from agricultural productivity and energy use to flood control, municipal and industrial water supply, and fish and wildlife management (Xu and Singh 2004). The tremendous importance of water in both society and nature underscores the necessity of understanding how a change in global climate could affect regional water supplies.

There is increasing interest in understanding how the rising concentration of greenhouse gases might affect climate (the mean and variability of temperature, precipitation, humidity, wind and other climate variables over several decades) at local and regional scales in the Nile basin. Being one of the very sensitive sectors, the changes and variability of climate can cause significant impacts on water resources by resulting changes in the hydrological cycle. The change on temperature and precipitation components of the cycle can have a direct consequence on the quantity of evapotranspiration component, and on both quality and quantity of the runoff component (Barnett et al. 2005). Fresh water resources, which are vital to all sectors and regions, are also prone to such direct effects. However, vulnerabilities to climate change vary regionally (Falloon and Betts 2006). East African countries like Ethiopia are particularly likely to experience adverse impacts from climate change due to their topographical settings and poor adaptation capacity (Kotir 2011). Any extreme change, either positively or negatively, on the hydrological and meteorological variables will have great potential to affect regional water resources. Hence, there is an imperious need to quantify the impact of climate change on water resources to support building adaption measures to mitigate the impacts of climate change. The Megech reservoir is one of the proposed projects in the Upper Blue Nile Basin which is included in Tana Sub-basin. The implementation of this reservoir will be expected to minimize the food scarcity from the surrounding area by irrigating 7311 ha of land and solving the domestic water demand of Gonder town. Generally, the introduction of irrigation will make farmers feel more secure about their basic food supply and enable them to diversify their crops based on local market demand and export opportunities. Even though there is huge uncertainty related to climate change in Ethiopia, there are few or no studies conducted on the existing and upcoming climate variability. Although the Megech reservoir is on construction phase, the impacts of climate change under future condition is not investigated yet. Therefore, understanding the impacts of climate change using different climate and hydrological model is more urgent than ever. For that reason, this study mainly deals with evaluating of the impact of climate change on the surface water availability of the Megech watershed. To evaluate the possible impact of climate change on the flow of Megeh river catchment for the future period as compared to the baseline period based on the downscaled climate scenario data using QSWAT model.

Methodology

The Megech River is about 75 km long, has a drainage area of about 850 km² and an average annual discharge of 11.1 m³/s. The catchment area at the dam site is 499.06 km² with a mean annual flow of 5.6 m³/s. The river, which flows generally in a southern direction and drains into Lake Tana, is one of the main streams flowing into Lake Tana from the North. Four major tributaries join the Megech River: two from the right bank and two from the left. The watershed is highly vulnerable to sheet, rill and gully erosions. In the upper catchment of Megech reservoir new gullies which directly run into the Megech river were formed as a result of the increased agricultural activities performed in the watershed. The geographical location of the river bed at the centre of the dam axis is 37°27′36″E and 12°30′10″N latitude with an elevation of 1877 m above sea level. The catchment is situated in the northern portion of the sub-basin with its geographical location between 12°45′45″ to 12°28′30″ latitude and 37°21′ 45″ to 37°37′45″ longitudes with catchment area of 499.06 Km² (Figure 1).

Figure 1. Location map of the study area.

Megech watershed is topographically characterized as a steep mountainous watershed with circular shape and more than 60% of the watershed found on more than 15% slope range. The elevation range in the watershed varies between 1863 and 2966 m a.s.l. The watershed highest elevation is 2,962m in its north-eastern part with mean elevation of 2338.194 m a.s.l (Figure 2).

Figure 2. Digital elevation model (DEM) of Megech watershed.

Megech watershed is covered by eight major soil groups according to FAO/UNESCO classifications; Rhodic Nitisols (27.82%), Eutric Leptosols (25.08%), Lithic Leptosols (17.63%) Haplic and Luvisols15.1%). The soil database is clipped and merged to the eight similar characteristics groups (Figure 3).

Figure 3. Megech watershed soil map.

The land use map of the study area was obtained from Ministry of Agriculture. Spatial distribution and specific land use parameters were required for modelling. QSWAT has predefined land uses identified by four letter codes and it uses to link land use maps to SWAT land use databases in the QGIS interface. The dominant land use of the watershed is rain-fed agriculture and cultivated land is in a various forms. More than three quarter of all land in the watershed has already been brought under cultivation, bush or shrub land, grazing land, forest/wood land and wetland/swap are other land cover types in the watershed (Figure 4).

Figure 4. Land use/land cover map of Megech watershed.

The climate of the Megech watershed is marked by a rainy season from May to October, with monthly rainfall varying from 67 mm in October to 306 mm in July. Mean annual rainfall is about 1,100mm in the upper part and about 1000 mm in the lower part. Rainfall over the Megech watershed is unimodal with nearly 79% of the annual rainfall occurring in the period from June to September. The dry season, from November to April has a total rainfall of about 8% of the mean annual rainfall. Dependable rainfall (85%) varies from less than 1.2 mm during the dry season to 88–225mm/month during the period of June to July/August, equivalent to 55–75% of the average values. Maximum temperatures vary from 23 °C in July to 30 °C in March, whereas minimum temperatures ranges from 11.5 °C in January to 15.6 °C in April & May.

Data collection and analysis

There are three ground based meteorological stations within and around the Megech watershed namely Gondar, Ambagiorgis and Maksegnit. The daily precipitation data for all three gauging stations were collected. These data later prepared in the.dbf format required by the model. The daily temperature were collected similar to the rainfall for the same period. Relative humidity, solar radiation and wind speed data were available only for principal stations or class-1. In this case SWAT can generate data for the rest stations by using weather generator. Precipitation is key meteorology variable that mostly affects the hydrological regime of the study area. Daily precipitation data of six stations for which data was obtained and have less than 30% missing data in the analysis period were used. Data for the missing period was filled using estimation technique. Arithmetic mean and normal ratio are the most commonly used methods for estimation of missing rainfall and temperature data sets.

Simple Arithmetic mean method was used where the mean monthly precipitation of all the index stations is within 10% of the station under consideration (station x) and calculated the missing data by Equation 1(1) $$\mathit{\text{Px}}=\frac{p1+P2+P3}{3}$$(1) . Whereas the mean monthly precipitation of one or more of the adjacent (index) stations differs from that of station x by more than 10% then the normal ratio method was used (Equation 2(2) $$\mathit{\text{Px}}=1/3[P1\frac{\mathit{\text{Nx}}}{N1}+P2\mathrm{ }\frac{\mathit{\text{Nx}}}{N1}+P3\mathrm{ }\frac{\mathit{\text{Nx}}}{3N3}]$$(2) ). (1) $$\mathit{\text{Px}}=\frac{p1+P2+P3}{3}$$(1) (2) $$\mathit{\text{Px}}=1/3[P1\frac{\mathit{\text{Nx}}}{N1}+P2\mathrm{ }\frac{\mathit{\text{Nx}}}{N1}+P3\mathrm{ }\frac{\mathit{\text{Nx}}}{3N3}]$$(2)

Where Px is the precipitation for the station with missed record, P1, P2, P3……P_N are the corresponding precipitation at the index stations N1, N2, N3……N_N are the long term mean monthly precipitation at the index stations and at station x under consideration respectively

Megech is one of the major rivers which contributes significant amount of inflow to Lake Tana. For this reason the Ethiopian ministry of water resources installed gauging station downstream of the river. The gauging station is located near the dam site and measures daily instantaneous flows of the river. Based on the long recorded flow data obtained from Ethiopian Ministry of Water Resources hydrology department, the average daily flow of the river is 11.1 m³/s.

Global Circulation Model (GCM) derived scenarios of climate change were used for predicting the future climates of the study area based on criteria proposed by the Intergovernmental Panel on Climate Change (IPCC). Among the wide range of GCM models HadCM3, (Hadley Centre for Climate Prediction and Research, England), ECHAM4 (Max Plank Institute, Hamburg, Germany), CGCM2 (Canadian Centre of Climate Modelling and Analysis), GFDL_R30 (Geophysical Fluid Dynamics Laboratory & NOAA), and CCSR/NIES (Centre for Climate Systems Research & Japanese National Institute for Environmental Studies) are commonly used. Use of average outputs of different GCMs can minimize the uncertainties associated with each GCMs and can result in plausible future climates for impact studies (Vidal and Wade 2008). However, for this study the HadGCM3 model was selected for the impact assessment investigations. Besides, the HadGCM3 GCM output is chosen since the model is widely used for climate change impact assessment and the results of HadCM3 can be easily downscaled using SDSM (Enyew and Hutjis 2015).

HadGCM3, a climate model developed by the Hadley Centre for Climate Prediction and Research, is commonly employed for hydrological modelling due to several reasons:

1. Global Climate Representation:

   • HadGCM3 provides a comprehensive representation of global climate patterns, including temperature, precipitation, and other climatic variables.

   • Its spatial resolution allows for capturing large-scale climate features that influence hydrological processes.

2. Coupling with Hydrological Models:

   • Researchers often couple HadGCM3 with hydrological models to simulate the impact of climate change on water resources.

   • By feeding HadGCM3’s climate projections into hydrological models, they can assess changes in river flow, groundwater recharge, and other hydrological components.

3. Scenario-Based Assessments:

   • HadGCM3 allows for exploring various climate scenarios (e.g. different greenhouse gas emissions pathways).

   • These scenarios help assess the potential hydrological impacts under different climate futures.

4. Historical Data Calibration:

   • Researchers calibrate hydrological models using historical climate data from HadGCM3.

   • This calibration ensures that the hydrological model accurately represents observed streamflow, evapotranspiration, and other hydrological processes.

HadGCM3 serves as a valuable tool for understanding the complex interactions between climate change and hydrology, enabling informed decision-making for sustainable water resource management

Statistical downscaling methods have several practical advantages over dynamical downscaling approaches. In situations where low–cost, rapid assessments of localized climate change impacts are required, statistical downscaling (currently) represents the more promising option (Chen et al. 2011). Thus, for this study SDSM 5.2-a decision support tool for the assessment of regional climate change impacts was used to downscale large scale predictors and it was freely downloaded from http://www.sdsm.org.uk. SDSM develops statistical relationships, based on multiple linear regression techniques, between large-scale (predictors) and local (predictand) i.e. rainfall and maximum and minimum temperature which were used as an input for hydrological modelling.

Large-scale predictor variable information are freely obtained from the Canadian climate impact scenario group of: http://www.cics.uvic.ca/scenarios/sdsm/select.cgi/. The National Center for Environmental Prediction (NCEP_1961–2001) reanalysis data and HadCM3 predictor variables are obtained on a grid by grid box basis. The watershed is found at a grid box of BOX_12.5N_37.5E. Hence the required predictor data that represents the watershed were downloaded from the nearest average location of the watershed.

Few meteorological stations have 100% complete and/or fully accurate data sets. Handling of missing and imperfect data is necessary for most practical situations. Simple Quality Control checks in SDSM enable the identification of gross data errors, specification of missing data codes and outliers prior to model calibration (Wilby and Dawson 2007). The rainfall and maximum and minimum temperature at Gonder station are comparatively consistent and lower missed values. After visual inspection and correlation analysis was done, Gonder station was used for downscaling climate variables due to their long-term and high-quality data.

Identifying the empirical relationships between gridded predictors and single site predictands (such as station rainfall, maximum and minimum temperature) is central to all statistical downscaling methods and is often the most time consuming step in the process (Wilby and Dawson 2007). The main purpose of the screen variables operation is to assist the user in the selection of appropriate downscaling predictor variables for model calibration. The decision process is also complicated by the fact that the explanatory power of individual predictor variables varies both spatially and temporally. Screening is identifying the downscaling predictors which have high correlation with the actual climate variables. The model correlates each predictands (observed rainfall, maximum and minimum temperature) at Gonder with National Center for Environmental Prediction (NCEP) downloaded predictor data. Annual analysis period was used which provides the predictor-predictand relationship all along the months of the year. As a QG1S plugin, QSWAT is discovered and made available by QGIS when the latter is started. QSWAT has the option of creating a new project or opening an existing project. A project consists of a project file (holding the state of QGIS, the loaded layers, user settings, etc.) and a project directory (where the project input and output files are stored). The main components of QSWAT are

• Watershed delineation,

• Hydrological Response Units (HRU) creation,

• Opening the SWAT Editor to complete input preparation and execute SWAT, and

• Visualization of results.

In order to utilize any predictive watershed model for estimating the effectiveness of future potential management practices one need to select values for the model parameters so that the model closely simulates the behaviour and performance characteristics of the study site. Refsgaard and Storm (1996) distinguished three types of calibration methods: the manual trial and-error method, automatic or numerical parameter optimization method; and a combination of both methods. According to Refsgaard and Knudsen (1996) the first method is the most common, and especially recommended for the application of more complicated models in which a good graphical representation is a prerequisite. For this research work, daily stream flows measured at Megech river gage station from 2001 to 2011 including the first three years warm up period as an input and manual calibration was adopted until the model objective functions reached a satisfactory level (i.e. D < 15%, R² >0.6 and ENS> 0.5).

Downscaling, or translation across scales, is a set of techniques that relate local and regional scale climate variables to the larger scale atmospheric forcing (Wilby et al. 2004). The downscaling approach was developed specifically to address present requirements in global environmental change research and the need for more detailed temporal and spatial information from GCM (Figure 5). Most impact models require information at a sub-grid scale featuring topography, clouds and land use-land cover (Halder et al. 2016). Downscaling bridges the gap between large and local scale climatic data. It tries to link what has been provided by the global climate modelers and what is needed by the decision makers. The translation across scales is based on the assumption that similar synoptic atmospheric patterns produce similar climatic conditions.

Figure 5. Conceptual representation of downscaling.

The Downscaling process plays a crucial role in driving impact assessment models such as drought analysis, water resources management, water demand availability, ecological impacts and risk and vulnerability assessments. A quantitative relationship between circulation and local climate in the form of the equation: (3) $$\mathrm{y}=\mathrm{f}\left(\mathrm{x}\right)\text{is established}$$(3) $$\mathit{\text{local climate response}}\left(y\right)=f(\mathit{\text{external}},\mathit{\text{larger scale}})$$

SDSM uses local surface variables (predictands) from observed meteorological stations and large scale atmospheric variables (predictors) from outputs of GCM runs both at the daily time-scale to develop quantitative relationships between them. These can be expressed in the most general form as (Chu et al. 2010): (4) $$\mathrm{\text{Rt}}=\mathrm{\text{FY}}({\mathrm{X}}_{\mathrm{T}},\mathrm{\theta });\text{For T}\le \mathrm{t}$$(4)

Where

θ = the parameter set,

F_Y = represents a multiple non-linear regression function, that in many cases reverts to a multiple linear regression model,

R_t = the local downscaled predictand variable at a site,

X_t = the large-scale predictor variable (from GCMs).

The procedures to be taken in the SDSM method are outlined in Figure 6 and can be summarized as follows:

Figure 6. Calibration results between observed and generated mean daily maximum temperature (°C) (a) and mean daily minimum temperature (°C) (b), in Gonder station.

• Prepare input data of predictands and GCM predictors at daily time scale

• Screen the most potential predictors.

• Fit the SDSM to reanalysis predictors and observed predictands.

• Drive this fitted SDSM with independent temporally reanalysis predictors.

• Compare the statistical properties of the results of the above steps with those of observed predictands. A good agreement implies that the SDSM can reconstruct the climatology of the observed local variables well, when driven by large scale observation.

• Drive the SDSM fitted to observed time series

• Compare the simulated time series in the above steps with observed time series. Good agreement implies that the GCM predictors can adequately simulate the local climate variables.

Table 1 lists the daily predictor variables held in the UKSDSM data archive. Ideally, candidate predictor variables should be physically and conceptually sensible with respect to the predictand, strongly and consistently correlated with the predictand, and realistically modelled by GCMs. For precipitation downscaling, it is also recommended that the predictor suite contain variables describing atmospheric circulation, thickness, stability and moisture content. In practice, the choice of predictor variables is often constrained by data availability from GCM archives.

Table 1. Model performance ratings based on the range of values for the RSE, NSE and PBIAS.

Performance rating	RSR	NSE	PBIAS (%)
Very good	0.00 ≤ RSR ≤ 0.50	0.75 <NSE ≤ 1.00	PBIAS < ±10
Good	0.50 < RSR ≤ 0.60	.65 <NSE≤ 0.75	±10 ≤ PBIAS < ±15
Satisfactory	0.60 < RSR ≤ 0.70	0.50 ≤ NSE ≤ 0.65	±15 ≤ PBIAS < ±25
Unsatisfactory	RSR > 0.7	NSE < 0.50	PBIAS ≥ ±25

The land phase of the hydrologic cycle is described by the transient water balance equation applied to water movement through the soil, namely: (5) $$\mathrm{\text{SWt}}=\mathrm{\text{Swo}}{\sum }_{i=1}^t(\mathrm{\text{Rday}},\mathrm{i}‐\mathrm{\text{Qsurf}},\mathrm{i}‐\mathrm{\text{Eact}},\mathrm{i}‐\mathrm{\text{Wseep}},\mathrm{i}‐\mathrm{\text{Qlat}},\mathrm{i})$$(5)

Where:

SW_t = the final soil water content after t days (mm water),

SW₀ = the initial soil water content (mm water),

R_day = the amount of precipitation on day I (mm water),

Q_surf = the amount of surface runoff on day i (mm water),

E_a = the amount of evapotranspiration on day i (mm water),

W_seep = amount of percolation and bypass flow exiting the soil profile bottom on day i

Q_gw = the amount of return flow on day i (mm water)t = the time (days).

SWAT provides two infiltration methods for estimating the surface runoff volume component from HRUs, namely, the SCS-curve number (CN) method or the Green & Ampt infiltration method (Kim and Lee 2008). The surface runoff in the SCS curve number method is computed by the following equation: (6) $$\mathrm{\text{Qsurf}}=\frac{\left(\mathit{\text{Rday}}-\mathit{\text{Ia}}\right){}^2}{(\mathit{\text{Rday}}-\mathit{\text{Ia}}+S)}$$(6)

Where:

Q_{surf =} the accumulated runoff or precipitation excess (mm H₂O),

R_day = the precipitation depth for the day (mm H₂O),

I_a = initial abstract which include surface storage, interception and infiltration

S = a retention parameter (mm H₂O). Retention parameter varies due to change in soil, land use, management, and slope and temporally due to change in soil water content.

The initial abstraction is commonly expressed as 0.2S, so that the above equation becomes (7) $$\mathrm{\text{Qsurf}}=\frac{\left(\mathit{\text{Rday}}-0.2S\right){}^2}{(\mathit{\text{Rday}}+0.8S)}$$(7)

Which means that Q_surf > 0, whenever R_day >0.2*S; The retention parameter S is defined by: (8) $$\mathrm{S}=25.4(\frac{1000}{\mathit{\text{CN}}}-10)$$(8)

Where

CN is the SCS-curve number, which ranges from 0 to 100, depending on the soil permeability, land use and the antecedent soil water conditions.

The peak discharge or the peak surface runoff rate q_peak is the maximum volume flow rate passing a particular location during a storm event. SWAT calculates the peak runoff rate q_peak based on a modified rational method as: (9) $$\mathit{\text{qpeak}}=\frac{\mathrm{\alpha }\mathit{\text{tc Qsurf Area}}}{3.6t\mathrm{ }\mathit{\text{conc}}}$$(9)

Where

q_peak = the peak runoff rate (m²s⁻¹),

α_tc = the fraction of daily precipitation that occurs during the time of concentration,

Q_surf = the surface runoff (mm),

Area = the sub-basin area (km²)

t_conc = the time of concentration (hr.), and 3.6 is a conversion factor.

The SWAT model has three options of potential evapotranspiration estimation. These are Hargreaves, Priestley-Taylor and Penman Monteith. The choice of a method among the three estimation methods is based on the data availability. For this study the Hargreaves PET method, which requires only temperatures – available from the downscaling simulations was adopted. The Hargreaves method is based on the following equation: (10) $$\mathrm{\lambda \lambda }\mathrm{E}0=0.0023\mathrm{H}0\left(\sqrt{\mathrm{\text{Tmax}}‐\mathrm{\text{Tmin}})}\right(\mathrm{\text{Tavg}}+17.8)$$(10)

Where

λ is the latent heat of vaporization (MJ/kg),

E_o is the potential evapotranspiration (mm/d),

T_max is the maximum air temperature for a given day (°C),

T_min is the minimum air temperature for a given day (°C), and

T_avg is the mean air temperature for a given day (°C).

H₀ is the extraterrestrial radiation (MJ/m²/d), which depends on the latitude and the day of the year and is computed using the solar constant and complex solar declination relationships as presented in the SWAT manual.

Once the potential evapotranspiration is determined, the actual evapotranspiration can be computed by considering transpiration from the plant canopy interception, transpiration from plants and evaporation from the soil. A simplified formula is given by Mekonen et al. (2009) for computation of actual evapotranspiration (11) $${\mathrm{E}}_{\mathrm{a},\mathrm{i}}=\mathrm{\text{min}}\left({\text{SW}}_{\mathrm{o}},{\text{LAIxE}}_{\mathrm{o}}\right)$$(11)

Where:

LAI = the leaf area index which can be taken as 0.6 during crop periods and 1 otherwise

The SWAT Groundwater flow contribution to a total streamflow is simulated by creating a shallow aquifer storage which is recharged by the water percolated from the bottom of the root zone (Arnold et al. 2000). The other component of groundwater is the deep, confined aquifer which doesn’t add to streamflow rather it contributes to flow outside the watershed. Water that percolates through the soil layers flows through the vadose zone before becoming the shallow or deep aquifer recharge. An exponential decay weighting function is employed to compute the recharge to both aquifers on a given day as: (12) $$\text{Wrchrg},\mathrm{i}=(1-\hspace{0.17em}\mathrm{\text{exp}}\hspace{0.17em}(-1/\mathrm{\delta }\text{gw}))*\text{Wseep}+\hspace{0.17em}\mathrm{\text{exp}}\hspace{0.17em}(-1/\mathrm{\delta }\text{gw})*\text{Wrchrg},\mathrm{i}-1$$(12)

WhereW_rchrg,i = the amount of recharge entering the aquifers on day i(mm)δ_gw = the delay time or drainage time of the overlying geologic formations (days)W_{rchrg, i−1} = the amount of recharge entering the aquifers on day i–−1 (mm).w_seep = the total amount of water exiting the bottom of the soil profile on day i (mm).

Measures the average difference between the simulated and measured values for a given quantity over a specified period (usually the entire calibration or validation period in the study). The percent difference can be calculated using the following equation. (13) $$D=100\left(\frac{\sum _{i=1}^n\mathit{\text{qsi}}-\sum _{i=1}^n\mathit{\text{qoi}}}{\sum _{i=1}^n\mathit{\text{qoi}}}\right)$$(13)

Where,

D = the percent of difference,

q_si = the simulated values in each model time step,

q_oi = the measured values in each model time step. Value close to 0% is best for D.

Regression coefficient (R²) is the square of the Pearson product–moment correlation coefficient and describes the proportion of the total variance in the observed data that can be explained by the model. The closer the value of R² to 1, the higher is the agreement between the simulated and the measured flows. It is calculated as using the following equation: (14) $${\mathrm{R}}^2={x\left\{\frac{n\left\{\sum _{i=1}^n(\mathit{\text{oi}}-o)(\mathit{\text{Si}}-S)\right\}\left(n-1\right)}{{⌈\sum _{i=1}^n{\left(\mathit{\text{oi}}-o\right)}^2⌉}^{0.5}\left|{⌈\sum _{i=1}^n{\left(\mathit{\text{Si}}-S\right)}^2⌉}^{0.5}\right|}\right\}}^2$$(14)

Where

n = the number of compared values

O_i = observed data, $\overline{\mathrm{O}}$ is observed mean

S_i: = simulated data, $\overline{\mathrm{S}}$ is simulated mean

A value of 1 for NSE indicates a perfect match between simulated and observed data values. The optimal value for PBIAS is 0. 0, negative for overestimation bias and positive for underestimation bias of a model (Moriasi et al. 2007). Performance of the model for calibration and validation can be evaluated using performance rating guide line (Table 1).

Results and discussions

The observed daily precipitation at Gonder station showed that the maximum amount of precipitation occurred within the major parts of the rainy season (June-September) and the maximum daily precipitation were recorded in May and October respectively. Mean annual maximum temperature shows positive difference is higher than the negative difference. The test conclude that it’s increasing through observed time period. Mean annual maximum temperature shows increasing and decreasing trend in observed time period. The maximum mean annual temperature observed in 1998 was 27.4 °C and the minimum observed in 1985 and 1989 is 25.8 °C with a positive trend value (0.86 °C) In the year 1989 the lowest (11.2 °C) and in 1995 the highest (14.4 °C) average monthly minimum temperature was observed. The trend line of minimum temperature showed significant increments by 2.5 °C. The results of comparisons of observed data among different stations for quality control, the calibration and validation of SDSM model, selection of predictor variables and scenario generations for base period and for the future time period were discussed. Up to six predictors were selected at a time and analyzed to investigate the percentage of variance explained by specific predictand–predictor pairs, then those predictors that have high explained variance was selected (Table 2).

Table 2. List of predictor variables that have better correlation with the predictands at Goder station with significant level of less than 0.05 (p < 0.05).

Predictands	Predictors(NCEP reanalysis data)	Notation	Partial r
Precipitation	surface divergence	ncepp_zhaf dat	0.073
	near surface relative humidity	nceprhumaf dat	0.068
	surface meridional velocity	ncepp_vaf.dat	0.047
	500hpa wind direction	ncepp5thaf dat	0.026
Maximum Temperature	mean Temperature at 2 M	nceptempaf. dat	0.468
	850 hpa meridional velocity	ncepp8_vaf.dat	0.446
	500pa zonal velocity	ncepp5_uaf.dat	0.230
	surface specific humidity	ncepshumaf. dat	0.135
Minimum temperature	mean Temperature at 2 M	nceptempaf. dat	0.405
	850hpa geopotential height	ncepr850af. dat	0.277
	surface air flow strength	ncepp_faf. dat	0.113
	surface wind direction	ncepp_thaf. dat	0.039
	850 hpa meridional velocity	Ncepp8_vaf.dat	0.026

The partial correlations indicate that on average surface divergence ncepp_zhaf. dat has the strongest association with local precipitation once the influence of all other predictors has been removed, whereas maximum and minimum temperatures are strongly correlated with mean temperature at 2 m, which shows their heavy dependence on regional temperatures.

Calibration and validation

It develops multiple regression equations between the metrological station data (predictand), regional scale and atmospheric (predictor) variables. Then calibration of climatic predictors from 1971 to 1990 periods followed, coefficient of determination R² that measures the percentage of the explained variance between modelled and observed variables is 45%, 83% and 62% for precipitation, maximum and minimum temperature. The weather generator generates twenty ensembles of synthetic daily weather series by using calibration. PAR file and NCEP-reanalysis data from the 1/1/1991 to 31/12/2000. Then mean of those 20 ensembles generated by weather generator for remaining 11 years used for validation gives 44%, 77% and 57% (Table 3).

Table 3. Calibration and validation R² values of SDSM 5.2 downscaled for precipitation, maximum and minimum temperature.

Predictands	R^2
	Calibration (1971–1990)	Validation (1991–2000)
Precipitation	0.45	0.44
Maximum temperature	0.83	0.77
Minimum temperature	0.62	0.57

Maximum and minimum temperature values give a better R² values, inferring that future projections would also be well replicated. But the result of precipitation is unsatisfactory because of its conditional nature. For maximum temperature, the model slightly over estimates in September and slightly under estimates in June. While in the remaining months the generated and observed values are more or less similar. But estimates of the minimum temperature between generated and observed values showed similar patterns (Figure 6).

The calibration results of precipitation in Figure 7, shown the model overestimates in some months However, the overall agreement was good.

Figure 7. Calibration results between observed and generated mean daily precipitation (mm).

Validation was done using an independent observed data for the period of ten year from 1991 to 2000. Here also twenty ensembles (runs) of daily values were generated and the averages of these ensembles were taken for the comparison. A good agreement was also found between the observed and simulated precipitation (R²= 0.44) though it is a conditional process and for minimum temperature (R² = 0.57) maximum temperature R² = 0.77) during validation. However, like the calibration, the validation results also found to be in good agreements between generated and observed precipitation, maximum and minimum temperature which proves the SDSM model ability to generate weather variables in Megech watershed under a given set of model parameters (Figure 8).

Figure 8. Validation results between observed and generated of Gonder station.

One of the criteria commonly used in evaluating the performance of any useful method is whether the historic (observed) condition can be replicated or not. For this study, the HadCM3A2a and HadCM3B2a were the two GCM output files used for the scenario generation. The regression weights produced during the calibration process were applied to the time series outputs of the GCM model based on the assumption that the predictor-predictand relationships under the current condition remain valid for future climate conditions. Twenty ensembles of synthetic daily time series data were produced for each of the two SRES scenarios for a period of 1961 to 2099 (139 years). Finally, the mean of the twenty ensembles for the specified period was produced for maximum and minimum temperature and precipitation. The climate scenario for the future period was developed from statistical down scaling using the GCM predictor variables for the two emission scenarios for 100 years based on the mean of 20 ensembles and the analysis was done based on 30 years period from 2011 to 2040, 2041 to 2070 and 2081 to 2099. The IPCC recommends 30 years of 1961–2000 as a climatological base period in impact assessment. Hence, for this research the period from 1971 to 2000 was taken as a base period with in which the comparison was made.

For the sake of comparison with the observed values, the generated values of the base period were averaged to a monthly time step. Maximum temperature the downscaled monthly, seasonal and annual average maximum temperature resulted in good agreement with the observed temperature for the baseline period in both A2a and B2a emission scenarios except slight variations but the overall agreement confirms that the SDSM is able to produce synthetic maximum temperature series given daily atmospheric predictor variables supplied by a GCM (either under present or future greenhouse gas forcing) (Figure 9).

Figure 9. Pattern of observed and downscaled mean monthly and annual maximum temperature.

The downscaled mean monthly minimum temperature for the emission scenario for the base period are shown in the Figure 11. The downscaled monthly annual mean minimum temperature resulted a good agreements with the observed mean monthly, and annual minimum temperature in both A2a and B2a emission scenarios with insignificant variations but the general trend shows similar patterns. As compared to the minimum temperature, a better agreement was found in maximum temperature (Figure 10).

Figure 10. Pattern of observed and downscaled mean monthly and annual minimum temperature.

As compared to the minimum and maximum temperature the precipitation could not able to replicate the observed data. This is due to complicated nature of precipitation processes and its distribution in space and time. This might be due to the conditional process (dependent on other intermediate processes like on the occurrence of humidity, cloud cover, and/or wet-days) and high variability in space and so the relatively coarse spatial resolution of the current generation of climate models is not adequate to fully capture that variability. As it can be seen from the Figure 11 the mean monthly generated value in both A2a and B2a scenario is over estimated. The largest monthly variation occurs from month of June to September which is main rainy season (June–September) the generated values for both scenarios are over estimated. These indicate that the SDSM model over estimates the precipitation in the major parts of rainy season. The annual mean values are also variations with the observed value.

Figure 11. Pattern of observed and downscaled mean monthly and annual precipitation.

The mean annual change in maximum temperature for the future period (2001–2099) for both A2a and B2a emission scenarios are shown. As it can be seen from Figure 4, 9a and b, the overall results (2001–2099) (Figure 12).

Figure 12. Change in average monthly maximum temperature in the future (2001–2099) for A2a.

Monthly and annual increase in maximum temperature for A2a scenario is greater than B2a scenario because A2a scenario represents medium to high scenario which produces more CO₂ concentration than B2a scenario which represents low to medium emission scenarios (Figure 13).

Figure 13. The mean annual change in maximum temperature for the future period (2001–2099).

Figure 14 Shows the mean monthly and annual changes of minimum temperature for the future period (2001–2099) of the two emission scenarios (A2a and B2a) as compared to the base period. The minimum temperature shows increasing trends for the future periods (2001–2099) for both A2a and B2a emission scenarios except in January and March which the temperature shows insignificant decreases from the period 2001–2070 with a range of 0.01–0.05. But in March (2071–2090) minimum temperature shows slight increment of 0.1^oc in A2a emission scenarios. The increment of minimum temperature varies from 0.01^oc in February (2020s) to 1.7^ocin March (2080s) of A2a emission scenarios. The overall results (2001–2099) for annual mean minimum temperature increases for both A2a and B2a emission scenarios with a significant increment of 0.6^oc on 2080s of A2a emission scenarios.

Figure 14. Future trend of minimum temperature change for the two scenarios a) HadCM3B2a. (b) HadCM3A2a.

As shown in the Figure 15 there are high fluctuation of minimum temperatures for both A2a and B2a scenarios. However, the overall average annual minimum temperature shows a slight increasing trend line for the future time period.

Figure 15. Trend of downscaled minimum temperature (1971–2099) for both B2a and A2a emission scenarios.

In the case of precipitation the results of this study was thus in line with the previous researches done on Tana basin, upper Blue Nile basin (Melke 2015). On monthly basis, the percentage change in precipitation is not systematic i.e. precipitation increases in some months and decreases in other ones. The percentage changes in precipitation increases in June on the period 2080s for A2a scenarios in which the precipitation increases by 5.6% and 5.2% and 5.1% in 2020s and 2080s respectively. The decrease in precipitation reaches to a maximum of 10.1% (2090s) for A2a scenario and 5% (2040s) for B2a scenario (Figure 16).

Figure 16. Percentage change monthly precipitations in the future (2001–2099) for A2a and B2a scenario from the base period.

The overall results of the period 2001 to 2099 for annual percentage change in precipitation showed an increasing trend for both scenarios (A2a and B2a). The increment of annual percentage change of precipitation 3% (2080s) for B2a scenario and for A2a scenario the increment ranges between 2.2% (2080s). The generated future scenarios for annual precipitation generally show an increasing trend with respect to the base period (Figure 17).

Figure 17. Trend of downscaled annual precipitation (1971–2099) for both B2a and A2a emission scenarios.

The calibration of stream flow was conducted depending on the sensitive parameters which were demonstrated as influential variables on the simulated water balance by the model. The parameters which has influence on the simulated flow were taken in to consideration from sensitivity analysis result. Using the most sensitive parameters, the 9 were used in calibration. Calibration of stream flow was carried out at the outlet of sub-basin 15 (Near the dam gauging station). Observed daily stream flows were adjusted on the monthly basis and Simulations runs were conducted on monthly basis to compare the modelling output with the measured daily discharge (Table 4).

Table 4. Flow sensitive parameters and their category of sensitivity.

Parameters	Code	Lower and upper bound	Initial value	Calibrated value
Initial SCS CN II value	CN2	±25%	default	+25%
Manning’s ‘n’ value for the main channel	CH_N2	0.01–0.3	0.1	0.15
Threshold water depth in the shallow aquifer for flow (mm)	GWQMN	0–5000	0	100
Groundwater delay (days)	GW_DELAY	0–500	30	100
Channel effective hydraulic conductivity (mm/hr)	CH_K2	0–150	50	100
Base flow alpha factor (days)	ALPHA_BF	0–1	0.48	0.1
Soil evaporation compensation factor	ESCO	0–1	0	1
Groundwater ‘revap’ coefficient	GW_REVAP	0.02–0.2	0.02	0.02
Threshold depth of water in the shallow aquifer for ‘revap’ to occur (mm).	REVAPMN	0–500	1	500

The calibration results in Figure 19 show that there is a good agreement between the simulated and measured monthly flows by the coefficient of determinations (R² = 0.77) and the Nash-Sutcliffe simulation efficiency (ENS = 0.76) and thus fulfilled the requirements suggested by Santhi et al. (2001) for R² > 0.6 and ENS > 0.5. There is good agreement between monthly observed and simulated flows at Megech River gauge station. The coefficient of determinations (R²) and Nash Sutcliffe simulation efficiency (ENS) between the observed and simulated monthly flow were found to be 0.77 and 0.73 respectively and these shows a very good correlation of the simulation results with the observed values. Furthermore, Figure 18 shows the hydrograph of monthly observed and simulated flows.

Figure 18. Time series of simulated and observed monthly Megech River flow for calibration period.

Thus, the validation check illustrates the accuracy of the model for simulating time-periods outside of the calibration period. The model performed as good in the validation period (2001–2005), as for the calibration period (1993–2000) at Megech gauge station as indicated in Table 4. Hence, the set of optimized parameters used during calibration process can be taken as the representative set of parameter to explain the hydrologic characteristic of the Megech watershed and further simulations using SWAT model can be carried out by using these parameters for any period of time (Figure 19).

Figure 19. Time series of simulated and observed monthly Megech River flow for the validation period.

As it can be seen from Table 5 the overall results (2001–2099) for average annual total flow showed a decreasing trend for both A2a and B2a scenarios as compared to the base period except B2a 2020'S scenario. The decline of total average annual flow volume ranges from 3.6 (A2a2080'S) to 0.3 (A2a2020'S) the increment ranges is only 0.8 in (B2a2020'S). SWAT model simulation result using HadCM3A2a and HadCM3B2a scenarios inputs indicated annual declination of flow by 4.4%, −4.1% and −18.9% 2050s and 2080s on HadCM3B2a −1.6%, −5.3% −20% (Figures 20 and 21).

Figure 20. Monthly change of Megech River flow HadCM3A2a scenarios comparing to the observed period.

Figure 21. Monthly change of Megech River flow HadCM3A2a scenarios comparing to the observed period.

Table 5. Monthly and annual average discharge.

Month	Base period	A2a			B2a
		2020'S	2050'S	2080'S	2020'S	2050'S	2080'S
January	2.9	4.9	6.9	4.8	8.4	7.0	10.6
February	3	3.7	4.0	8.2	6.6	7.2	7.2
March	3.8	4.9	8.3	6.4	7.4	10.5	6.5
April	4.3	3.7	3.9	4.0	6.2	6.5	7.0
May	8.6	13.6	12.0	8.6	17.4	10.0	8.9
June	20.7	26.4	25.1	19.5	24.1	16.5	14.5
July	52	31.3	28.7	25.4	27.3	23.7	21.4
August	60.8	30.0	30.5	29.0	27.4	27.7	27.1
September	33.6	30.4	28.2	23.1	29.2	28.5	25.4
October	17.3	28.4	27.1	21.4	27.9	27.8	22.2
November	6.2	20.5	18.7	14.3	25.7	24.1	14.5
December	3.6	15.5	12.0	8.6	18.7	18.3	10.5
Annual	18.1	17.8	17.1	14.5	18.9	17.3	14.7

Conclusion

Human beings, at the larger scale, are vulnerable to climate change as recurrent floods and droughts continue to bring misery and keep on claiming many lives all over the world. Regarding to the integrated water resources management in basins development approach, it is required to make indicative predictions of the impacts of climate change on the spatial and temporal variability of precipitation and temperature and their influences on stream flows. Understanding the problem is key for the solution and predicting the level of climate change impact on water resources as a prerequisite for planners and decision makers to reduce, prevent and/or to find the possible adaptation measures. Hence, the impact of climate change on Megech River was carried out to address part of the global problem by showing the possible indicative predictions of climate changes. The study confirmed that the Statistical downscaling Model (SDSM) is able to simulate climatic events. The calibration and validation results of SDSM showed that the model is able to simulate the climatic variables (precipitation and temperature) which follow the same trend with the observed one. Even though, the precipitation is a conditional process and high special variability, the overall result from SDSM is well correlated with the observed precipitations. Hence, SDSM can predict the future climatic events under changing conditions based on the assumption that the predictor-predictand relationships under the current condition remain valid for future climate conditions. The hydrological model results of calibration and validation indicate that SWAT model is able to accurately explain the hydrological characteristic of Megech watershed. Flow results show that average monthly and annual inflow volume changes mainly corresponding to the change in precipitation except in the month of September and Bega season in which the flow volume increases even though the precipitation decreases.

Authors’ contributions

SA has made substantial contributions in conception design, acquisition of data, interpretation of results and leading the overall activities of the research. He has given also the final approval of the version to be published. AT and HL contributed in designing, data collection and analysis of this research. All authors read and approved the final manuscript.

Competing interests

The authors declare that they have no competing interests.

Consent for publication

All authors read the manuscript and agree to publication.

Ethics approval and consent to participate

The authors hereby declare that, this manuscript is not published or considered for publication elsewhere.

Acknowledgements

This study would never be completed without the contribution of many people to whom we would like to express our gratitude. The administrative kebele’s development agents, district agricultural officials, local youths, in each of the sampling sites were indispensable for the successful completion of the field work and ground truth.

Availability of data and materials

The dataset supporting the conclusions of this article is included within the article.

Disclosure statement

No potential competing interest was reported by the authors.

Additional information

Funding

Funded by Bahir Dar University and Abay Basin Authority.