The HEC-RAS hydrological model

This tool enables hydraulic modelling of a water’s permanent and temporary flow patterns in artificial canals and natural channels, including rivers (^{US Army Corps Engineers & Hydrologic Engineering Center, 2006}). This software’s hydraulic simulation is based on deterministic differential equations that enable the prediction of water level dynamics that occur during events with high rainfall that cause flooding. Flood levels are defined by cross-sectional profiles. The dynamics of the water and channel behaviour are simulated, including: cross-sections having any type of geometry along a channel, different depths of water and variable flow along a channel in sub-critical or super-critical flow conditions, having hydraulic effects due to natural or artificial transverse obstacles in the channel (^{Sarhadi, Soltani, & Modarres, 2012}; ^{Mohammadi, Nazariha, & Mehrdadi, 2014}).

HEC-RAS software (for the simulation model used in this research) uses a continuity equation (^{US Army Corps Engineers & Hydrologic Engineering Center, 2006}) that describes the conservation of mass for a one-dimensional system, as well as calculates storage terms:

where x = distance along the channel, t = time, Q = flow, A = cross- section area, S = storage and q ₁ = lateral input per unit of distance.

The following calibration parameters were used in the modelling:

Hydrographs: graphs enabling the flow rate or flow to be observed at a given point on the current (^{Chow et al., 1994}).

Flow curves or calibration curves: graphical representations of the relationship between the water level and its respective flow (^{Salazar & Chaparron, 1990}).

Cross-sections: these define a river’s shape and geometric characteristics and must be topographically connected so that they define the longitudinal profile.

where K = cross-section, n = Manning coefficient of roughness for the section, A = the section’s flow area and R = the section’s hydraulic ratio (area/wetted perimeter).

Manning’s roughness coefficient: also called the roughness coefficient, which enables a channel’s runoff resistance to be estimated (^{Ruberto, Carreras, & Depettris, 2003}). When there are several Manning coefficients (nc) for a channel’s roughness. The main channel is divided into N parts, each having a wetted perimeter Pi and a roughness coefficient n _i :

where n _c = composite roughness coefficient; P = main channel wetted perimeter; P _i = wetted perimeter of section I, and n _i = roughness coefficient per section.

Research using this software usually evaluates how well a hydraulic model predicts floods, in order to identify vulnerable areas, critical infrastructure and the affected land use value (^{Sarhadi et al., 2012}; ^{Zazo et al., 2015}). Studies have demonstrated that HEC-RAS modelling enables different scenarios to be evaluated for forecasting areas of flooding (^{Guida et al., 2015}). This software is also used for optimising the geometric visualization of areas prone to flooding, which can be subsequently visualised using a geographic information system (GIS) ^{(Sarhadi et al., 2012}). ^{Mohammadi et al. (2014}) simulated flooding and hydraulic conditions in flood areas for different return periods (recurrence intervals) using HEC-RAS, HEC-GEORAS and GIS models in a case study, presenting their results as risk analysis and flood damage.

Artificial neural network (ANN) model

A standard ANN structure (Figure 1) consists of a set of neurons organised into a hierarchy of layers (input, hidden and output) constituting an autonomous functional system (^{Chen, Chen, Chou, & Yang, 2010}). The following elements can be identified with this type of intelligent system: input and output variables and synaptic weights (the intensity of interaction between neurons and propagation, activation and output functions) (^{Komatsu et al., 2014}). The amount of layers and neurons represents one of the most important parameters in ANN since that determines a system’s efficiency.

Figure 1 Backpropagation in an artificial neural network (ANN) structure (^{Chen et al., 2010}). 

One of ANNs’ advantages is that they are useful tools for modelling when the ratio of entry data to output data is unknown (which is why this type of model is called a black box) (^{Chau, Wu, & Li, 2005}; ^{Wang, Wang, Lei, Jiang, & Song, 2011}), enabling complex systems to be modelled based on their mathematical composition (i.e. hydrological processes) (^{Dawson et al., 2006}). Another benefit is that ANNs can produce output from a specific combination of inputs, and their response capability concerning managing non-linear data (^{Santillán, Fraile-Ardanuy, & Toledo, 2014}; ^{Cervantes-Osornio, Arteaga-Ramírez, Vázquez-Peña, Ojeda-Bustamante, & Quevedo-Nolasco, 2013}).

Available data

The Bogota River basin is located in Colombia’s Cundinamarca department. It has a 5 891 km² surface area, representing around 32% of the department’s total surface area. The Bogota River is the basin’s main river; it runs for a total of 308 kilometres from an altitude of 3 300 metres above sea level (masl) in the municipality of Villapinzón to its outlet into the Magdalena River at 280 masl in the municipality of Girardot (^{CAR, 2006}).

The Bogota River is divided into 3 sub-basins: upper, middle and lower. The stretch being studied was in the upper basin (Figure 2) between the Puente Florencia hydrological station in the municipality of Gachancipá (upstream) and the Tocancipá station in the municipality of Tocancipá (downstream). This stretch is characterised by having hourly flow frequency records available and stability regarding such records (i.e. there is no reservoir or other large water body nearby that significantly alters the basin’s hydrological pattern).

Figure 2 The study segment (authors, with cartography supplied by the CAR and IGAC). Regional Autonomous Corporation (Corporación Autónoma Regional - CAR), Agustín Codazzi Geographical Institute (Instituto Geográfico Agustín Codazzi, IGAC). 

Twelve small sub-basins were identified as contributing flows to the section being studied, which added together, considerably increase the river’s level during periods of heavy rainfall. The rational formula method (estimating peak runoff rate at a specific location as a function of drainage network area, runoff coefficient, and mean rainfall intensity) was used for calculating each sub-basin’s contribution to flow, using the water level curves corresponding to each sub-basin and flow order (^{Horton, 1945}).

Physical model produced using HEC-RAS software

The physical model implemented in HEC-RAS required establishing the calibration parameters, with which the flows were simulated at the output of the model. The parameters used for simulating the model’s output flows were as follows: hydrographs from the Puente Florencia and Tocancipá stations, flow curves or calibration curves, cross- sections and Manning roughness coefficient.

Data from April, May, October and November 2011 and 2013 were taken for simulating the (branch-network) flow model as flooding occurred during these dates due to heavy rainfall.

Concerning the flow curves, upstream at the Puente Florencia station, a maximum of roughly 60 m³/s can occur at a 5-meter water level height. And downstream at the Tocancipá hydrological station a maximum of roughly 50 m³/s at a 5-meter water level height can occur.

Regarding the model’s cross-sections, bathymetry was used for taking measurements at points in the field; 151 sections from the study segment were used, with distances varying from 100 to 800 meters in length, depending on the shape of the channel (i.e. measurements were made at shorter distances in areas having very tight curves).

Manning’s roughness coefficient was calibrated based on constant friction with the surface and the surface with the least friction on the sides of the channel.

Calibrating the one-dimensional hydraulic model simulated in HEC-RAS began by identifying a simple hydrograph showing a wave without distortions, during a period of time in which average flows occurred in the study segment, using Manning coefficients ranging from 0.021 - 0.04 for all cross-sections (^{Santos, Cubillos, & Vargas, 2008}; ^{Cook & Merwade, 2009}). With the aforementioned characteristics in mind, the period from the 12^th to the 23^rd of July 2010 was chosen, based on which the Puente Florencia station’s hydrographs and the calculated flows from the 12 sub-basins were entered in the (HEC-RAS) database.

Three scenarios (January, April to June, and October) were simulated for 2011 and 2013 after calibrating Manning coefficients. The parameter concerning the last cross-section (i.e. the Tocancipá station or the model’s output) was configured with a normal 0.0001 depth value, this being suitable for situations where flow approaches uniform rate (^{US Army Corps Engineers & Hydrologic Engineering Center, 2006}).

Artificial neural network (ANN) model

MATLAB Neural Network Toolbox (2013) software was used for ANN simulation, with flow data from the Puente Florencia station and from the 12 sub-basins along the study segment as input, and flows at the Tocancipá station as the model’s output.

The data was normalised between -1 to 1 (^{Matworks, 2013}), giving an input entry for the model using the Puente Florencia station’s flows and those from the 12 sub-basins along the study segment. The model’s output was a vector from the ANN flow data calculated for the Tocancipá station.

Accurately training the ANN and its forecasts involved dividing the data into two parts. Data from September 2009 to December 2012, including February, March, July, August and September 2013, were used for training (70%), whereas data (the remaining 30%) from January, April, May, June and October 2011 and 2013 were used for forecasting.

MATLAB’S Neural Network Toolbox was configured using backpropagation training (^{Kia et al., 2011}; ^{Chen et al., 2010}). The Levenberg-Marquardt backpropagation algorithm was used for the learning function (trainlm), this being the fastest algorithm for this type of training with large networks. It adjustment function performs better for recognising a target system’s patterns (^{Matworks, 2013}). The following parameters were configured in the Toolbox in order to run the model: a maximum 2 000 iterations (repetitions), 1e-05 minimum gradient and a maximum of 6 validation reviews for evaluating the model’s quality.

A multilayer structure was used for training every scenario (^{Kia et al., 2011}; ^{Siou, Johannet, Borrell, & Pistre, 2011}), modifying the amount of layers (2 to 20) and neurons (2 to 50). Altogether, 168 scenarios were trained as input layers and divided into hidden layers according to their propagation function: 85 had a sigmoid-sigmoid configuration and 83 a sigmoid-linear configuration. A forecast was simulated for every scenario and a MATLAB programme used the results for ascertaining the model’s efficiency.

The models’ statistical evaluation criteria

After the simulations had been made, resulting in the Tocancipá station’s output hydrographs for each period, these were compared to real data for the same periods of time. The following six statistical indexes were used for data analysis, which have been used in most articles consulted and as a method to evaluate the performance of simulation models) (^{Dawson, Abrahart, & See, 2007}): mean absolute error (MAE) (^{Singhal & Swarup, 2011}), mean squared error (MSE) (^{Gomes & Ludermir, 2013}), mean absolute percentage error (MAPE) (^{Lewis, 1982}), root-mean-squared error (RMSE) (^{Singhal & Swarup, 2011}), Pearson’s correlation coefficient (CC) (^{Lin, Hedayat, Bikas, & Yang, 2002}), and concordance correlation coefficient (ρc) (^{Lin, 2011}).

The results of this research were also compared to HydroTest Statistical Assessment of Hydrological Forecasts, which evaluated 20 statistical measures reported by hydrological modelling studies (^{Dawson et al., 2007}). Four HydroTest metrics were used for evaluating real data and modelled data (HEC-RAS, ANN sig-lin and ANN sig-sig) (i.e. 30% of the data selected for validation).

Comparing the physical model to artificial intelligence models

A literature search revealed investigations that compared mathematical models (such as linear regression or multiple regression) to intelligent artificial systems, concluding that intelligent systems had a higher real data to simulated data ratio (^{Aqil et al., 2007}; ^{Firat & Güngör, 2007}; ^{Kisi, Shiri, & Nikoofar, 2012}; ^{Karimi, Kisi, Shiri, & Makarynskyy, 2013}), thereby providing a better forecast than mathematical models. However, the search did not find evidence of hydraulics or traditional hydrology studies that made a comparison with a physical model, which is why this comparison was made. It was found that traditional statistical criteria (CC, MAE, MAPE, RSME) were used in such research but none involved analysis using the concordance correlation coefficient (ρc), which indicates the relationship between a model’s precision and its accuracy (^{Firat & Gürgör, 2007}).

The results were used for comparing the HEC-RAS simulated physical model to the best two ANN MATLAB models. Table 4 gives the results for the three best models. The data suggests that the models had very similar forecasts; the sigmoid-sigmoid ANN, the sigmoid-linear ANN and HEC-RAS models are shown in order of effectiveness.

Figure 3 compares the best three models used in the research (HEC-RAS, sigmoid-sigmoid ANN and sigmoid-linear ANN) and the real data to simulated data ratio using a reference line. It should be noted that 30% of the total flow data was used to validate the models in this simulation.

Figure 3 Real to simulated data ratio for each model. 

For the Tocancipá station, real output flow data compared to simulated data was found above and below the reference line in this figure. If data were above it then the flow forecast would be overestimated and if the data were below the flow forecast it would have been underestimated.

It can be seen that the HEC-RAS model shown in Figure 3a underestimated output flow, since most of the data were below the reference line, meaning that it could not predict flooding levels. Whereas data predicted by the ANN models (sig-sig Figure 3b and sig-lin 3c) had a more homogeneous dispersion, with data above and below the reference line, indicating that it would have a greater possibility of predicting flooding levels corresponding to high flows at the model’s output.

Table 4 shows that the sigmoid-sigmoid configuration resulted in the best ANN model. Hydrographs were then drawn after selecting the best intelligent model in order to compare real flow data to simulated flow data for a period of heavy rainfall, as in April 2011.

As can be seen, the simulated values in the HEC-RAS model’s hydrograph shown in Figure 4a were found to be lower than the real ones, meaning that the physical model did not properly predict the real flows that occurred during that period. The model would thus not be reliable for predicting future flooding events.

Figure 4b shows that the selected ANN model better predicted real flows during the same period as that shown in Figure 4a; however, it fluctuated around the real data (as seen in the scatter plot / hydrograph in Figure 3b).

Figure 4 Comparing real hydrographs to simulated hydrographs from the best two models. 

Validating the model

Table 5 compares the statistics calculated with MATLAB® to those calculated with HydroTest, along with three additional statistical criteria: the R-squared coefficient of determination (RSQR) (^{Pearson, 1896}), Willmott’s index of agreement (IoAd) and the coefficient of efficiency (CE) (^{Ablan, Marquez, Rivas, Molina, & Querales, 2011}).

The CC, MAE and RMSE values were very close (varying by tenths or hundredths), indicating that they were correctly found in the analysis performed by this research. The results of RSQR ranged from 0.801 to 0.847, indicating satisfactory models, since this was close to 1.0 (^{Pearson, 1896}); it should be noted that the ANN model’s sigmoid-sigmoid configuration was very close to being a good forecast model. Willmott’s index of agreement (IoAd) showed that the results for the best two simulated models were good, given values over 0.9, and generally very similar, with a range of 0.93 to 0.95. The coefficient of efficiency (CE) revealed a large difference between both prediction models, rejecting the HEC-RAS physical model (0.7837), which resulted in a value less than 0.8. The other two artificial intelligence models were found to be satisfactory, with values ranging from 0.80 to 0.84 (the latter corresponding to the sigmoid-sigmoid ANN model) (^{Dawson et al., 2007}).