Effect of climate change on streamflow in the Gelana watershed, Rift valley basin, Ethiopia

Water problems have occurred for ages and nowadays water managers and planners are facing substantial uncertainties in the future demand and availability of water (Kuma et al. 2021). Population growth coupled with climate changes and their probable hydrological effects are increasingly contributing to these problems. In this case, serious action is essential for the safety of water resources. Climate change is one of the biggest problems, ranking as the second-highest risk to livelihoods as stated in the Global Risk Report (Khoi et al. 2021). It is believed to have led to the changes in global patterns of the water cycle and caused the redistribution of water resources in time and space (Shiferaw et al. 2018). Its problem is global stress, for instance, in the study conducted by Chou et al. (2014); the warming starts in central and southeastern Brazil and progresses strongly toward the north. According to Zhang et al. (2018), the effects of the future climate change on streamflow dynamics of the Grand and Thames rivers in Canada will result in heterogeneity in streamflow. These problems, coupled with the high dependence of African economies on agriculture and the direct consumption of natural resources, will have dramatic negative consequences for the economy (IPCC 2014). Currently, one-fourth of the population of Africa is facing high water stress, the magnitude of which is expected to increase in the next 40 years, and the adaptation mechanisms developed by farmers are not sufficient to cope with current and future climate variability (IPCC 2013; Wagesho et al. 2013).

Ethiopia, which is situated in the Horn of Africa, has a large population with a fast growth rate and an extremely diverse climate ranging from an equatorial rainforest in the south to the desert-like conditions in the northeast. The economy is heavily reliant on agriculture, which supports 41% of the national income and 80% of the workforce, undoubtedly making it a victim of this global challenge (Tarekegn et al. 2021).

Ethiopia is often called the water tower of East Africa. However, according to the Global Climate Risk Index (2015), it is one of the 10 African countries most affected in the period 1996–2015. However, Ethiopia's contribution to global greenhouse gas emission is limited (Ghebrezgabher et al. 2016). The increment in temperature may increase the water demand that affects the water resource availability (Geleta & Gobosho 2018; Dibaba et al. 2020). The western, northern, southern, and central parts of the country are characterized by an annual decline (Wagesho et al. 2013). Accordingly, crop failure and malfunction of existing water infrastructures are the most common problems in Ethiopia.

The previous studies show that many parts of the Ethiopian Rift valley basin are sensitive to climatic changes. In the recent study conducted by Shiferaw et al. (2018) in the Ilala watershed, northern Ethiopia, the minimum and maximum temperature increases and high evapotranspiration loss reduced the surface runoff and effect on groundwater recharge in the watershed for the future periods. A similar conclusion was made by Legesse et al. (2015) in the hydrology regime of the Didessa catchment. The climate variation results in changes in streamflow in Lake Ziway catchments (Gadissa et al. 2019). The increment in temperature may increase the water demand that affects water resource availability in the Finchaa watershed (Geleta & Gobosho 2018). Moreover, the study carried out by Dibaba et al. (2020) supported these reports. Similarly, the climate change impact was assessed in the Bilate watershed and it was seen that the average total annual flow at the outlet of the watershed will decrease up to 2041, whereas for the 2041–2070 periods it will increase (Wagesho et al. 2013).

The soil and water assessment tool (SWAT) model is computationally effective for assessing the relationships between the climate change and the watershed hydrological process (Shiferaw et al. 2018; Tarekegn et al. 2021).

The principle of parameter regionalization is to provide effective hydrological information from gauged to ungauged catchments (Guo et al. 2021). The principal component analysis (PCA) is the best method to get a better correlation and group the optimized parameters in physically significant components. The multiple regression techniques can apply in modeling the hydrological responses (surface runoff) from the watersheds (Sharma et al. 2015).

The hydro-meteorological time-series trend analysis for the future is necessary. The Mann–Kendall (MK) test analysis is more applicable and employed by many researchers. It is distribution-free and tolerates missing values (Belihu et al. 2018; Orke & Li 2021).

The Lake Abaya water level depth is shallow compared to other Rift valley basin lakes with a maximum depth of 13.1 m, due to the fact that it is highly sensitive to change in volume. The lake water level fluctuation results from water exchange characteristics within its tributary watersheds in response to climatic and hydrological factors within natural amplitudes. The Gelana river is a tributary of the Lake Abaya sub-basin, Rift valley basin. Its watershed area experiences drastic changes over the year due to deforestation of vegetation and being replaced by cultivated land. Most of the time there is a sudden occurrence of flooding in the Gelana watershed area, resulting in loss of life and damage to agricultural products and properties of the community. It is also characterized by fast population growth rates; irrigated agriculture is expanding from time to time, and demand for water supply is increasing. These factors, in combination with the future climate change, impact the result in water stress within and around the Gelana watershed. Understanding the general trends of the future climate variables, such as precipitation and maximum and minimum temperatures, is necessary to minimize the adverse risks. However, the effect of climate change on hydrology and streamflow was not investigated in the Gelana watershed. Therefore, the main objective of this study is to evaluate and quantify the effects of climate change on streamflow in the Gelana watershed by using the regional climate model (RCM) and the SWAT under two representative concentration pathways (RCPs). The study plays an important role in allowing planners, decision-makers, and other concerned parties to integrate their duties regarding climate change within and around the watershed.

The Gelana river is located in the Lake Abaya sub-basin, Rift valley basin in the Southeastern part of Ethiopia within the Oromia Regional State (Gelana district) and South Nation Nationality Peoples Regional States (Amaro special district). The watershed is located at 5° 25′–6° 18′ North latitude and 37° 50′–38° 20′ East longitude about 450 km south of Addis Ababa as shown in Figure 1. The Gelana River originates from the Yirga Chefe area, Gedeo zone that drains into and eventually meets the southeastern side of Lake Abaya with a drainage area of about 3,364.6 km².

Terrain data of the Gelana watershed digital elevation model (DEM) (12.5 × 12.5 m resolution) were downloaded from the Alaska satellite facility (https://asf.alaska.edu/) and used to obtain spatial information, to delineate the watershed, sub-basins, and to calculate its parameters. In addition, it was used for slope classification. The SWAT model requires the climate data for streamflow simulation. The observed climate data (maximum and minimum temperatures, rainfall, relative humidity, wind speed, and solar radiation) from 1987 to 2019 were collected from the Ethiopian National Meteorological Service Agency. The outcome of data analysis depends on the quality and completeness of data. However, missing data are a common problem in hydrological studies; so, climate data filling and reconstruction are vital (Hamzah et al. 2021). Therefore, the average and normal ratio methods were used to fill the missing meteorological data records because of their computational simplicity and significant precision. The double-mass curve analysis was used for checking the consistency of recorded data. The homogeneity of the selected station's rainfall records was checked by the non-dimensional parameterization method. Streamflow data are necessary for calibration and validation of the SWAT model and used to evaluate the climate change effects. The Gelana streamflow gauged at Tore and Yirga Chefe stations from 1980 to 2015 was collected from the Ministry of Water, Irrigation, and Energy (MoWIE). Similarly, the Gidabo, Kulfo, and Hare river streamflow data were collected and used for regionalization purposes.

Land use land cover (LULC) and soil map were used to divide the basin into sub-basins based on the hydrological response unit (HRU) and also used to classify vegetation types and soil types that have an impact on the streamflow and water balance components of the study area. LULC of 2020 was classified based on supervised classification by using Google Earth images and field survey by checking with many ground control points using GPS considering each land-use type. The classification process and analysis of the different LULC classes were done using single Landsat satellite images covering the Landsat 8 OLI/TIRS (path 168, rows 56) with a 30×30 m grid cell. The Landsat images were downloaded from the United States Geological Survey (USGS) website using earth explorer (https://earthexplorer.usgs.gov/). The selection of the Landsat satellite image dates was influenced by the quality of the image especially for those with limited or low cloud cover. Thus, the images were downloaded with less than 10% cloud cover.

The Landsat image was then processed in ERDAS IMAGINE 2014 software. According to the LULC classification in 2020, nine major land uses and land cover types were identified in the Gelana watershed. The classes are coded by four letters according to the SWAT database as Agricultural land (AGRL), Forest-Evergreen (FRSE), Shrubs land (RNGB), Bare land (BARR), Forest-mixed (FRST), Settlement (URBN), Grass land (PAST), Wetland (WETL), and Water bodies (WATR). However, the major part of the watershed was covered by agricultural land which covered about 1,336.23 km² (39.70%) of the watershed area, and the lowest part of the watershed was covered by water body which accounts for about 0.19 km² (0.01%) of the watershed from the whole study area.

After the image was classified, the accuracy assessment of classifications by generating a set of random 400 points was done in ArcGIS 10.4. Then, the value of each random point was identified from the Google Earth image with the help of some of the field survey points. The Kappa coefficient can be used as a measure of agreement between classification by model and reality. If Κ≤0.5 shows rare agreement, 0.5≤Κ≤0.75 shows a medium level of agreement, 0.75≤Κ≤1 shows a high level of agreement, and Κ=1 shows perfect agreement (Dibaba et al. 2020). The result of accuracy shows that the total (overall) accuracy of land use and land cover is 89.25% and the Kappa coefficient (Κ) is 87.8% which shows a high level of agreement, so it is acceptable.

The SWAT model requires a basic physical–chemical property of the soil types such as the texture, chemical composition, physical properties, and available moisture content, hydraulic conductivity, bulk density, and organic carbon content of different layers (Dibaba et al. 2020). Its simulation depends on the soil data to determine the variety of hydrological characteristics found in each sub-basin within the watershed. The soil map of the study area was obtained from the MoWIE that was extracted by using the Gelana watershed shape file. Based on the dominant characteristics, the soil of the study area was classified into four major groups: Humic Nitisols (45.02%), Chromic Luvisols (25.64%), Eutric Fluvisols (18.52%), and Eutric Vertisols (10.82%). From the identified soil, Humic Nitisols is the dominant soil type covering an area of 1,514.90 km², whereas the Eutric Vertisols soil type covered the lowest area (363.16 km²) of the watershed.

The collected data were effectively analyzed by using software and materials listed as follows: Arc-GIS was used to obtain the physical parameters and spatial information of the watershed, to generate the climate data from CORDEX-Africa to the watershed, and to assess the LULC accuracy. Arc-SWAT 2012 was used to delineate the Gelana watershed, assess the hydrology/water balance components, and simulate the streamflow for present and future periods. SWAT-CUP 2012 was used for sensitivity analysis, calibration, validation, and uncertainty analysis of the SWAT model. Google Earth was used to handle coordinate and elevation data, and provide viewing and conversion to ERDAS software. ERDAS IMAGINE 2014 was used to stack the satellite images and to classify the 2020 land use/land cover of the Gelana watershed. CMhyd software was used to extract and correct the bias of climate data (rainfall and maximum and minimum temperatures) obtained from two RCMs. RAINBOW software was used to check the homogeneity of the hydrological data (the Gelana streamflow gauged at Tore and Yirga Chefe stations). IBM SPSS statistics software was used to check the correlations and covariances of the physical characteristics of the catchments, to develop the PCA, and stepwise regression equation. PCPSTAT and DEW02 were used to prepare and arrange the weather data generator algorithm in the SWAT model.

RCMs were used due to high resolution and resolved information and their ability to model atmospheric processes and land cover changes explicitly, and many variables are available with better representation of some weather extremes than in global climate models (GCMs). The downscaled RCM data were obtained from https://esgf-data.dkrz.de/search/cordex-dkrz. The CORDEX-Africa model data at longitude 0.44° and latitude 0.44° horizontal resolution and a multi-model ensemble of RCMs with their driving GCMs provide the boundary conditions (Kuma et al. 2021).

Africa was nominated as the target area of CORDEX for three major reasons. These are the high vulnerability of this region in many sectors resulting from climate variability, the relatively low adaptive capacity of its economies, and significant changes in temperature and precipitation patterns (Giorgi et al. 2009; Dibaba et al. 2019). The climate variables were extracted from the output of CORDEX-Africa RCMs by using Arc-GIS to the study area by using the meteorological data gauging station's latitude and longitude coordinates. In addition, based on the vintage, resolution, validity, and representativeness, among the several RCMs, the RACMO22T and RCA4 models were selected as recommended by Geleta & Gobosho (2018) and Dibaba et al. (2019).

The driving forces for future expansion of emissions of substances are technological innovation, energy choice, socioeconomic and demographic growth, and their relationships. Therefore, the Gelana watershed is characterized by the deforestation of vegetation and replaced by cultivated land or agricultural activities due to fast population growth rates, the lower incomes in our country and technology enhancements such as industry development and burning of fossil fuel. Similarly, in other countries that will increase the carbon dioxide concentration, because it is earth climate system, the climate change and warming of temperature is a global circulation and its effect is not bounded in the unique region. By considering these pre-requisites according to Vuuren et al.’s (2011) explanation, the RCP4.5 and RCP8.5 scenarios were used for the study.

The RCM data of rainfall and minimum and maximum temperatures were extracted by ArcGIS 10.4.1 software through a multidimensional tool and a NetCDF table view. The grid point data were extracted, which are the nearest for observed meteorological stations by using the latitude and longitude of observed stations.

The performance of the selected two RCMs (RACMO22T and RCA4) to simulate the daily observed precipitation and minimum and maximum temperatures of the Gelana watershed was statistically evaluated in the base period. Similarly, work is done for all bias correction methods. The four statistical approaches are Pearson correlation coefficient (R), root mean square error (RMSE), Nash–Sutcliffe efficiency (NSE), and percent bias (PBIAS). The RMSE is a measure of the absolute error of the climate model in simulating certain climate variables. The smaller the RMSE, the better the model and vice versa. The Pearson correlation coefficient is a measure of the strength of the relationship between model simulations and observations and has the limits of 1 and −1. The Pearson correlation coefficient of 1 indicates a perfect positive correlation between model simulation and observed data, while −1 indicates a perfect negative correlation between the two. Typically, the value of greater than 0.5 is taken as satisfactory. The NSE indicates how well the simulation matches the observation and it ranges between and 1, with NSE = 1 meaning a perfect fit. The higher this value, the more reliable the model is in comparison to the mean (Fang et al. 2015). The PBIAS measures the average tendency of the simulated data to their observed counterparts. Positive bias values indicate overestimation, while negative bias values indicate underestimation by the climate model. The optimal value of the PBIAS is 0.0, with low-magnitude values indicating accurate model simulations (Fang et al. 2015; Geleta & Gobosho 2018; Hamzah et al. 2021).

The results indicate that the RACMO22T climate model performed better than the RCA4 model in reproducing daily precipitation and the maximum and minimum temperatures in the Gelana watershed. The related result was seen by Dibaba et al. (2019). However, climate modeling by using multi-RCMs helps to minimize uncertainties comparatively to the use of single RCMs (Dibaba et al. 2020; Kuma et al. 2021). Using an ensemble mean of climate models was preferred in the analysis of rainfall and temperature (Kuma et al. 2021). Therefore, both the RACMO22T and RCA4 climate models were selected for this study to model climate variables and streamflow.

The bias corrections of RCM (RACMO22T and RCA4) simulation raw data were carried out by the linear scaling (LS), power transformation (PT), and distribution mapping (DM) of precipitation methods that were used for adjusting the precipitation. Whereas LS, variance scaling (VARI) of temperature, and DM of temperature methods were used for adjusting the maximum and minimum temperatures.

Based on the comparison of statistical performance indicator results, all the selected bias correction techniques performed well for climate model simulations for all climate variables. However, based on their range of variability and their ability to bring the raw RCM simulations closer to observations, the PT method performed better and was selected for precipitation, while the DM method performed better and was selected for maximum and minimum temperature correction.

A similar study was carried out on hydrologic responses to climate and land use/land cover changes in the Bilate catchment, Southern Ethiopia (Kuma et al. 2021). A similar conclusion was also drawn by Teutschbein & Seibert (2012), and it is supported by Geleta & Gobosho (2018) in the Finchaa watershed in Ethiopia. Moreover, there are differences between the correction methods, but the DM reproduces very well for temperature bias correction in the Finchaa watershed (Dibaba et al. 2020).

While LS accounts for a bias in the mean, it does not allow differences in the variance to be corrected. Therefore, PT can be used to specifically adjust the variance statistics of a precipitation time series (Teutschbein & Seibert 2012). The PT method utilizes a nonlinear approach in an exponential form to correct the mean and variance of the precipitation series (Luo et al. 2018; Zhang et al. 2018).

The local intensity scaling (LOCI) method for the precipitation method aims to simultaneously correct the precipitation intensity and frequency. Initially, the rainfall intensity threshold for each month is confirmed. Accordingly, the number of wet days in RCM precipitation that exceed this threshold matches the number of days for which observed precipitation was determined. This approach can effectively eliminate the drizzle effect because too many drizzly days are often included in original RCM outputs (Luo et al. 2018).

The DM of the temperature method is applied to correct the distribution function of the RCM outputs and to align them with the observed distribution function. This can be done by creating a transfer function to shift the occurrence distributions of temperature. It is based on the assumption that both the RCM-simulated and -observed climatic variables obey a specific frequency distribution (Teutschbein & Seibert 2012). Concerning temperature, the Gaussian distribution (normal distribution), with location parameter and scale parameter, is often assumed to agree with the optimal temperature distribution (Luo et al. 2018). The DM method uses the transferring of function to adjust the cumulative distribution of estimated data to the cumulative distribution of observed data and it reproduces temperature very well (Dibaba et al. 2020).

The SWAT is a physically based semi-distributed hydrological model developed in the early 1990s by the Agriculture Research Service of the United States Department of Agriculture and using a process-based, basin-scale and continuous-time model (Neitsch et al. 2002b). The SWAT model has important applications for assessing the relationships between climate change and the watershed hydrological process.

The SWAT input files were organized and the model was set to run; at the end, it simulates the streamflow and evapotranspiration to evaluate the climate change effect on hydrology in the Gelana watershed.

Daily weather data, rainfall, temperatures (maximum and minimum temperatures), solar radiation, and relative humidity were required for the SWAT modeling. The 33 years including a 2-year warm-up period of the seven meteorological stations from January 1st 1987 to December 31st 2019 were used for SWAT simulation depending on data availability.

The SUFI-2 (sequential uncertainty fitting version 2) embedded in the SWAT-CUP was used for the sensitivity analysis, calibration, and validation of the hydrological model in this study. This algorithm was used for parameter optimization and for time-consuming large-scale models, and it was found to be quite efficient (Abbaspour et al. 2015). The sensitivity analysis was used to identify the most sensitive hydrological parameters for the simulated streamflow of the Gelana watershed. The entire number of the model parameters involved in the hydrological process was reduced into the most effective and sensitive parameters with a different degree of sensitivity based on its low P-value and a high absolute value of t-Stat (Dibaba et al. 2020). Calibration of the hydrological model is the process of estimating model parameters by comparing the model prediction with the observed data for the same condition (Arnold et al. 2012; Dibaba et al. 2020). It was done alongside identifying the most sensitive parameters by comparing model-simulated streamflow with observed streamflow data for the period of January 1989–December 2007. Validation is used to test the calibrated model without further parameter adjustments with an independent dataset, and the results are compared with the remaining observational data to evaluate the model prediction (Dibaba et al. 2020). It involves running a model using parameters that were determined during the calibration process and comparing the predictions to observed data not used in the calibration (Arnold et al. 2012). It is the last step of the modeling, proving the performance of the model for simulated flows in periods different from the calibration periods of January 2008–December 2015 and finally comparing the simulated flow with observed streamflow data during the validation period.

The goodness of fit of the model simulation with the observed streamflow was expressed by statistics techniques such as the coefficient of determination, NSE, and PBIAS (Abbaspour et al. 2015; Dibaba et al. 2020).

where ⁠, 97.5% and ⁠, 2.5% represent the upper and lower boundaries of the 95PPU, and σ_obs is the standard deviation of the measured data.

The regionalization method was used for transforming hydrological information from gauged catchments to ungauged catchments. The principle of parameter regionalization is to provide effective hydrological information from gauged to ungauged catchments (Guo et al. 2021). The observed streamflow data are not gauged at the outlet of the Gelana watershed, due to the fact that the hydrological model optimized parameters were transferred from gauged donor watersheds to the ungauged watershed.

PCA is a multivariate technique in which several related variables are transformed into a smaller set of uncorrelated variables. Large datasets are increasingly common and every so often problematic to interpret. The PCA is a technique for reducing the dimensionality of such datasets, increasing interpretability, and minimizing the loss of information. It is developed by forming new uncorrelated variables that continually maximize variance. Finding such new variables, the principal components, reduces to solving an eigenvalue or eigenvector problem, and the new variables are defined by the dataset at hand, not a priori, hence making PCA an adaptive data analysis technique (Wuttichaikitcharoen & Babel 2014).

The principal component loading is the best method to get a better correlation and group the optimized parameters into physically significant components. According to Sharma et al. (2015), some of the parameters are strongly correlated with the components. It has been screening out the parameters or variables of least significance and regrouping the remaining variables into the physically significant factors. The multiple regression techniques can apply in modeling the hydrological responses such as surface runoff from the watersheds. Furthermore, to evaluate and analyze the variability of physical variables contributing to the hydrological cycle, it showed proper cyclicity at seasonal to annual timescales (Syed et al. 2004).

To develop a regionalization method, the optimized parameters during calibration in the SWAT model at the gauged watersheds in the upstream of the Gelana watershed and its nearby watersheds were assessed. The watersheds used for regionalization were the Gelana river gauged at Tore, the Gelana river gauged at Yirga Chaffe, the Gidabo river gauged at Aposto, the Hare river gauged at Arba Minch, and the Kulfo river at Sekala. The Gelana river gauged in Tore and Yirga Cheffe, the Gidabo river, and the Hare river was a donor, while the Kulfo river is the neighboring watershed and is used for the validation of regionalized parameters. The targeted location was the Gelana watershed outlet near the Lake Abaya.

The neighboring (or nearby) watersheds must share physical attributes. Eighteen physical catchment characteristics were selected and determined for correlation and to make an equation with optimized SWAT parameters. These characteristics are categorized as two climate descriptors such as mean annual rainfall (MAR) and mean annual PET, three soil descriptors such as saturated hydraulic conductivity (Ksat), available water capacity of the soil layer (Swc), and bulk density moist (BDM), seven LULC descriptors such as Forest-Evergreen (% FRSE), Agricultural (% AGRL), Shrub land (% RNGB), Forest-Mixed (% FRST), Bare land (% BARR), Grass land (% PAST), and Settlement (% URBN), and six topographical descriptors such as drainage area, mean elevation (ME), length of the longest flow path (LLP), topographic wetness index (TWI), aspect, and flow accumulation (FA), which were prepared for developing the regionalization method.

The hydro-meteorological time-series trend analysis for future predicted climate variables and streamflow is necessary. It can be performed by using parametric and non-parametric tests. The non-parametric statistical trend analyses are not affected by the normal distribution of data, and they are less sensitive to outliers and better suited to normally undistributed data (Orke & Li 2021). These advantages make a non-parametric test for the trend analysis of hydro-climatic variables over a parametric analysis.

MK rank-based non-parametric statistical tests are applied in many hydro-climatological studies, which are used for detecting monotonic trends in time-series data (Belihu et al. 2018). The bigger the absolute magnitude of trend, the more powerful are the tests; as the sample size increases, the tests become more powerful; and as the amount of variation increases with time series, the power of the test decreases. The presence of significant trends in the variables using the MK test is more applicable and employed by many researchers. It is distribution-free and tolerates missing values (Belihu et al. 2018). Therefore, the MK test was used for analysis of hydro-climatic time series in this study. It is a statistical test mostly used for the analysis of trends in climatic and hydrologic data time series. According to this test, the null hypothesis (H₀) designates that there is no trend in the series, while the alternative hypothesis (H₁) describes that there is a trend (Orke & Li 2021). The null hypothesis is that a data series is serially independent and identically distributed with no trend.

The computational procedure for the MK test considers the time series of n data points and x_i and x_j as two subsets of data where i = 1, 2, 3, …, n − 1 and j = i + 1, i + 2, i + 3, …, n. The data values are evaluated as an ordered time series. Each data value is compared with all subsequent data values. If a data value from a later period is higher than a data value from an earlier period, the statistic S is incremented by 1. On the other hand, if the data value from a later period is lower than a data value sampled earlier, S is decremented by 1. The net result of all such increments and decrements yields the final value of S.

Test statistic Z_s evaluates the statistical significance of the trend. If Z_s is positive, it indicates an upward trend, while negative denotes decreasing trends. For a given level of significance α, the null hypothesis of no trend in the time series will be rejected if > ⁠. More specifically, > indicates the existence of an increasing trend (Belihu et al. 2018; Orke & Li 2021).

The summary of methodology to achieve the objective in the current study is shown in Figure 2.

The sensitive parameters that affect the SWAT model output with its rank were assessed by using SWAT CUP (SUFI-2). The parameters with small sensitivity values do not significantly affect the SWAT model output. However, medium, high, and very highly sensitive values that affect the model output significantly were used to calibrate the hydrological model.

The results of SWAT sensitivity analysis, ranges, and its fitted value in the Gelana watershed at the Tore gauging station are presented in Table 1. The result indicated that ALPHA_BF.gw (baseflow alpha factor in days), RCHRG_DP.gw (deep aquifer percolation fraction), CH_K2.rte (effective hydraulic conductivity in main channel alluvium), CN2.mgt (Soil Conservation Service (SCS) runoff curve number), GWQMN.gw (threshold depth of water in the shallow aquifer required for return flow to occur in mm), SOL_K(..). sol (saturated hydraulic conductivity), SLSUBBSN.hru (average slope length), and HRU_SLP.hru (average slope steepness) were found to be the eight topmost sensitive parameters. Similarly, the remaining parameters were listed to their sensitivity order and are also sensitive parameters in the Tore gauging station in the Gelana watershed.

The selected 19 sensitive parameters were used for calibration and validation of the hydrological model with SWAT CUP (SUFI-2) on monthly time steps from 1989 to 2015. Two-thirds of the total continuous time-series flow gauging period from 1989 to 2007 were selected as the calibration period, and the remaining one-third period, i.e., from 2008 to 2015, was used for validation of the hydrological model gauged in the Tore station in the Gelana watershed. Moreover, the SWAT model calibration and validation are shown in Figures 3 and 4.

The calibration results on mean monthly flow show that the SWAT model can capture the observed streamflow with R², NSE, and PBIAS of 0.74, 0.74, and 4.8, respectively. The average monthly flow validation indicates R², NSE, and PBIAS of 0.71, 0.67, and 7.5, respectively. Furthermore, the p-factor, which is the percentage of observations bracketed by the 95PPU, brackets 78%, and 66% of the observation, and r-factor equals 0.89 and 0.76 during calibration and validation periods, respectively. Generally, the hydrological model showed a reasonably good agreement in the Tore gauging station. Therefore, according to descriptions by Moriasi et al. (2007), the SWAT model has an acceptable performance in the Gelana watershed.

The observed streamflow data are not gauged at the outlet of the Gelana watershed, due to the fact that the 19 SWAT sensitive parameters were transferred from gauged rivers to an ungauged outlet to estimate the climate change effect in the streamflow at the watershed outlet.

The calibration and validation were carried out to check the efficiency of the SWAT model performance of donors and validation rivers. The results indicate that the parameter donors and the regionalized model validation rivers show a good agreement between the observed and simulated flow. Then, the 19 optimized SWAT model parameter values at gauged neighbor catchments were assessed.

The physical catchment characteristics that relate to the topography of the watersheds were extracted from the DEM. LULC and soil were obtained from the LULC map and soil map, respectively. Similarly, the physical catchment characteristics under climate descriptors were obtained from the meteorological data as made available by the National Meteorological Agency.

The inter-correlations of 18 physical catchment characteristics were calculated. Depending on the inter-correlation values, all physical catchment characteristics were in the range of moderate to strong correlation. The correlation matrix of physical catchment characteristics reveals strong correlations (correlation coefficient greater than ±0.9), good correlation (correlation coefficient greater than ±0.75), and moderate correlation (correlation coefficient greater than ±0.60). Therefore, it is very difficult at this stage to group the parameters into components and attach any physical significance. Furthermore, the correlation matrix is subjected to the PCA.

By using IBM SPSS statistics software, the PCA applying the Varimax rotation indicated that three principal components with an eigenvalue greater than 1.00, as shown in Table 2, accounted for the total cumulative variance of 95.4% as per their eigenvalues. The first component (PC1) has explained about 55.0% of the variance in the catchment characteristics, whereas the second component (PC2) explained 33.4%, and the third component (PC3) explained 7.0%. Therefore, 95.4% of the variance is explained by only the three components.

The bold values in Table 3, by the extraction method (PCA, rotation method: Varimax with Kaiser normalization), indicate highly correlated variables and factor loadings in the PCs.

The first component (PC1) is strongly correlated (loading of more than ±0.90) with Swc, BDM, and PET, and moderately correlated (loading of more than ±0.60) with FRST, Aspect, Ksat, AGRL, FRSE, and MAR based on higher loading factors which may be termed as soil descriptor component. The second component (PC2) is strongly correlated with Area, LLP, and FA, and moderately correlated with BARR, RNGB, and ME, which may be termed as a topographical descriptor component. The third component (PC3) is correlated as good (loading of more than ±0.75) with URBN and PAST, and moderately correlated with AGRL, FRSE, and MAR, which may be termed as the LULC descriptor component.

To select prominent variables for subsequent regression analyses, the first three variables with the highest factor loadings and greater than ±0.60 were selected as representative variables of each of the PCs. A threshold of 0.60 was used for identifying a reliable factor in this study according to Wuttichaikitcharoen & Babel (2014). The useful PCA has been screening out the parameters of least significance and regrouping the remaining variables into the physically significant factors (Sharma et al. 2015).

These are evident from the results that nine of the parameters are highly correlated with some of the components. However, other parameters show a moderate correlation. To screen out parameters having less significance in explaining the component variance, FRST, Aspect, Ksat, FRSE, MAR, BARR, RNGB, and ME were screened out from the analysis.

Then correlation matrix and principal component matrix are obtained for nine parameters as follows: for PC1: Swc, BDM, and PET were selected, for PC2: Area, LLP, and FA were used, and URBN, PAST, and AGRL were considered for PC3. All nine factors were assumed to be the forcing factors of streamflow with positive and negative effects, which can be used subsequently as predictor variables in regression analysis.

A stepwise multiple linear regression between the predictor physical catchment characteristics and SWAT model-sensitive parameters was assessed. After the stepwise regression by IBM SPSS Statistics software, the correlation coefficients were identified and the physical relevancy of each index to each parameter was checked. The catchment indices that showed better correlation and those that are hydrologically relevant to each SWAT model parameter regarding the catchment response were selected and regressed over each model parameter.

The regression equation was developed after testing its statistical significance through R², t-test, and P-value of the regression statistics. Finally, the regression equations for each SWAT parameter with a function of the physical characteristics and principal components were developed. Then, by substituting the physical characteristics and principal components of the Gelana river at the outlet near the Lake Abaya, the most sensitive parameters were transposed to the ungauged outlet of the Gelana watershed as shown in Table 4.

Before using the transposed model parameters in the ungauged catchment at Gelana outlet, the performance of the regionalization method was calibrated and validated in the neighbor catchment. However, the catchment that was used for validation was not used to derive the regional equation for the ungauged watershed. Therefore, the Kulfo catchment was not used to derive the regionalized equation for transposed parameters, and it is the nearest watershed to the ungauged part of the Gelana watershed outlet, located in a similar climate condition (semi-arid zone) to the ungauged part of the study watershed.

The overall results show the correlation coefficient (R² = 0.62) and the Nash–Sutcliffe simulation efficiency (NSE = 0.51) during the calibration, and R² = 0.66 and NSE = 0.52 during validation of the Kulfo watershed. This result demonstrates a satisfactory agreement between observed and simulated values. In general, the regionalization model has acceptable performance in a nearby catchment; therefore, the transposed parameters can determine the Gelana watershed outlet (ungauged part) simulation.

Climate change projection was performed by using the ensemble mean of two RCMs for the near future (2031–2050), middle future (2051–2070), and last future (2071–2090) under two emission scenarios relative to the baseline (1996–2015) over the Gelana watershed.

The ensemble mean annual of the two RCMs suggested decreasing precipitation in the three future periods under both scenarios as shown in Table 5. The changes show rainfall decline by 11.28% under RCP4.5 and 12.52% under RCP8.5 in.the near future period. For the middle future, it declines by 15.12% under RCP4.5 and 10.08% under the RCP8.5 scenario. Finally, it is expected to reduce by 7.21% under RCP4.5 and 4.85% under RCP8.5 in the last future period. A large reduction in rainfall was also confirmed by Dibaba et al. (2020) using four different RCMs in the Finchaa sub-basin, Ethiopia. There are variations among both the RCMs, and the RCA4 climate model has shown the maximum reduction under RCP4.5 and RCP8.5 scenarios in all study years.

There are four seasons, namely winter (Bega) including December, January, and February; spring (Belg) including March, April, and May; summer (Kiremt) including June, July, and August; and autumn (Tseday) including September, October, and November. The more profound changes in projected rainfall were seen in seasonal bases compared to the annual bases in the watershed. It shows a declining trend in all seasons, except in winter and spring (from the period 2071 to 2090 for the RCP4.5 scenario). However, it is expected that there will be a higher decline in rainfall changes in the wet season (higher rainfall) than in the dry season, as shown in Figure 5. In particular, in the spring season, it will decline by 16.68% in RCP4.5 from 2051 to 2070, and in the autumn season, it reduces by 29.50% in RCP4.5 from 2051 to 2070. In addition, in the summer season, it is expected to decrease by 23.52% from 2071 to 2090 in the RCP8.5 scenario compared with the base period. Meanwhile, in winter (dry season), it is expected to face a higher incremental rainfall change of 73.43% from 2051 to 2070 under the RCP4.5 scenario.

Therefore, for wet seasons, all projections indicate a higher decreasing signal in the Gelana watershed. This condition may affect the economy of the farmers because of the subsistence of rain-fed agricultural activity in the watershed. The supportive result was reported in the Bilate watershed by Kuma et al. (2021) who stated that it is expected to decrease during the short rainy season (February–May) and the long rainy season (June–September). Corresponding results state that the warm phases of the El Niño Southern Oscillation have been associated with reduced rainfall in the main wet season (July, August, and September) in Ethiopia, causing severe drought, famine, and desertification, but with enhanced rainfall in early February which mainly affects Southern Ethiopia (Legesse et al. 2015). In addition, the decreasing pattern of precipitation intensity in southern and southwestern regions of Ethiopia was reported by Dibaba et al. (2020). The consistent study carried out in the Central Rift Valley Basin using the ensemble mean of RCMs shows that the average projected rainfall will decrease by 7.97 and 2.55% under RCP4.5 and RCP8.5, respectively, for the future scenario period (Gadissa et al. 2019). The climate trend analysis from 2010 to 2039 exposed to the rainfall decline ranges from 50 to 150 mm and is related to lower harvest and poor pastoral rangelands across the southern and eastern parts of Ethiopia. Similarly, the annual rainfall will vary in the range of +10 to −40% by 2080, in which time it will rise outside the growing season and decline during the growing seasons (Kuma et al. 2021).

Regarding the ensemble mean of the RCMs, the maximum temperature will increase on average by 1.18, 2.01, and 2.12 °C under medium emission scenarios (RCP4.5) for the near, middle, and last future, respectively. Likewise, it will increase by 1.52, 2.51, and 3.48 °C for the near, middle, and last future, respectively, under high emission scenarios (RCP8.5). The overall results show in Table 6 that maximum temperatures increase under both RCPs throughout the study years considered showing the warming trends in the watershed. RCA4 has shown the highest increase in maximum temperature.

It indicates that the RCP8.5 scenario is warmer than the RCP4.5 scenario in the Gelana watershed. A high increase was seen in the wet season (autumn) and dry season (winter) and a large increase was also seen in summer and spring seasons. The highest seasonal rise of maximum temperature is estimated in winter (dry season) by 15.88% under the RCP8.5 scenario during 2071–2090, while under the RCP4.5 scenario the lowest rise (+2.40%) is expected in spring (wet season) during 2031–2050.

The highest seasonal rise in the maximum temperature is expected by the high emission scenario (RCP8.5) than the medium emission scenario (RCP4.5). Likewise, the highest rise will be expected in the last future and the middle future than the near future.

A related study was done in the Bilate watershed; in the mid-century from 2051 to 2080, the highest temperature changes are likely to become warmer than those in the present century from 2021 to 2050 under RCP4.5 and RCP8.5 scenarios (Kuma et al. 2021). In addition, the highest temperature change in the RCP8.5 scenario is also confirmed by Dibaba et al. (2020).

The mean annual minimum temperature in the Gelana watershed will increase, on average, by 1.30, 2.36, and 2.61 °C under medium emission scenarios for the near, middle, and last future, respectively. Likewise, it will increase by 1.56, 2.75, and 4.19 °C for the near, middle, and last future, respectively, under high emission scenarios as shown in Table 7. It also varies seasonally. However, the seasonal variations are higher for the minimum temperature than the maximum temperature. The seasonal changes were expected to rise in wet and dry seasons. The maximum seasonal rise is predicted in the summer season by 35.53% under the RCP8.5 scenario during 2071–2090, while under the RCP4.5 scenario it is estimated as the lowest rise (+7.97%) in winter (dry season) during 2031–2050.

Generally, minimum temperatures increase under both RCPs throughout the study years considered showing the warming trends in the watershed. Moreover, the magnitude of changes by the mean ensemble is maximum for higher emission scenarios (RCP8.5) compared to the RCP4.5 scenario. Likewise, it is higher in the last future than the middle and recent future under both emission scenarios. For all future periods under scenarios, the mean monthly variations will be greater compared to mean annual variations. In addition, the monthly and annual variations are higher for the minimum temperature than the maximum temperature.

The related study was reported on the Rift valley basin. A study by Kuma et al. (2021), in the Bilate watershed, stated that the minimum temperatures in the mid-century from 2051 to 2080 are warmer than from 2021 to 2050 under RCP4.5 and RCP8.5 scenarios. Similarly, the highest minimum temperature change in RCP8.5 is also confirmed by Dibaba et al. (2020).

Therefore, the projections of maximum and minimum temperatures in the three future time horizons under RCP4.5 and RCP8.5 scenarios are within the range predicted by IPCC (2014) and by some other researchers. Several climate change reports in different regions in Ethiopia have shown that the temperature is expected to increase. But, the magnitude of change varies with the procedures of downscaling and climate model types. According to Dibaba et al.’s (2020) report, the study on climate risk and adaptation country profile states that the mean annual temperature is expected to increase by 1.1–3.1 °C by the 2060s and 1.5–5.1 °C by 2090s. A study in the Kulfo watershed by Demmissie et al. (2018) clearly stated that the mean annual temperature is expected to increase by 0.5 °C from 2.5 to 3 °C and from 4.5 to 5 °C in the 2050s and 2080s, respectively.

According to the ensemble mean of the RCMs, the PET (Table 8) will increase on average by 9.19, 12.30, and 12.90% under the medium emission scenario (RCP4.5) for the near, middle, and last future, respectively. Likewise, it will increase by 10.28, 13.91, and 17.85% for the near, middle, and last future, respectively, under high emission scenarios (RCP8.5). The seasonal variations were also evaluated by the ensemble mean of the two RCMs. So, it will be expected to increase in wet and dry seasons with respect to the baseline period. The maximum increases were seen in the summer, autumn, and winter seasons compared to the spring season. However, the maximum seasonal rise of PET was predicted in summer by 22.05% under the RCP8.5 scenario during 2071–2090, while under the RCP4.5 scenario it shows the minimum rise (+8.25%) in winter (dry season) during 2031–2050.

Likewise, a higher seasonal rise in PET is expected under the high emission scenario (RCP8.5) compared to the medium emission scenario (RCP4.5). In this manner, higher increases will be expected in the last future (2071–2090) and in the middle future (2051–2070) than in the near future (2031–2050). Therefore, it is expected there will be a higher increase in all seasons for the future periods with respect to a baseline period in the Gelana watershed.

The changes in temperature and PET are correlated positively in the Gelana watershed. Accordingly, the rise in temperature resulted in an increase in the PET.

The corresponding study was carried out by another researcher in the Rift valley basin watersheds. Kuma et al. (2021) clearly stated that in the Bilate watershed PET was expected to increase under RCP4.5 and RCP8.5 scenarios in future periods from 2021 to 2080. These changes follow the direction of temperature changes. A 2.0 °C increase in the mean annual temperature resulted in about a 19% increase in PET under RCP8.5 in the mid-century from 2051 to 2080 in the watershed. This result also agrees with a study carried out in the Finchaa watershed by Dibaba et al. (2020). The consistent study done in the Central Rift Valley Basin demonstrates that the results from the ensemble mean of RCMs show that the average temperature will increase by 1.9 and 2.7 °C under RCP4.5 and RCP8.5, respectively. This corresponds to a 4.89 and 6.59% increase in PET under RCP4.5 and RCP8.5, respectively (Gadissa et al. 2019).

The changes in rainfall and temperature were used to predict the effect of climate change on the water balance components of the watershed. The projected MAR and temperature for future periods from 2031 to 2090 under both RCP scenarios cause significant changes in water balance components in the Gelana watershed. The sensitivity of water balance components to climate change is varied in future periods, and the effect is explained with the integration of rainfall, temperature, and evapotranspiration. The hydrological responses of the watershed to climate changes were investigated with the ensemble mean (Table 9) of the two models for three future periods under RCP4.5 and RCP8.5 scenarios.

The hydrological responses of the watershed were considered in terms of the process that contributes to the water balance components of annual surface runoff, groundwater, total water yield, and actual evapotranspiration. The results show that the mean annual surface water runoff, groundwater, and total water yield will decrease, whereas actual evapotranspiration will be expected to increase in all future periods in the watershed.

The highest percentage of decreasing water balance components was forecasted using an ensemble mean (Table 10) of two model results, surface runoff by 22.23%, groundwater by 42.54%, and total water yield by 35.89% in middle future under the RCP4.5 scenario and, similarly, surface runoff decline by 11.80% in near future, groundwater by 37.0% in the middle future, and total water yield by 30.16% in near future under the RCP8.5 scenario corresponding to the base period.

The annual actual evapotranspiration in the Gelana watershed will be expected to increase in all future study periods under RCP4.5 and RCP8.5 scenarios. The highest percentage of variation was forecast in the higher emission scenario and in the last future period. The maximum increasing value was estimated in the middle future by a rate of 8.07% under the RCP4.5 scenario and also rose by a rate of 14.94% in the last future under the RCP8.5 scenario with respect to the base period.

Generally, the SWAT model simulations for the near future, middle future, and last future showed that the decline of rainfall and the increase of maximum and minimum temperatures will lead to reduced surface runoff, groundwater, and total water yield, whereas actual evapotranspiration and PET will be expected to increase in the Gelana watershed. Consequently, the increase in temperature resulted in increased PET and actual evaporation, which could be a critical factor for the reduction of total water yield. As the forms of water are exposed to losses, owing to the changes in temperature, evaporation is also a factor for the future period reduction of the surface runoff and groundwater in the Gelana watershed. This is expected considering the warming trends of temperature due to climate changes. Therefore, the Gelana river is a tributary for the Lake Abaya sub-basin in a Rift valley basin, and thus, the decline of water yield may also result in reducing the volume of water in Lake Abaya in the future period from 2031 to 2090. This result agrees with a study on the impacts of the climate changes on the Finchaa watershed, which states decreasing precipitation results in reduced surface runoff, groundwater, and total water yield due to increasing temperatures amplified by the increase in evapotranspiration (Dibaba et al. 2020). Moreover, the result is consistent with other research reports in the Bilate watershed by Kuma et al. (2021) and concluded that a 15.39% reduction in the annual rainfall caused about 9.22, 11.11, and 10.25% reductions in water yield, groundwater, and surface runoff, respectively, under the RCP4.5 scenario from 2051 to 2080. Similarly, a 2.0 °C increase in the mean annual temperature resulted in about 19 and 13.9% increases in PET and ET, respectively, under RCP8.5 in.the mid-century (2051–2080).

The projected changes in rainfall and temperature, under two RCP scenarios, cause a significant variation of the future streamflow of the Gelana watershed. However, the other weather variables, such as solar radiation, relative humidity, and wind speed observed in the baseline period were considered in the future period projection without making any change, because changing these variables may not have a significant impact in modeling the climate change scenarios on local hydrology (Dibaba et al. 2020).

The annual streamflow in the Gelana river will be expected to decline by 21.09% under RCP4.5 and 32.30% under RCP8.5 in the near future period. Under the middle future, it will decline by 44.14% under RCP4.5 and 24.27% under the RCP8.5 scenario. Finally, streamflow will be reduced by 22.10% under RCP4.5 and 15.61% under RCP8.5 in the last future period.

The reasons may be the decreasing rainfall and the increasing maximum and minimum temperatures, together with PET in the Gelana watershed in all study years under both RCP scenarios.

The results from this study agreed with the findings of similar studies in different areas. For example, the study by Gadissa et al. (2019) shows that a 10% increase in PET gives a decrease in streamflow of about 15.1% for the Lake Ziway basin. A study carried out by Demmissie et al. (2018) states that the average annual streamflow of the Kulfo river is projected to decrease by 8.2% in the 2080s, and increasing temperature by 0.5 °C decreases streamflow rates by 2.99% in the 2050s.

The maximum and minimum seasonal changes in flow projections for future periods were studied. The predicted streamflow demonstrates a declining trend in all seasons, except in winter for the middle and last future under the higher emission scenario, and in spring for the last future period under the RCP8.5 scenario, and for near future in the RCP4.5 scenario, and similarly, in summer for near and last future under the medium emission scenario. A higher reduction of flow changes is expected in the wet season than the dry season as shown in Figure 6. Moreover, the maximum decrease of flow variation is expected in the autumn season. Specifically, it declines by a rate of 75.85% under RCP4.5 from 2051 to 2070 in months including September, October, and November. Similarly, in all future periods for both RCP scenarios, the streamflow shows a large decline in the autumn season.

In addition, the reduction of flow will be expected in the summer season by 43.25% from 2071 to 2090 under the RCP8.5 scenario. It is also expected there will be a higher decrease compared with the base period in spring and winter seasons, whereas, in winter (December, January, and February), it will be expected to look at a higher incremental flow by 38.75% from 2051 to 2070 and by 36.92% from 2071 to 2090 under the RCP8.5 scenario. Similarly, it is expected there will be a maximum rise in the summer season of 34.17% in the last future under the medium emission scenario.

In general, as with rainfall trends, all streamflow projections show a higher reduction trend in wet seasons than dry seasons in the Gelana watershed. Similar results were studied by Shanka (2017) in the Gidabo watershed, who states that during autumn and winter seasons, the total average seasonal runoff showed a reduction in all future horizons.

The trends of the climate variables and streamflow were tested in terms of annual time series using the non-parametric MK trend test at a significance level of 5% for both climate models under RCP scenarios. According to the MK trend test analysis, the positive and negative values of the S statistic indicate upward and downward trends, respectively. In addition, the significance is rejected (insignificant) when the absolute value of Z is less than 1.96, whereas the significance should be accepted (significant) when the absolute value of Z is greater than 1.96 at a significance level of 5%.

The projected annual rainfall over the Gelana watershed shows a slightly increasing trend for the future 60 years from 2031 to 2090 as shown in Table 11. However, compared with the base period (observed), it will decrease with different rates for various scenarios and climate models. The projected annual rainfall shows an increasing trend for the future years from 2031 to 2090 at a rate of 1.53 mm/year (statistically insignificant) and 3.79 mm/year (statistically significant) in RCP4.5 and RCP8.5, respectively, in the RACMO22T climate model as shown in Figure 7. Similarly, the projected annual rainfall shows the increasing trend at a rate of 0.74 mm/year in RCP8.5, whereas the statistically insignificant decreasing trend at a rate of 0.44 mm/year RCP4.5 in the RCA4 climate model is shown in Table 12 and Figure 8. This variation of increment and decrement trends was due to the complexity in the nature of rainfall and its reliance on topographic and physical factors (Dibaba et al. 2020). Related study of the non-significant trend in annual rainfall was reported in the middle Gidabo catchment (Belihu et al. 2018; Belihu et al. 2020).

The projected annual average temperature clearly shows an increasing trend (statistically significant), as shown in Tables 11 and 12 and Figures 9 and 10. An increasing trend will be 0.02 and 0.05 °C/year for RCP4.5 and RCP8.5, respectively, for the RACMO22T model and 0.04 and 0.06 °C/year for RCP4.5 and RCP8.5, respectively, for the RCA4 model. The related result of the trend analysis indicates that the annual temperatures increased significantly in the Bilate watershed (Orke & Li 2021).

The projected annual PET shows an increasing trend (statistically significant), as shown in Tables 11 and 12 and Figures 11 and 12, in both climate models under RCP4.5 and RCP8.5 scenarios. It shows an increasing trend for future periods at a rate of 0.83 and 1.81 mm/year under both RCPs, respectively, for the RACMO22T model. Similarly, an increasing trend at a rate of 1.28 and 2.09 mm/year under both RCPs, respectively, for the RCA4 climate model. This result indicates that due to climate changes, the increment in maximum and minimum temperatures results in an increase in PET in the watershed. In addition, the maximum increase seems to be in a higher emission scenario than the medium emission scenario. This shows the direct relationship between temperature and evapotranspiration. The increase in PET suggests an increased crop water requirement in future period crop production. Thus, the design of irrigation infrastructure should take this into account.

The streamflow is mainly dependent on the amount of rainfall falling on its watershed area and the evapotranspiration amount released into the atmosphere from the waterbody and land. Hence, the variations in rainfall and temperature can significantly influence streamflow patterns. The projected annual streamflow of the Gelana river for the future 60 years from 2031 also shows a slightly increasing trend as presented in Tables 11 and 12 and Figures 13 and 14 at a rate of 0.06 m³/s/year (statistically insignificant) and 0.09 m³/s/year (statistically significant) in the RACMO22T model under both RCP scenarios, respectively, and 0.03 m³/s/year in the RCA4 model under RCP8.5 scenarios (which is statistically insignificant). However, its slight decreasing trend (which is statistically insignificant) is shown in Tables 11 and 12 and Figure 14 at a rate of 0.05 m³/s in the RCA4 model under RCP4.5 scenarios. A corresponding result was found in the Rift valley basin watershed. The annual streamflow showed a decreasing trend which is not significant in the Gidabo watershed (Belihu et al. 2018).

In general, the projected future streamflow comparatively with the baseline period resulted in a decreasing trend in the river. This reduction is due to the increase of the maximum and minimum temperatures coupled with the increased PET in the Gelana watershed in both climate models under RCP4.5 and RCP8.5 scenarios.

However, for only future period trend analysis, particularly after the 2050, a little rise in the rainfall results in a slight increase in the streamflow under RCP4.5 and RCP8.5 scenarios in the RACMO22T model and under the RCP8.5 scenario in the RCA4 model. This indicates that the streamflow is more dependent on the amount of rainfall falling on its watershed area than the evapotranspiration amount released to the atmosphere. It may also be due to the watershed being located in two types of climatic conditions: the upstream part is in a humid region (low mean temperature, high rainfall) and the downstream part is in a semi-arid region (high mean temperature, low rainfall). The results show that the future streamflow in the Gelana watershed will be more sensitive to change in rainfall than a change in temperature and PET. This indicates the direct relationship between rainfall and streamflow in the watershed. This was in agreement with the results in the Kulfo watershed (Demmissie et al. 2018; Shah et al. 2021). In addition to these climatic factors, some non-climatic, i.e., anthropogenic actions (like LULC change), may also affect the streamflow of the Gelana river.

The population growth together with climate changes increases uncertainty in the future demand and availability of water. This study aims to investigate the climate change effect on streamflow in the Gelana watershed using the SWAT for three consecutive periods of 2031–2050, 2051–2070, and 2071–2090. Climate variables were downscaled from two different RCMs (RACMO22T and RCA4) under RCP4.5 and RCP8.5 scenarios from CORDEX-Africa. The performance of the climate models and bias correction methods was evaluated using the Pearson correlation coefficient (R), RMSE, NSE, and PBIAS. PT and DM performed the best for precipitation and temperatures, respectively, and were used for further analysis. The SWAT model was calibrated (and validated) for the 1989–2007 (2008–2015) period at the Tore station in the Gelana watershed; the model showed good performance. The general results indicate that an ensemble mean of two RCMs shows that the highest changes in annual maximum and minimum temperatures and PET are expected to increase by 2.12 °C, 2.61 °C, and 12.90%, respectively, in the last future period. Whereas the rainfall and streamflow are expected to decrease by 15.12 and 44.14%, respectively, in the middle future period under the RCP4.5 scenario. In the same way, maximum and minimum temperatures and PET are expected to increase by 3.48 °C, 4.19 °C, and 17.85%, respectively, in the last future period. However, the decline in the rainfall and streamflow by 12.52 and 32.30%, respectively, are predicted in the near future period under the RCP8.5 scenario. Correspondingly, increases in maximum and minimum temperatures results in an increase in evapotranspiration, which may have also been the cause of the decreased surface runoff by 22.23%, groundwater by 42.54%, and total water yield by 35.89% in the middle future period under the RCP4.5 scenario. Therefore, the Gelana river is a tributary for the Lake Abaya sub-basin in a Rift valley basin; thus, the decline of total water yield may also result in a reduction in the volume of water in the Lake Abaya for the future periods. Moreover, the projected rainfall and streamflow are expected to face a higher decline in wet seasons, although the highest seasonal rise of maximum temperature was estimated in dry seasons. Additionally, the RCP8.5 scenario projection is warmer than the RCP4.5 scenario. Likewise, higher increases will be expected in the last future than the near future period in the Gelana watershed. These conditions explain the integration of precipitation, temperature, and evapotranspiration to water balance and their sensitivity to climate change. Furthermore, this suggests that the farmers, governmental institutions, and other concerned organizations in and around the Gelana watershed should look for additional small-scale water sources and adopt drought-resistant crops, integrate their duties with climate change, and select the probable adaptation measures. Trends of the rainfall, average temperature, PET, and streamflow at the watershed outlet on an annual basis for the future periods under both RCP scenarios in the RACMO22T and RCA4 climate models have been estimated by using the MK test analysis.

This study considered only precipitation, maximum and minimum temperatures, and assumed all other parameters to be constant for future periods. However, change in land use/land cover scenarios and other weather variables may have some impacts on streamflow and water balance components. Therefore, this study should be extended by considering these parameters. Moreover, the results were based on two RCMs. However, the green house gas (GHG) emission scenarios are constructed on technological growth, socioeconomic development, and demographic growth, which add uncertainty into climate change impact evaluations. So, to diminish such uncertainties, it is often recommended to use multiple RCMs. Thus, future studies should consider multiple RCMs with different emission scenarios for the climate change analysis.

All relevant data are included in the paper or its Supplementary Information.