Application of HEC-HMS Model to Estimate Daily Streamflow in Badddegama Watershed of Gin Ganga Basin Sri Lanka

Rapid urbanization causes stress in limited water resources. Water managers need quantified streamflow as high, medium and low flow regions to assess and manage water resources. This requires calibrated and validated models together with a rationalized method for the selection of a proper model, components, data, temporal resolution, objective function, and performance criteria. The current preference is to use process-based models for streamflow modelling. The widely used process-based HEC-HMS 4.2.1 version model was chosen to simulate daily streamflow at Baddegama watershed (749 km2) while using 2007 to 2012 period for calibration, and 2013 to 2017 for verification. Simple canopy, initial deficit and constant loss, SCS direct runoff, and recession baseflow were selected as model components. Semi-automatic optimization was done with RMSE objective function. Model was calibrated and validated with RMSE of 3.0 mm/day and 3.5 mm/day. Model with RMSE value of 6.2 mm/day for high flows and 2.3 mm/day for medium flows indicate better capability of model on flooding management and water resource. The HEC model showed very good matching of daily flows in the high and intermediate regions but reflected a poor matching in the low flow region. Hence the use of HEC-HMS model for the considered regions can be recommended for flood and water resources management. The use of HECHMS for low flow management must be with caution.

Existing water resources and future water demands should be assessed in order to evaluate the impacts on society and the environment because of any planned water infrastructure development [1]. Hence, such assessments performed in a reliable manner are vital for sustainable water resources management.
Hydrological processes in a watershed are complex. Representation of processes in a watershed using scientific methods is needed to simulate these hydrological processes accurately. A sufficient collection of a continuous series of streamflow (SF) data is necessary to perform sustainable water resource management using scientific watershed models [2]. However, lack of past observations and advances in technology have caused problems when using hydrological models to quantify streamflow [3]. Loucks [4] had mentioned that hydrological models are versatile tools to manage watersheds by evaluating the upcoming effects of proposed water management plans. Along with the advantage of scientific representation, a process-based model enable the incorporation of hydrological processes with space-time variability [5].
Comprehensive physically-based distributed models require information such as land use, meteorological data and properties of the watershed for their parameterization and adjustment [6]. Though there is a preference to use elaborate models, the data constraints must be carefully considered when an appropriate model is selected for a particular catchment to ensure optimum watershed management [7], [8]. However, main problem is the lack of guidance to select a suitable rainfall-runoff model. A suitable hydrologic model can be selected by using several selection criteria.
Criteria of selection mainly depend on the specification of a particular project [9].
Mimicking the physical processes in a watershed by using a process-based model requires the inclusion of fundamental governing equations to represent each representative process of the hydrologic system [6]. The most widely accepted physically-based models which are lately developed or regularly updated can be named as GSSHA, HSPF, KINEROS2, HEC HMS, AnnAGNPS, SWAT, GSSHA, WinSRM, PRMS, HYPE, WetSpa, and MIKE-SHE [10], [11]. With compared to recently developed and updated models, HEC-HMS is providing more options to simulate [11].
HEC HMS is a free watershed model widely applied in the world and to a certain extent in Sri Lanka [12]- [15]. HEC HMS model is relevant to different kinds of geographic areas for solving watershed issues [16]. In case of Sri Lanka, Siriwardana & Wijesekera indicated that HEC HMS model would be a better option for planning irrigation reservoirs in ungauged watersheds [14], Jayadeera & Wijesekera concluded that capability of model on flood and water resource management at Kalu Ganga basin by evaluation of potential water management [13]. Present work is carried out to estimate and assess daily streamflow at Baddegama watershed by application of HEC-HMS 4.2.1 version model for sustainable water resource management.

Study Area
Gin ganga, the fifth largest river in Sri Lanka, is located approximately between longitudes 80°08" E and 80°40" E, and latitudes 6°04" N and 6°30" N. Gin Ganga basin is entirely located in the wet zone and faces frequent flooding. Sandy clay loam is the main soil type in the watershed and temperature varies from 24° C to 32° C. The watershed has a considerable rainforest cover in its upper watershed [17]. The river gauging stations and rainfall stations with monthly variation of rainfall are shown in Figure 1. The studied watershed area covers approximately 749 km 2 .

Data and Data Checking
Daily rainfall (RF) data of six rainfall stations at Thawalama, Anningkanda, Deniyaya, Neluwa, Hiniduma, Baddegama, monthly evaporation data at Kottawa (nearly 18 km away from Baddegama river station) and daily streamflow at Baddegama were collected for the period from 2007 to 2017. Kottawa is the only nearest evaporation station which is having high percentage of available data for the selected period. The latest data period from 2007 to 2017 has been selected because higher percentage of availability of data and lesser percentage of data missing.
Land use, topographic and contour data in 1:50000 scale were collected from the Department of Survey. Annual water balance (AWB) was performed to check the water budget at each water year. Visual data checks were done to identify anomalies of the pattern. Daily rainfall vs streamflow variations at each station were plotted and compared. Double mass curves were plotted for all the stations to check consistency of data and no data inconsistency could be observed.

3.2.2
Basin Model A lumped model was developed for Baddegama watershed because lumped modelling approach is a good alternative to the complex physically-based models when the main focus is only streamflow prediction [20]. The model components, objective function, flow thresholds for evaluation, model evaluation criteria, calibration and verification requirements were selected for this task by comparing other options with its relevance to the purpose of the current study. As a model structure, only one sub basin was selected from the basin model to represent a lumped model. There was no reach and junction components since river routing was not considered for lumped model.

Canopy Model
Model's user manual clearly states that canopy interception should be incorporated with an initial deficit and constant method in continuous modelling [21], [16]. In order to represent interception and evapotranspiration, user can select the canopy component from the three methods that are available, as dynamic, simple and gridded simple canopy [16]. Simple canopy is preferable for this study because of the compatibility with the loss method and applicability in continuous process [22]. This value may be very small because of the high impervious percentage of 40% at Baddegama. Impervious percentage was determined by the percentage of homestead and rock land use percentage from land use distribution at Baddegama. Initial canopy storage was estimated by optimization while considering both dry and wet seasons to determine evapotranspiration.

Precipitation Loss Model
Various methods are available to simulate losses. Event modelling consists of several options as initial constant, SCS CN (Soil Conservation Service Curve Number), Exponential, Green Ampt, and Smith Parlange. For primary continuous modeling, usually Deficit and constant loss method was used. This study also used Deficit and constant loss method because it considers the regain of initial loss after a protracted duration of no rainfall. Also it is most suitable for continuous simulation and having a smaller number of parameters [23].
Initial deficit, maximum deficit and constant loss rate are the parameters of selected loss method and those were estimated by optimization. Initial values for the parameters were estimated as follows.
Initial value of constant loss rate was taken as 1.27 mm/hr which is the minimum loss rate given for hydrological soil type C [21]. The weighted CN value of 79.9 was computed by averaging CN values for each land use type which is from the standard tables for AMC (antecedent moisture condition) II and type C hydrological soil group [24]. SCS abstraction method formula was used to determine maximum potential retention which was taken as the maximum deficit. The initial deficit parameter was computed by the formula which provided initial abstraction as equal to maximum retention times 0.2 [24]. Table 3 shows this initial parameter values of the Baddegama watershed.

Transform Model
HEC HMS provides seven options to convert excess rainfall into direct runoff [16]. SCS method was used because i) SCS method adaptation in various environments and generate better results; ii) only two types of parameters make easy calculation; and iii) reliable and excellent results can be obtained as same as complex models [25]. Lag Time is the only parameter in SCS UH method and it was calculated by considering the relationship of lag time with time of concentration (TC). TC was calculated using the Kirpich formula. Lag time for watershed was calculated as 3220 min.

Baseflow Model
To represent baseflow, five methods are available in model [16]. Baseflow was calculated by the exponential recession method because its applicability in continuous process [26], lesser number of parameters and the compatibility with loss model. Sri Lankan studies used baseflow recession because its ability to automatically reset after each storm event [13]- [15]. Recession baseflow model consists of three parameters such as initial discharge, recession constant and ratio to peak. In this, ratio to peak was taken as 0.2 and it was justified by comparing the behaviour with the observed hydrographs. Groundwater flow recession constant was taken as 0.95 by selecting from typical values proposed for watershed area range 300 km 2 -1600 km 2 [21]. Initial discharge at the beginning was taken as the initial baseflow value [22].

Model Evaluation
Model evaluation was carried out using six criteria. They are 1) annual water balance; 2) hydrograph matching; 3) flow duration curve (FDC) matching; 4) FDC matching for high flows; 5) FDC matching for medium flows; and 6) FDC for low flow matching.

Flow Threshold Separation
According to the order of magnitude of probability of exceedance and variation of corresponding FDC gradient, 15% and 70% probability of exceedance, were selected as high and low flow thresholds (Figure 2) for Baddegama watershed [27].

3.3.2
Objective Function In this study, RMSE was selected as the objective function because it is most widely used to compare observed and estimated streamflow in HEC HMS [28]. HEC HMS model has two automatic optimization search algorithms, namely Univariate Gradient (UG), and Nelder and Mead (NM) methods [29]. The NM search algorithm is used to improve search direction [30]. The NM method is used more frequently than the UG method. This is mainly because the NM method uses a downhill simplex to assess all parameters simultaneously [20]. Considering these reasons, NM was used for this study. Fully automatic calibration results were not good and a semi-automatic calibration approach was used to optimize the parameters.

3.3.3
Simulation Controls HEC HMS guide recommends that simulation time interval should not exceed 0.29 times lag time for a sub basin [21]. Accordingly, 6 hours simulation time interval was selected. Beginning and end of the date from 2007 to 2012 period was taken as calibration period. A five time repetition of calibration data set was used to warmup the model to reduce the effect of initial conditions [31]. Warmup period was adequate to stabilize the initial conditions for water year. Data from 2013 to 2017 was used for model verification.

Watershed Status
Average annual rainfall, streamflow and evaporation over the study period are 3857 mm, 2802 mm and 913 mm, respectively. Runoff coefficient during the period 2007 to 2017 varied between 0.7 and 0.8. The highest runoff coefficient value was indicated in the 2015/16 water year while the 2013/14 water year indicating the lowest value for the runoff coefficient. The lowest and highest evaporation values can be seen in 2014/15 and 2008/09 water years, respectively. Baddegama daily streamflow indicates high low flows during less rainfall period throughout the calibration and verification period. This may be due to the closeness of river gauge to the sea and flat terrain. However, data collection station's authorities confirmed that there is no backwater effects at the Baddegama gauge even though there is a saltwater barrier. These were the observations about the streamflow at Baddegama.  overall matching in hydrograph (Table 1). In Figures 3 and 5, most of the high peaks were captured by the model but there were shifts in the magnitude of peak flow occurrence in small peaks.  Figure 4) also showed that the model does respond well during low flow periods but verification ( Figure 6) period indicates poor matching. In both calibration and verification periods, model does not reach the best fitted of highest flow values. The annual water balance in calibration and verification periods were 3 mm and 187 mm/year (Table 2), respectively. The optimised parameters are in Table 3.

Model Component Selection
The crucial decision of this study was the selection of appropriate model components according to the study objective. Canopy, precipitation loss, transform and baseflow components need most relevant method to calculate streamflow accurately. According to several criteria and justification by literature, sub models were selected as simple canopy, initial deficit and constant, SCS UH for transform, and exponential recession base flow. Hydrograph (Figure 3 and Figure 5) for the entire period showed that similar pattern and magnitude of simulated streamflow with compared to observed. Therefore, selected components are capable to simulate most accurate streamflow at Baddegama.

5.2
Model Performance Evaluation criteria and objective function are key factors to assess the model. In this study, a total six number of criteria were considered to evaluate the model. Table 2 indicates that lesser annual water balance percentage error for calibration with compared to verification. However, hydrograph error matching is not showing that much error difference between calibration and verification. FDC for high flow error shows that higher deference of error between calibration and verification periods. Hence, one evaluation criteria is not sufficient to justify the model performance. Therefore, evaluation of all these six categories together causes to indicate a better picture of developed model performance.
RMSE objective function was taken to optimized Baddegama model. In this study mm/day is the unit for RMSE. RMSE statistic considers an average of deference between observed and simulated. Therefore, when optimization running, it is more favourable to high flows. However, just focusing only on low flow region affects the matching of medium and high flow regions during optimisation. In optimisation, RMSE objective function was able to give best fit to high and moderate flows with more considerable matching to low flow regions. Hence, objective function is key factor to develop a reliable model.

5.3
Model Calibration Satisfactory RMSE value of 3.08 and 3.55 mm/day indicated in hydrograph matching in calibration and verification, respectively. According to model evaluation guideline, ratio of RMSE to the standard deviation of measured data is considered to evaluate model performance. The ratios are 0.46 and 0.53 for calibration and validation, respectively.
It indicates very good model performance rating for calibration, and good performance rating for verification ([32], [33]). Hydrographs in both periods reveal that matching of flow pattern is not satisfactory, especially in low flows, while RMSE reflects acceptable less error values. Hence, hydrograph matching in low flows needs more improvement. However, similar pattern fittings of peaks in hydrographs reveal Baddegama model's capability to use in flood management.
During the wet season all regions of flow show a considerable level of matching at both periods. However, during dry season, models were not capable to respond to lesser rainfall and high rainfall after dry or no rainfall period. This may be due to spatial variabilities of rainfall over the watersheds and may be some observe data errors in some years. These reasons may lead to a high objective function value during calibration and verification.
FDC plots are illustrated in Figure 4 and Figure  6 for both periods. It is observed that the model does not respond well during low flow period in verification compared to calibration period. In Baddegama, 15% is a high flow margin and 70% is a low flow margin. During calibration period, FDC shows a 6.28 mm/day value of RMSE for high flow and it is slightly underestimated with better matching. Medium flow is showing good matching and slight overestimation with 2.32 mm/day value of RMSE. Low flow matching shows perfect matching with 1.21 mm/day value of RMSE. In verification period, FDC shows 7.45 mm/day value of RMSE for high flow matching and it is underestimated with better matching. The medium flow is matching perfectly with 2.58 mm/day value of RMSE. Low flow matching at Baddegama shows underestimation with moderate matching of 1.21 mm/day RMSE value. Watershed's flow duration curve fitting for medium flow is very good and this may lead to a model capable of using in water resource management better than in drought management.
Average annual water balance error is -0.3 % for -3 mm in calibration while it is -16.7 % for -187 mm for verification. As per Table 2