## Abstract

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.

## HIGHLIGHTS

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

## INTRODUCTION

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).

In the above equations, *S _{t}*,

*I*, and

_{t}*Q*are storage, inflow, and outflow magnitudes, respectively, at time

_{t}*t*.

*a*is a storage time constant and

*χ*is a weighting factor which is defined between 0 and 0.5. In Equation (1), Δ

*S*is assumed as the time rate of change of storage volume (Karahan

_{t}*et al.*2013; Moghaddam

*et al.*2016).

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.

## ROUTING PROCESS OF THE NL-MUSK MODEL

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 Δ*S _{t}* 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 (

*Q*) is calculated by using

_{t}*step 3*expressed in Figure 1 that is made from Equation (3), in which the unit time step is used.

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 *Q _{t}* and observed outflows

*Q*

_{obs}_{,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.

## DATA

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 (*R ^{2}_{inflow–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.

## METHODS

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.

*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).

*W*(

*t*) is the inertia weight, which is linearly decreasing coefficient based on the iteration number

*t*.

*c*

_{1}and

*c*

_{2}are positive acceleration coefficients controlling

*pbest*and

*gbest*components of Equation (6), while

*r*

_{1}and

*r*

_{2}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 (

*iter*).Considering the recommendation of Eberhart & Shi (2000),

_{max}*θ*was set to 4.1 (

*c*

_{1}=

*c*

_{2}= 2.05) and the constant multiplier

*cfac*is thus 0.7298.

*W*

_{min}and

*W*

_{max}were taken as 0.001 and 1.0, respectively. The position and velocity update operations in the PSO are shown schematically in Figure 3.

### 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).

*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,..,

*iter*) of LM, the parameter updating process is performed by the following equation:

_{max}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

*CR*), which was previously tested by He & Lin (2016), was used (Equation (11)).where

*fitness*is equal to zero for the targeted

_{desired}*SSE*value,

*fitness*

_{0}is the best fitness value among those initially produced,

*fitness*is the stable fitness value that is assumed to not be significantly changed after the

_{h}*h*th step.

## RESULTS AND DISCUSSION

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 *iter _{max}* 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

*N*. For the sake of practicality, it was taken as 20. This population size was also recommended by Ouyang

_{p}*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 cms^{2} for the Brutsaert, Viessman and Lewis, Sutculer, and River Wye data, respectively.

Algorithm . | Data type . | fitness (cms_{mean}^{2})
. | fitness (cms_{Std}^{2})
. | CR
. _{mean} | h
. _{mean} |
---|---|---|---|---|---|

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 |

Algorithm . | Data type . | fitness (cms_{mean}^{2})
. | fitness (cms_{Std}^{2})
. | CR
. _{mean} | h
. _{mean} |
---|---|---|---|---|---|

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 |

*fitness _{mean}* and

*fitness*are the mean and standard deviation statistics of the fitness values obtained in the last iteration of all runs, respectively.

_{std}*CR*and

_{mean}*h*are the mean convergence rate and mean stable iteration achieved for 30-run average fitness values, respectively.

_{mean}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 *CR _{mean}* 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 (

*R*

^{2}

*), 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*

_{inflow–outflow}*R*

^{2}

*. 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*

_{inflow–outflow}*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.

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 *R ^{2}*, 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.

Data type . | Calibrated parameters . | Statistical Performances . | |||||||
---|---|---|---|---|---|---|---|---|---|

χ
. | a
. | m
. | SSE (cms^{2})
. | RMSE (cms)
. | RSR
. | NS
. | R^{2}
. | DPO (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 type . | Calibrated parameters . | Statistical Performances . | |||||||
---|---|---|---|---|---|---|---|---|---|

χ
. | a
. | m
. | SSE (cms^{2})
. | RMSE (cms)
. | RSR
. | NS
. | R^{2}
. | DPO (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 |

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.

## CONCLUSIONS

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.

## CONFLICT OF INTEREST

The authors declare that they have no conflict of interest.

## ACKNOWLEDGEMENTS

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.