Hydrological Response of Land Use Land Cover Change on Stream Flow of Anger Watershed, Blue Nile Basin, Ethiopia

Stream flow of Anger watershed (7,982 km2), located in Western Oromia Regional State, Ethiopia was simulated by using soil and water assessment tool (SWAT) model. Assessment of layer staking, mosaicking, image classification and accuracy was carried by ArcGIS. The main target of this study was to determine the SWAT Model applicability and impact of past land-use land-cover change on stream flow of the Anger watershed. The time series data from 1990 to 1999 was used for model parameter calibration, and data from 2000 to 2004 were used to validate the model using the input parameter set. The performance of the SWAT model was evaluated by both R2 time series plots and ENS the statistical measures. The results of the model during calibration for monthly stream flow was R2= 0.90, ENS= 0.76 and during validation R2= 0.84, ENS= 0.63 respectively. Validation of image classification is achieved through Mosaicking image and Google Earth image. The performance of the Arc GIS model was evaluated by using the overall accuracy and kapa coefficient. The agreement between classified land use image and ground truth of the year 1986, 2000 and 2010 were showed a good agreement with result of overall accuracy=88%, 85% and 87% and kapa coefficient=81%, 79% and 82% respectively. Simulation of steam flow of the watershed under the impacts of past land-use land-cover changes can be determined by using SWAT model as a tool for water resources planning and management in this watershed.


Introduction
Water resources are currently under severe pressure because of impacts of climate change and human activities, which include land-use land-cover change, increasing population growth, and economic development [1]. Assessing water resources becomes a complex task that must consider many aspects, of which climate and land-use land-cover change and variability are identified as key factors [2]. Whereas land-use land-cover change can cause changes in hydrological components in a sub watershed, such as evapotranspiration, stream flow and groundwater. Thus, separating the effects of climate and land-use land-cover changes on stream flow is important to produce accurate predictions of stream flow simulation and to provide useful information to land-use planning and water resources management [3]. Semi-distributed hydrological models, which use input parameters directly representing land surface characteristics, have been applied to assess the impact of land-use land-cover change on stream flow in water resource management areas [4].
Some research findings identified that, land-use land-cover change was the main driving factors for the decrease or increase in stream flow simulation. In Anger sub watershed, there is no studies considering land-use land-cover change have been performed in the sub watershed, therefore, this research focus on stream flow simulation under the impact of land-use land-cover change in Anger sub watershed. Hence, Anger sub watershed stream flow simulation under impact of land-use land-cover change was particularly important in improving water management efficiency and benefiting various water use needs such as irrigation and small-scale hydro power on Anger sub watershed. The overall objective of this study is to model stream flow under the impact of land-use land-cover change by using SWAT model. The specific objectives are: (1) to calibrate and validate SWAT model for study area (2) to quantify the effect of past landuse land-cover change on stream flow of sub watershed through scenario simulation. This study provides important information that decision-makers will need in order to assist with water resources management in the Anger sub watershed under changing land-use/ land-cover conditions.

Description of the Study Area
The Anger sub watershed is located in the East Wollega Zone, Oromia Regional State, Ethiopia, between latitudes 1 0 . 95° N to 11.80° N and longitudes 36.70° E to 37.40° E (Figure 1). The sub watershed has an area of 7,982 km 2 and covers parts of fifteen woredas, namely: -Wayu Tika, Bila Seyo, Sibu Sire, Guto Gida, Sasiga, Gudaya Bila, Bilo Jegonfoy, Abe Dongoro, Horo, Gida Kiremu, Yaso, Limu, Haro, Jarti Jadeja, and Amuru. The climate of the A n g e r watershed is 'tropical highland monsoon' with an average annual rainfall of 1,604 mm. Most of the rain falls during the months of June to September with peaks occurring during July to August and it is virtually dry from November through to April. As the watershed is located in a high rainfall area, it receives frequent torrential showers and frequent flash floods during the rainy season. The mean monthly maximum and minimum temperature of the area varies from 22 to 30°C and 11 to 15°C respectively ( Figure 2).
The major landform of the watershed includes flat to gently sloping, undulating plains, hills and mountains. The western part of the watershed is characterized by highly rugged, mountainous and rolled topography with steep slopes and the lower part is characterized by a valley floor with flat to gentle slopes. Elevation in the watershed varies from 868 to 3,144 m above mean sea level. The major portion of the watershed is under intensive cultivation and maize, barley and wheat are the major crops grown in the watershed. Shrub land, grazing land, forest, woodland and wetland/swamp are other land cover types in the watershed.
The catchment has a wide range of soil types mainly dominated by Haplic Alisosls, Haplic Acirisols, Rahodic Nitisols and Dystiric Leptosols. The largest portion of the watershed is characterized by Haplic Alisosls. And also; Haplic Nitisols, Eutric Vertisols and Eutiric Regosols were part of watershed.

SWAT Model Description
The Soil and Water Assessment Tool (SWAT) is a physical process-based model to simulate continuous-time landscape processes at a catchment scale [5]. The catchment is divided into hydrological response units (HRUs) based on soil type, land use and slope classes that allows a high level of spatial detail simulation.
The major model components include hydrology, weather, soil erosion, nutrients, soil temperature, crop growth, pesticides agricultural management and stream routing [5]. The model predicts the hydrology at each HRU using the water balance equation, which includes daily precipitation, [5].
Where: Q surf is the accumulated runoff or rainfall excess (mm) R day is rainfall depth for the day (mm) S is retention parameter. SWAT model provided three options for potential evapotranspiration calculation. Penman Monteith [10], Priestley Taylor [11], and Hargreaves [12] methods. The methods were selected based on data availability.
SWAT model partitions groundwater in to two aquifer systems: a shallow, unconfined aquifer, which contributes return flow to stream within the watershed, and deep, confined aquifer which contributes return flow to streams outside the watershed [13]. SWAT model uses the following equation to simulate ground water for both confined and un confined aquifer.
where: -Aqshi: the amount of water stored in the shallow aquifer on a day I (mm), aqsh i-1: the amount of water stored in the shallow aquifer on day i-1 (mm), Wrchrg: the amount of recharge entering the aquifer on day I (mm), Qgw: The ground water flow or base flow, or return flow, in to the main channel Wrevap: amount of water moving in to soil zone in response to water deficiencies Wdeep: the amount of water percolating from the shallow aquifer in to the deep aquifer on day I (mm), and Wpumpsh: the amount of water removed from the shallow aquifer by pumping SWAT model have both variable storage and Muskingum methods of flow routing to route the flow of the sub watershed.
Where: ∆V stored is a volume of during the time step (m³) V in is a volume of outflow during the time step (m³) V out is changed in the volume of storage during the time step (m³)

Geographical Information System (GIS)
Geographical information system (GIS) is the power full model that used for image processing and image classification. To identify change in land use land cover distribution in the Anger watershed over 25-year period from 1986-2010 land sat imageries of seven band were used. during the period of 1986, 2000 and 2010 land sat TM and ETM+ were used respectively. The image data files were downloaded in zipped files from the Global land cover facility (GLCF) website and extracted to Tiff format files. And also, composed the band using the RGB color composition and clip the image covering the complete study area.   Image classification is the process of assigning of pixels of continuous raster image to the predefined land-use landcover classes. The result of the classification is mostly affected by various factors such as classification methods, algorithms, collecting of training sites etc.
In remote sensing there are various image classification methods. Those are supervised, unsupervised and hybrid. The supervised classification type was applied for this study and the land cover map was produced based on the pixel based supervised classification thought the step such as: selecting of the training sites which are typically representative for the land-use land-cover classes, image enhancement and image. Using these approaches same training sites were collected as from each image (1986, 2000 and 2010) and perform the classification using the Maximum Likelihood classifier.

Filling Missing Data
Some precipitation stations may have short breaks in the records because of the absence of the observer or because of instrumental failures. It is often necessary to estimate or fill in this missing record. The missing precipitation of a station is estimated from the observations of precipitation at some other stations as close to and as evenly spaced around the station with the missing record as possible. There are methods to fill in missing data. These are the arithmetic mean method, normal ratio method, and inverse distance weighting method. Arithmetic mean method can be used to fill in missing data when normal annual precipitation is within 10% of the gauge/station for which data are being reconstructed. The normal ratio method is used when the normal annual precipitation at any of the index station differs from that of the precipitation station by more than 10%. In the absence of Normal annual rainfall for the stations, the inverse distance weighting method can be used to fill the missing data. According to [7] the two formulas are described below.

i) Arithmetic Mean Method
Where Px is the missing precipitation, Pi is precipitation at ? ℎ station and n is the number of nearby stations. (6) P x is the missing precipitation for any storm at the interpolation station x, P i is the precipitation for the same period for the same storm at the i th station of a group of index stations. N x is the normal annual precipitation for station x, and N i is the normal annual precipitation value for the i th station.

Data Homogeneity Test
Homogeneity analysis was used to separate a change in the statistical properties of the time series data. The causes can be either natural or man-made. These include alterations to land use and relocation of the observation gauging station. Therefore, in order to select the representative meteorological station for the analysis of areal rainfall estimation, checking homogeneity of group stations is essential, the homogeneity of the selected gauging stations daily rainfall records was carried out by Non-Dimensional Equation: [8].

Data Consistency Test
The double mass analysis will use for checking the consistency of a hydrological or meteorological record and is considered to be an essential tool before taking it for analysis purpose.
The accumulated totals of the gauges are compared with the corresponding totals for a representative group of the nearby gauge. The double mass curve is based on the principle that when each recorded data comes from the same parent population, they are consistent [9] If the cumulative plot of DMC shows different, gradient it should be corrected as:

Accuracy Assessment
Accuracy assessment is part of the image classification process and its objective was to evaluate the total number of correctly classified pixels divided by the total number of ground truth pixels. The most widely used classification accuracy is in the form of error matrix which can be used to derive a series of descriptive and analytical statistics [15]. User's accuracy and producer accuracy measured the correctness of each category with respect to errors matrix. The users' accuracy is defined as the probability that a reference pixel has been correctly classified as well as the producer accuracy is defined as the probability that a pixel classified on the map represents that class on the ground.
The accuracy of thematic maps was determined by the constructed matrices along kappa statistics in order to test whether any difference exists in the interpretation. Kappa statistics (Eqn.) considers a measure of overall accuracy of image classification and individual category accuracy as a means of actual agreement between classification and observation. [15] suggested the use of subjective kappa value as < 40% as poor, 40-55% fair, 55-70% good, 70-85% very good and > 85% as excellent. The generally the estimate of Kappa was computed as follows: Where, Po= proportion of observed agreements Pc= proportion of agreement expected by chance.

Model Performance Evaluation by Using SWAT-CUP
SWAT-CUP is a public domain program linking the SUFI-2 procedure to SWAT. SWAT-CUP provides a decision making frame work that incorporates a semi-distributed approach using both manual and automated calibration incorporating sensitivity and uncertainty analysis. The Sequential Uncertainty Fitting, version 2 (SUFI-2) is one of the uncertainty analyses programs that is incorporated in an independent program called SWAT calibration and Uncertainty Program (SWAT-CUP) [6], that perform uncertainty analysis due to both parameter and model uncertainties. Its main function is to calibrate SWAT and perform validation, sensitivity and uncertainty analysis for a watershed model created by SWAT. SUFI-2 is developed for a combined calibration and uncertainty analysis to find parameter uncertainties while calculating smallest possible prediction uncertainty band. In SUFI-2, parameter uncertainty accounts for all sources of uncertainties such as uncertainty in driving variables (e.g., rainfall), parameters, conceptual model, and measured data (e.g., observed flow, sediment) [6].
The performance SWAT was evaluated using the statistical measures to determine the quality and reliability of simulate when compare to observed values. The performance of model in simulating stream flow and sediment yield is evaluated using SWAT-CUP. The selection and use of specific efficiency criteria and the interpretation of the results can be a challenge for even the most experienced hydrologist since each criterion may place different emphasis on different types of simulated and observed behaviors. For this study coefficient of determination (R 2 ) and Nash-Sutcliffe Efficiency (NSE) was choose for the criteria of model performance.
i. Coefficient of Determination (R 2 ): is an indicator of strength of relationship between observed and simulation value. Coefficient determination ranges from 0 to 1, with higher values indicating less error variance, and typically values greater than 0.5 are considered acceptable [14]. The following formula used to calculate R 2 : Where: R 2 is coefficient of determination, Oi is observed hydrologic variable, Om is mean observed hydrologic variable, Si is model simulated output, Sm is mean of the model simulations and N is total number of observations. ii. Nash-Sutcliffe Efficiency (NSE): is determines the relative magnitude of the residual variance compared to the measured data variance is calculated as follows.
Where: S is model simulated output, O is observed hydrologic variable, Om is mean of the observations that NSE uses as a benchmark against which performance of the hydrologic model is compared and n is total number of observations. NSE indicates how well the plot of observed versus simulated data fits the one ratio to one line. According to [14] NSE ranges from negative infinity to 1.0, where NSE is 1.0 shows that a perfect model. NSE 0.0, implies the observed mean is as good a predictor as the model. Generally, the value between 0.0 and 1.0 acceptable model performance where below 0.0 is not acceptable performance and the model simulation can be judged a satisfactory if NSE >0.50

Root Mean Square Error Observation Standard Deviation Ratio (RSR)
RSR varies from the optimal value of "0", which indicate zero roots mean square error (RMSE) or residual variation and therefore perfect model simulation, to a large positive value. Generally, if the value of RSR ≤ 0.70 the model simulation can be considered as satisfactory [15]. Where: Si is model simulated output; Oi is observed hydrologic variable and Sm is mean of the model simulations.

Percent Bias (PBIAS)
PBIAS is the deviation of data being evaluated, expressed as a percentage. If PBIAS ± 25% for streamflow and PBIAS ± 55% for sediment, the model simulation can be judged as satisfactory. Where: Si is model simulated output, Oi is observed hydrologic variable.

Over All Methodology
The following framework illustrates the general workflow of the study Figure 7. General framework of the study.

Accuracy Assessment
The accuracy assessment is used to identify the precision of the classified image. It was evaluated by using both the original mosaic image and the Google Earth image as a reference. Therefore, the following randomly selected points 130, 131 and 141 were used for validation during the year of 1986, 2000 and 2010 images respectively. Tables 4, 5 and 6 show a confusion matrix for the three Landsat images.

Overall Accuracy
The overall accuracy was achieved by dividing the total number of accurate pixels by the total number of pixels in the confusion matrix. For the maps of 1986, 2000 and 2010 the result was 88%, 87% and 92% respectively. According to [15], the minimum accuracy value for reliable land cover classification is 85%.

Kappa Coefficient
The kappa coefficient gives the agreement between classified image and reference or ground truth. A kappa coefficient equal to 1 means perfect agreement where as a value close to zero means that the agreement is no better than would be expected by chance. The results indicated that the kappa coefficient for the maps of 1986, 2000 and 2010 were 81%, 79.5% and 83.2% respectively. Therefore, based on tables 4, 5 and 6 the classification carried out in this study produces a kappa coefficient with substantial agreement for the year 2000 and almost prefect agreement for year 1986 and 2010.

Producer's Accuracy
The producer's accuracy was found by dividing the number of properly classified pixels in the class by the total number of pixels of the class in the reference data. Its overall result ranges from 58%% to 100%. The lowest values were achieved due to misclassified of different land-use land-cover classes.

User's Accuracy
The User's accuracy was the ratio between the total number of pixels properly fitting to a class and the total number of pixels allocated to the same class. its overall accuracy ranges from 78% to 100%. The lowest value was achieved due to misclassified because of the similarity spectral properties of forest and water body.

Land Use Land Cover Change Analysis
The maps of the year 1986, 2000 and 2010 were generated from the land sat TM and ETM+ imaginary. The classified land-use land-cover of the watershed shown that there is an increasing of agricultural land and built-up area. And also, decreasing of forest, shrub, grassland and water bodies. In general, during 25-year period the agricultural land was increased by 39.3% whereas the Shrub and grass land decreased by 32.9% and also forest area decreased by 6.6%.

Sensitivity Analysis
The simulated stream flow for the watershed was evaluated by using the monthly observed flow data for identifying the most sensitive parameter and for further calibration of the SWAT model. Twenty-six flow parameters were checked for sensitivity and seven of them were found to be highly sensitive. For, detail information presented in (Table 8) below. From those parameters Groundwater delay (GW_DELAY), SCS runoff curve number (CN2), Average slope steepness (HRU_SLP), Soil evaporation compensation factor (ESCO), Base flow alpha factor for bank storage (ALPHA_BNK), Channel effective hydraulic conductivity (CH_K2) and Manning's value for main channel (CH_N2) were sensitive and ranked from 1 up to 7 respectively.

Flow Calibration
The simulated and observed stream flow had a good agreement with the default value of parameters in Anger watershed. After sensitivity analysis has been done, the calibration of stream flow was done manually and Auto calibration respectively for the period of 10 years from 1990-1999. The result of agreement between observed and simulated stream flow (Figure 10) shown by statical determinant with value of NSE, (R 2 ) 0.76 and 0.90 respectively.

Uncertainty
It is necessary to consider uncertainties in predicting stream flow of the watershed. There are different sources of uncertainties. Those Sources of uncertainty are: -Input and Output data.

Model Validation
After calibration, validation was carried out for 6 years i.e from 2000-2005 were performed. The model validation also showed a good agreement between simulated and observed monthly flow as shown ( Figure 11) with the NSE, R 2 of 0.63 0.84 respectively.

Land-Use Land-Cover Change Impact Analysis
Effects of land-use land-cover change on average annual surface runoff of Anger sub watershed shown that surface runoff volume increased due to land use land cover change

Conclusion
In Anger sub watershed, the SWAT model was used to simulate the effect of past land use land cover change on runoff volume change and Arc GIS was used to process preliminary data, extract, layer stake, mosaic satellite images, and image classification. Performance of the model was confirmed, with a NSE and R 2 .
Due to the effect of urbanization, deforestation and agricultural expansion the past land use land cover changes increased runoff volume of the sub-basin by 12.54 Mm 3 , 25.29 Mm 3 during 1986-2000 and 2000-2010 respectively.
As research findings shown that, much of the original forest in the sub watershed has already converted to agricultural lands, for instance agricultural land increased by 39.3% and built-up area increased by 0.5% and also, shrubland, forest area, and water body decreased by 32.9%, 6.6%, and 0.3% respectively. These changes would have aggravated problem related to water scarcity in dry period and hill-slope erosion during wet period The accuracy of LULC Classification was achieved by:overall accuracy, user's accuracy, and producer's accuracy and kappa coefficient. The kappa coefficient shows good agreement between classified image and ground truth with 81%, 79.5% and 83.2% during 1986, 2000 and 2010 respectively.