A framework for event-based flood scaling analysis by hydrological modeling in data-scarce regions

Flood scaling, which is extensively applied to flood prediction in ungauged and poorly gauged river basins of different sizes, is usually expressed by an equation that relates peak discharge to catchment descriptors such as the drainage area, slope of the watershed, land use, and rainfall characteristics; among these descriptors, the drainage area is the most widely used explanatory variable to predict flood peak discharge (Jothityangkoon & Sivapalan 2001; Menabde & Sivapalan 2001; Galster et al. 2006; Al-Rawas & Valeo 2010; Ishak et al. 2011; Han et al. 2012; Farmer et al. 2015; Furey et al. 2016; Lee & Huang 2016). The power law function recommended worldwide is Q = αA^θ, where Q is the peak discharge for a given drainage area A, α is the scaling intercept, and θ is the scaling exponent (Gupta et al. 1996; Medhi & Tripathi 2015). If a reference river basin with drainage area A₁ has a flood peak discharge Q₁ and θ is determined, then the flood peak Q₂ in an ungauged basin with drainage area A₂ can be estimated by the following equation: Q₁/Q₂ = (A₁/A₂)^θ.

Most previous studies focused on the relation between the quantile-based peak discharge and drainage area. The peak discharge shows simple scaling or multiscaling for snowmelt-generated or rainfall-generated floods, which means that the scaling parameters (especially the scaling exponent) remain invariant or vary, respectively, with the flood return period (Gupta & Dawdy 1995). This hypothesis has been tested in many river basins of different sizes from small-scale to mesoscale basins (Ogden & Dawdy 2003; Gupta et al. 2007; Yue & Gan 2009; Kroll et al. 2014; Ayalew et al. 2015). Generally, the scaling exponent is less than 1 in quantile-based flood scaling but might exceed 1 in event-based flood scaling. Moreover, flood scaling analyses are often conducted in homogeneous and nested watersheds, in which the flood generation mechanisms are similar (Eash 2001; Ogden & Dawdy 2003; Furey & Gupta 2005, 2007). Accordingly, the function between the flood quantile and watershed area can be established, and the scaling parameters can be obtained. However, the fitted function may not have satisfactory precision; in this situation, other controlling factors, such as the corresponding quantile-based rainfall intensity, should be taken into consideration in addition to the drainage area.

Although the quantile-based scaling exponent is widely used for future flood prediction, the physical mechanism cannot be interpreted because the peak discharges of floods with the same return period in different catchments may not be generated by the same rainfall events (Liu et al. 2019). Therefore, Gupta et al. (2015) proposed performing flood scaling research by using event-based floods for statistical analysis. As such, observed floods that occur simultaneously in different catchments can be selected, and the spatial and temporal distributions of rainfall can be considered. Furey & Gupta (2005) presented a new index to reflect the temporal variability of excess rainfall, and the results showed that the scaling parameters depended mainly on the duration and amount of excess rainfall.

Together with statistical analysis, hydrological modeling can provide intuitive results regarding how meteorological and land surface factors impact flood scaling parameters (Yang et al. 2020). For simulated flood events, other controlling factors of the scaling parameters, such as the channel slope, flow velocity, and rainfall distribution, can be identified by changing the hydrological model parameters (Mantilla 2007; Venkata 2009). Therefore, it is better to use physically based distributed hydrological models that contain all the important parameters reflecting flood generation processes. For example, Ayalew et al. (2014) applied the CUENCAS rainfall–runoff model to three river basins in Iowa to investigate how the rainfall intensity, duration, hillslope overland flow velocity, and channel flow velocity influenced the scaling parameters; they found that the key factors affecting the scaling parameters were the rainfall duration and overland flow velocity.

The main aims of this paper are to (1) perform a simulated flood scaling analysis and evaluate the effects of the rainfall characteristics and land surface and (2) propose a framework for estimating the scaling exponent by hydrological modeling. Compared with previous studies, the novelty of this study is that a framework for determining the flood scaling exponent is proposed by hydrological modeling in data-scarce regions, especially for multiscaling. Hence, the proposed method can be used in other data-scarce regions on a global scale.

The flood scaling analysis was carried out in the Zijingguan catchment, situated in the northeastern part of China. The catchment has a drainage area of 1,760 km². The long-term average annual precipitation is approximately 650 mm, 70–80% of which is concentrated in the flood season (July, August, and September). Flood events are mostly generated by rainstorms with a high intensity and short duration.

Hourly observed flood data at the Zijingguan hydrological station were used for flood modeling. The corresponding hourly rainfall data at seven rainfall gauging stations located in the Zijingguan catchment were obtained from the Hydrology and Water Resource Survey Bureau of Hebei Province. The hourly data range from 1956 to 2018. Since no large floods occurred after 2000, some flood events before 2000 were selected for simulation and flood scaling analysis.

To establish a hydrological model, a digital elevation model (DEM) with a resolution of 30 m was downloaded from the Geospatial Data Cloud (http://www.gscloud.cn/). Remotely sensed land use data from 1980 to 2000 were collected from the Institute of Geographic Sciences and Resources of the Chinese Academy of Sciences (http://www.resdc.cn/) and were classified into forestland, grassland, agricultural land, water area, urban area, and unused land. Soil data were downloaded from the Environmental and Ecological Science Data Center for West China (http://westdc.westgis.ac.cn/).

The Hydrologic Engineering Center-Hydrologic Modeling System (HEC-HMS) model developed by the US Army Corps of Engineers is widely used for watershed hydrological simulations. This model has been verified to be applicable in semi-arid regions throughout China (Wang & Sun 2019). In addition, the Soil Conservation Service curve number (SCS-CN) method is often used to calculate the amount of generated runoff, and the SCS unit hydrograph method is used for flow routing. Therefore, we selected these two modules for the model. The CN value is determined based on the spatial distributions of the soil water content and vegetation types (A et al. 2019) and will vary considerably under wet and dry antecedent conditions (Han et al. 2018).

This study aims to analyze the flood scaling exponent in a data-scarce region. Therefore, the Zijingguan catchment with a drainage area of 1,760 km² (Figure 1), which has only one hydrological station (Zijingguan station), is selected. Scaling analysis is conducted by using the simulated flood peak discharges from 11 subcatchments. The HEC-HMS model, which is applicable in semi-humid and semi-arid regions, is employed to simulate flood events due to its data availability and structure (Alfy 2016). The simulated peak discharges and corresponding subcatchment areas are fitted with power functions, and then the scaling exponent is obtained for each flood event.

Based on the simulated flood events and scaling analysis, a framework that links event-based scaling exponents to quantile-based scaling exponents can be proposed. There are two main assumptions: (1) the input rainfall at each rain gauge is spatially uniform and (2) the rainfall events have the same frequency as the generated flood events. The idea and procedures are as follows:

1. A typical rainfall event is selected from the observed data. A typical rainfall event should have a large rainfall amount and high rainfall intensity. Generally, the largest rainfall event is selected for its rare occurrence as the low-frequency design flood.

2. The typical rainfall event is amplified. The typical rainfall event is amplified by multiple constants, such as 1.5, 2, 2.5, … , to form hypothetical rainfall events with the same temporal distribution as the typical rainfall event.

3. The amplified rainfall events are used as inputs to drive the established HEC-HMS model, and the corresponding flood events at the outlets of the subcatchments and Zijingguan catchment are output with various CN values ranging from 10 to 100.

4. The flood scaling is analyzed for all the simulated flood events generated by the hypothetical rainfall events. Then, using the relationship between the scaling exponent and rainfall intensity (total amount), the CN value is obtained, and the existence of simple scaling or multiscaling is concluded.

5. If the flood shows simple scaling, then the scaling exponent can be used in quantile-based flood scaling.

6. If the flood shows multiscaling, then the rainfall and the corresponding flood events are assumed to have the same frequency. The rainfall with a specific frequency is obtained from rainfall frequency analysis, and the corresponding scaling exponent is calculated according to the relationship established in step 4.

The framework is illustrated in Figure 2.

The Zijingguan catchment was divided into 11 subcatchments, and 17 outlets were selected with drainage areas ranging from 3 to 1,760 km². Fourteen observed flood events were selected for HEC-HMS calibration and validation, including large floods (>5-year return period with a peak discharge greater than 500 m³/s) and small floods (<5-year return period with a peak discharge smaller than 500 m³/s). The model parameters were adjusted artificially to fit the observed flood hydrographs, and comparisons between the simulated and observed flood events are listed in Table 1 to estimate the model performance. Some simulated and observed flood events are shown in Figure 3.

The RE values of the peak discharges for all 14 flood events are less than 20%, ranging from −1.98 to 19.82%, and the RE values of the flood depths range from −13.62 to 17.93%. The NSE coefficients of the 14 flood events are between 0.51 and 0.95 with an average of 0.77. Therefore, the calibrated model is capable of capturing the peak discharges, and the RE and NSE coefficient have acceptable accuracies. The simulated flood hydrographs fit the observations very well during both the calibration and the validation periods, which implies that the performance of the HEC-HMS model is sufficient to simulate flood events in the Zijingguan catchment and can be used in the following sections to analyze the flood scaling in this area.

According to the simulated flood events in each subcatchment, the scaling of the peak discharge for each flood event was analyzed with the subcatchment area as the only explanatory variable. The scaling parameters for each flood event are listed in Table 2. The power function fit the plots very well with acceptable coefficients of determination except for the flood event that occurred on 30 July 1988, verifying that flood scaling exists in this catchment.

The scaling exponent θ varies for different flood events in the range from 0.5396 to 0.9922 with an average of 0.7829, which is less than 1. These results reflect multiscaling, and large flood events show greater scaling exponents and coefficients of determination (Figures 4 and 5). The values of θ for large floods range from 0.7621 to 0.9234 with an average of 0.8634, while those of small floods range from 0.5396 to 0.9922 with an average of 0.7225. Although the values of θ for two small flood events (occurring on 13 August 1973 (730,813) and 6 August 1975 (750,806)) are relatively high, the R² is approximately 0.66, which is lower than that of large floods. These results imply that there must be other factors, such as rainfall and land surface properties, influencing the statistics (Fang et al. 2018).

Due to the large drainage area of the Zijingguan catchment, rainfall is not spatially uniform. Figure 6 shows the spatial distributions of the two rainfall events. Each rainfall event was distributed nonuniformly over the whole catchment, especially the small rainfall events. This may be the reason for the small values of the fitted R² in the above section.

To research the effects of the spatial rainfall distribution on the flood scaling properties, the areally averaged rainfall calculated by the Thiessen polygon method was input at each of the seven rain gauges. Then, the simulated flood peak discharge at each subcatchment was subjected to flood scaling analysis, as shown in Figure 7 and Table 3.

The peak discharge simulated by uniform rainfall shows multiscaling as well. Compared with the fitted power function determined by the observed rainfall, all coefficients of determination increase, demonstrating stronger flood scaling (Figure 8) with the lowest R² of 0.6149. The scaling exponent also varies with an average of 0.7829, which is the same as the value under the condition of nonuniform spatial rainfall. For small flood events, subcatchment W150 has a relatively higher peak discharge, as shown in Figure 7 because of its higher calibrated CN value (CN = 65). W150 is located at the lower part of the catchment, and thus, the initial soil moisture content may be higher than that in the other subcatchments, leading to a higher runoff generation capacity and higher peak discharge. However, for large flood events, the peak discharge of W150 shows obvious scaling properties. This may be because the subcatchments are unsaturated during small flood events but become saturated for large flood events.

On the basis of the above 14 scaling exponents θ simulated by the areally averaged rainfall, the relationships between the maximum 1 h rainfall (P₁) and θ and between the total rainfall (P) and θ were analyzed (Figure 9). The scaling exponent has no obvious relation with P₁ or P, which is uncommon. This may be because the rainfall intensity and total rainfall amount from different rainfall events have different temporal distributions. Therefore, the effects of the rainfall intensity on flood scaling need to be studied with the same spatial and temporal distributions.

The rainfall event that occurred on 6 August 1963 was selected as the typical event, based on which some hypothetical rainfall events were obtained by amplifying the typical rainfall event to drive the HEC-HMS model. The flood hydrographs were simulated under different initial CN values. Then, a scaling analysis was conducted, as shown in Figure 10. The larger the maximum 1 h rainfall amount and CN value are, the larger the scaling exponent. When the CN value is sufficiently large (CN ≥ 70), the exponent becomes constant regardless of the rainfall intensity. For CN values less than 70, the scaling exponent increases with increasing maximum 1 h rainfall and finally approaches a constant. From this figure, if we determine the CN value and the design maximum 1 h rainfall, we can obtain the scaling exponent and make flood predictions.

Based on the above flood scaling analysis, flood predictions can be conducted by implementing a hydrological simulation according to the following steps. (1) A hydrological model is established first and then calibrated and validated. (2) A typical rainfall event with a high rainfall intensity and a large total amount is selected, and some hypothetical rainfall events are formed by amplifying the typical rainfall event to keep the temporal distribution of the hypothetical rainfall events the same. (3) The hypothetical rainfall events are areally averaged and used to drive the hydrological model to simulate the peak discharge for each subcatchment. (4) Flood scaling analysis is conducted by using the simulated peak discharge and the area of each subcatchment, and the relations among the scaling exponent, rainfall and CN values are built, as shown in Figure 10. (5) According to the assumption that the rainfall and the corresponding flood have the same frequency, from Figure 10, the scaling exponent of a specific quantile can be obtained by interpolation under any rainfall intensity and the CN value.

Flood scalings in ungauged basins should be analyzed by analyzing homogeneous regions with sufficient flood data. In watersheds possessing only one or a few hydrological stations, flood scaling analysis cannot be conducted by fitting the observed flood peak and watershed attributes. Therefore, hydrological modeling is an approach for performing scaling research by simulating the flood processes in each subwatershed.

During the flood scaling analysis, the spatial rainfall distribution was found to be a controlling factor for the scaling exponent. Spatially uniform rainfall resulted in better fitting of the simulated peak discharge and watershed area and resulted in a different scaling exponent. Hence, the difference between the scaling exponent obtained from fitting the observed rainfall and that obtained from fitting the areally averaged rainfall cannot be neglected. In practice, we use the scaling exponent obtained from fitting the observed peak discharge. However, this may lead to a significant error.

Similarly, land surface conditions significantly influenced the scaling exponent. However, for different rainfall events, the relationship between θ and the maximum 1 h or total rainfall is not good. Since all 14 flood events were simulated by using the same CN value, the temporal rainfall distribution must have a significant influence on θ. Furey & Gupta (2005) presented a variable to represent the effective rainfall duration and analyzed the relationship with the scaling exponent. Therefore, the influence of the temporal rainfall distribution needs further research.

In addition, the simulated peak discharge of each subcatchment was used for flood scaling analysis. Thus, the result must be affected by the peak discharge simulation error. In this catchment, the peak discharge simulation errors ranged from −1.98 to 19.82%, with the largest simulation error obtained for the event that occurred on 6 August 1963. However, this flood event was selected to illustrate the multiscaling analysis framework. If the result is applied to a data-scarce region, the simulation accuracy should be improved further.

Gupta & Dawdy (1995) noted that floods caused by rainstorms exhibited multiscaling, and the results of our study further verify their statements. However, our results were obtained in semi-humid and semi-arid regions with multiscaling, whereas in humid regions, the initial CN value is large, and simple scaling may be inferred from the findings in this study.

In this study, we proposed a method to conduct flood scaling analysis by hydrological modeling in data-scarce regions. The scaling exponents were very different between using the observed rainfall and the areally averaged rainfall as model inputs, and the coefficient of determination was larger for spatially uniform rainfall. Therefore, the spatial rainfall distribution has a significant impact on flood scaling.

Based on how the spatial rainfall distribution and CN values influence the scaling exponent, a framework was proposed to estimate the scaling exponent in a data-scarce region. We recommend that the areally averaged rainfall of an event be used as the input to the hydrological model; accordingly, the correlation curves relating the scaling exponent, rainfall intensity, and CN values can be drawn. As long as the flood frequency and CN value are determined, the scaling exponent can be obtained from these correlation curves.

This research is supported by the National Key Research and Development Program of China (2018YFC0407902) and the National Natural Science Foundation of China (No. 51779165). We extend sincere thanks to the Hydrology and Water Resource Survey Bureau of Hebei Province for providing the research data.

Data cannot be made publicly available; readers should contact the corresponding author for details.