Kinematic wave hydrologic model of the Turbio River basin, Guanajuato, Mexico

Modelo hidrológico de onda cinemática de la cuenca del río Turbio, Guanajuato, México

A model is a representation of a real system, natural or built by man (Dingman, 2002), which is simpler than a real system, so it represents only some of its characteristics (Doodge, 1986). In hydrology, in general, models are usually divided into two categories: a) physical, that is, laboratory-scale, and b) mathematical, which describe the system and the phenomena that occur in it through equations. In this sense, several authors recognize that these models enable us to better understand how a certain system works and predict its behavior, which can improve its control and operation (Dingman, 2002; Doodge, 1986; Ponce, 1994).

According to Ponce (1994), many hydrologic models, especially those related to surface water, are developed in the frame of reference of a basin. The model of a basin has become a set of mathematical abstractions that describe relevant phases of the hydrological cycle in order to simulate the conversion of precipitation into runoff. Dingman (2002) states that most current research in hydrology is focused on improving the prediction of water levels in aquifers, runoff variability, water quality (from the time it comes down from the slopes and crosses the basin until it empties into the sea), and the effects of land-use and climate change on the hydrological balance.

In Mexico, the use of a specialized program for hydrologic modeling is recent. Studies have been reported with the Soil and Water Assessment Tool (SWAT), Hydrologic Engineering Center-Hydrologic Modeling System (HEC-HMS) and Hydrologic Engineering Center-River Analysis System (HEC-RAS) (Vargas-Castañeda, Ibáñez-Castillo, & Arteaga-Ramírez, 2015), which may be freely downloaded, although the methodologies of Mexican federal institutions encourage the use of HEC family programs for the case of defining urban areas prone to flooding. The Centro Nacional de Prevención de Desastres (CENAPRED, 2011) recommends using HEC-RAS and hydrologic methodologies corresponding to the runoff curve number (CN) and the unit hydrograph (UH) of the Soil Conservation Service (SCS), both included in the HEC-HMS. Both HEC-HMS and HEC-RAS (rainfall-runoff models) are more recommended for flood prevention and warning through simulations by increasing events and with increases over time. On the contrary, SWAT is more recommended for watershed management with simulations of long periods (months, years, etc.) and with one-day increments; that is, the SWAT specialty is not flood prevention and warning (Vargas-Castañeda et al., 2015).

Based on the experience gained in various studies, the series of steps to be executed in a typical rainfall-runoff model, in the HEC-HMS program, can be: 1) convert the rainfall depth into runoff depth with the CN, 2) convert the runoff depth into a hydrograph by means of the UH and 3) route the flood hydrograph through a channel with the Muskingum method and through a dam with the continuity equation, that is, through a volume balance (Juárez-Méndez, Ibáñez-Castillo, Pérez-Nieto, & Arellano-Monterrosas, 2009; Miranda-Aragón, Ibáñez-Castillo, Valdez-Lazalde, & Hernández-de la Rosa, 2009; United States of America Corps of Engineers [USACE], 2000).

A theoretical alternative with a greater physical basis to that already described, also contained in the HEC-HMS program, is the kinematic wave (KW) approach. Since its development, by Lighthill and Witham (1955a, 1955b), it has been widely used in environmental and water sciences. Examples of its application include surface and current flow, baseflow, unsaturated flow, macropore flow, furrow and bank flow, river hydraulics, glacier movement, erosion, and microbial sediment, solute and chromatographic transport, among others. These applications are used in various branches of hydrology, including surface hydrology, vadose zone hydrology, fluvial and coastal hydrology, irrigation, subsurface hydrology and water quality (Singh, 2001). In the last two decades, KW equations have been used in several cases and models (Ajami, Gupta, Wagener, & Sorooshian, 2004; Cheah, Ball, & Cox, 2008; Chua, Wong, & Sriramula, 2008; Chua & Wong, 2010; Huang & Lee, 2013; Jinkang, Shunping, Youpeng, Xu, & Singh, 2007; Kazezyilmaz-Alhan & Medina, 2007; Mejia & Reed, 2011; Rai, Upadhyay, & Singh, 2010).

The aim of the present study was to apply the HEC-HMS methodology based on KW theory in the Turbio River basin, which has a history of damage due to the occurrence of extreme rainfall. It is hypothesized that this methodology can increase the model’s accuracy in predicting storm runoff, which is an important element in water use and the protection of the population from flood disasters.

The Turbio River basin (Figure 1) is located at the western end of the state of Guanajuato, on the border with Jalisco, Mexico. Its border to the North, which in turn is the highest part, is the Sierra de Comanja and to the South the Lerma River, which in turn is the lowest elevation point. The Guanajuato Civil Protection Agency reports that 67 % of the basin area is located in this state and only 33 % in Jalisco.

Figure 1
Turbio River basin to the "Las Adjuntas" hydrometric station.

Table 1 shows the most important characteristics of the basin under study. The rains and their floods that have occurred since the 1970s have resulted in flooding in the municipalities of southwestern Guanajuato and caused damage in aspects such as health, property damage, interruption of economic activities and agricultural losses. The major floods so far of this millennium occurred in 2003 and 2007; therefore, in this work the rainfall-runoff behavior of the rainfall event that occurred from July 24 to July 31, 2007 was simulated, since detailed data of storms recorded by automatic weather stations (AWSs) are required to calibrate the model (Table 2). The floods of 2003, although they were more severe than those of 2007, were not modeled because this technology was not yet available in the basin.

Table 1
Characteristics of the Turbio River basin, Guanajuato, Mexico.

Table 2
Rainfall depth in the Turbio River basin from July 24 to 27, 2007.

¹AWS: automatic weather station.

The various rainfall-runoff model options developed in this work were carried out in the HEC-GeoHMS (USACE, 2013) and HEC-HMS (USACE, 2000, 2010a, 2010b) programs. The basin model and its physical parameters were obtained with the HEC-GeoHMS geotechnical extension designed for ArcGIS 10.1 (Environmental Systems Research Institute [ESRI], 2012), and were executed with the rainfall information represented through 13 hyetographs spatially distributed in the Turbio River basin with the Thiessen polygons method. These hyetographs fed the model at 1-h intervals. On average, the rainfall event of the evaluated period was 114 mm, the minimum value was 72 mm and the maximum 152 mm (Table 2).

For the event studied, three rainfall-runoff model alternatives were executed. In all three, the conversion of rainfall depth to direct runoff depth was carried out with the CN methodology. The methods that varied among the models were: a) the transformation of runoff depth to a direct runoff hydrograph (unit hydrograph or kinematic wave) and b) flood routing in channels (Muskingum or kinematic wave). In order to compare the predictive capacity of the models, three different combinations were made: 1) unit hydrograph-Muskingum (UH-Musk), 2) kinematic wave-Muskingum (KW-Musk) and 3) kinematic wave-kinematic wave (KW-KW). The three model alternatives executed the flood routing in the dams located in the Turbio River basin. These variants were developed in HEC-HMS and generated a hydrograph at the exit of the "Las Adjuntas" station.

Model calibration

The quantitative measure of the goodness of fit between the simulation results and observed discharge is called the objective function and is equal to zero if these are identical. The objective function of minimum errors is reached when the values of the parameters that best reproduce the measured hydrograph are found (USACE, 2010b).

Among the objective functions offered by the HEC-HMS are those that minimize the root mean square error (RMSE). However, when it comes to minimizing the prediction error in peak flows, it is recommended to use the peak-weighted RMSE objective function, which is a modification of the standard RMSE function, which gives greater weight to flows above average and less weight to those below average (Cunderlik & Simonovic, 2004; USACE, 2010b). Therefore, in this work, when calibrating the following error function is minimized:

$PWRMSE=\sqrt{\frac{\sum _{t=1}^N{\left(Q_{obs}^t-Q_{sim}^t\right)}^2\frac{Q_{obs}^t+{\stackrel{-}{Q}}_{obs}}{2{\stackrel{-}{Q}}_{obs}}}{N}}$(6)

where PWRMSE is the objective function to minimize peak-weighted RMSE, Q^t_sim and Q^t_obs are the flows simulated and observed at time t, N is the number of observations in the period considered and ${\stackrel{-}{Q}}_{obs}$ is the arithmetic mean of the flows observed. The calibration uses measured flow data; however, the program proposes value range intervals that are possible for the USACE (2000) parameters to have, although the modeler can limit these values more.

Once the reported simulations were carried out, the necessary parameters were calibrated to obtain a better fit, taking as reference the hydrograph measured in "Las Adjuntas." To calibrate, both the CN and some flood routing parameters (Muskingum or KW) were used. After different tests, it was more useful to calibrate the models by fitting the CN for each sub-basin in a range of ± 10 %, a common practice among modelers (Arnold et al., 2012). The HMS model of the Turbio River presented in this work has 43 sub-basins.

Finally, at the moment that the HEC-HMS minimized the errors, between the simulated and observed values, only the curve numbers of 17 sub-basins were varied in UH-Muskingum and KW-Muskingum. In the UH-Musk model, the largest CN fit was -10 %, but on average they varied -6 %, and in the KW-Musk, the largest fit was -13 %. In the KW-KW model, the best CN fit achieved was -3 % in eight sub-basins, since when the fit was made in the same sub-basins as in the previous ones the results worsened.

An important guide to calibrate the models was the sensitivity analysis provided by the model. For example, when calibrating the UH-Muskingum model, it was observed that when a parameter does not provide improvement in the objective function it reports a value of 0.0. In Table 3, it is observed that when the CN of sub-basins W840 and W860 is varied, an improvement in the objective function is not reported, while in sub-basins W580 and W490 the CN improves the objective function.

Table 3
Sensitivity analysis of the runoff curve parameter in several sub-basins of the model in the unit hydrograph-Muskingum simulation.

¹CN: runoff curve number.

According to Moriasi et al. (2007), the quality of a hydrologic model can be evaluated based on the value taken by the Nash-Sutcliffe efficiency index (NSE) (Table 4).

Table 4
Evaluation of hydrologic models according to Moriasi et al. (2007).

¹NSE: Nash-Sutcliffe efficiency index.

$NSE=1-\frac{\sum _{t=1}^N{\left(Q_t^{obs}-Q_t^{sim}\right)}^2}{\sum _{t=1}^N{\left(Q_t^{obs}-{\stackrel{-}{Q}}^{obs}\right)}^2}$(7)

Additionally, the percent bias (PBIAS), which measures the average tendency of the simulated data to be greater or less than the observed values, was calculated (Moriasi et al., 2007):

$PBIAS=\frac{\sum _{t=1}^N\left(Q_t^{obs}-Q_t^{sim}\right)}{\sum _{t=1}^{obs}{Q}_t^{obs}}\times 100$(8)

A PBIAS of 0.0 indicates a very accurate model; if the value is positive it indicates that the model underestimates and if it is negative that the model overestimates. According to Moriasi et al. (2007), a PBIAS of less than 15 % is desirable when calibrating flows. However, Cunderlik and Simonovic (2004) mention that at the time of calibrating an hourly-type hydrograph, the lowest flows affect the PBIAS, and it is more desirable that there is no such bias in the largest flow values when it comes to a flood warning.

Figure 2 shows the hydrographs resulting from the three simulated model alternatives; the one that presents the values obtained from the "La Adjuntas" station shows a peak flow of 51.7 m³·s^-1 on July 28 at 4:00 p.m. This hydrograph is the benchmark for comparison with the three simulated models. In this sense, the KW-Muskingum method gives better results; however, the KW-KW model is the best at predicting the time at which the peak flow occurs, although it has the greatest error in volumes.

Figure 2
Hydrograph of simulated and observed flows from July 24 to 30, 2017 at the "Las Adjuntas" station of the Turbio River, Guanajuato, Mexico.

Table 5 shows the observed and simulated peak flows for each of the methodologies used. In this work, the UH-KW combination is not presented because the NSE reported values of -0.20, which is not acceptable.

Table 5
Peak flows of the simulations.

Table 6 shows the fits of the surface hydrologic model for the three methodological combinations, in addition to the statistics to choose the best model. According to the NSE, all three models are satisfactory (Moriasi et al., 2007). The best combination for the Turbio River rainfall-runoff model is the KW-Musk one, and the worst is the UH-Musk method. In the best methodological combination, the conversion of runoff depth to a hydrograph (flow) was made with the KW approach and the flood routing in channels with the Muskingum method, although it can be seen that the most used methodology in Mexico (UH-Muskingum) also shows a good fit, although of lower quality. In reality, the UH-Muskingum combination is preferred in Mexico because it is easier to have data available to feed the hydrologic model.

Table 6
Efficiency indices obtained in the three models.

¹NSE: Nash-Sutcliffe efficiency index; r²: coefficient of determination; PBIAS: percent bias; RMSE: root mean square error. *In all three methods, the hourly rainfall depth is converted into runoff depth with the runoff curve number.

Figure 3 shows the dispersion behavior of the simulated flow predictions. In the KW-Muskingum and KW-KW methodologies, the simulated values present more error when it comes to predicting small flows, so that their r² is lower than that of the UH-Muskingum method since the error occurs in small flows.

Figure 3
Dispersion of simulated versus observed flows in the Turbio River basin, Guanajuato, Mexico.

The best simulation of the rainfall-runoff process in the Turbio River basin is given by the kinematic wave-Muskingum combination; that is, the runoff depth becomes a hydrograph with the kinematic wave model and flood routing in channels is carried out with the Muskingum method. The Nash-Sutcliffe efficiency index obtained with this model was 0.709, which is considered good. The second-best simulation was with the kinematic wave-kinematic wave model, and the third with the traditional unit hydrograph-Muskingum.

According to the Nash-Sutcliffe coefficient, the three methods give satisfactory results, but the kinematic wave-Muskingum combination gives better results in terms of errors when estimating total flows, peak flow and runoff volume.

The authors thank the Guanajuato State Water Commission (CEAG) for providing the rainfall data from its automatic weather stations.