Journal of Earth Science & Climatic Change
Open Access

Land and water resources are crucial physical assets in the subsistence agricultural economy [1]. As a result, any changes in land use and land cover have a direct impact on food production and alternative economic activities. Land-use changes studied at different scales ranging from small scale to large scale indifferent parts of Ethiopia indicates that, agriculture has gradually expanded from gentle slope land to steeper slopes [2]. Land use and land cover (LULC ) could lead to significant changes in evapo-transpiration, soil moisture content, infiltration rates, water flow regimes, surface runoff, and soil erosion through interactions with vegetation, topography, soils, geology and climate processes [3-5]. The impact of land use and land cover changes on surface hydrology, surface energy balance, and surface roughness is not linear function or straightforward but rather complex to warrant any generalization as it is dependent on the environmental processes in general and scale of the watershed, seasons, climate, and soil conditions in particular [6]. Since the knowledge about the impact of land use and land cover changes was still limited at community level, insights into the consequences of land use and land cover change on hydrology and surface properties have been explained at small spatial, observable scales [3]. The authors such as Tesfaye and Bogale, Assefa, Yilma and Takala et al., conducted the land use land cover change analysis and the results mainly show the progressive expansion of agricultural land in Gibe III watershed. Amsalu (2010) also justified that, Gibe III watershed is one amongst such land resources which are subjected to the land use and land cover dynamics. Accordingly, the condition of land under little vegetative cover is expected to be subjected to high surface runoff amounts, low infiltration rate, and reduced groundwater recharge, eventually leading to the lowering of water tables and intermittence of once-perennial streams (Woldeamlak and Geert). The selected watershed was contributing as a source of water resource for hydropower generation, irrigation, sustainable biodiversity conservation and eco-tourism services for about five National Parks [7]. The conversion of natural forest, grazing land, shrubs and woodlands to agricultural lands in this catchment area was prominent during the last 20 years (Paper I). These changes were primarily due to anthropogenic activities, which indicate expansion of agricultural land including change of uncultivable, woodlands to cultivated lands as well as overgrazing, leading to accelerated soil erosion.

However, land use and land cover change at the area is loading stress on the sustainability of aforementioned benefits [8]. In the current world of increasing population and limited resource distribution, the complex relationship between human development and the environment is dynamic and seems to be major cause for land degradation. On the other hand, the human being is also responsible to rehabilitate the environment as a whole and watershed specifically based on the priority area of development action plan. In view of this case, this research seeks to investigate the LULC change effects on hydrological responses like stream flow and sediment yield condition at Gibe III watershed, Omo river basin that could show directions for policy makers to perform preventive action.

Therefore this study was planned to know the impact of land use and land cover change on sediment yield and stream flow in the selected watershed based on LULC map of 1994, 2006 and 2018. The finding of this research is expected to contribute to understanding of the impact of land-use dynamics in the stream flow and sediment load of the area.

Description of the study area

The study area is located around 450 km south of Addis Ababa having latitude 6.6o–9.4oN and longitude of 35.78°–38.42°E. The catchment area is about 32154.16 km² with a long term mean flow estimated to be 438.2 m³s^-1 [7]. It has monthly mean stream flow ranging from 60 m³/s in March to 1,500 m³/s in August [9]. The Gibe III Hydropower scheme comprises a 243 meter high dam which will create a reservoir with a surface area of some 200 km² and live storage of some 11,750 million m3 and a surface power house equipped with ten power generating units and a switchyard (EEPCO, 2009).

Climate

The climate of the Omo-Gibe River Basin varies from a hot arid climate in the southern part of the floodplain to a tropical humid one in the highlands that include the extreme north and north-western part of the Basin. Intermediate between these extremes and for the greatest part of the basin the climate is tropical sub-humid. The amount of mean annual rainfall decreases with a decrease in elevation throughout the Omo-Gibe catchments and ranges from 1,200 mm-1,900 mm.

Moreover, the rainfall regime is uni-modal for the northern and central parts of the basin and bimodal for the south. The average annual rainfall calculated over the whole Gibe III dam watershed where the dam is located is 1,426 mm and 75- 80% of the annual rainfall distribution occurs during five months from May to September. The mean annual temperature in the basin varies from 16oC in the highlands of the north to over 29°C in the lowlands of the south [7].

Topography

The topography of the Omo-Gibe basin is characterized by its physical variation. The northern two-thirds of the basin has mountainous to hilly terrain cut by deeply incised gorges of the Omo, Gojeb, and Gilgel-Gibe rivers (Figure 1), while the southern one-third of the basin is a flat alluvial plain punctuated by hilly areas. The basin in general lies at an altitude range of 333-3570 m.a.s.l while the study site (watershed) has an altitude range of 681‑3570 m.a.s.l and the plains of the Lower Omo lies between 400-500m.a.s.l (EEPCO, 2009).

Figure 1: Location of the Gibe-III watershed.

Land Use

The major land-use types in the study area are forest lands, urban areas, rangeland, agricultural land, water body, and other built-up areas [10]. In a very broad term, most of the northern catchments of the Omo-Gibe basin are, under extensive cultivation with increased land pressure, i.e., the expansion of cultivated areas into marginal lands at the expense of woodlands. The flatter and poorly drained bottomlands of the northern catchments are usually not cultivated but are used for dry season grazing and eucalyptus tree plantations. The main gorges of the basin are relatively unpopulated and support open woodland and bush land with grasses. The eastern part of the basin has some of the most densely populated and intensively farmed areas. The south of the basin is more sparsely populated with coverage of natural vegetation through deforestation, which is increasing at an alarming rate (EEPCo, 2009).

SWAT model Setup

In this study, the Arc SWAT2012 version model was applied to estimate impact of land use land cover change on stream flow and sediment yield. The model was selected after thorough evaluation of hydrological models based on predefined criteria for achieving the objectives of the study. SWAT was selected because its data requirements are not complex (eg. SWAT need DEM, Land Use/ Land Cover, Soil and Weather data as input for setup); it is simple to perform model sensitivity and uncertainty analysis; its applicability for complex catchment; it is useful for spatial data analysis; it applied as the interface with Geographic Information System (GIS) and its potential for continuous review and improvements [11].

SWAT model is more preferable than others due to its use of SWAT-CUP program interface for sensitivity analysis, calibration, validation and uncertainty analysis procedures that is very crucial for this study [12]. As a physically based model, SWAT use Hydrologic Response Units (HRUs) to describe spatial heterogeneity in terms of land cover, soil type and slope within a catchment. The SWAT model uses two steps for the simulation of hydrology: the land phase and routing phase. The land phase controls the amount of sediment, nutrient and pesticides loading to the main channel in each sub-basin. Routing phase defines the movement of water, sediments, and nutrients through the channel network of the catchment to the outlet. The land phase of the hydrologic processes is simulated by the model based on the water balance equation in [13] defined as:

(1)

Where;

SW_t = final soil water content (mm),

SW₀ = initial soil water content on day i (mm),

t = time (days),

R_day= amount of precipitation on day i (mm),

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

E_a= amount of evapo-transpiration on day i (mm),

W_seep= amount of water entering the vadose zone from the soil profile on day i (mm) and

Q_gw= amount of return flow on day i (mm).

SWAT offers two methods for estimating surface runoff: the Soil Conservation Service (SCS) curve number (CN) procedure (SCS, 1972) and the Green and Ampt infiltration method (Green and Ampt). Using daily rainfall amounts, SWAT simulates surface runoff volumes and peak runoff rates for each HRU. Since Green-Ampt method need subdaily data, SCS curve number method will be used to estimate surface runoff volumes as is less data intensive [14]. The SCS curve number surface runoff equation (SCS, 1972) is:

(2)

Where;

Q_surf = daily accumulated surface runoff or rainfall excess (mm),

R_day = rainfall depth for the day (mm),

Ia= initial abstractions which includes surface storage, interception and infiltration prior to runoff (mm), and

S = retention parameter (mm).

The retention parameter varies spatially due to changes in soils, land use, management and slope and temporally due to changes in soil water content and is defined as:

(3)

Where;

CN is the curve number for the day. Runoff will only occur when R_day>I_a (=0.2S). The hydrological model component estimates the runoff volume and peak runoff rate that are in turn used to calculate the runoff erosive energy variable.

SWAT calculates the peak runoff rate using a modified rational method (Neitsch et al., 2005) as;

(4)

Where;

= peak runoff rate (m³s^-1)

= fraction of daily rainfall that occurs during time of concentration

= surface runoff volume (mm H₂O)

A= basin area in km²

= time of concentration for the basin (hr) and 3.6 is unit conversion factor

Additional information about runoff calculation can be found in SWAT2005 theoretical documentation (Neitsch et al., 2005).

The SWAT model calculates the surface erosion caused by rainfall and runoff within each HRUs using the Modified Universal Soil Loss Equation (Arnold and Fohrer, 2005) as;

Where;

Sed is the sediment yield on a given day (metric tons),

Qsurf is the surface runoff volume (mm ha^-1),

qpeak is the peak runoff rate (m³s^-1),

Area hru is the area of the HRU (ha),

K USLE is the USLE soil erodibility factor (0.013 metric ton hr2/ (m³-metric ton cm)),

C USLE is the USLE cover and management factor,

P USLE is the USLE support practice factor,

LS USLE is the USLE topographic factor and

CFRG is the coarse fragment factor

The sediment routing model [15] that simulates the sediment transport in the channel network consists of two components operating simultaneously: deposition and degradation. The details of the USLE factors and the descriptions of the different model components can be found in SWAT theoretical documentation (Neitsch et al.,).

Inputs of SWAT Model: The data required for the SWAT model were determined following the information given in Neitsch et al. The spatial data include Digital Elevation Model (DEM), land use land cover, and soil data while non-spatial data include observed daily meteorological and hydrological data.

Digital elevation model (DEM): DEM stands for the elevation of any point in a watershed at a specific spatial resolution. Ii was used as input for SWAT hydrologic model to delineate the watershed, extract information about the topography or slope/ elevation of the watershed and to analyze the drainage patterns of the land surface terrain. Subbasin parameters such as slope gradient, slope length of the terrain, and the stream network characteristics were also derived from the DEM.

For this study 30 m by 30 m DEM of Omo basin of Ethiopia was collected from NASA EARTHDATA. The obtained DEM was made to be mosaic to one raster dataset and clipped to area of interest by using ArcGIS software v10.1.

Land Use Land Cover data (LULC): Study of the effect of land use and land cover change on hydrological processes specially stream flow and sediment yield requires LULC data and the study used the image downloaded and processed from USGS Earth Explore website for three years (1994, 2006 and 2018) as shown in (Figure 2) below for modeling the impact on LULC change on stream flow and sediment yield. The detail sources and preprocessing was explained under paper I. SWAT has predefined land uses identified by four letter codes and it uses these codes to link land use maps to SWAT land use databases in the GIS interface. The land use and land covers in the watershed were with the SWAT code is indicated in (Table 1) as Urban and exposed rocks, (URMD), Water body (WATR), Shrubs and grasslands (RNGE), Woodlands (RNGB), Dense forest (FRSE) and Agricultural lands (AGRL). The LULC definition and corresponding SWAT code.

Figure 2: LULC of 1994, 2006 and 2018 (After image processing in paper I).

Table 1: SWAT land use/cover definition, code and coverage in study area.

Soil Data

Important soil data required by the SWAT model was obtained from FAO-UNESCO Soil map database and from the Omo Gibe River Basin master plan document which was provided by MoWIE. The major soil classification map of the watershed area was made according to the FAO-UNESCO soil classification system. The identified name of soil class was taken from the user soil of the soil map attribute table from which look up (Table 2,3) was prepared and linked to the SWAT data base. The soil types distinguished and used as SWAT input were Chromic Luvisols (LVx), Dystric Vertisol (VRd), Eutric Vertisols (VRe), Humic Alisol (ALu), Humic Nitisols (NTu), and Lithic Leptosol (LPq) with coverage of 13%, 4.48%, 14.78%, 27.24%, 32.04% and 8.46% respectively with the dominant soil class of Humic Nitosols and Humic alisols as shown in (Table 4).

Table 2: Area coverage of major LULC in Gibe III watershed in 1994, 2006 and 2018.

Table 3: Land use and land cover description in the map.

Table 4: Area coverage of soil types of the study watershed used in SWAT data base.

Meteorological Data: Daily Meteorological data (rainfall, maximum and minimum temperature) for station in and around the study watershed required for the task under this study was obtained from Ethiopian National Metrological Agency while hydrological data (stream flow and sediment concentration data) as well as spatial data including soil data and qualitative data describing the watershed was obtained from Ministry of Water, Irrigation and Energy (MoWIE). Since SWAT require relative humidity, sunshine hour and wind speed data as input for station, these data were obtained for synoptic station, Jima, and SWAT weather generator was used to obtain the data for the remaining sites (Figure 3).

Figure 3: Soil types of the watershed from FAO-UNESCO Soil map database

Hydrological Data: Hydrological data including discharge and sediment data was obtained from MoWIE so as to calibrate and validate SWAT model at Abelti and Gojeb sub-watersheds. The calibration and validation of SWAT model was done with 21 years discharge data (1990-2010) while the sediment yield was derived by using rating curve (Figure 4) developed from observed discharge data and small number of observed sediment data obtained from MoWIE. Developing rating curve equation for fulfilling sediment data shortage was also applied by different scholars [16-20] for different watershed during the study of hydrological characteristics of the watershed and to estimate sediment load in the structures. The detail of data coverage was shown on (Table 5) below.

Figure 4: Sediment concentration for Abelti and Gojeb sub-watersheds.

Table 5: The discharge data used for the study obtained from MoWIE.

For Abelti sub-watershed (6)

For Gojeb sub-watershed (7)

Where;

is sediment concentration (mg/l)

Q is discharge at the point of gauge (m³/s)

is the sediment amount in (t/day)

Data quality assessment and preprocessing

Missing data estimation

The observation from rainfall station may have a break in the records because of instrument failure or absence of the observer [21] which requires estimation of the missing records from the neighboring station. The arithmetic mean and normal ratio method were among the testing approaches. Arithmetic mean method can be used to fill in missing data when normal annual precipitation is within 10% of the gauging station for which data are being reconstructed. The normal ratio method is used when the normal annual precipitation at any of the normal index station differs from that of the precipitation by more than 10%. Subramanya, [22] suggested to apply inverse distance weighing when normal annual rainfall for the station is unavailable. In this study the arithmetic and normal ration methods were applied as needed (Eq. 9).

Arithmetic mean method (9)

Py = Amount of precipitation estimated for the time t, from Xa and Xb values of precipitation at time t and TN is total number of years

Normal ration method (10)

Where;

px is annual precipitation at station X that to be estimated at M stations, 1, 2, 3 …m, with the annual precipitation values are p1, p2, p3 ……pm, respectively. The normal annual precipitation N1, N2, N3… Ni at each of the above (m+1) stations including the station x is known.

The quality control can be done by visual inspection, filling of missing data if there is any, by normal ratio method and double mass curve. This can help identify if there are any gaps or unphysical peaks in data series and correct them before the data is used or input to the model.

In this study, Missing hydrological daily data were filled by developing correlation between the station with missed data and any of the nearby stations, having best relation for common data period of specified season (Table 6).

Table 6: The sediment data obtained from MoWIE.

Consistency analysis and homogeneity test of rainfall data

Double mass-curve method that use the annual cumulative total rainfall of the station under study as ordinate and the average annual cumulative total of neighboring stations (base stations) as abscissa was employed to check whether the data is consistent or not. The homogeneity test was done for data of each rainfall stations using XLSTA which tests the homogeneity of station data comparing with long term mean at 5% error level and 95% confidence. Null hypothesis (Ho) says ‘the data is homogenous’ while the Alternative hypothesis (H_a) against Ho says, ‘the data is not homogeneous’. If the P-value is greater than 0.05, the data is homogenous (Table 7). Accordingly, except Dedo and Butajira station data, all others were homogenous and useful for further analysis. Based on the length of year of data and spatial representation of the study area, about eight stations such as Wolyta, Hosaina, welkite, Asendabo, Woliso, Jima, Shebe and Bonga were selected. Multiple linear regression analysis was found to be the best fitting method to fill the missed data based on the available daily climate data (Figure 5,6).

Table 7: Homogeneity test for annual rainfall data (mm).

Figure 5: Consistency test for annual rainfall using double mass curve method for stations.

Figure 6: Meteorological station used for the climate data input.

Preparation of Data for SWAT model

Representative eight weather stations were selected for the watershed modeling process. These stations were having full length quality daily data for precipitation, maximum and minimum temperature. However only one station has relative humidity, sunshine hours and wind speed data and therefore, SWAT weather generator was used to derive these data for remaining stations. To prepare weather generator data, the precipitation and temperature data statistics were required and produced by running pcpSTAT.exe and dew02.exe software for rainfall and temperature respectively).

General approaches of the study

The study structural skeleton of the approach (input/output relationships) for this study was shown in Figure 7. It shows the schematic representations of the steps followed and the parameters interlinked to each other for estimating the impact of LULC change on stream and sediment yield.

Figure 7: Flow chart for SWAT input, processes and output for the study.

Evaluation of the effectiveness of SWAT model: Evaluation of the model performance using observed data was done using the performance measures such as; coefficient of determination (R2), Nash and Sutcliffe simulation efficiency (Ens) and percent bias (PBIAS) [23-25].

The coefficient of determination (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:

(11)

Where, R2 is the regression coefficient,

q_si is the simulated values,

q_{si_av} is the average simulated value,

q_oi is the measured value,

q_{oi_av} is the average measured value,

n is the number of computed values.

Nash and Sutcliffe simulation efficiency (ENs) indicates the degree of fitness of observed and simulated data. The value of ENS ranges from 1.0 (best) to negative infinity. ENS value of 0.0

(12)

Where,

E_NS is the Nash and Sutcliffe simulation efficiency, q_si is the simulated values, qoi is the measured value, q_oiav is the average measured value, n is the number of computed values. As Bekele (2009), SWAT developers assumed an acceptable calibration for hydrology at, R²>0.6 and E_NS>0.5. Goodness-of-fit is quantifiable using percent bias (PBIAS). It assesses the average tendency of simulated data to exhibit under or overestimate (positive or negative BIAS) value respectively.

(13)

Where;

PBIAS is percent of model bias,

q_sim and q_obs are the simulated and observed values and

n is the number of computed value.

Sensitivity analysis, Calibration and Validation of the Model Sensitivity Analysis

Soil and Water Assessment Tool (SWAT) have a large number of parameters that can lead to over-parameterization problem. SWAT input parameters are process based and must be held within a realistic uncertainty range. The determination of the most sensitive parameters for a given watershed or sub-watershed is the prior activity in the calibration and validation process in SWAT model [26]. Sensitivity analysis is a process of testing and identifying model parameters that affects most the output from the model when changed (VanGriensven et al.,). In other words, sensitivity analysis is the process of determining the rate of change in model output when changes in model input parameters occur. Accordingly, a parameter is said to be sensitive if the change in that parameter causes large change on model output. In general, identifying sensitive parameters prior to model calibration helps to allow the possible reduction in the number of parameters that must be calibrated there by reducing the computational time required for modeling. Once the sensitivity analysis is done calibration can be performed for limited number of influential parameters. The current version of SWAT model, SWAT2012, provides the algorithmic techniques for sensitivity analysis. Two types of sensitivity analysis, global Sensitivity and one-at-a-time sensitivity analysis, are allowed when using SUFI2 (Sequential Uncertainty Fitting version 2).

Local (one-at-a-time) sensitivity analysis

The one-at-a-time (OAT) sensitivity analysis is performed for one parameter at a time only by keeping the value of other parameters constant. OAT sensitivity analysis shows the sensitivity of a variable to changes in a parameter if all other parameters are kept constant at some reasonable value. This constant value can be the value of parameters from the best simulation (simulation with the best objective function) of the last iteration. According to Abbaspour (2013), OAT sensitivity analysis has a drawback that the correct values of other parameters that are fixed are never known.

Global sensitivity analysis

Global sensitivity analysis performs the sensitivity of one parameter while the value of other related parameters are also changing. Global sensitivity analysis uses t-test and p-values to determine the sensitivity of each parameter. The t-stat provides a measure of the sensitivity (larger in absolute values are more sensitive) and the p-values determine the significance of the sensitivity. A p-value close to zero has more significance. The drawback related to global sensitivity analysis is, it needs a large number of simulations to get most sensitive parameters [26]. The performance of the SWAT model was evaluated through sensitivity analysis, calibration, and validation.

Calibration and validation

To define calibration, it is test of a model with known input and output information that is used to adjust or estimate factors for which data are not available. Model calibration and validation is a mandatory procedure when using physically based models like SWAT for a watershed hydrological process and sediment yield analysis. Once the models are calibrated, the result can be taken as a representative value and can be used for further analysis [12,26]. The reliability of the results of such models depends on the quality of the hydro-climatic and hydraulic input data. The tiresome and lengthy calibration and validation work of SWAT model needs a wide range of professional expertise and reliable data sources since calibration is tuning of model parameters based on checking simulated output results against observations to ensure the similar response over time. Calibration is comparing the model results (discharge or sediment) obtained with the use of historic model input data to recorded discharge or sediment data to get approaching (similar) results. In this process, model parameters varied until recorded flow patterns are accurately simulated. It is important to select values for the model parameters so that the model closely simulates the behavior and performance characteristics of the study site in order to utilize the predictive watershed model for estimating the effectiveness of future potential management practices. Validation is comparison of the model outputs with an independent data set without making further adjustments. Model validation confirms the applicability of the watershed-based hydrologic parameters derived during the calibration process (Abbaspour, 2013).

The sequential uncertainty fitting (SUFI-2) found in SWATCUP was used to calibrate and validate the SWAT model as it accounts for possible sources of uncertainties such as uncertainty in variables, in model, in parameters and in measured data (Abbaspour, 2007). SUFI-2 algorithm, in particular, is also suggested by Yang et al., as a suitable for calibration and validation of the model because it represents uncertainties of all sources. The degree to which all uncertainties are accounted for is quantified by a measure of P-factor which is the percentage of measured data bracketed by the 95% prediction uncertainitiey-95PPU. The 95PPU is calculated at 2.5% and 97.5% levels of the cumulative distribution of an output variable obtained through Latin Hypercube sampling.

For this study, the Gibe III watershed was partitioned into four sub watersheds namely, Gojeb sub-watershed, Abelti sub-watershed, Wabi sub-watershed and Un-gauged sub-watershed (Figure 8) and sensitivity analysis, Calibration and Validation of the Model was done at Abelti and Gojeb a sub-watersheds only and the model best parameter after calibration and validation was transferred by regionalization to Wabi and un-gauged sub-watersheds.

Figure 8: Gibe III sub-watersheds.

Making the process of calibration at the sub watershed level could increase confidence to capture the variability in the predominant processes for each of the sub watersheds instead of determining at global (watershed-wide) processes [12].

The sensitivity analysis was done for flow and sediment separately since some parameters are sensitive to flow and sediment, some sensitive to flow only and others sensitive to sediment only. Therefore, the sensitivity of the parameters for flow and sediment was done separately to know the most important parameters using discharge and sediment data for a period of 21 years (1990 to 2010). The data of three years (1990-1992), eleven years (1993-2003) and seven years (2004-2010) were used for warm up period, calibration and validation respectively. Since the sediment data was generated from rating curve, the period of calibration and validation were similar to flow [27-29]. Arnold et al., also applied SWAT to estimate average annual sediment loads for five major Texas river basins sediment yield from rating curves and observed that the SWAT predicted sediment yield reasonably well with estimated sediment yields obtained from rating curves (Figure 9). After sensitivity analysis of parameters at two sub-watersheds, the SWAT model was calibrated and validated for each period by dividing the stream flow and sediment yield data in to two (calibration data and validation data) [30,31].

Figure 9: Calibration and validation of monthly mean flow at Abelti sub watershed.

The performance of SWAT model was evaluated during calibration and validation processes using coefficient of determination (R2), Nash and Sutcliffe efficiency (Ens) [23] and Percent difference between simulated and observed data and percent bias (PBIAS) as recommended by Moriasi et al. (2007) [32-35]. The general performance of NSE in SWAT according to (Moriasi, et al.,) is NSE > 0.65 is very good, NSE between 0.5 and 0.65 is adequate, NSE >0.5 is satisfactory and NSE < 0.5 is unsatisfactory both for calibration and validation. After routine procedures of parameter selection, sensitivity analysis, calibration and validation, SWAT model was used for estimation of stream flow and sediment yield from three land use types (lulc1994, lulc2006 and lulc2018). The calibration and validation performance was evaluated based on [23,36,37] as shown in the (Table 8) below.

Table 8: General performance ratings for the recommended statistics.

Regionalization of parameters to un-gauged sub watersheds

Watershed modeling requires basically the discharge and sediment data for the calibration of the selected model. However, the challenge of insufficient or absence of discharge measurements forces the modeling work to find methods to transfer hydrological information from gauged to un-gauged catchments by regionalization principle. Regionalization method helps to identify similar or proxy catchments, and transfer a model parameter set calibrated on a gauged donor catchment to the target catchment [38-41]. Catchments having apparently similar physical characteristics are assumed to have a similar hydrological behavior (Figure 10). Acreman and Sinclair, [42] and Nathan and McMahon, (1990) took part for application of catchment similarity principle based identifying similar basin. In this study, two subwatersheds such as Wabi and un-gauged sub-watersheds have no observed data for calibration and validation of SWAT model and regionalization was used to transfer best parameters from donor catchment to them. Before transferring the best parameters obtained at donor catchment, physical similarity approach was used to know if two adjacent catchments were similar [43-45] using the equation 14. Physical similarity measure was applied based on comparison of catchment characteristics including catchment topography, land cover types, and soil type. These characteristics are assumed to be major drivers of the hydrological processes and catchment runoff response. The physical similarity among catchments was measured by means of a weighted Euclidean distance:

Figure 10: Calibration and validation of monthly mean sediment yield at Abelti sub watershed.

(14)

Where,

S = similarity index of catchment a to catchment b

= the Euclidean distance between catchment a and b,

J = catchment descriptor at the a^th catchment and b^th catchment respectively

= the weight attributed to the j^th catchment descriptor

Application of equation 14 involves measures generally having different units and scales, and therefore requires a standardization of the descriptors. The standardization was carried out by dividing each descriptor by the maximum of the descriptor as follows;

(15)

Where, = the value of the catchment descriptor at the k^th catchment before standardization. Weights were given by:

(16)

Where;

Is the difference among the j^th descriptor of the catchments

The Gibe III watershed has shown land use changes between consecutive three years (1994, 2006 and 2018) as discussed in chapter I which could affect the hydrological processes specially, stream flow and sediment yield of the watershed. The analysis result of modeling the impact using SWAT model was discussed as follows.

SWAT model application at Abelti sub-watershed

Sensitivity analysis for the model parameters: After tedious and repeated efforts of parameter selection, procedural steps in a SWATCUP setup, the most sensitive parameters were obtained for Abelti gauge station for stream flow and sediment yield as shown at Table 9 below.

Table 9: Rank of Sensitive parameters for flow and sediment modeling at Abelti station.

The sensitivity analysis result denotes that, the SCS_CN for moisture condition II (Cn2), Lateral flow travel time (LAT_TTIME), groundwater ″revap″ coefficient (GW_REVAP.gw), Average slope steepness (HRU_SLP), soil bulk density (SOL_BD.sol), Manning's "n" value for overland flow (OV_N), Effective hydraulic conductivity in main channel alluvium (CH_K2, Average slope length (SLSUBBSN), Soil conductivity (SOL_K) and soil available water capacity (Sol_ AWC) were found to be top ten sensitive parameters that could affect the hydrological process in the area among all tested twentytwo parameters. Similarly, for the case of sediment, about eight most sensitive parameters were observed and believed to have high effect on the simulated values. Accordingly the sensitive parameters for sediment yield in the watershed were identified as USLE cover or management factor (USLE_C), the SCS_CN for moisture condition II (Cn2), exponential factor for channel sediment routing (SPEXP), USLE support practice factor (USLE_ P), Manning's "n" value for the main channel (CH_N2), Average slope steepness (HRU_SLP), linear factor for channel sediment routing (SPCON and Effective hydraulic conductivity in main channel alluvium (CH_K2) were found to be highly sensitive parameters [46-48]. Although these parameters are very important, considering the limitations of parameterization with respect to physical meaning of parameter in the area (parameters expected to affect stream flow and sediment yield but not top sensitive), additional parameters were considered for calibration and validation for both variables. After intensive simulation and iteration, about eleven parameters for flow and seven (not include common parameters of flow & sediment) parameters were selected (Table 10).

Table 10: Calibrated parameters of flow and sediment yield at Abelti sub watershed.

Calibration and Validation of flow and sediment: Flow calibration was performed with great caution to prevent over-parameterization, loss of meaning of parameters and loss of reasonability of the output. It was done with eleven parameters for flow while seven parameters of sediment. The selected sensitive parameters, parameter value ranges and fitted value result was shown in Table 10 below.

Evaluation of flow and sediment yield calibration and validation result: The performance of the model to simulate the stream flow during the calibration and the validation periods was evaluated based on the computed results of the indicators and the suggested model performance rating standards. At Abelti station, the computed statistical indicators shown in (Table 11) for stream flow resulted to R2, NSE and BIAS of 0.90, 0.87 and -5.8% for calibration and 0.82, 0.77 and 14.8% for validation period respectively. According to Moriasi et al. recommendations, the statistical indicators show very good model performance for calibration both calibration and validation periods. In terms of sediment yield result, the model showed better simulation statistics during calibration and the validation periods. The computed statistical indicators shown in Table 11 for sediment yield at Abelti resulted for 0.87, 0.86 and -6.0% for R², NSE and PBIAS at the at the calibration period while 0.75, 0.73 and 10.9% for R², NSE and PBIAS for validation period indicating that the sediment yield is well reproduced by the model at the station. The match between monthly measured and simulated flows in the calibration and validation period were demonstrated by the acceptable results of R², NSE and PBIAS as well as the mean monthly values obtained.

Table 11:  SWAT model performance evaluation statistic for flow and sediment calibration and validation at Abelti sub-watershed.

The mean monthly observed and simulated flow for calibration was 218.75 m³/s and 231.79 m³/s respectively while that of validation resulted to observed amount of 226.37 m³/s and simulated amount of 192.92 m³/s. The flow result was simulated well approaching to the observed discharge. This confirms that the model could be applicable for the use in the watershed. On the other hand, the monthly sediment yield at the Abelti catchment was found to be 1.05 Mton and 1.11 Mton for sediment obtained from rating curve and simulation respectively while 1.01 Mton and 0.9Mton for observed and rating curve result for validation respectively. The good agreement between observed and simulated values indicates the better performance of the model.

The model performance was very good for stream flow and good for sediment yield modeling. For the gauging stations considered for model calibration, the stream flow was measured continuously and the data has been in a good quality as compared to the suspended sediment yield data. The sediment data has been generated from the sediment rating curves and its result cannot be depended better than flow. However, the model performance for sediment yield in Abelti is in an acceptable range and therefore it can be a basis for further analysis or for sediment management at watershed.

SWAT model application at Gojeb sub-watershed

Sensitivity analysis for the model parameters

As it was done for Abelti sub watershed, the same procedure was followed to identify sensitive parameters for both flow and sediment and the result were shown in (Table 12) below. According to Narsimlu et al., [49], the ranks of sensitivity, sensitive parameters can be selected depending on global sensitive analyses p-value and t-statistic. The larger t-stat value in absolute values or P-value close to zero has more significance or more sensitive taking prior rank in sensitivity analysis. At table 12, the rank for each parameter was assigned using superior value of t-stat and smaller value of p-value since they indicate most sensitive parameters (Abbaspour et al.,) From twenty-two flow parameters, the top eight parameters including; groundwater delay (GW_DELAY), Base flow alpha factor for bank storage (ALPHA_BNK), the SCS curve number (CN2), the Maximum canopy storage (CANMX), effective hydraulic conductivity in the main channel (CH_K2), Groundwater "revap" coefficient (GW_REVAP), base flow alpha factor (ALPHA_ BF), and threshold depth of water in the shallow aquifer required for return flow to occur (GWQMN) were identified as top eight sensitive while top parameter such as Average slope steepness (HRU_SLP), USLE cover or management factor (USLE_C), Channel cover factor (CH_COV2, Average slope length (SLSUBBSN) and Initial residue cover (kg/ha)(RSDIN) were found to be top five parameters affecting sediment yield (Table 12). During sensitivity analysis, parameters that were expected to affect sediment process such as slope run off and channel parameters were taken in to consideration as sediment yield parameters for calibration and validation parameters identified in sensitivity analysis that influence predicted outputs are often used to calibrate a model (VanGriensven et.al.,) (Figure 11).

Table 12:  Sensitivity analysis for stream flow and sediment yield at Gojeb sub watershed.

Figure 11: Representative photomicrograph of histopathological features in pulmonary necropsies.

Calibration and validation of flow and sediment yield

Overall statistical results at Gojeb catchment showed better performance of SWAT model to simulate the stream flow and sediment yield during the calibration and validation periods. Accordingly the computed values of performance rating evaluation statistics showed 0.81, 0.80 and 0.0% for R², NSE and PBIAS at the at the calibration period and 0.78, 0.76 and 9.7% for R², NSE and PBIAS for validation period for stream flow respectively. According to Moriasi et al., the statistical indicators show very good model performance at calibration and validation period. The figure below showed that the model underpredicted the peak flow at about ten events. The (Figure 12) below showed that over-estimation was observed in low flows while high flows were slightly underestimated. The scholars such as Chu and Shirmohammadi, Spruill et al., [50,51] and Kanishka and Eldho, (2020) also agreed that the SWAT model mostly under-predicts peak flows. The similar situation of the model was also reported by. Approaching results of mean observed and simulated discharge values were obtained for both flow calibration (62.02 m³/s & 61.99 m³/s) and validation periods (66.37 m³/s & 59.92 m³/s) respectively indicating the potential applicability of the model in the area (Table 13).

Figure 12: Area rainfall (inverted) and observed discharge (upright) for comparison of linkage.

Table 13: Calibrated parameters of flow and sediment yield at Gojeb sub watershed.

Evaluation of flow and sediment yield calibration and validation result

The figure 13 indicates the direct linkage of monthly area rainfall characteristics with discharge in the watershed which confirms the quality of flow data at the outlet. It was also observed that the time of concentration can be seen from the chart that the peak flow occurs after peak rainfall in all peaks of both events. The sediment yield result at Gojeb was not reproduced better like Abelti station. However the model showed good simulation statistics during calibration and the validation periods with R², NSE and PBIAS of 0.73, 0.73 and 5.4% for calibration and 0.60, 0.60 and 2.1% for validation period respectively (Table 14) indicating that the sediment yield is well reproduced by the model in order to apply the model. The mean observed and simulated sediment yield were obtained for both calibration (56×10³ ton & 53×10³ ton) and validation periods (54.6×10³ ton & 53.5×10³ ton) respectively. The result shows that the model can be considered robust. Therefore, it could be said that, the model performance was very acceptable for stream flow and sediment yield modeling in the sub-watershed. Since SWAT used the observed discharge from rating curve, it showed discrepancy in terms of sediment yield simulation result especially at peak values. To determine the sediment yield from the watershed, the model under estimated the sediment yield during high flows where the simulated runoff was less than the observed runoff. However, the model performance for sediment yield in is in better performance range and the model result can be a basis for sediment management at watershed scale [52].

Figure 13: Calibration and validation of monthly mean flow at Gojeb sub watershed.

Table 14:  SWAT model performance evaluation statistic for Stream flow and sediment.

Uncertainty analysis for Abelti and Gojeb sub watersheds

The degree to which all uncertainties are accounted for is quantified by a measure of the p-factor, which is the percentage of measured data bracketed by the 95% prediction uncertainty (95PPU) (Abbaspour et al.,). Similarly, r-factor is used as the measure quantifying the strength of a calibration or uncertainty analysis which is the average thickness of the 95PPU band divided by the standard deviation of the measured data (Abbaspour). Theoretically, a P-factor of 1 (all observations bracketed by the prediction uncertainty) and R-factor of 0 (achievement of rather small uncertainty band) indicate that the simulation exactly corresponds to the measured data but p-factor ranges between 0 and 100%, while that of r-factor ranges between 0 and infinity. The degree to which they are away from these numbers can be used to judge the strength of our calibration. When acceptable values of r-factor and p-factor are reached, then the parameter uncertainties are the desired parameter ranges (Moriasi et al.,) [53-55].

The result of flow calibration and validation at Gojeb gave p-factor of 0.70 and 0.67 respectively while r-factor 0.69 and 0.66 respectively. Calibration and validation of sediment indicated the p-factor of 0.55 and 0.25 respectively while r-factor 0.73 and 0.38 respectively. The result indicates that 70% of measured data was bracketed by 95 ppu range of uncertainty in the output for flow calibration while 67% during validation. The result is inacceptable range. During the calibration and validation, the thickness of uncertainty band was obtained as 95ppu of the r-factor 0.69 & 0.66 for flow and r-factor of 0.73 & 0.38 for sediment respectively. Setegn et al., [13] justified that low p-factor and large r-factor in uncertainty simulation could be due to the error in the rainfall and temperature input data which may be the same reason in this study during validation periods to get 38% r-factor. The area has complex topography and highly converting land use land cover condition which could also increase uncertainty more than model parameter errors [56].

Parameter regionalization using catchment physical similarity analysis

The physical parameters such as LULC result, soil property coverage and slope level of the sub-watersheds indicated in (Table 15) below were used for physical similarity index determination based on the equation 14 to equation 16 explained above and the result was shown in Table 15 below.

Table 15: Parameters applied for physical similarity index determination.

The index value more than 0.6 shows high physical similarity of watersheds [43]. According to the result at (Table 16), Abelti and Wabi are physically similar with index value of 0.74. In addition Abelti and un-gauged sub-watershed are also similar with index value of 0.63.

Table 16: The result of physical similarity of sub-watersheds to each other.

More similarity was observed between Un-gauged sub-watershed and Wabi than un-gauged catchment with Abelti with index value of 0.75 and 0.63 respectively. However, due to limitation of data at Wabi, calibration and validation of the SWAT model was not done at Wabi so that parameter transfer to Un-gauged catchment was carried out from Abelti which was the donor for both sub-watersheds. On the contrary, Gojeb sub-watershed is not similar to any of three catchments [57-59].

Effect of LULC Change on stream flow and sediment yield

Since there was high expansion of agricultural lands in the expenses of other lands in the area, simulation of the impacts of land use and land cover change on stream flow and sediment yield was considered to be the most important. The mean annual stream flow and sediment yield of the study watershed was simulated at three different land use and land cover conditions since highly changing land use could affect stream flow timely (Figure 14).

Figure 14: Model equation for observed and simulated flow at calibration and validation.

Twenty-five years (1993-2017) climate data was used with each LULC map (1994, 2006 and 2018) as input to calibrated SWAT model for simulation and the output was shown in the table 17. Form each sub-watershed, high flow rate was obtained at Un-gauged (285, 168.6 and 241.6 m³/s) while low amount of flow was obtained at Wabi subwatershed (27.4, 20.82 and 20.73.6 m³/s) at 1994, 2006 and 2018 LULC respectively.

Table 17: Mean monthly stream flow and sediment yield at different land use and land cover by sub-watersheds and Gibe III.

The result in Table 17 shows that Gibe III has mean annual stream flow of 614.29 m³/s, 446.09 m³/s and 515.93 m³/s was obtained under each LULC map (1994, 2006 and 2018) respectively. The study conducted by Teshome and Koch (2013) also obtained 521, 552 and 530 m³/s for the 2001-2010, 2011-2020, 2021-2030 periods respectively using SWAT model at Gibe III which confirms the current result. The mean annual flow of 438 m³/s was reported by Tesfaye and Bogale (2014) at Gibe III which was less than the current (2018 LULC) simulated amount by 17% while equal to that of 2006 LULC. It is clearly observed that the flow is decreasing with significant decreasing change was obtained between 1994 and 2006 while less difference observed between 2006 and 2018. Although there was agricultural land coverage increase from 1994 to 2006, the reduced flow may be due to agricultural practice that can increase infiltration during 2006 than 1994 because broad definition of agricultural land in Table 3 above could hide the role of different practices with agriculture (Figure 15).

Figure 15: Calibration and validation of mean monthly sediment yield at Gojeb sub watershed.

In each LULC the sediment yield was higher at Abelti sub-watershed (1.0, 4.9 and 4.7 Mton during 1994, 2006 & 2018 respectively) and at Un-gauged sub-watershed (1.4, 1.1 and 1.2Mton during 1994, 2006 & 2018 respectively) which was related to the higher sub-watershed area coverage (Figure 8). In terms of sediment load, the increasing trend was significant between 1994 & 2006 than between 2006 & 2018 land use changes. However the mean monthly amount of sediment load obtained at three different land use land cover were 3.1 Mtons, 7.2 Mtons and 7.3 Mtons from 1994, 2006 and 2018 land use land cover conditions which indicate high risk of sedimentation. This may be attributed to degradation of forest and shrub land being changing to agricultural land during lulc of 1994 to 2006 (Figure16).

Figure 16: Model equation for observed and simulated sediment at calibration and validation.

The total annual sediment load was also observed to increase from 1994 to 2018. It was shown that the annual sedimentation load during 2018 at Abelti, Gojeb, Un-gauged, Wabi and Gibe III were 56.7, 8.5, 14.5, 8.3 and 88.2 Mton respectively (Table 17). The study conducted by Devi et al. indicated that annual load of sedimentation at Abelti was 45Mton during the study period which is 26% lower than current amount at Abelti. Since the study was before 10 years, the current result is consistent to it. According to Das, the dry bulk densities of soil materials ranges between 1.0 to 1.6 g/cm³, but, soils with large amounts of organic material could have dry bulk density less than 1.0 g/cm³. Taking the average dry bulk density of 1.3 g/cm³, the annual total volume of reservoir that could be occupied by 88.2 Mtone sediment will be 67.8 Mm³. The result indicates that under this rate of sedimentation 50% of the reservoir storage volume 11,750 Mm³ of Gibe III will be occupied after 86 years or 80% will be filled by sediment after 138 years. Anonymous (2008) said that a high rate of sedimentation is anticipated in the Gibe III reservoir where one-third of its space is reserved for sediment to accumulate over time [60].

The progressing increase of sedimentation in reservoir increases the risk of reduction in storage capacity of dam. Overall increase in sedimentation at Gibe III watershed in 2018 relative to 1994 LULC was 133% with average annual increase by about 44%. It indicates that, making other drivers of sedimentation set constant, the contribution of LULC for sedimentation in the watershed is very high [61].

The Monthly distribution of stream flow and sediment yield at Gibe III for the selected LULC was shown at (Table 18) above. Thus the flow was highly reduced from 1994 to 2018 in all months while in each month effect of LULC was very high from 1994 to 2006 than between 2006 and 2018. The monthly distribution of the sediment load shows the maximum sediment load occurs during August month and minimum in the month of February.

Table 18: Monthly mean stream flow and sediment yield at different LULC at Gibe III watershed.

The non-significant change in stream flow simulated between 2006 and 2018 was related to the minimum change in agricultural land between the respective years (Figure 17).

Figure 17: Change of major LULC in each of sub-watersheds and Gibe III watershed for the years 1994, 2006 and 2018.

Since Gibe III watershed has different sub-watersheds that contribute flow and sediment to it, it was better to know which subwatershed contribute the load and what was governing factors in generation of runoff. Accordingly, figure 17 shows that, Abelti, Wabi and un-gauged sub-watersheds were affected by intensive agriculture which covers more than 80% of their area. However Gojeb subwatershed is dominated by forest and shrub as well as woodlands than agriculture in all the selected periods. The watershed of Gibe III that contributes 80% of the Omo basin flow (WWDSE and SDCSE, 2013) has agricultural land coverage more than 70% (Abdella) that could play major role in sedimentation.

Land and water resources were basically interlinked to human livelihood improvement and ecosystem conservation. The dynamics of land use and land cover changes play great role in managing water resource potential in general and hydroelectric power service specially. The impact of LULC changes on stream flow and sediment yield in Gibe III watershed, Omo-Gibe basin Ethiopia under the LULC conditions of 1994, 2006 and 2018 was assessed using SWAT model. The SWAT model was calibrated and validated at Abelti sub-watershed and Gojeb sub-watershed using SWATCUP SUFi2 algorism using observed discharge data of 21 years (1990-2010) and sediment estimation from rating curve from discharge data.

Method of physical similarity index for regionalization was used to transfer best parameter of SWAT to other two sub-watersheds (Wabi and Ungauged). After selecting flow and sediment sensitive parameters, calibration and validation of the SWAT at Abelti and Gojeb subwatersheds, regionalization technique was employed to transfer best parameter to un-gauged catchments from Abelti which was tested to be to similar catchments for regionalization and finally flow and sediment load at Gibe III was estimated. Performance of SWAT was evaluated using the evaluation statistics such as R², NSE and PBIAS for accuracy assessment while P-factor and r-factor for uncertainty analysis. Accordingly, SWAT satisfactorily estimated flow and sediment yield at the Abelti sub-watershed with R², NSE and BIAS of 0.90, 0.87 and -5.8% during calibration and 0.82, 0.77 and 14.8% during validation of stream flow while 0.87, 0.86 and -6.0% during calibration and 0.75, 0.73 and 10.9% during validation of sediment yield indicating that the flow and sediment yield is well simulated by the model in the study watershed. The mean monthly observed and simulated flow amount of 218.75 m³/s and 231.79 m³/s at calibration while 226.37 m³/s and 192.92 m³/s at validation respectively. In the case of Gojeb sub-watershed, R², NSE and PBIAS of 0.81, 0.80 and 0.0% during calibration period and 0.78, 0.76 and 9.7% during validation for stream flow respectively. The observed and simulated values for flow calibration (62.02 m³/s & 61.99 m³/s) and validation (66.37 m³/s & 59.92 m³/s) respectively were obtained confirming that the model could be used for simulation of hydrological processes in the watershed. Sediment calibration resulted to R², NSE and PBIAS of 0.73, 0.73 and 5.4% during calibration and 0.60, 0.60 and 2.1% during validation period respectively. The mean observed and simulated sediment yield was (56×10³ tone & 53×10³ tone) during calibration and (54.6×10³ tone & 53.5×10³ tone) during validation respectively. The mean annual stream flow at Gibe III watershed was obtained as 614.29 m³/s, 446.09 m³/s and 515.93 m³/s was obtained under each LULC map (1994, 2006 and 2018) respectively which indicate high reduction in flow from 1994 to 2018.

Similarly, mean monthly amount of sediment load obtained at three different land use land cover were 3.1 Mtons, 7.2 Mtons and 7.3 Mtons from 1994, 2006 and 2018 LULC conditions which indicate high risk of sedimentation. It was shown that the annual total sediment load during 2018 at Abelti, Gojeb, Un-gauged, Wabi and Gibe III were 56.7, 8.5, 14.5, 8.3 and 88.2 Mton respectively. About 67.8 Mm³ storage volume of reservoir will be occupied by sediment per year indicating that 50% of the volume of reservoir storage (11750 Mm³) would be filled by sediment after 86 years of its initial period. Institutionalized and well organized community based watershed management could reduce the sedimentation of the water body in general and the dam in particular. Furthermore, it is recommended that the gibe III watershed more focus on proper management of Abelti sub-watershed since it is major contributor for sedimentation. In general, from this study, the following recommendations could be forwarded;

Since the study did not included the exploration of major driving forces fot expansion of agricultural land, father study is required on socio-economic factors that governs the LULC change in the area.

Hydrologic model like SWAT could be applied to predict the potential impacts of land use change on the stream flow and sediment yield for improving watershed ecosystem service in general that could help stakeholder and decision makers to do best option during land and water resource planning and development.

Understanding the future consequence of climate change on water resource potential in the future using climate models and scenarios is very important to prevent risk of flow reduction and sediment load increment in the future.

On the other hand, sediment sources in the watershed are spatially variable so that it should be identified and prioritized based on degradation level in order to implement soil and water conservation action.