Regionalization of a Rainfall-Runo ﬀ Model: Limitations and Potentials

: Regionalized lumped rainfall-runo ﬀ (RR) models have been widely employed as a means of predicting the streamﬂow of an ungauged watershed because of their simple yet e ﬀ ective simulation strategies. Parameter regionalization techniques relate the parameter values of a model calibrated to the observations of gauged watersheds to the geohydrological characteristics of the watersheds. Thus, the accuracy of regionalized models is dependent on the calibration processes, as well as the structure of the model used and the quality of the measurements. In this study, we have discussed the potentials and limitations of hydrological model parameter regionalization to provide practical guidance for hydrological modeling of ungauged watersheds. This study used a Tank model as an example model and calibrated its parameters to streamﬂow observed at the outlets of 39 gauged watersheds. Multiple regression analysis identiﬁed the statistical relationships between calibrated parameter values and nine watershed characteristics. The newly developed regional models provided acceptable accuracy in predicting streamﬂow, demonstrating the potential of the parameter regionalization method. However, uncertainty associated with parameter calibration processes was found to be large enough to a ﬀ ect the accuracy of regionalization. This study demonstrated the importance of objective function selection of the RR model regionalization.


Introduction
Estimating the continuous long-term runoff of ungauged watersheds is a task frequently required in hydrological analyses and design [1][2][3][4][5]. The RR models have been used as a tool to predict streamflow in gauged and ungauged watersheds; the simple structure of a lumped RR model makes it popular, especially for investigations that focus on the overall responses (rather than detailed internal transport processes) of a watershed [6][7][8][9]. The parameters of a lumped RR model represent the spatially aggregated hydrological characteristics of a watershed, and parameter values tend to vary nonlinearly with the spatial scales of RR modeling. Furthermore, for the same reasons, the parameter values cannot be directly measured in the field; instead, the parameters of an RR model have to be calibrated to observations so that its prediction can achieve reliability. Studies have shown that the parameter values of a lumped RR model are associated with the geohydrological characteristics of watersheds [8,[10][11][12]. The statistical relationships between the parameter values and watershed characteristics that are derived from gauged watersheds have been employed to predict the streamflow of ungauged watersheds, which is called parameter or model regionalization [2,5,13,14].
The regression approach is one of the most widely used regionalization methods [4,5,[15][16][17]. In this approach, the parameters of an RR model are calibrated using observations from gauged watersheds that can then be used to quantify the characteristics of the watersheds. The regression approach statistically links the calibrated parameter values to quantified watershed characteristics. The relationship is then applied to ungauged watersheds that have geohydrological features statistically similar to those of the gauged ones [4,5,13,18,19]. The quality of parameter regionalization is dependent on many factors, including model structure, the quality of observations used for regionalization, hydrological variables of interest, calibration techniques, and regression methods. However, large uncertainty and the consequential limited reliability prevent model regionalization from being widely employed in hydrological analyses and design [16,[20][21][22].
There are many sources that introduce uncertainty in the regionalization process whose effects have been extensively discussed in literature [20][21][22][23][24]. However, the extent of uncertainty remains unclear even though many methods have been proposed to estimate the uncertainty in hydrological modeling [2,21,22,25,26]. Calibration can be a major underlying uncertainty source because calibrated parameter values are directly related to the selected watershed characteristics in a regionalization approach [8,19,20,27,28]. Researchers have argued that parameter values cannot be specifically determined [27,29,30], and significantly different parameter sets may be obtained, depending on the objective functions used in the calibration processes [8,[31][32][33][34][35]. However, the conventional practice of model regionalization does not explicitly consider such aspects.
Tank models are a type of lumped hydrological models, and they have been applied to study a wide range of watersheds owing to their computational and conceptual simplicity [8,13,36,37]. Tank models tend to have more parameters as compared to other parsimonious models because they are designed to represent non-linear responses of a watershed using multiple linear equations [8]. Thus, parameter calibration can be challenging, especially when streamflow observations are limited [36,37]. The models have been widely used in East Asia, including Korea, Japan, and Taiwan, and studies have demonstrated their accuracy and performance in applications for humid and mountainous watersheds [8,9,17,38,39].
Various regional Tank models have been developed for water resources planning and management for ungauged watersheds [4,13,14,17,18,40,41]. However, there is substantial room for improvement of the regional models as methods and techniques for calibration and hydrological analyses are improving along with advanced computing resources. In the previous regionalization practices, for instance, model parameters were calibrated with only a single objective function, such as root mean square error (RMSE), mean squared error (MSE), and Nash and Sutcliffe efficiency (NSE) [42], which are over-sensitive to high values (and large outliers) but relatively unresponsive to low flow [8,43,44]. It is known that Kling-Gupta efficiency (KGE or KGE ) [45,46] can simultaneously take into account the multiple aspects of model evaluation, including correlation, bias, variance, and variability [47]; thus, they can serve as effective alternative objective functions for parameter regionalization.
This study explored ways to improve the accuracy of regionalizing an RR model; we also investigated the effect of objective function selection on regionalization results. A Tank model with three layers was used as a lumped hydrological model to be regionalized in this study, and the model was applied to simulate the streamflow hydrographs of 49 watersheds located in Korea. Two objective functions, RMSE and KGE sqrt , were employed in parameter calibration to investigate the impacts of objective function selection. Regionalization results obtained by using previous methods were compared with the results of the newly proposed method to demonstrate the efficacy of the latter.

Tank Model
There are many different possible configurations of a Tank model depending on the number of storage (tanks or layers) and the locations of outlets. The 3-Tank model is one of the most common forms of a Tank model; it represents the hydrological processes of a watershed with the help of three tanks (or storages) that are vertically connected with outlets located at the bottom and side of each tank ( Figure 1) [8,9,36]. The height of the side-outlets and the size of each tank are calibrated to observed streamflow so that the model can predict the overall watershed responses to rainfall events. In the model, water running out of the side outlets of the first tank (or the top layer), second tank, and third tank (or the bottom layer) represents surface runoff, intermediate runoff, and baseflow, respectively [8,9,13]. Water passing through the bottom outlets means infiltration (from the first tank to the second) or percolation (from the second tank to the third). The amount of runoff or infiltrated/percolated water through the outlets is assumed to be linearly proportional to the storage of a Tank as follows (Equations (1) and (2)): where i is the tank order, j is the side-outlet order, n is the number of tanks, m is the number of side outlets for each tank, q ij is the runoff (mm) from the jth side outlet in the ith tank, a ij is the side-outlet coefficient (dimensionless) for the jth side-outlet in the ith tank, ST i,t is the storage of the ith tank (mm), h ij is the height of the side outlet for the jth side outlet in the ith tank (mm), I i,t is the amount of water infiltrated in the ith tank (mm), and b i is the bottom-outlet coefficient (dimensionless) for the ith tank. The ST i,t+1 in the tanks, the storage for the next time step t + 1, can be expressed as follows (Equations (3) and (4)): where P t is the precipitation at time t (mm), and ET i,t is the actual evapotranspiration of the ith tank at time t (mm). The ET i,t is constrained by the volume of storage in each tank (Equations (5) and (6)). It is calculated by subtracting the evapotranspiration in the upper tanks from the total actual evapotranspiration (ET a,t ): The ET a was estimated using the guidelines of the Food and Agriculture Organization of the United Nations (FAO) as follows ( [48]; Equation (7)): where K c,t is the single crop coefficient (dimensionless), K s,t is the soil water stress coefficient (K s,t ) (dimensionless), and ET 0 is the potential evapotranspiration (mm), which was estimated using the FAO Penman-Monteith (PM) approach with gauged solar radiation, air temperature, wind speed, and relative humidity [48]. The K c,t is determined with the area-weighted averages of crop coefficient values for different land uses [36,48]. The K s was calculated from the simulated total watershed storage (ST t ) and a modified power function as represented in Equation (8) [36,40]: where SECP is a soil evaporation compensation parameter, which is positively correlated to K s . The SECP parameter has been reported to vary from 0.001 to 0.1, depending on watershed soil texture [8,40,49]. We adopted the feasible ranges of the Tank model parameters that were set considering the hydrological characteristics of Korean watersheds [8,36] (Table 1).
The was estimated using the guidelines of the Food and Agriculture Organization of the United Nations (FAO) as follows ( [48]; Equation (7): where , is the single crop coefficient (dimensionless), , is the soil water stress coefficient ( , ) (dimensionless), and is the potential evapotranspiration (mm), which was estimated using the FAO Penman-Monteith (PM) approach with gauged solar radiation, air temperature, wind speed, and relative humidity [48]. The , is determined with the area-weighted averages of crop coefficient values for different land uses [36,48]. The was calculated from the simulated total watershed storage ( ) and a modified power function as represented in Equation (8) [36,40]: where is a soil evaporation compensation parameter, which is positively correlated to . The parameter has been reported to vary from 0.001 to 0.1, depending on watershed soil texture [8,40,49]. We adopted the feasible ranges of the Tank model parameters that were set considering the hydrological characteristics of Korean watersheds [8,36] (Table 1). Table 1. Value ranges of 3-Tank model parameters [8,36].

Parameter
Description Min. Max.
Side-outlet coefficient for the first side outlet in the first tank (dimensionless) 0.08 0.5 Side-outlet coefficient for the second side outlet in the first tank (dimensionless) 0.08 0.5 ℎ Height of side outlet for the first side outlet in the first tank (mm) 5 60 ℎ Height of side outlet for the second side outlet in the first tank (mm) 20 110 Bottom-outlet coefficient for the first tank (dimensionless) 0.1 0.5 Figure 1. Schematics of the 3-Tank model structures [8,9,36].

Regionalization of the Tank Models
The Tank models have been widely employed as a tool to predict the streamflow of ungauged watersheds in Korea, Japan, and Germany [4,13,14,17,18,40,41] (Table 2). The model parameters are estimated from the relationship between parameter values and watershed characteristics, which are developed using observations made in gauged watersheds. The Tank model was originally developed with four layers (4-Tank) [37], and the parameter values of the 4-Tank models were generally derived from selected watershed features in Japan and Germany [4,17]. MRE, MSE, and RMSE refer to mean relative error, mean square error, and root mean square error, respectively. NA refers to information not available. W_Slope (%), Length (km), and Form (-) refer to the watershed slope, flow length, and form factor, respectively. Upland is an area where crops are cultivated under non-ponded conditions, and paddy is an area that can maintain a ponded water depth to grow rice. Forest (%), Upland (%), and Paddy (%) refer to the percentages of the forest, upland, and paddy area, respectively.
Yokoo et al. [17] regionalized the 12 parameters of the 4-Tank model using flow measurements made in 12 watersheds that had drainage areas with a range from 100 to 805 km 2 in Japan; the watershed characteristics were derived from topography, soil type, geology, and land-use. In their study, a multiple linear regression model successfully identified the statistical relationship between the model parameters and the watershed characteristics. Amiri et al. [4] examined whether changes in landscape metrics, including shape index, perimeter-area ratio, patch size, and patch density of land-use, could affect the calibrated values of the 4-Tank model parameters from 30 catchments (53-737 km 2 ) located in Germany. They found that multiple regression models could successfully explain the relationship between calibrated parameter values and a highly varying landscape, and landscape metrics should be included in the regionalization of conceptual rainfall-runoff models.
In Korea, a three-layer Tank model, modified with two side outlets on the top layer, has been widely employed, especially to study small-or medium-sized watersheds where flow travel time is short, and the recession limb of a streamflow hydrograph is steep [9,36]. Various regionalization approaches have been carried out in Korea to estimate ungauged streamflow [13,18,40,41]. Kim and Park [13] developed regression equations to estimate reservoir inflow in ungauged watersheds using only the information of drainage area and land-use composition. Similarly, Huh et al. [40], Kim et al. [41], and An et al. [18] considered the topographical characteristics, such as stream length, slope, and form factor when developing the regression relationship.

Objective Functions
Two objective functions were employed in the parameter calibration to see how the selection of objective functions could affect the regionalization of the Tank model: RMSE and modified KGE (KGE sqrt ) [45,46]. The RMSE (mm), which is relatively responsive to high flow or peak flow, has been widely used for calibrating hydrological models (Equation (9)): where O and S represent the observed and simulated discharge (mm), respectively, n is the number of time steps at a time step i. The KGE was proposed as an alternative performance statistics [45]. The modified version, KGE , was introduced to ensure that the bias and variability ratios are not cross-correlated to take into account the multiple aspects of model evaluation, including correlation, bias, and variability simultaneously [46]: where r is the Pearson correlation coefficient between simulated and observed streamflow (dimensionless), β is the bias ratio (dimensionless), α is the standard deviation ratio (dimensionless), γ is the variability ratio (dimensionless), cov is the covariance between observation and simulation, µ is the mean runoff (mm or cms), and σ is the standard deviation (mm or cms). As the modified KGE tends to be sensitive to large values, a square root transformation (or the Box-Cox power transformation with a lambda value of 0.5) (KGE sqrt ) was applied to the original flow before computing KGE to reduce the degree of expected bias toward high flow [8,47].

The Automatic Parameter Calibration Algorithm
Automatic calibration with selected objective functions was performed by the shuffled complex evolution algorithm (SCE), which is a sampling-based heuristic search strategy [50,51]. Previous studies have proved its applicability to the calibration of hydrological models [8,45,46,52], and the details of SCE are well described in literature [50,53]. The SCE sampling of this study employed 15 complexes and 21 points (or populations) per complex. The sampling continued until the differences between the objective function values sampled in the last 10 points were less than 0.1%.

Model Evaluation Statistics
The calibrated Tank models were evaluated using four performance statistics commonly used in hydrological modeling: (1) NSE [42], (2) a log-transformed NSE (NSE ln ) [34], (3) percent bias (PBIAS) [51], and (4) flow duration curve (FDC) index (NSEFDC) [8,54]. The NSE is known as a sensitive index to high flow, and it was used to evaluate high flow accuracy [54,55]. The NSE ln has an increased influence of low values as compared to the original NSE because of the logarithm, and thus it was employed to access low flow simulation [8,34,54]. The NSEFDC is a type of Nash-Sutcliffe efficiency designed to measure the similarity between FDCs [8,54,56], indicating the flow variability index. The values of NSE, NSEln, and NSEFDC are close to 1 when there is a complete agreement between simulated and observed streamflow, but they can become large negative values (−∞) when the discrepancy between them is wide. The measures the overall tendency or bias that the simulated data have compared to the observed [57]. More detailed equations of NSE, NSE ln , NSEFDC, and PBIAS can be seen in [8].

Relating Watershed Characteristics to Model Parameter Values
Multiple linear regression was carried out to find out the regression relationship between the calibrated Tank model parameter values (dependent variables) and watershed characteristics (independent variables): where y i is the calibrated parameter of the Tank model, x 1 , x 2 , . . ., x n−1 are the watershed characteristics, including surface geology and land use types, β 0 , β 1 , . . . , β n−1 are coefficients, and ε i is a model constant [58]. The logarithmic transformation for independent variables was considered to explore the best regression model. Three strategies of variable selection methods, forward, backward, and stepwise selections, were employed to determine the optimum number of independent variables; the strategy that provided the highest adjusted R 2 was finally selected. No more than five independent variables were used in each regression.

Study Watersheds and Their Characteristics
Forty-nine watersheds draining 3.8 to 2990.7 km 2 throughout South Korea were selected on the basis of the variability of their locations, land uses, and topography in this study [8,59] (Figure 2). The study watersheds represent a wide range of hydrological conditions, such as drainage areas (Area), watershed and channel slopes (W_Slope and S_Slope), drainage densities (Density) [60], flow lengths (Length), form factors (Form), and the percentages of forest (Forest), rice paddy fields (Paddy), and uplands (Upland) (Figure 2). These watersheds were randomly divided into two groups. The first group included 39 watersheds to be used in the development of the regional equations. The second group consisted of 10 watersheds to be used to verify the developed regression equations.
The precipitation records of weather stations associated with the study watersheds were obtained from the Korea Meteorological Administration (KMA), and the areal average precipitation was determined using the Thiessen polygon method [61]. Other daily weather variables, including temperature, relative humidity, mean wind velocity, and solar radiation, were obtained from the KMA and used to calculate the potential evapotranspiration (PET) by the FAO-PM method [48]. Topographic characteristics, including drainage area, watershed mean slope, channel slope, drainage density, flow length, and form factor, were calculated using 30-m digital elevation models (DEMs) provided by the National Geographic Information Institute (NGII). The percentages of forest, upland, and paddy areas were calculated from a land-use map obtained from the Ministry of Environment (MOE).

Parameter Calibration
We defined the and ′ as the cases of using and ′ as an objective function, respectively. The performance statistics were compared in terms of high and low flow, FDC, and runoff volume. The statistical significance of differences between the model performance statistics was investigated using a paired t-test at a significance level of 0.05.
Overall, the parameter calibration provided satisfactory performance in most study watersheds, Observed daily streamflow data of the study watersheds were compiled from the Korea Ministry of Land, Infrastructure, and Transport (MOLTM) and Seoul National University [59,62]. The length of available flow records varies from one watershed to another, and all study watersheds included at least five years of data. A split sample test scheme was employed to calibrate and validate the models; at least three and two years of streamflow records were used for calibration and validation, respectively. In addition, the calibration periods were set to include wet, average, and dry years [63,64]. The first two years of weather records were used to stabilize the hydrological variables of the model so that calibration results would be minimally affected by arbitrary and rough assumptions made for the initial conditions. In this study, the parameters of 3-Tank models were calibrated to streamflow observations made at the outlets of the 39 study watersheds with two different objective functions, RMSE and KGE sqrt (Table 1 and Figure 2).

Parameter Calibration
We defined the optRMSE and optKGE sqrt as the cases of using RMSE and KGE sqrt as an objective function, respectively. The performance statistics were compared in terms of high and low flow, FDC, and runoff volume. The statistical significance of differences between the model performance statistics was investigated using a paired t-test at a significance level of 0.05.
Overall, the parameter calibration provided satisfactory performance in most study watersheds, especially in the case of the optKGE sqrt ( Table 3). The optRMSE yielded slightly higher NSE values (or more accurate peak flow prediction) than optKGE sqrt in the calibration of the 3-Tank model parameters (p < 0.001) ( Figure 3). However, the differences between the NSE values were not statistically significant in the validation (p > 0.05). The use of KGE sqrt provided significantly better accuracy in predicting low flow (NSEln) and flow variability (NSEFDC) as compared to that of optRMSE (p < 0.001). In terms of PBIAS (or the overall water balance), the two objective functions provided no significant difference for both evaluation periods (p > 0.05). Such results imply that optKGE sqrt can provide more balanced evaluation aspects than optRMSE in model parameter calibration, as KGE is more responsive to low values than RMSE [8,45,47]. In this study, we included study watersheds where the corresponding models provided "satisfactory" performance (0.50 < NSE and PBIAS < ±15) so that we could exclude the influence of models that do not apply to the watersheds on parameter regionalization.   intermediate runoff; thus the three parameters were critical to the shapes of recession and baseflow parts of the streamflow hydrographs. Such findings indicate that the selection of an objective function in model calibration can substantially influence the hydrological analyses, including regionalization, by providing different mathematical representations for a watershed under consideration and by creating parameter uncertainty. In the following sections, we have demonstrated how the objective function selection can influence the regionalization of model parameters and their uncertainty.

Regionalization of the 3-Tank Model
The correlation structure between calibrated parameter values and selected hydrologic features of the watersheds was investigated to identify key watershed characteristics ( Figure 5). At least six of the parameters were correlated (| | ≥ 0.4, < 0.05) to one or more of the watershed features. The five parameters that were associated with the first layer (or the top tank) of the 3-Tank model, including , , ℎ , ℎ , and , turned out to be correlated to topographic factors, such as W_Slope, Length, S_Slope, and Form. Such a correlation structure was expected, as the first layer of a Tank model is usually introduced to simulate hydrological processes happening on the ground surface, such as direct runoff generation and routing [8,9,13].
In the case of using the objective function of , close correlation structures (0.7 > | | > 0.5, < 0.01) were found in between Form and the heights of two side outlets of the first layer, ℎ and ℎ ( Figure 5). When ′ was employed as the objective function, however, ℎ showed a close relationship with W_Slope ( = 0.58, < 0.001), and ℎ was associated with Form ( = −0.42, < 0.05) and Upland ( = −0.43, < 0.05). The outlet heights, ℎ and ℎ , of the first layer control the quick runoff response and high (and peak) flow of a watershed, respectively, and they control the surface storage capacity at the beginning of an event (ℎ ) and the total surface storage capacity (ℎ ) of the upper layer [8,17,37,65]. Thus, the findings imply that direct runoff of the study watersheds is relatively heavily controlled by Form, W_Slope, and Upland than the other watershed features, which is corroborated by our understanding and previous studies [8,[66][67][68].

Regionalization of the 3-Tank Model
The correlation structure between calibrated parameter values and selected hydrologic features of the watersheds was investigated to identify key watershed characteristics ( Figure 5). At least six of the parameters were correlated (|r| ≥ 0.4, p < 0.05) to one or more of the watershed features. The five parameters that were associated with the first layer (or the top tank) of the 3-Tank model, including a 11 , a 12 , h 11 , h 12 , and b 1 , turned out to be correlated to topographic factors, such as W_Slope, Length, S_Slope, and Form. Such a correlation structure was expected, as the first layer of a Tank model is usually introduced to simulate hydrological processes happening on the ground surface, such as direct runoff generation and routing [8,9,13].
In the case of using the objective function of RMSE, close correlation structures (0.7 > |r| > 0.5, p < 0.01) were found in between Form and the heights of two side outlets of the first layer, h 11 and h 12 ( Figure 5). When KGE sqrt was employed as the objective function, however, h 11 showed a close relationship with W_Slope (r = 0.58, p < 0.001), and h 12 . was associated with Form (r = −0.42, p < 0.05) and Upland (r = −0.43, p < 0.05). The outlet heights, h 11 and h 12 , of the first layer control the quick runoff response and high (and peak) flow of a watershed, respectively, and they control the surface storage capacity at the beginning of an event (h 11 ) and the total surface storage capacity (h 12 ) of the upper layer [8,17,37,65]. Thus, the findings imply that direct runoff of the study watersheds is relatively heavily controlled by Form, W_Slope, and Upland than the other watershed features, which is corroborated by our understanding and previous studies [8,[66][67][68]. The parameters of the second tank ( , ℎ , and ) and third tank ( ) determine the shapes of recession and baseflow-only parts of a streamflow hydrograph [8]. The calibrated values of are relatively strongly correlated to Paddy, regardless of the types of objective functions. In the case of calibrating with ′ , ℎ was found to be correlated to W_Slope ( = −0.51, < 0.01) or Forest ( = −0.50, < 0.01), and was also correlated to Upland ( = −0.51, < 0.01). The calibrated values of were associated with Density ( = −0.38, < 0.05) when ′ was used as the objective function; however, the parameter did now show any statistically significant relationship with the watershed characteristics in the case of the objective function. The parameters of the second tank (a 2 , h 2 , and b 2 ) and third tank (a 3 ) determine the shapes of recession and baseflow-only parts of a streamflow hydrograph [8]. The calibrated values of a 2 are relatively strongly correlated to Paddy, regardless of the types of objective functions. In the case of calibrating with KGE sqrt , h 2 was found to be correlated to W_Slope (r = −0.51, p < 0.01) or Forest (r = −0.50, p < 0.01), and b 2 was also correlated to Upland (r = −0.51, p < 0.01). The calibrated values of a 3 were associated with Density (r = −0.38, p < 0.05) when KGE sqrt was used as the objective function; however, the parameter did now show any statistically significant relationship with the watershed characteristics in the case of the RMSE objective function.

Performance of the Regionalized 3-Tank Models
The two regional models (regRMSE and regKGE sqrt ) were developed by relating the calibrated values of parameters to the watershed characteristics (Tables4 and 5). Regression equations were assigned to parameters that were at least "moderately" (|r| ≥ 0.35) correlated to any of the watershed features. When the correlation structure between the parameter values and watershed features was weak (|r| < 0.35), the median of the parameter values that had been calibrated to individual study watersheds was used to represent the overall average value of the parameter [15]. When KGE sqrt was used as the objective function in the calibration, stronger correlation structures were found between the parameter values and watershed characteristics as compared to the optRMSE ( Figure 5, Tables4 and 5).  The accuracy of the regionalized 3-Tank models was evaluated by comparing the model performance statistics provided by the calibrated and regionalized models for the 39 study watersheds employed in the regionalization (Figure 6). A one-way analysis of variance (ANOVA) was carried out to determine the statistical significance of differences between the performance statistics of the four groups (optRMSE, optKGE sqrt , regRMSE, and regKGE sqrt ). Subsequently, a post-hoc Tukey honest significant difference (HSD) test was performed to facilitate a pairwise comparison of the performance statistics provided by the models [8,69]. Overall, the two regionalized models ( and ′ ) provided similar accuracy to that of the calibrated models ( and ′ ). There was no statistically significant difference between the values achieved by the four models ( > 0.05), which implied that regionalization can predict high (or peak) flow at the level of accuracy similar to those of the calibrated models. In terms of and , however, ′ yielded better performance than the , presumably because ′ yielded better accuracy than ( < 0.01). However, ′ and ′ provided similar accuracy ( > 0.05). The slightly underestimated the overall runoff volume (e.g., positive ) as compared to ′ ; this may be attributed to ′ providing more balanced views on model performance than the [8,45,47]. The performance of the regionalized models ( and ′ ) was further investigated by applying them to the other 10 study watersheds that were not used in the regionalization processes (Figure 7). The sample size was small ( = 10 ); therefore, the nonparametric Wilcoxon signed-rank test was conducted to test the significance of any differences between the performance statistics provided by the two regionalized models at a significance level of 5%. The ′ model yielded , , and significantly better than those of Overall, the two regionalized models (regRMSE and regKGE sqrt ) provided similar accuracy to that of the calibrated models (optRMSE and optKGE sqrt ). There was no statistically significant difference between the NSE values achieved by the four models (p > 0.05), which implied that regionalization can predict high (or peak) flow at the level of accuracy similar to those of the calibrated models. In terms of NSE ln and NSEFDC, however, regKGE sqrt yielded better performance than the regRMSE, presumably because optKGE sqrt yielded better accuracy than regRMSE (p < 0.01). However, optKGE sqrt and regKGE sqrt provided similar accuracy (p > 0.05). The regRMSE slightly underestimated the overall runoff volume (e.g., positive PBIAS) as compared to regKGE sqrt ; this may be attributed to KGE sqrt providing more balanced views on model performance than the RMSE [8,45,47].
The performance of the regionalized models (regRMSE and regKGE sqrt ) was further investigated by applying them to the other 10 study watersheds that were not used in the regionalization processes ( Figure 7). The sample size was small (n = 10); therefore, the non-parametric Wilcoxon signed-rank test was conducted to test the significance of any differences between the performance statistics provided by the two regionalized models at a significance level of 5%. The regKGE sqrt model yielded NSE ln , NSEFDC, and PBIAS significantly better than those of regRMSE while they provided statistically similar NSE values. Such results highlight the potential of regKGE sqrt as a strategy for RR model regionalization.  The performance of the two models was compared with that of other regionalized models using the 10 validation watersheds ( Table 6). As seen in the comparison, the two models regionalized in this study outperformed the others in terms of and . The 3-Tank models regionalized by Kim et al. [41] provided a performance that was comparable to that of this study. However, An et al. [18] reported relatively poor performance as compared to others. We applied one-way ANOVA and HSD tests to see if there were significant differences between the performance statistics provided by the 3-Tank models regionalized for the 10 validation watersheds (Figure 8).  The performance of the two models was compared with that of other regionalized models using the 10 validation watersheds ( Table 6). As seen in the comparison, the two models regionalized in this study outperformed the others in terms of NSE and PBIAS. The 3-Tank models regionalized by Kim et al. [41] provided a performance that was comparable to that of this study. However, An et al. [18] reported relatively poor performance as compared to others. We applied one-way ANOVA and HSD tests to see if there were significant differences between the performance statistics provided by the 3-Tank models regionalized for the 10 validation watersheds (Figure 8). Table 6. Comparison of the performance of the two regional models (regRMSE and regKGE sqrt ) and others prepared for the 10 validation watersheds. The five regionalized models employed different regression equations that explained the relationships between parameter values and watershed characteristics, but there was no statistically significant difference between the NSE values provided by them (p > 0.05). In terms of NSEln and NSEFDC, however, regKGE sqrt provided significantly better performance as compared to ones that had been developed in Kim and Park [13] and Kim et al. [41]. The 3-Tank model regionalized by An et al. [18], yielded NSEln and NSEFDC values similar to those by regKGE sqrt (p > 0.05), but the model significantly underestimated runoff volume (positive PBIAS) compared to the other models (p < 0.05). The regionalized model of Huh et al. [40] provided a level of efficiency similar to that of regKGE sqrt , but regKGE sqrt outperformed the model in terms of the model evaluation criteria that are commonly employed in hydrological modeling practices [57] (Table 6). , and ) provided by and other regional models for 10 verification watersheds.

Model
The five regionalized models employed different regression equations that explained the relationships between parameter values and watershed characteristics, but there was no statistically significant difference between the values provided by them ( > 0.05). In terms of and , however, ′ provided significantly better performance as compared to ones that had been developed in Kim and Park [13] and Kim et al. [41]. The 3-Tank model regionalized by An et al. [18], yielded and values similar to those by ′ ( > 0.05), but the model significantly underestimated runoff volume (positive ) compared to the other models ( < 0.05). The regionalized model of Huh et al. [40] provided a level of efficiency similar to that of ′ , but ′ outperformed the model in terms of the model evaluation criteria that are commonly employed in hydrological modeling practices [57] (Table 6).
The relatively worse performances of the models reported by Kim and Park [13] and Kim et al. [41], reproducing low flow, might be related to the models being developed on manually calibrated parameters; thus, the local optima could be used for regionalization (Table 2). In addition, the studies did not calibrate models by considering low flow efficiencies in the development process. Huh et al. [40] and An et al. [18] employed automatic optimizations (Table 2), and their models provided better low flow performances. However, the models produced relatively worse performances as compared The relatively worse performances of the models reported by Kim and Park [13] and Kim et al. [41], reproducing low flow, might be related to the models being developed on manually calibrated parameters; thus, the local optima could be used for regionalization (Table 2). In addition, the studies did not calibrate models by considering low flow efficiencies in the development process. Huh et al. [40] and An et al. [18] employed automatic optimizations (Table 2), and their models provided better low flow performances. However, the models produced relatively worse performances as compared to regKGE sqrt as these studies used RMSE as an objective function (Table 2), which may yield less balanced results than the case when KGE sqrt is used.
In this study, we compared the prediction performance of two regionalized 3-Tank models that were calibrated with two different objective functions, optRMSE and optKGE sqrt . From the comparison, we found that the models could predict high flow and water balance of the 49 study watersheds at acceptable levels of accuracy (Table 6, Figures 6-8). We also saw that the use of RMSE that has been widely employed as an objective function in the regionalization studies provided relatively poor performance in reproducing low flow and FDC compared to the use of KGE sqrt . Such a result is not surprising because KGE sqrt considers flow variability more explicitly when evaluating model accuracy while NSE sqrt is still sensitive to peak or high flow [47]. It is worth noting that optKGE sqrt gained significantly higher NSEln and NSEFDC values but slightly lower NSE efficiency as compared to optRMSE (Figure 3). Such a finding suggests that when low flow is one of the modeling outputs of interest, we should not solely rely on statistics including NSE and R 2 that have been commonly used but also know that the statistics are very sensitive to high flow. Instead, combining alternative statistics, such as NSEln and NSEFDC, could provide a balanced view point to model evaluation [8,9,55,70].
The KGE and its variants are kinds of the Euclidian distance (ED) between the optimal and ideal points for correlation, bias, variance, and variability (Equations (10)-(15)). Pfannerstill et al. [71] proposed another form of the ED (ED DM , Equation (17)) as an objective function. The two ED measures, KGEs and ED DM , are similar to each other because both consider pairwise differences (r and NSE), water balance (β and PBIAS), and flow variability (α, γ, RSRQ0Q5, RSRQ5Q20, and RSRQ20Q70): where RSRQ0Q5, RSRQ5Q20, and RSRQ20Q70 are RMSE-observations standard deviation ratio (RSR) [64] for the FDC segments of 0% to 5%, 5% to 20%, and 20% to 70%, respectively. We found that the values of ED DM are highly correlated with those of KGE sqrt in the study dataset (R 2 = 0.79; Figure 9). Such a finding implies that the ED-based statistics including KGEs and ED DM could serve as an objective function to efficiently count for the multiple evaluation aspects in an RR model regionalization.
comparison, we found that the models could predict high flow and water balance of the 49 study watersheds at acceptable levels of accuracy (Table 6, Figures 6-8). We also saw that the use of that has been widely employed as an objective function in the regionalization studies provided relatively poor performance in reproducing low flow and FDC compared to the use of ′ . Such a result is not surprising because ′ considers flow variability more explicitly when evaluating model accuracy while is still sensitive to peak or high flow [47]. It is worth noting that ′ gained significantly higher and values but slightly lower efficiency as compared to (Figure 3). Such a finding suggests that when low flow is one of the modeling outputs of interest, we should not solely rely on statistics including and that have been commonly used but also know that the statistics are very sensitive to high flow. Instead, combining alternative statistics, such as and , could provide a balanced view point to model evaluation [8,9,55,70].
The and its variants are kinds of the Euclidian distance (ED) between the optimal and ideal points for correlation, bias, variance, and variability (Equations (10)-(15)). Pfannerstill et al. [71] proposed another form of the ED ( , Equation (17)) as an objective function. The two ED measures, s and , are similar to each other because both consider pairwise differences ( and ), water balance ( and ), and flow variability ( , , where 0 5, 5 20, and 20 70 are RMSE-observations standard deviation ratio ( ) [64] for the FDC segments of 0% to 5%, 5% to 20%, and 20% to 70%, respectively. We found that the values of are highly correlated with those of in the study dataset ( = 0.79; Figure 9). Such a finding implies that the ED-based statistics including s and could serve as an objective function to efficiently count for the multiple evaluation aspects in an RR model regionalization.

Summary and Conclusions
This study explored ways to accurately regionalize an RR model, i.e., the 3-Tank model, and evaluated the performance of models regionalized on the basis of the calibration results made with different objective functions. We also demonstrated the impacts of objective function selection on RR

Summary and Conclusions
This study explored ways to accurately regionalize an RR model, i.e., the 3-Tank model, and evaluated the performance of models regionalized on the basis of the calibration results made with different objective functions. We also demonstrated the impacts of objective function selection on RR model regionalization. Results showed that there was no significant difference between the performance of optRMSE and optKGE sqrt when predicting high flow and water balance. However, the optRMSE provided poorer performance in reproducing low flow and FDC than optKGE sqrt . The regKGE sqrt provided high flow modeling performance similar to that of the calibrated models (optRMSE and optKGE sqrt ). In addition, regKGE sqrt was superior to the regRMSE when predicting low flow, water balance, and FDC in the study watersheds. Such evaluation results suggested that regKGE sqrt can serve as an effective objective function when regionalizing a daily rainfall-runoff simulation model for ungauged watersheds in Korea. The regionalization methods proposed in this study should be applicable to other watersheds, even though the regression equations developed using the method might not be suitable for them.