In the literature about the parameter estimation of the nonlinear Muskingum (NL-MUSK) model, benchmark hydrographs have been subjected to various metaheuristics, and in these studies the minor improvements of the algorithms on objective functions are imposed as ‘state-of-the-art’. With the metaheuristics involving more control variables, the attempt to search global results in a restricted solution space is not actually practical. Although metaheuristics provide reasonable results compared with many derivative methods, they cannot guarantee the same global solution when they run under different initial conditions. In this study, one of the most practical of metaheuristics, the particle swarm optimization (PSO) algorithm, was chosen, and the aim was to develop its local search capability. In this context, the hybrid use of the PSO with the Levenberg–Marquardt (LM) algorithm was considered. It was detected that the hybrid PSO–LM gave stable global solutions as a result of each random experiment in the application for four different flood data. The PSO–LM, which stands out with its stable aspect, also achieved rapid convergence compared with the PSO and another hybrid variant called mutated PSO.

  • The routing process of different floods was done with the nonlinear Muskingum model.

  • The parameters of the model were calibrated by variants of PSO.

  • A new variant of PSO was obtained by integrating the Levenberg–Marquardt algorithm.

  • The new type of optimization algorithm showed better results.

  • It was seen that the hybrid model converged to optimum results faster.

Graphical Abstract

Graphical Abstract
Graphical Abstract

The determination of hydrograph characteristics is an essential issue for water resources planners, especially in making various flood protection studies and designing hydraulic structures (Kar et al. 2014, 2015). If the flows lead to flood events, the design procedure is implemented as flood routing, and it can be established as the change of the flood wave moving in the river bed, channel, or storage reservoirs (Gill 1978; Karahan et al. 2013). Flood routing is such a procedure, that is used to predict the variation in flood wave features including magnitude, flow rate, and hydrograph shape. In this concept, routing by the lumped system approach is generally called hydrologic routing, and routing by the distributed system is referred to as hydraulic routing (Karahan et al. 2013). The hydraulic methods (e.g., dynamic, kinematic, and diffusion wave models), in which flow routing is computed as a function of space and time throughout the system, are fundamentally based on the use of Saint-Venant equations (Kaya & Ulke 2012). On the other hand, the hydrologic methods are more simple and favorable owing to their lower demand for input and ability to enable the estimation of the outflow discharge as the function of the inflow discharge of the selected river channel (Karahan et al. 2013). For natural rivers and channels, the most frequently used hydrologic method is the Muskingum model developed by McCarthy (1938).

The classical Muskingum model is based upon two main equations that demonstrate the continuity equation (Equation (1)) and storage equation (Equation (2)) as follows (Tung 1985):
(1)
(2)

In the above equations, St, It, and Qt are storage, inflow, and outflow magnitudes, respectively, at time t. a is a storage time constant and χ is a weighting factor which is defined between 0 and 0.5. In Equation (1), ΔSt is assumed as the time rate of change of storage volume (Karahan et al. 2013; Moghaddam et al. 2016).

Since the first developed Muskingum method was grounded on a supposition that the storage volume has a linear relation with the weighted value of inflow and outflow, and that cannot be feasible for all cases, Gill (1978) proposed a three-parameter nonlinear Muskingum (NL-MUSK) model, in which updated storage function was defined as the following equation:
(3)

This model includes an additional power parameter m to enhance the relation between accumulated storage and weighted flow. It is possible to consider the nonlinear model as a linear routing model when m is equal to one. The unknown parameters of the NL-MUSK are a, χ, and m to be calibrated.

In the literature, to make the parameter estimation of the NL-MUSK, the numerous efforts based on nonmetaheuristic techniques can be cited. With a historical review, the segmented least-squares method (S-LSM) (Gill 1978), the Hooke–Jeeves pattern search with the Davidon–Fletcher–Powell algorithm (HJ + DFP) (Tung 1985), the nonlinear least-squares method (NLLSM) (Yoon & Padmanabhan 1993), the Lagrange-multiplier method (LAMM) (Das 2004), Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm (Geem 2006), and the Nelder–Mead simplex (NMS) algorithm (Barati 2011) have been used for the parameter estimation step. However, since none of the techniques is approved to reach an optimal solution of NL-MUSK globally, the recent studies have been generally focused on the population-based metaheuristic ones. The metaheuristic algorithms such as genetic algorithm (GA) (Mohan 1997), harmony search (HS) (Kim et al. 2001), particle swarm optimization (PSO) (Chu & Chang 2009; Moghaddam et al. 2016), immune clonal selection algorithm (ICSA) (Luo & Xie 2010), parameter-setting-free harmony search (PSF-HS) algorithm (Geem 2011), differential evolution (DE) algorithm (Xu et al. 2011), shuffled frog leaping algorithm (SFLA) (Orouji et al. 2013), honey bee mating optimization (HBMO) algorithm (Niazkar & Afzali 2015), and backtracking search algorithm (BSA) (Yuan et al. 2016) have been used so as to find optimal parameters defined in the NL-MUSK. In addition to the above-cited methods, Kirdemir & Okkan (2019) performed five modern metaheuristic algorithms, namely vortex search algorithm (Dogan & Olmez 2015), gases Brownian motion optimization (Abdechiri et al. 2013), water cycle algorithm (Eskandar et al. 2012), flower pollination algorithm (Yang 2012), and colliding bodies optimization (Kaveh & Mahdavi 2014). They compared the performance of these methods in the parameter estimation of NL-MUSK, and it was concluded that all of the methods yielded close objective function values but varied convergence speed. Moreover, it is also possible to encounter several NL-MUSK modeling studies in which the hybrid optimization techniques are employed. For instance, Haddad et al. (2015) used two different hybrid algorithms in order to optimize the same model. They integrated SFLA with NMS and GA with the generalized reduced gradient method. Ouyang et al. (2014) also presented a hybrid scheme to estimate the NL-MUSK parameters by employing PSO hybridized with NMS. In another hybrid approach developed by Karahan et al. (2013), the BFGS was used as a local search algorithm to strengthen the HS for NL-MUSK modeling.

According to all the literature given above, the NL-MUSK has been handled with many algorithms, and many investigations have been completed through some flood examples evaluated as benchmark data (for example, Wilson data). However, it is also clear that usage of metaheuristic methods in NL-MUSK modeling do not make any significant difference to the objective function which is generally assigned as the sum of squared error (SSE) function. In our opinion, examining a three-parameter flood-routing model with intensive algorithms in terms of complex and intrinsic control variables is both difficult and unfavorable, because they generally make the modeling process cumbersome or impractical (e.g., Haddad et al. 2015; Niazkar & Afzali 2015). Moreover, it is essential to measure the convergence capability of algorithms instead of praising them with minor improvements in SSE. Unfortunately, many studies which focus on NL-MUSK modeling, do not handle the issues simultaneously referring to algorithm parsimony, rapid convergence, and robustness except for a few (e.g., Ouyang et al. 2014). At this stage, it will be more remarkable for engineering practice to make a decision and also to equip an available practical algorithm with less mathematical effort (Tigkas et al. 2016). Therefore, in this study, we only focus on the PSO algorithm, which is one of the most practical methods among the other metaheuristics due to its fewer control parameters. Besides the PSO's inherent capability, the Levenberg–Marquardt (LM) algorithm is also employed to improve the local searching abilities of PSO (Okkan & Kirdemir 2020). According to our recent literature review of NL-MUSK modeling, there is no application in the context of the hybrid approach based on the combined use of PSO and LM (termed as PSO–LM).

The remainder of the presented study is arranged as follows. We give the background of NL-MUSK in the second section. Then, we explain the employed datasets in the third section. The methods involving the fundamental principles of the PSO, hybrid PSO–LM, and the convergence examination exposed to algorithms are detailed in the fourth section. The fifth section consists of results and discussion. The final section provides a brief conclusion to the article.

After assigning candidate values for the parameters of NL-MUSK, the initial storage is calculated (step 1 expressed in Figure 1). Initially, the outflow can be assumed to be equal to the inflow (Karahan et al. 2013). In steps after t = 1, the time rate of change of storage volume ΔSt is computed. Then, the amount of storage to be used in the next step is estimated (step 2 expressed in Figure 1). Finally, the next outflow (Qt) is calculated by using step 3 expressed in Figure 1 that is made from Equation (3), in which the unit time step is used.

Figure 1

The elementary computing steps of NL-MUSK.

Figure 1

The elementary computing steps of NL-MUSK.

Close modal

Once the above operations denoted in Figure 1 are repeated for all time, the outflow values produced by NL-MUSK should be compared with measured outflows. In this study, the objective (fitness) function to be minimized within the employed algorithms is the sum of the squared errors between Qt and observed outflows Qobs,t.

In some cases, infeasible parameter values of the NL-MUSK may result in negative outflows and storage values. To deal with this problem, a penalty function (Pen) is integrated to the routing process by referencing Karahan et al. (2013, 2015). According to this approach, if the negative outflow or storage values are produced, their absolute values are multiplied by a relatively large value of Pen. In this study, the use of a value of 1,000 for Pen as a result of various trials was found to be consistent.

In the study, four different datasets were used as numerical examples in the parameter calibration of NL-MUSK. Among them, the data compiled from Brutsaert (2005) and Viessman & Lewis (2003) are both hypothetical and conceivably do not relate to any real-life measurements. These data can be accessed from the studies performed by Yuan et al. (2016) and Bagatur & Onen (2018), respectively. Besides, the real data obtained from the Sutculer flood event in Turkey, which occurred on 4 November 1995, were also exposed to the modeling strategy in order to understand the response of the proposed approach in different character data. In the Sutculer flood event, 9–88 streamflow gauging station data are used as inflow hydrograph, while 9–89 station data are used as outflow hydrograph. For the locations of these stations, readers can view the study performed by Kaya & Ulke (2012), while the related data can be accessed from Lee et al. (2018). Another example is the data obtained from the 1960 flood event at the River Wye in the UK. As a nearly 70 km branch of the river has minimal lateral inflow contribution, it is considered ideal data for the flood-routing exercises (Karahan et al. 2013). The inflow and outflow hydrographs of the examples are given in Figure 2. Here, the first and the last data have a typical hydrograph shape with a single peak, while the others have relatively fluctuated waves including two major peaks. Also, the existing determination coefficient between the inflow and outflow hydrographs (R2inflow–outflow), which were specified in the right column in Figure 2, provides an idea of the degree of nonlinearity between these hydrographs. For instance, it can be understood that the relationships between inflow and outflow rates are closer to linear for Sutculer, while they are more complex for other cases.

Figure 2

Flood hydrographs used in the study (left column) and the scattering between the inflow and outflow data (right column).

Figure 2

Flood hydrographs used in the study (left column) and the scattering between the inflow and outflow data (right column).

Close modal

In this section, we first briefly explain the PSO algorithm. Next, the hybrid algorithm termed as PSO–LM is introduced. The last subheading of this section describes a relatively new measure of the convergence examination called the average convergence rate.

The main principle of the PSO

The PSO is a population-based metaheuristic algorithm brought forward by Kennedy & Eberhart (1995). It is typically based on the social behavior of flocks of birds. The algorithm was frequently practiced in the modeling and optimization stages of hydro-environmental problems such as reservoir operation (Nagesh & Janga 2007; Reddy & Kumar 2007), water quality modeling (Chau 2005; Afshar et al. 2011), analysis of pollutant levels (Lu et al. 2002), water resources planning topics (Shourian et al. 2008), and so on.

The PSO procedure starts with a random distribution of a flock of birds into the food zone. A matrix composed of Np particles and d parameters can be written in matrix form as follows:
(4)
In the above matrix, the particle in row i is expressed as follows:
(5)
Bird swarms flocks (particles) try to sense how far they are from the food resources they are looking for. For this purpose, they calculate the fitness function values for the coordinates in which they are located. Since the particles do not already know the locations of the foods, they interact with two different terms and act iteratively in space. The first of these terms, pbest, represents the best location information that any particle has ever reached, while the second term, gbest, denotes the best position obtained from all the particles in the population. In the algorithm, there is also a vector v which gives the iterative velocity variation of the particles (Equation (6)). Then, with Equation (7), the particles are allowed to move to their new position in space according to the calculated velocities (Kang & Zhang 2016).
(6)
(7)
Here, W (t) is the inertia weight, which is linearly decreasing coefficient based on the iteration number t. c1 and c2 are positive acceleration coefficients controlling pbest and gbest components of Equation (6), while r1 and r2 refer to random numbers that are uniformly distributed between 0 and 1. Besides, cfac is the constriction factor to be used to control the changes in velocity and enhance the convergence capability of the PSO, as shown in Equation (8). The above operations are repeated between t = 0 and the maximum number of iterations (itermax).
(8)
Considering the recommendation of Eberhart & Shi (2000), θ was set to 4.1 (c1 = c2 = 2.05) and the constant multiplier cfac is thus 0.7298.
Shi & Eberhart (1998) have integrated the variable inertia weight given in Equation (9) into Equation (6). In the study, Wmin and Wmax were taken as 0.001 and 1.0, respectively. The position and velocity update operations in the PSO are shown schematically in Figure 3.
(9)
Figure 3

Schematic representation of the standard PSO procedure (left corner), the Jacobian matrix calculation required for LM (right corner), and the PSO–LM hybrid algorithm procedure (bottom).

Figure 3

Schematic representation of the standard PSO procedure (left corner), the Jacobian matrix calculation required for LM (right corner), and the PSO–LM hybrid algorithm procedure (bottom).

Close modal

Hybrid PSO–LM algorithm

In recent years, the researchers interested in the parameter estimation of NL-MUSK have focused on the hybrid optimization techniques that combine the global search abilities of metaheuristic algorithms with local fine-tuning properties of the nonmetaheuristic ones such as BFGS (Karahan et al. 2013) and NMS (Ouyang et al. 2014). Even if all these efforts increase the possibility of approaching the global solution, the combinations intensified the algorithmic structure in terms of the control variables. So how logical is it that the complex algorithms are made more effective by using local search additions? As a matter of fact, considering that almost all metaheuristic themed algorithms are designed with inspiration from GA and/or PSO, it is argued that concentrating on one of the simplest structured algorithms, namely the PSO, is more consistent. Although PSO is considered to be a popular and fast adaptable method to any problem, its local search ability is relatively weak as probably observed in other metaheuristics, and what's more, it also does not guarantee the same global solution for each run (Zhang et al. 2007; Okkan & Kirdemir 2020).

On the other side, Newtonian approaches, which are considered as local search methods, are hypersensitive to the initial solution. Therefore, their integration into a population-based algorithm can be more effective (Okkan & Kirdemir 2020). In the case of using the algorithms like BFGS, the inversion of the Hessian matrix is an obstacle during the numerical solution. This issue can be tackled by the LM algorithm, which uses an approximate solution of the Hessian matrix (Okkan & Kirdemir 2018). The LM algorithm is operated with the Jacobian matrix (J) which is composed of first-order partial derivatives according to the related parameters. After the finite difference approach (e.g., forward type) is evaluated for calculation of J matrix, in each iteration step (k = 1,2,..,itermax) of LM, the parameter updating process is performed by the following equation:
(10)

Here, e represents the error vector calculated between observations (target values) and model predictions, while λ denotes the Marquardt parameter that is adaptively updated.

In the LM, λ is multiplied by a certain decay rate β when SSE function decreases in a certain iteration step, and it is divided by β when SSE increases in another step (Okkan & Serbes 2012). In this study, as suggested by Hagan & Menhaj (1994), values of 0.01 and 0.1 were assigned for λ0 (initial λ) and β, respectively.

In the study, the gbest component of the PSO was exposed to the LM algorithm to improve the local search behavior of the process. In fact, a similar procedure was previously implemented by Zhang et al. (2007) in the training of feed-forward neural networks (FFNWs). In another study performed by Okkan & Kirdemir (2020), PSO (without constriction factor) was coupled with LM for the robust calibration of conceptual rainfall–runoff models. Zhang et al. (2007) have stated that PSO converges quickly during the initial stages of a global search, but around global optimum, the search process slows down. In this concept, a hybrid algorithm combining PSO with a gradient back-propagation algorithm has been shown to provide convergent speed and convergent accuracy in the training step of FFNW. That's why we maintain in our study that this parallel operation will directly contribute to obtaining a solution with the less objective function call and increase the convergence performance of NL-MUSK calibration. With this hybrid approach, the advantageous aspects of the two methods (only PSO and only LM) are revealed, while their shortcomings are disabled (Okkan & Kirdemir 2020). The conceptual demonstration of the hybrid approach, also referred to as PSO–LM algorithm, is given in Figure 3.

Convergence examination

In the population-based evolutionary algorithms, it is also important to make out how rapidly they converge to the desired optimal point along with the number of function calls (Okkan & Kirdemir 2020). In this study, a geometric convergence rate (CR), which was previously tested by He & Lin (2016), was used (Equation (11)).
(11)
where fitnessdesired is equal to zero for the targeted SSE value, fitness0 is the best fitness value among those initially produced, fitnessh is the stable fitness value that is assumed to not be significantly changed after the hth step.
When a threshold value ɛ is taken, the required step h can be determined by using the following equations:
(12a)
(12b)

To examine the credibility of the suggested hybrid PSO–LM technique and to compare it with PSO, four flood-routing examples, which were mentioned in the data section of this paper, were considered. Furthermore, it was provided to compare PSO–LM with another hybrid type, the elitist-mutated version of PSO (PSO-Mutated) proposed by Kang & Zhang (2016) in order to present more detailed analyses to readers with two different hybrid schemes. All the algorithms were coded in the MATLAB environment, in which the SSE function was minimized. The parameter ranges of NL-MUSK were selected as χ= 0.01–0.5, a= 0.01–5.0, and m= 0.01–5.0. Here, a relatively large range was chosen to compel the algorithms in a wide solution space where the penalty constant may be needed to prevent generating negative outflows or storage values. Though such algorithms were run with a huge number of iterations in the literature, a feasible number of generation were exerted in this study similar to Niazkar & Afzali (2015), and itermax was taken as 1,000 for PSO. Since both PSO–LM and PSO-Mutated use extra local tuning operation at each iterative step, the assessment should be done based on the number of function evaluations instead of iteration numbers. So, the maximum generation number was set to 500 for them. Furthermore, there is no straightforward criterion when choosing population size Np. For the sake of practicality, it was taken as 20. This population size was also recommended by Ouyang et al. (2014).

Due to stochastic features of PSO variants, 30 independent runs were operated to check over the robustness of the performances. After operating the hybrid approach depending on the scheme stated in Figure 3, the variations concerning fitness values obtained for each run were extracted throughout the objective function calls. Then, the fitness values obtained in each step during the 30-run are averaged. These results were stored for four different groups of flood data and given in a comparative manner including the PSO and PSO-Mutated variant (see Figure 4). Standard PSO did not produce qualified solutions because the convergence lines of PSO crossed the y-axis higher in almost all data types. On the other hand, it can be clearly seen from Figure 4 that the PSO–LM and PSO-Mutated rapidly converge to the solution during function calls. Even if the responses indicated that the hybrid variants succeeded in comparison to the PSO for all datasets, looking at the performance of the algorithms in the ultimate step would be more descriptive. Thus, the arithmetic mean and the standard deviation statistics of SSE values detected in the last function call for 30 runs were also computed (see Table 1). These descriptive statistics regarding all flood examples indicate that the PSO may get trapped in a local minimum. In contrast to the PSO and PSO-Mutated, the PSO–LM has reached almost the same SSE fitness value in all of the simulations because the standard deviation around the mean SSE was calculated as 9.44E-09, 4.29E-08, 1.41E-06, and 8.12E-09 cms2 for the Brutsaert, Viessman and Lewis, Sutculer, and River Wye data, respectively.

Table 1

Performance scores of the algorithms employed

AlgorithmData typefitnessmean (cms2)fitnessStd (cms2)CRmeanhmean
PSO Brutsaert 16,175.936 20,140.533 0.0255 322 
Viessman and Lewis 76,875.047 14,279.914 0.0068 588 
Sutculer 779.988 1,198.632 0.0065 620 
River Wye 45,044.086 30,444.387 0.0049 570 
PSO-Mutated Brutsaert 12,394.369 9.496 0.0103 596 
Viessman and Lewis 73,436.667 90.668 0.0048 792 
Sutculer 561.052 0.269 0.0176 252 
River Wye 37,955.327 37.460 0.0040 746 
PSO–LM Brutsaert 12,391.819 9.44 × 10−9 0.0565 112 
Viessman and Lewis 73,399.293 4.29 × 10−8 0.0193 190 
Sutculer 560.913 1.41 × 10−6 0.0179 200 
River Wye 37,944.145 8.12 × 10−9 0.0152 198 
AlgorithmData typefitnessmean (cms2)fitnessStd (cms2)CRmeanhmean
PSO Brutsaert 16,175.936 20,140.533 0.0255 322 
Viessman and Lewis 76,875.047 14,279.914 0.0068 588 
Sutculer 779.988 1,198.632 0.0065 620 
River Wye 45,044.086 30,444.387 0.0049 570 
PSO-Mutated Brutsaert 12,394.369 9.496 0.0103 596 
Viessman and Lewis 73,436.667 90.668 0.0048 792 
Sutculer 561.052 0.269 0.0176 252 
River Wye 37,955.327 37.460 0.0040 746 
PSO–LM Brutsaert 12,391.819 9.44 × 10−9 0.0565 112 
Viessman and Lewis 73,399.293 4.29 × 10−8 0.0193 190 
Sutculer 560.913 1.41 × 10−6 0.0179 200 
River Wye 37,944.145 8.12 × 10−9 0.0152 198 

fitnessmean and fitnessstd are the mean and standard deviation statistics of the fitness values obtained in the last iteration of all runs, respectively. CRmean and hmean are the mean convergence rate and mean stable iteration achieved for 30-run average fitness values, respectively.

Figure 4

The change in SSE values through the number of objective function calls used for standard PSO, PSO-Mutated, and PSO–LM. (Each variant was run 30 times; 30-run average values are shown.)

Figure 4

The change in SSE values through the number of objective function calls used for standard PSO, PSO-Mutated, and PSO–LM. (Each variant was run 30 times; 30-run average values are shown.)

Close modal

In addition to rendition of the fitness performances, the convergence of the algorithms was also investigated. An algorithm's CR value is naturally expected to be close to one. For mathematical benchmark functions with a global minimum of zero (e.g., Rastrigin, Griewank functions), this assessment is definitely much more coherent. However, in the case of evaluating flood-routing or other hydrological models, derived CR values may seem relatively fewer. Since the expected value for fitness is zero (we had to choose zero because we do not know the global minimum solution), CR values can be small or large depending on the data characteristics. As seen in Table 1, the complexity of the relationship between inflow and outflow data, and the fluctuation in the series may have caused the variability on CRs. To address this issue, empirical CRmean relations were obtained by taking into account both the number of major peaks in the hydrographs and the existing determination coefficient between the inflow and outflow hydrographs (R2inflow–outflow), which were given in the right column in Figure 2. While Brutsaert and River Wye data have one major peak, in other cases this number can be considered as two. The regression relations and the efficiency of empirical estimates from these relations have been denoted in Figure 5. It can be derived from the evinced scattering that consistent relationships are achieved for both PSO and PSO–LM. It is understood from the relations given in Figure 5 that the convergence rate of both algorithms decreases depending on the number of major peaks, while it increases due to R2inflow–outflow. Of course, it may not be reliable to use a limited number of data to generalize such a relationship. As we focused only on a hydrological method, extra physical parameters may be required to explain those uncertainties. So, in our concept, it will be more appropriate to inspect the CRs relatively. In other words, it will be consistent to query the calibration capabilities of the algorithms, not the data case. In this regard, the PSO–LM gave relatively qualified results. Approximately 2–3 times faster convergence was achieved in all data types compared with PSO (see Table 1). In this determination, combining PSO with LM based on the usage of Jacobian matrix has a big role. It has been also proved that by using PSO–LM appropriate results can be obtained with two to three times less function calls compared with the PSO. When it is compared with the other NL-MUSK studies in point of iteration settings (e.g., Ouyang et al. 2014; Yuan et al. 2016; Lee et al. 2018), the proposed PSO–LM is absolutely applicable. However, it should be noted that the selection of the number of function calls (or stopping criteria) is problem specific.

Figure 5

Empirical CR estimates based on the number of major peaks and the determination coefficient between the observed inflow and the observed outflow.

Figure 5

Empirical CR estimates based on the number of major peaks and the determination coefficient between the observed inflow and the observed outflow.

Close modal

The calibrated NL-MUSK model was also questioned with the additional measures after PSO–LM proved its success in both precise calibration of parameters and convergence. When referencing Niazkar & Afzali (2015), the DPO criterion, which is expressed as the absolute value of the deviation in peak of observed and routed outflows, can be used for the magnitude of peak consideration. Moreover, the closeness of form and size of the hydrograph can be examined assessing the criterion of the variance explained (Moghaddam et al. 2016). In this context, the coefficient NS, which is a modification of the determination coefficient R2, was evaluated. The formula of NS coefficient can be found in Nash & Sutcliffe (1970) and Moriasi et al. (2007). The ultimate parameter estimations procured from PSO–LM and the measures which were calculated are given in Table 2. On account of the fact that the performance limits NS > 0.75 and RSR < 0.5 have been met in accordance with Moriasi et al. (2007), the computed values obtained for all flood events point out a high-performance modeling (here, RSR standardizes root-mean-square error (RMSE) using the standard deviation of outflow observations). As shown in Figure 6, it can be concluded again that both observed and predicted outflows are soundly matched for all data. Additionally, the DPO criteria, which are evaluated in the peak value estimation of the outflow hydrograph, also appear to be relatively small compared with the total flow volume for all flood examples.

Table 2

The NL-MUSK parameters estimated by PSO–LM and the related statistical performance measures

Data typeCalibrated parameters
Statistical Performances
χamSSE (cms2)RMSE (cms)RSRNSR2DPO (cms)
Brutsaert 0.12533 1.14046 1.07666 12,391.819 19.6785 0.0302 0.9991 0.9993 49.969 
Viessman and Lewis 0.16719 0.07614 1.44595 73,399.293 55.3019 0.1272 0.9831 0.9832 48.87 
Sutculer 0.010 0.99802 1.01024 560.913 4.3240 0.0936 0.9909 0.9954 7.501 
River Wye 0.40924 0.07924 1.58148 37,944.145 33.4066 0.1492 0.9771 0.982 97.808 
Data typeCalibrated parameters
Statistical Performances
χamSSE (cms2)RMSE (cms)RSRNSR2DPO (cms)
Brutsaert 0.12533 1.14046 1.07666 12,391.819 19.6785 0.0302 0.9991 0.9993 49.969 
Viessman and Lewis 0.16719 0.07614 1.44595 73,399.293 55.3019 0.1272 0.9831 0.9832 48.87 
Sutculer 0.010 0.99802 1.01024 560.913 4.3240 0.0936 0.9909 0.9954 7.501 
River Wye 0.40924 0.07924 1.58148 37,944.145 33.4066 0.1492 0.9771 0.982 97.808 
Figure 6

Comparison of routed outflow hydrographs obtained by PSO–LM operations for all cases used.

Figure 6

Comparison of routed outflow hydrographs obtained by PSO–LM operations for all cases used.

Close modal

The study also aimed to emphasize which parameter of the NL-MUSK model is more sensitive. As the exponent parameter m represents the nonlinearity of the NL-MUSK model, this parameter is undoubtedly the most sensitive one. It is shown in Figure 7 through the example of Brutsaert data that the m parameter is relatively more sensitive than the other two parameters. This again shows that the calibration of the model is sensitive and needs to be optimized with a qualified algorithm.

Figure 7

Relative sensitivity curves of the NL-MUSK parameters. (As an example, an application is shown only through Brutsaert data.)

Figure 7

Relative sensitivity curves of the NL-MUSK parameters. (As an example, an application is shown only through Brutsaert data.)

Close modal

The parameter estimation is the most important process of the NL-MUSK, which is a popular hydrologic flood-routing method. In the literature about parameter calibration of the NL-MUSK, it is possible to encounter a number of applications varying from the derivative-based algorithms to metaheuristics. In fact, metaheuristic algorithms can be said to be more popular thanks to their global search capabilities and their inherent structures mimicking individuals or phenomena in nature (Kirdemir & Okkan 2019). However, in many cases, with the metaheuristic algorithms that require more control parameters than the number of parameters in the NL-MUSK, the performed applications within a limited solution space have begun to move away from practicality and robustness. Reaching accurate results with less time and less effort is a generally approved principle in engineering, and the developed methods are expected to be in harmony as well. In this regard, Okkan & Kirdemir (2020) have shown that the parsimonious hybridized techniques contain minimal uncertainty in terms of their control parameters and this will therefore be the reason for their preference in such calibration experiments.

In the study presented, to test the suitability of hybrid algorithms in NL-MUSK modeling, four flood data were exposed to different PSO variants. The comparisons between the PSO, PSO-Mutated, and PSO–LM revealed that the PSO–LM satisfied both the stability in each random run and the fast convergence throughout the objective function calls. On the grounds that the proficiency about the statistical criterions including NS, RSR, and DPO was achieved, it was concluded that both observed and predicted outflows were thoroughly matched. Moreover, it is anticipated that the proposed PSO–LM algorithm can be adjusted to other flood-routing models incorporating lateral flow (Karahan et al. 2015; Moghaddam et al. 2016), large-scale water quality models (Afshar et al. 2011), multi-purpose reservoir operation studies (Nagesh & Janga 2007), as well as the other engineering issues subject to continuous optimization.

The authors declare that they have no conflict of interest.

The modeling strategy applied in this study was developed in the scientific research project titled ‘Evaluation of global optimization algorithms in conceptual hydrological model calibration (Grant No. 2017/134)’ managed by the first author. So, he would like to thank ‘Scientific Research Projects Unit of Balikesir University/Turkey’ for their financial support. The authors are also grateful to two anonymous reviewers for providing valuable comments, which we considered to improve the quality of this paper. As mentioned in the Data section, we also thank the researchers who shared the flood data.

Abdechiri
M.
Meybodi
M. R.
Bahrami
H.
2013
Gases Brownian motion optimization: an algorithm for optimization (GBMO)
.
Appl. Soft. Comput. J.
(
13
),
2932
2946
.
doi: 10.1016/j.asoc.2012.03.068
.
Afshar
A.
Kazemi
H.
Saadatpour
M.
2011
Particle swarm optimization for automatic calibration of large scale water quality model (CE-QUAL-W2): application to Karkheh reservoir
.
Water Resour. Manage.
25
(
10
),
2613
2632
.
https://doi.org/10.1007/S11269-011-9829-7
.
Bagatur
T.
Onen
F.
2018
Development of predictive model for flood routing using genetic expression programming
.
J. Flood Risk Manage.
(
11
).
doi: 10.1111/jfr3.12232
.
Barati
R.
2011
Parameter estimation of nonlinear Muskingum models using Nelder-Mead simplex algorithm
.
J. Hydrol. Eng.
16
(
11
),
946
954
.
https://doi.org/10.1061/(asce)he.1943-5584.0000379
.
Brutsaert
W.
2005
Hydrology: An Introduction
.
Cambridge University Press
,
New York
.
Chau
K.
2005
A split-step PSO algorithm in prediction of water quality pollution
. In:
Advances in Neural Networks – ISNN 2005. Lecture Notes in Computer Science
, Vol.
3498
(
Wang
J.
Liao
X. F.
Yi
Z.
, eds.).
Springer
,
Berlin, Heidelberg
.
https://doi.org/10.1007/11427469_164
.
Chu
H. J.
Chang
L. C.
2009
Applying particle swarm optimization to parameter estimation of the nonlinear Muskingum model
.
J. Hydrol. Eng.
14
(
9
),
1024
1027
.
https://doi.org/10.1061/(asce)he.1943-5584.0000070
.
Das
A.
2004
Parameter estimation of Muskingum models
.
J. Irrig. Drain. Eng.
130
(
2
),
140
147
.
https://doi.org/10.1061/(asce)0733-9437(2004)
.
Dogan
B.
Olmez
T.
2015
A new metaheuristic for numerical function optimization: vortex search algorithm
.
Inf. Sci.
(
293
),
125
145
.
doi: 10.1016/j.ins.2014.08.053
.
Eberhart
R.
Shi
Y.
2000
Comparing nonlinear inertia weights and constriction factors in particle swarm optimization
. In:
Proceedings of the 2000 Congress on Evolutionary Computation, CEC00
. pp.
84
88
.
Eskandar
H.
Sadollah
A.
Bahreininejad
A.
Hamdi
M.
2012
Water cycle algorithm – a novel metaheuristic optimization method for solving constrained engineering optimization problems
.
Comput. Struct.
(
110–111
),
151
166
.
doi:10.1016/j.compstruc.2012.07.010
.
Geem
Z. W.
2006
Parameter estimation for the nonlinear Muskingum model using the BFGS technique
.
J. Irrig. Drain. Eng.
132
(
5
),
474
478
.
https://doi.org/10.1061/(asce)0733-9437(2006)132:5(474)
.
Geem
Z. W.
2011
Parameter estimation of the nonlinear Muskingum model using parameter-setting-free harmony search algorithm
.
J. Hydrol. Eng.
16
(
8
),
684
688
.
https://doi.org/10.1061/(asce)he.1943-5584.0000352
.
Gill
M. A.
1978
Flood routing by Muskingum method
.
J. Hydrol.
36
(
3–4
),
353
363
.
https://doi.org/10.1016/0022-1694(78)90153-1
.
Haddad
O. B.
Hamedi
F.
Fallah-Mehdipour
E.
Orouji
H.
2015
Application of a hybrid optimization method in Muskingum parameter estimation
.
J. Irrig. Drain. Eng.
141
(
12
).
https://doi.org/10.1061/(ASCE)IR.1943-4774.0000929
.
Hagan
M. T.
Menhaj
M. B.
1994
Training feed forward techniques with the Marquardt algorithm
.
IEEE Trans. Neural Netw.
5
(
6
),
989
993
.
https://doi.org/10.1109/72.329697
.
He
J.
Lin
G.
2016
Average convergence rate of evolutionary algorithms
.
IEEE Trans. Evol. Comput.
20
(
2
),
316
321
.
https://doi.org/10.1109/tevc.2015.2444793
.
Kar
K. K.
Yang
S. K.
Jung
W. Y.
2014
Determination of surface runoff for recent years rainfall over Hancheon Basin in Jeju Island
. In:
The 23rd Korean Environmental Sciences Society Conference
,
Daegu, South Korea
, pp.
6
8
.
Kar
K. K.
Yang
S. K.
Lee
J. H.
2015
Assessing unit hydrograph parameters and peak runoff responses from storm rainfall events: a case study in Hancheon Basin of Jeju Island
.
J. Environ. Sci. Int.
(
24
),
437
447
.
doi:10.5322/jesi.2015.24.4.437
.
Karahan
H.
Gurarslan
G.
Geem
Z. W.
2013
Parameter estimation of the nonlinear Muskingum flood-routing model using a hybrid harmony search algorithm
.
J. Hydrol. Eng.
18
(
3
),
352
360
.
https://doi.org/10.1061/(ASCE)HE.1943-5584.0000608
.
Karahan
H.
Gurarslan
G.
Geem
Z. W.
2015
A new nonlinear Muskingum flood routing model incorporating lateral flow
.
Eng. Optim.
47
(
6
),
737
749
.
https://doi.org/10.1080/0305215X.2014.918115
.
Kaveh
A.
Mahdavi
V. R.
2014
Colliding bodies optimization : a novel meta-heuristic method
.
Comput. Struct.
(
139
),
18
27
.
doi:10.1016/j.compstruc.2014.04.005
.
Kaya
B.
Ulke
A.
2012
Solution of kinematic wave model using DQM and Sutculer flood example
.
Teknik. Dergi.
23
(
2
),
5869
5884
(in Turkish)
.
Kennedy
J.
Eberhart
R.
1995
Particle swarm optimization
.
Proc. IEEE Int. Conf. Neural Netw.
(
4
),
1942
1948
.
https://doi.org/10.1109/ICNN.1995.488968
.
Kim
J. H.
Geem
Z. W.
Kim
E. S.
2001
Parameter estimation of the nonlinear Muskingum model using harmony search
.
J. Am. Water Resour. Assoc.
37
(
5
),
1131
1138
.
https://doi.org/10.1111/j.1752-1688.2001.tb03627.x
.
Kirdemir
U.
Okkan
U.
2019
General review of calibration process of nonlinear Muskingum model and its optimization by up-to-date methods
. In:
Decision Support Methods for Assessing Flood Risk and Vulnerability
(
Karmaoui
A.
, ed.).
IGI Global
,
Hershey, PA
,
USA
, pp.
61
83
.
Lee
E.
Lee
H.
Kim
J.
2018
Development and application of advanced Muskingum flood routing model considering continuous flow
.
Water
(
10
),
1
21
.
https://doi.org/10.3390/w10060760
.
Luo
J.
Xie
J.
2010
Parameter estimation for the nonlinear Muskingum model based on immune clonal selection algorithm
.
J. Hydrol. Eng.
15
(
10
),
844
851
.
https://doi.org/10.1061/(asce)he.1943-5584.0000244
.
McCarthy
G. T.
1938
The Unit Hydrograph and Flood Routing
.
Army Engineer District
,
Providence, RI
.
Moghaddam
A.
Behmanesh
J.
Farsijani
A.
2016
Parameters estimation for the new four-parameter nonlinear Muskingum model using the particle swarm optimization
.
Water Resour. Manage.
(
30
),
2143
2160
.
https://doi.org/10.1007/S11269-016-1278-X
.
Mohan
S.
1997
Parameter estimation of nonlinear Muskingum models using genetic algorithm
.
J. Hydraul. Eng.
123
(
2
),
137
142
.
https://doi.org/10.1061/(ASCE)0733-9429(1997)123:2(137)
.
Moriasi
D. N.
Arnold
J. G.
Liew
M. W. V.
Bingner
R. L.
Harmel
R. D.
Veith
T. L.
2007
Model evaluation guidelines for systematic quantification of accuracy in watershed simulations
.
Trans. ASABE
50
(
3
),
885
900
.
http://dx.doi.org/10.13031/2013.23153
.
Nagesh
K. D.
Janga
R. M.
2007
Multipurpose reservoir operation using particle swarm optimization
.
J. Water Resour. Plann. Manage.
133
(
2
),
192
201
.
https://doi.org/10.1061/(asce)0733-9496(2007)133:3(192)
.
Nash
J. E.
Sutcliffe
J. V.
1970
River flow forecasting through conceptual models part I—a discussion of principles
.
J. Hydrol.
10
(
3
),
282
290
.
https://doi.org/10.1016/0022-1694(70)90255-6
.
Niazkar
M.
Afzali
S. H.
2015
Assessment of modified honey bee mating optimization for parameter estimation of nonlinear Muskingum models
.
J. Hydrol. Eng.
20
(
4
).
https://doi.org/10.1061/(asce)he.1943-5584.0001028
.
Okkan
U.
Kirdemir
U.
2018
Investigation of the behavior of an agricultural-operated dam reservoir under RCP scenarios of AR5-IPCC
.
Water Resour. Manage.
32
(
8
),
2847
2866
.
https://doi.org/10.1007/s11269-018-1962-0
.
Okkan
U.
Kirdemir
U.
2020
Towards a hybrid algorithm for robust calibration of rainfall-runoff models
.
J. Hydroinf.
.
Okkan
U.
Serbes
Z. A.
2012
Rainfall–runoff modeling using least squares support vector machines
.
Environmetrics
23
(
6
),
549
564
.
https://doi.org/10.1002/env.2154
.
Orouji
H.
Haddad
O. B.
Fallah-Mehdipour
E.
Mariño
M. A.
2013
Estimation of Muskingum parameter by metaheuristic algorithms
.
Proc. Inst. Civil Eng: Water Manage.
166
(
6
),
315
324
.
https://doi.org/10.1680/Wama.11.00068
.
Ouyang
A.
Li
K.
Truong
T. K.
Sallam
A.
Edwin
H. M. S.
2014
Hybrid particle swarm optimization for parameter estimation of Muskingum model
.
Neural Comput. Appl.
(
25
),
1785
1799
.
https://doi.org/10.1007/S00521-014-1669-Y
.
Reddy
M. J.
Kumar
D. N.
2007
Optimal reservoir operation for irrigation of multiple crops using elitist-mutated particle swarm optimization
.
Hydrol. Sci. J.
52
(
4
),
686
701
.
https://doi.org/10.1623/hysj.52.4.686
.
Shi
Y.
Eberhart
R. C.
1998
Parameter selection in particle swarm optimization
. In:
International Conference on Evolutionary Programming EP 1998 Evolutionary Programming VII
. pp.
591
600
.
Shourian
M.
Mousavi
S. J.
Tahershamsi
A.
2008
Basin-wide water resources planning by integrating PSO algorithm and MODSIM
.
Water Resour. Manage.
22
(
10
),
1347
1366
.
https://doi.org/10.1007/s11269-007-9229-1
.
Tigkas
D.
Christelis
V.
Tsakiris
G.
2016
Comparative study of evolutionary algorithms for the automatic calibration of the Medbasin-D conceptual hydrological model
.
Environ. Proc.
3
(
3
),
629
644
.
https://doi.org/10.1007/s40710-016-0147-1
.
Tung
Y. K.
1985
River flood routing by nonlinear Muskingum method
.
J. Hydraul. Eng.
111
(
12
),
1447
1460
.
https://doi.org/10.1061/(asce)0733-9429(1985)111:12(1447)
.
Viessman
W. J.
Lewis
G. L.
2003
Introduction to Hydrology
.
Pearson Education
,
New York
.
Xu
D. M.
Oiu
L.
Chen
S. Y.
2011
Estimation of nonlinear Muskingum model parameter using differential evolution
.
J. Hydrol. Eng.
17
(
2
),
348
353
.
https://doi.org/10.1061/(asce)he.1943-5584.0000432
.
Yang
X.
2012
Flower pollination algorithm for global optimization, in: unconventional computation and natural computation 2012
.
Lect. Notes Comput. Sci.
7445
,
240
249
.
Yoon
J. W.
Padmanabhan
G.
1993
Parameter-estimation of linear and nonlinear Muskingum models
.
J. Water Resour. Plann. Manage.
119
(
5
),
600
610
.
https://doi.org/10.1061/(asce)0733-9496(1993)119:5(600)
.
Yuan
X.
Wu
X.
Tian
H.
Yuan
Y.
Adnan
R. M.
2016
Parameter identification of nonlinear Muskingum model with backtracking search algorithm
.
Water Resour. Manage.
30
(
8
),
2767
2783
.
https://doi:10.1007/s11269-016-1321-y
.
Zhang
J.-R.
Zhang
J.
Lok
T.-M.
Lyu
M. R.
2007
A hybrid particle swarm optimization–back-propagation algorithm for feed forward neural network training
.
Appl. Math. Comput.
185
(
2
),
1026
1037
.
https://doi:10.1016/j.amc.2006.07.025
.