Land-Use and Land Cover Is Driving Factor of Runoff Yield: Evidence from A Remote Sensing-Based Runoff Generation Simulation

Abstract

The spatial distribution of water storage capacity has always been the critical content of the study of saturation-excess runoff. Xin’anjiang model uses the water storage capacity curve (WSCC) to characterize the distribution of water storage capacity for runoff yield calculation. However, the mathematical and physical foundations of WSCC are unclear, which is impossible to simulate runoff generation with complex basins accurately. To fill this gap, we considered the dominant role of basin physical characteristics in water storage capacity and developed a new integrated approach to solve the spatial distribution of water storage capacity (L-WSCC) to account for the spatiotemporal dynamics of their impact on runoff generation. The main contribution of L-WSCC was to confer WSCC more physical meaning and the spatial distribution of water storage capacity was explicitly represented more accurately, so as to better express the runoff generation and provide a new approach for runoff yield calculation in non-data basin. L-WSCC was applied to Misai basin in China and promising results had been achieved, which verified the rationality of the method (the mean Nash–Sutcliffe efficiency (NSE):0.86 and 0.82 in daily and hourly scale, respectively). Compared with WSCC, the performance of L-WSCC was improved (mean NSE: 0.82 > 0.78, mean absolute value of flood peak error (PE): 12.74% < 21.66%). Moreover, the results of local sensitivity analyses demonstrated that land-use and land cover was the major driving factor of runoff yield (the change of mean absolute error (ΔMAE): 131.38%). This work was significant for understanding the mechanisms of runoff generation, which can be used for hydrological environmental management and land-use planning.

1. Introduction

Floods are the most common and highly destructive natural disaster in the world [1,2,3,4]. According to relevant research, the annual casualties caused by floods over the world exceed 50% of the total natural disasters, and the economic losses exceed 30% [5,6]. With the population surging and the acceleration of urbanization, the frequency of extreme rainfall events has become higher and higher, and the frequency, intensity, loss and the ranges of impact of flood disasters have shown a significant upward trend [7,8]. China is one of the countries with the most serious flood disasters in the world, and the flooding issues have become a comprehensive problem that cannot be ignored in the sustainable development [9]. Therefore, it is urgent to implement reasonable and scientific flood management measures to minimize the damage caused by floods.

Accurate runoff simulation is essential for flood management and plays a vital role in sustainable water resources management in reducing losses caused by flood disasters [2,10,11,12,13]. Water storage capacity curve (WSCC) is the key to the calculation of saturation-excess runoff and is applied to describe the distribution of basin water storage capacity for runoff simulation [14]. WSCC is represented by empirical Generalized Pareto Distribution, which has been extensively applied worldwide for a wide range of watershed scales, hydrology environmental conditions, flood simulation, and provides a scientific way to analyze the process of rainfall-runoff [15,16,17,18].

Although the WSCC method performs well in humid and subhumid regions for nature basins, it cannot characterize its connection to the internal characteristics of the basin due to its empirical origins [19,20,21,22]. This limitation causes that it performs only fairly even poorly in complex basins: the findings of Jayawardena et al. suggested that WSCC is less effective for runoff simulation in wet and dry transition seasons [23]; Zhuo et al. found that the soil moisture schemes of WSCC are generally not realistic in comparison with the reality [24]. Bartlett et al. pointed out that WSCC does not describe the spatial and temporal variability of runoff, which limited its application [25,26]. In addition, studies also showed that WSCC may fail to capture the vertical partial runoff [27]. These poor performances may be due to the WSCC method assuming a uniform basin with spatially uniform runoff processes [28]. Moreover, the WSCC method also lacks a description of the spatial variability of runoff and the physical meaning of WSCC is not clear enough, which is a generalized model [29,30]. This generalization ignores the influence of underlying surface changes on runoff yield [31,32,33]. Hence, when WSCC is applied to complex basins, significant errors can and do occur when calculating runoff yield.

To move beyond the limitations of the WSCC method, the related research have been done to study the relationship between the WSCC and related physical factors based on the sufficient remote sensing data [19,27]. Among them, water storage capacity and land-use types were strongly correlated because landscape had significant impacts on soil physical properties [34,35]. In vegetation region, many studies had been undertaken and the results showed that vegetation structures such as leaf area index, vegetation rooting depth, root depth, vegetation coverage and canopy density would profoundly affect water storage capacity [36,37,38,39,40]. Additionally, the effects of topography and terrain on water storage capacity were of great concern [25,30,41,42,43,44]. A significant number of studies had tried to calculate the water storage capacity using topographic index [22,45,46], slope [19,47,48], and elevation [30]. Furthermore, many studies evaluated the soil characteristics such as soil type [25], soil texture [30], soil wetness index [49], soil layer depth [29], soil hydraulic characteristics [19], soil moisture [50], soil humidity, and soil wetness [51], on water storage capacity to represent the soil water availability [30,52,53,54]. The effect of watershed area on water storage capacity was also studied [19,55].

As can be seen from the above studies, the WSCC is most closely associated with land-use, topography, soil, and vegetation. Consequently, many scholars have described the spatial distribution of water storage capacity based on these factors to produce a more reasonable runoff generation response: the results of Yao et al. showed that the water storage capacity of a basin can be derived by the soil texture and vegetation type [30]. Wang et al. solved the spatial distribution of water storage capacity through the relationship between the soil wetting ratio and the storage index, and achieved good outcome [56]. Guo et al. studied the water storage capacity based on topographic index and incomplete beta distribution and the results showed the new method produced a more reasonable runoff simulation [57]. Moreover, Gao et al.’s research in 404 watersheds showed that the spatial distribution of water storage capacity is the result of the combined effect of remote sensing factors such as slope, land-use, soil, and topography [19].

Although there have been many studies on the spatial distribution of water storage capacity, these studies have only conducted research from a certain aspect. No method has been proposed to solve the spatial distribution of water storage capacity by integrating various remote sensing factors to date. Moreover, no report quantified the influence of each factor on the runoff yield. To fill this gap, a comprehensive method considering various factors for spatial distribution of water storage capacity based on abundant remote sensed data was proposed in our study, called L-WSCC. The L-WSCC has physical basis in a spatial distribution of water storage capacity that is derived from distributions of spatially heterogeneous variables such as land-use, terrain geomorphology, and soil. The object of the paper is to construct a physics-based water storage capacity method through various remote sensing factors to improve the physical basis of the original empirical WSCC. L-WSCC can provide a new representation of the runoff generation, so as to extend the traditional WSCC with a spatial characterization of runoff variability. In addition, this paper also quantitatively assessed which factors dominated the impact on runoff yield and how much it accounted for. It was expected that the proposed method could enhance the physical meaning of WSCC to improve its runoff simulation capability.

2. Materials and Methods

2.1. Study Area and Data Collection

In this study, Misai watershed was employed as study area. It is a complex forest watershed located in Quzhou city, Zhejiang Province, China, with an area of 797 km² (29.1°–30.0° N, 118.0°–119.0° E) (Figure 1a). Its altitude varies between 105 m and 1253 m, and the overall terrain decreases from northwest to southeast (Figure 1b). The watershed has one major tributary, Majinxi, one of the sources of the Qiantang River (Figure 1b). Misai watershed belongs to the mid-subtropical zone and has frequent monsoon activities. This region has a mean annual temperature of 17 °C, with average minimum and maximum temperatures of 12.9 °C and 21.3 °C, respectively [58]. The climate is humid with an average annual rainfall from 1500 mm to 2000 mm. 60% of the rainfall occurs during short heavy rainfall in summer. The humic acrisols (FAO soil classification, ACu) soil is the main soil type and the soil moisture within 200 cm is about 352.80 kg/m² (https://ldas.gsfc.nasa.gov/gldas/, accessed on 9 November 2021). The watershed is rich in forest resources, accounting for 86.3% of the area, including woodland, shrub forest, and sparse woodland.

Flood disasters have occurred frequently in this region. According to statistics, flood disasters of varying degrees occurred in the region (big floods in 1955, 1969, and 1999) and the flood disasters of 2006 and 2017 were highly devastating, which caused heavy loss on property and society [59,60]. It is, therefore, of great importance to investigate the hydrological regime in this region. The study area is chosen because of the threat of flooding and its compliance with the saturation-excess runoff.

The hydrological data used in this paper mainly come from 6 hydrological observation stations (included 6 rainfall stations, 1 runoff observation station and 1 evaporation station) (Figure 1b), including rainfall, evaporation, and runoff data (daily data) and 17 floods data (hourly data) (Table S1).

The remote sensing data used in this paper mainly include DEM data, Landsat data (land cover data), soil data, and vegetation type map. The DEM data set is provided by Geospatial Data Cloud site, Computer Network Information Center, Chinese Academy of Sciences (http://www.gscloud.cn, accessed on 9 November 2021). Five different resolutions are obtained through resampling: 30 m, 100 m, 200 m, and 500 m to study the change of sub-grid so as to weaken the scale influence as much as possible (see Section 2.3). The corresponding terrain information (relief, slope, and aspect) are shown in Figure S1. Landsat satellite data are taken from NASA website (https://earthexplorer.usgs.gov, accessed on 9 November 2021). Seven phases of TM/MSS satellite images with less than 10% cloud cover were obtained. The synthesis time of the image is about June of 1982–1988. Based on remote sensing images, the land-use classification of the study area is obtained through supervised classification. The land-use map is shown in Supplementary Materials (Figure S2). The results are shown in Table S2. The data of soil types comes from the global soil data released by FAO, which is in accordance with the FAO/UNESCO System of Soil Classification. The distribution of soil types is shown in Figure S3a. The vegetation data are from the vegetation map of the people’s Republic of China (1:1,000,000). The vegetation area of the study area accounts for a large proportion, mainly conifer and broadleaved mixed forests, which are mainly located in the northern mountainous area (Figure S3b).

2.2. WSCC and L-WSCC

The WSCC is the main form of calculating runoff yield in saturation-excess runoff. It is the cumulative frequency curve of part of the area varying with the storage capacity and represents the surface distribution of water storage capacity of the basin [14]. Rainfall-runoff relationships can be obtained through WSCC. According to experience and experimental data, the WSCC is parabolic usually. Please see Text S1 for more details.

In this paper, the distribution of water storage capacity is calculated based on different land-use combined with relevant soil, terrain, vegetation, and other data, aiming to generate a WSCC with obvious physical meaning in entire watershed. On the basis of land-use and combined with the similarity of water storage capacity of different land-use types, the watershed is divided into five sub-regions: vegetation region (Figure 2a), bare land region (Figure 2b), reservoir region (Figure 2b), impervious surface region (hardened surface and buildings) (Figure 2d), and another region (Figure 2e).

Vegetation processes are often ignored in conceptual runoff [61]; however, for the watershed with good vegetation, the vegetation has a great influence on the runoff. Relevant studies have shown that the water storage capacity of vegetation mainly includes the canopy interception layer, the litter interception layer, and the water storage of soil and roots [24,27,62,63]. Therefore, the water storage capacity of the vegetation region is calculated as follows (more details in the Text S2 of Supplementary Material):

$${\mathrm{W}}_{\mathrm{i}}={\mathrm{f}}_{\mathrm{i}}\times {\mathrm{LAI}}_{\mathrm{i}}\times {\mathrm{h}}_{\mathrm{i}}+{\mathrm{l}}_{\mathrm{i}}+{\mathrm{N}}_{\mathrm{i}}\times \frac{{\mathrm{D}}_{\mathrm{i}}-{\mathrm{H}}_{\mathrm{i}}}{100},$$

(1)

where ${\mathrm{W}}_{\mathrm{i}}$ is the storage capacity of the i-th grid located in the vegetation area (mm); ${\mathrm{f}}_{\mathrm{i}}$is the vegetation coverage of the i-th grid (0~1)$;{ \mathrm{LAI}}_{\mathrm{i}}$is the vegetation leaf area index of the i-th grid; ${\mathrm{h}}_{\mathrm{i}}$ is the maximum water holding capacity of vegetation leaves in the i-th grid (mm); ${\mathrm{l}}_{\mathrm{i}}$ is litter interception (mm), obtained by remote sensing techniques or related research [64,65,66,67]; ${\mathrm{N}}_{\mathrm{i}}$ is non-capillary porosity (%), obtained by soil type; ${\mathrm{D}}_{\mathrm{i}}$ is thickness of aeration zone (mm), will be described later; ${\mathrm{H}}_{\mathrm{i}}$ is the depth of litter layer, generally is 50 mm. ${\mathrm{f}}_{\mathrm{i}}$, ${\mathrm{LAI}}_{\mathrm{i}}$, and ${\mathrm{h}}_{\mathrm{i}}$ can be determined by remote sensing data, and their physical meaning is clear.

For bare land, soil water storage capacity is the storage capacity between the field capacity and the wilting coefficient [14], which can be solved by soil and topographic index [68,69]:

$${\mathrm{W}}_{\mathrm{i}}=({\mathsf{\theta }}_{\mathrm{f},\mathrm{i}}-{\mathsf{\theta }}_{\mathrm{w},\mathrm{i}})\times {\mathrm{D}}_{\mathrm{i}},$$

(2)

where ${\mathrm{W}}_{\mathrm{i}}$ is the tension water storage capacity of the i-th grid in bare land area, ${\mathsf{\theta }}_{\mathrm{f},\mathrm{i}}$ is the field capacity (%), ${\mathsf{\theta }}_{\mathrm{w},\mathrm{i}}$ is the wilting coefficient (%), ${\mathrm{D}}_{\mathrm{i}}$ is thickness of aeration zone (mm). ${\mathsf{\theta }}_{\mathrm{f},\mathrm{i}}$ and ${\mathsf{\theta }}_{\mathrm{w},\mathrm{i}}$ can be determined according to the soil type and can be obtained by soil water characteristics in soil-plant-atmosphere-water [70].

However, there is a difference between the WSCC and the real full storage state of the basin [24,71]. The modified calculation way of the water storage capacity in bare land area is as follows, which should be added:

$${\mathrm{W}}_{\mathrm{mod}-\mathrm{i}}={\mathrm{N}}_{\mathrm{i}}\times \frac{{\mathrm{D}}_{\mathrm{i}}}{100}\times \frac{1}{10},$$

(3)

where ${\mathrm{W}}_{\mathrm{mod}-\mathrm{i}}$ is the modified water storage capacity of the i-th grid in bare land area (mm).

For reservoir region and impervious surface region, if many small reservoirs or dams exist in the basin, the storage needs to be considered separately [72,73]. The hardened surface and other buildings have no or extremely small water storage capacity because the thickness of the aeration zone is 0 [74], so the water storage capacity is 0, ${\mathrm{W}}_{\mathrm{i}}=0$.

Finally, for other regions, such as wetlands, lakes, and swamps, their water storage needs to be considered additionally if there are a large number of lakes, swamps, and wetlands in the basin. When calculating the water storage capacity of the basin, it can be adjusted according to the water storage data of swamps and lakes (interannual, monthly and seasonal). The change of wetland water storage capacity is small for a rainfall secondary flood, but its annual runoff changes greatly, which needs to be analyzed based on the wetland data over the years.

The thickness of aeration zone is critical for solving water storage capacity in Equations (1)–(3). At present, the thickness of aeration zone calculation is mostly based on the topographic index [68]:

$${\mathrm{D}}_{\mathrm{i}}=\overline{\mathrm{D}}-{\mathrm{S}}_{\mathrm{zm}}\times \left[{\mathrm{TI}}_{\mathrm{i}}-\mathsf{\lambda }\right],$$

(4)

where,${\mathrm{D}}_{\mathrm{i}}$ is the thickness of aeration zone on the i-th grid, $\mathrm{T}{\mathrm{I}}_{\mathrm{i}}$ is the topographic index on the i-th grid, $\mathsf{\lambda }$ is the mean value of topographic index, $\overline{\mathrm{D}}$ is the average saturated groundwater surface depth at the initial time, and ${\mathrm{S}}_{\mathrm{z}\mathrm{m}}$ is the maximum storage depth in unsaturated area.

$\overline{\mathrm{D}}$ and ${\mathrm{S}}_{\mathrm{z}\mathrm{m}}$ have a relationship:

$$\overline{\mathrm{D}}=-{\mathrm{S}}_{\mathrm{zm}}\mathrm{Ln}\frac{{\mathrm{Q}}_{\mathrm{b}}}{{\mathrm{Q}}_0},$$

(5)

${\mathrm{Q}}_{\mathrm{b}}$ is the base flow at the outlet of the basin at the time of drought, ${\mathrm{Q}}_0$ is the base flow at the outlet of the basin when the soil reaches saturation.${\mathrm{S}}_{\mathrm{z}\mathrm{m}}$ can be obtained by base flow (digital filtering method) or the mass curve technique (MCT) [63,75,76].

Thickness of aeration zone can be obtained by Equation (4), whereas the thickness of impervious surface and water is 0, accordingly it needs to be adjusted by land-use. The scale problem of topographic index in Equation (4) has been discussed [77], and fitting method and multiscale method have been proposed [57,78,79]. A combination method is proposed to minimize the effect of scale. First, the multi-resolution method of topographic index is expressed as:

$${\mathrm{TIM}}_{\mathrm{i}}=\frac{1}{\mathrm{n}}{{\displaystyle \sum }}_{\mathrm{j}=1}^{\mathrm{n}}\mathrm{Ln}\left(\frac{{\mathsf{\alpha }}_{\mathrm{j},\mathrm{i}}}{{\tan \mathsf{\beta }}_{\mathrm{j},\mathrm{i}}}\right),$$

(6)

$$\mathrm{M}=\frac{1}{\mathrm{k}}{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{k}}{\mathrm{TIM}}_{\mathrm{i}} ,$$

(7)

where $\mathrm{T}\mathrm{I}{\mathrm{M}}_{\mathrm{i}}$ is the topographic index under multi-resolution;${\mathsf{\alpha }}_{\mathrm{j},\mathrm{i}}$ is the single width catchment area of point i under the j-th scale, ${\mathsf{\beta }}_{\mathrm{j},\mathrm{i}}$ is the slope angle of point i under the j-th scale, $\mathsf{\lambda }\mathrm{M}$ is the mean value of multi-resolution topographic index, $\mathrm{n}$ is the scale dimension, which is generally between 3 and 5 depending on the detail of data and watershed conditions, and $\mathrm{k}$ is the number of watershed grids.

Although the multi-resolution topographic index obtained above reduces the scale dependence to a certain extent, yet it depends on the parameter $\mathrm{n}$. This paper selects a variety of fitting functions (Figure S4) and finally applied the generalized beta distribution of the second kind (best fits) to represent the watershed topographic index, which can be written as:

$$\mathrm{f}(\mathrm{TIM})=\frac{{\mathrm{r}}_3}{\mathrm{bB}\left({\mathrm{r}}_1,{\mathrm{r}}_2\right)}{\left(\frac{\mathrm{TIM}}{\mathrm{b}}\right)}^{{\mathrm{r}}_1{\mathrm{r}}_3-1}{\left(1+{\left(\frac{\mathrm{TIM}}{\mathrm{b}}\right)}^{{\mathrm{r}}_3}\right)}^{-\left({\mathrm{r}}_1+{\mathrm{r}}_2\right)},$$

(8)

where, $\mathrm{B}\left({\mathrm{r}}_1,{\mathrm{r}}_2\right)$ is Beta function, ${\mathrm{r}}_1$,${\mathrm{r}}_2$, and ${\mathrm{r}}_3$ are shape parameters and${\mathrm{r}}_1$ > 0, ${\mathrm{r}}_2$> 0, ${\mathrm{r}}_3$ > 0. The parameters ${\mathrm{r}}_1$ and ${\mathrm{r}}_2$ jointly determine the skewness of the distribution, and ${\mathrm{r}}_1$ controls the shape of the left tail, ${\mathrm{r}}_2$ controls the shape of the right tail, ${\mathrm{r}}_3$ determines the global shape and affects the peak value of probability density function, b is the position parameter and its value depends on the $\mathrm{TIM}$.

2.3. Data Preprocessing and Model Setting

According to the requirements of Equations (1)–(8), the data preprocessing for L-WSCC mainly includes the following steps: (1) Get the topographic index of the watershed; (2) obtain the thickness of aeration zone; (3) with soil parameters and vegetation parameters, the final distribution of water storage capacity under different land-use is calculated (Figure S5). Please see Text S3 for detailed processes.

To compare the performance of the WSCC and L-WSCC, two sets of simulations were implemented on hourly and daily scale. The PEST (parameter estimation) automatic calibration algorithm [80] was adopted for the calibration of the WSCC and other modules (evapotranspiration, runoff sources partition, and flow concentration) parameters based on Xinanjiang model. PEST is a model-independent nonlinear parameter estimation and uncertainty analysis program, which is a flexible and generic calibration tool [80]. It has the advantages of fast optimization speed and good robustness [81], which can be used in any of these hydrologic models and has been widely applied [82,83]. For the L-WSCC, the parameters were obtained by remote sensing data (mentioned in Section 2.3 and Section 2.4), the other modules parameter values were directly adopted from the calibrated WSCC. The parameters are listed in Table S3.

2.4. Evaluation Criteria

To evaluate the performance of WSCC and L-WSCC on daily scale, Nash–Sutcliffe efficiency (NSE), Kling–Gupta Efficiency (KGE), and percent bias error (BE) are calculated for the calibration (1982–1985) and validation (1986–1988) period. Daily NSE and BE values are compared with the values for basin-scale models recommended by related study, and model performance was classified as “very good”, “good”, “satisfactory”, and “unsatisfactory” [84] (Table S4). In addition, the appearance time of peak ($\mathrm{ARPT}$) and flood peak error (PE) are also used to evaluate model performance on the hourly scale. In order to study the magnitude of impact of different factors on runoff yield, the change of mean absolute error (ΔMAE) is also used. Please see Text S4 for more details.

3. Results and Discussion

3.1. Results on Land-Use Change

In order to fully demonstrate the degree of land type transfer of Misai watershed from 1982 to 1988, the rate of land-use change is used to conduct the research [85] and the results are shown in

Table 1. Rate of land-use change of the study area from 1982 to 1988.

	Year	1982	1983 (%)	1984 (%)	1985 (%)	1986 (%)	1987 (%)	1988 (%)
Land Use
Cultivated land		-	−18.41	26.00	2.53	−42.42	23.12	12.83
Forest land		-	2.31	−8.03	−0.01	16.00	−8.04	7.05
Grass land		-	753.00	−14.39	1.48	26.52	51.92	−99.42
Water		-	1.67	−0.96	−2.57	1.88	2.42	−0.73
Urban land		-	−38.15	110.17	−45.90	14.29	14.17	−45.83

. Due to the forest land and cultivated land accounting for about 95% of the watershed (Table S2), hence they play an essential role in hydrological processes, especially in runoff generation processes at Misai watershed. Therefore, we focus on changes in these two land-use types. Forest and cultivated land changed dramatically in 1986 (Figure S6 and Table 1), the rate of land-use change is 16.00% and −42.42%, respectively, which will have the greatest impact on the hydrological modelling (mentioned in Section 3.2, Section 3.3 and Section 3.4).

3.2. Calibration and Validation on Different Time-Scales

Comparing simulated and observed discharge at the outlet of the watershed, a close fit was observed with an R² value above 0.7 whether calibration or verification period (Figure 3a,b). Additionally, the simulation effects of WSCC and L-WSCC can well meet the simulation requirements depending on the NSE and KGE (Figure 3c,e) and simulated streamflow hydrographs matched the observed data well (Figure S7). No significant model over or underestimation was found, since the obtained parameter values are in the recommended range. Nevertheless, compared with WSCC, L-WSCC significantly improved the flood volumes: L-WSCC improved the BE by 3.10% and 12.02% during calibration and validation period, respectively (Table S5). Even worse, there was a sudden change in the simulation performance of WSCC in terms of BE in validation period (Figure 3d), especially, WSCC gets “unsatisfactory” result in 1988, which may be due to the fact that WSCC uses the initial calibration parameters whereas the land-use has changed strongly at this time (Table 1). However, the L-WSCC performs equally well indicating that integrated approach considers the impact of land-use on runoff yield, so it is more flexible (than WSCC) to handle differences in the rainfall and runoff in validation period. This corresponds to the purpose of our research.

In addition to the daily scale, 17 representative floods from 1982 to 1988 are chosen to evaluate the simulation performance of L-WSCC on rainstorm process (hourly scale).

The hydrographs of these floods are shown in Figure S8 and Figure 4a. The performance statistics are shown in Figure 4b. WSCC and L-WSCC resulted in at least satisfactory performance (NSE > 0.5) during the calibration and validation period (Figure 4b). The mean value of NSE were 0.78 and 0.82 of WSCC and L-WSCC considering the entire period, indicating that the WSCC and L-WSCC can reasonably predict the general shapes of the storm hydrographs in Misai watershed.

L-WSCC achieves better results compared to WSCC on the hourly scale with the NSE of L-WSCC increased from 0.76 to 0.85 during the validation period, PE decreased from values of approximately 21.56% to 10.10% and from 21.79% to 16.50% during the calibration and validation period, respectively (Table S6 and Figure 4b). When the BE is less than 20%, we think it is qualified. The BE qualified rate of WSCC is 70% in calibration period, and all qualified in the validation period, 80% and 86% in L-WSCC, respectively. The BE decreased from approximately 12.96% to 11.50% compared to WSCC during the calibration period and from 11.07% to 10.33% during the validation period (Figure 4b).

Despite the good performance in reproducing flood volumes, the PE was unsatisfactory (PE > 20%) for some floods (40.00% of the floods during the calibration period and 57.14% during the validation period of WSCC, while only 10.00% of the floods during the calibration period and 28.57% during the validation period of L-WSCC). This result indicated that WSCC has markedly poorer performance for the validation period, which experiences higher variability in runoff than the calibration period, suggesting that WSCC cannot be accepted in the validation period. Hence, the parameters of WSCC and L-WSCC calibrated were able to reproduce runoff generation process but failed to capture the behavior of confluence process at the basin scale. These bad simulations occurred after 1985 (64% in all bad events). This is because of the fact that the land-use in the watershed has changed dramatically in 1986 and affected the confluence of the watershed, resulting in the unreliable simulation results.

Figure 4a demonstrates the visual contrast of observed and simulated discharge time series to investigate the similarity of discharge in the watershed. The overall shape of the discharge hydrograph matches well with the observed discharge. However, the peak discharge is degraded during both runs. This degradation is mainly due to frequent negative updates (reduction) in the surface runoff, which is more affected by land-use change while L-WSCC weakens the impact of land-use on it. Hence, the simulated runoff series from the L-WSCC corresponds fairly well with the observed data (compared with WSCC), even though minor discrepancies still exist for the peak discharge simulation.

3.3. Seasonal Variation of Simulation

The observed and simulated discharge of every month have been analyzed and the mean absolute error (MAE) are shown in Figure 5. L-WSCC runoff simulation results are better than WSCC in the rainy season (May–July). The main reasons are: most of the remote sensing data used are from around June (Section 2.1). The water storage capacity calculation for runoff generation from May to July is more accurate. For other months, such as from September to April, the vegetation growth status is less than that in rainy season, which will influence the simulation results of L-WSCC. Therefore, the improvement effect of L-WSCC during this period is not as obvious as that in the rainy season (Figure 5). For WSCC, the runoff simulation results are better for the rainy season because WSCC is suitable for the simulation of large floods. Therefore, how to improve the simulation effect in dry season needs further research.

3.4. Simulation in Different Flood Types

To understand the simulation performance of WSCC and L-WSCC in different flood types, 21 floods are divided into three typical flood processes: big (B: peak discharge > 1000 m³/s), medium (M: 500 m³/s < peak discharge ≤ 1000 m³/s), and small (S: peak discharge ≤ 500 m³/s) floods according to the observed peak discharge. The simulated and observed values of WSCC and L-WSCC are shown in Figure S9. The discharges relationship in big floods is very good, but not so good in small and medium floods, especially small floods. This is noteworthy as the initial soil deficit has a greater impact on small floods and this is also consistent with the phenomenon that small floods are more easily disturbed by human activities. The NSE/KGE of the L-WSCC is greatly improved compared with the WSCC in small and medium floods, especially for small floods (Table S7 and Figure 6). However, there is very little difference between the WSCC and L-WSCC for big floods (Figure 6). This is expected because the watershed is close to saturation-excess state, and all rainfall produces runoff yield in big floods. Although there is a certain improvement in terms of BE, the flood volume of some events is obviously less (lower 20%), which may be caused by artificial interception. In summary, L-WSCC has a significant improvement of small and medium floods, but almost no improvement for large floods. This corresponds to the impact of land-use on saturation-excess runoff.

3.5. The Local Sensitivity Analysis of Different Driving Factors

To investigate the magnitude of impact of different factors on runoff yield, the sensitivity of model input parameters was evaluated through local sensitivity analysis by adjusting different factor parameters by ±10% to construct different scenarios: (a) land-use scenarios. The area of bare land and vegetation is increased by 10% (according to Section 2.3), and the scenarios are land-use scenario 1 and land-use scenario 2, respectively. (b) Soil scenarios. The soil parameters are increased and decreased by 10%, and the scenarios are soil scenario 1 and soil scenario 2 respectively. (c) Topography scenarios. The topography parameters are increased and decreased by 10%, and the scenarios are topography scenario 1 and topography scenario 2 respectively. Figure 7a shows the relationship between rainfall and runoff under different scenarios. An upward or downward shift in the relationship indicated that the relationship between rainfall and runoff has been altered by land-use, soil, and topography, correspondingly. However, the original relationship between rainfall and runoff is the same. As the parameter changes, the rainfall-runoff relationship changes in different directions (Figure 7a). Remarkably, the bare land (land-use scenario 1) has an enhancing effect on the runoff, while the vegetation (land-use scenario 2) has a significant inhibitory effect. Similarly, changes of soil parameters and topography parameters also cause changes in different trends of runoff (above and below the L-WSCC line, respectively). In addition, it should be noted that the degrees of influence for different factors on the runoff are different; land-use has the greatest, followed by topography, with the lowest for soil. Moreover, vegetation changes have the deepest impact on runoff (the farthest distance to the line of the L-WSCC). Moreover, the higher the rainfall, the greater the increase in the runoff difference, suggesting that the relationship between rainfall and runoff is more sensitive to different factors with increasing rainfall. Figure 7b shows the changes of MAE in different scenarios and we can also find that different land-use has a great influence on the runoff yield. When the bare land region is changed, the ΔMAE is about 20%, but when the vegetation region is changed, the ΔMAE changes to 113.97–173.14% with the average is 131.38%. For the soil scenarios, the ΔMAE is relatively small (about 20%), indicating that the runoff yield is less sensitive to soil factors. For the topography scenarios, although ΔMAE has changed, the variation range is between the land-use and soil factors (about 60%), demonstrating that the sensitivity of topography factors is between the land-use and soil factors. These results suggest that the land-use is the largest of the explanatory variable, revealing that land-use is the most significant driver for runoff yield.

3.6. Strength and Recommendations for Future Studies

This study systematically evaluated the integrated approach proposed by remote sensing and WSCC in the Misai watershed. Compared to WSCC, L-WSCC obtained more satisfactory results (Section 3.2, Section 3.3 and Section 3.4), revealing its potential in hydrological modeling considering human activities at watershed scale. These results were related to the difference in the runoff generation structure of the WSCC and L-WSCC, in which L-WSCC considered the underlying surface change on runoff yield. The results of local sensitivity analysis (Section 3.5) suggested that land-use plays an important role in runoff yield and is consistent with our study. The information and knowledge obtained from remote sensing data and WSCC were coupled together, facilitated addressing the problem and predicting the runoff yield, and provided a calculation method in ungauged basins [86,87]. In addition, runoff generation parameters of L-WSCC were obtained from remote sensing data without calibration, thus, L-WSCC can eliminate a source of parameter uncertainties. Taking this approach further, it is also applicable to the runoff generation part of other rainfall-runoff models based on saturation-excess runoff such as VIC [88], GR4 J [89], PDM [90], ARNO [91], HBV [92], G2G [54], unifying the distribution of WSCC for all these models.

Although study objectives generally are achieved, this study has limitations that should be outlined and perhaps targeted in forthcoming studies. First, the influence of climatic factors on runoff yield is not considered in L-WSCC. Despite this will not have an impact on the results of this study, it cannot be ignored if long-term runoff predictions are carried [93]. Second, although multi-resolution and fitting method are used to minimize the effect of topography scale, the accuracy of land-use, soil, vegetation, and other data has not been considered. The solution of the problem will have maximized the accuracy of our model. Furthermore, seasonal variation of simulation results indicates that the acquisition time of remote sensing data has a certain influence on the simulation results. Higher temporal resolution data for the improvement of simulation needs further study. Moreover, the study period is quite long ago and the basin hydrology may be substantially altered. Therefore, obtaining more recent data for comparative research is the direction of our further research. Although this study has obtained good results in the short-term hydrological simulation (1982–1988), the long-term hydrological process of the basin will inevitably be affected by climate change, so it is worth exploring how to take the climate factors into account for the long-term hydrological simulation. Notwithstanding these limitations, this study outlines a comprehensive runoff generation process for flood simulation based on remote sense and WSCC. L-WSCC considers the impact of human activities on runoff yield, enhances the physical meaning of WSCC, and can better understand the hydrological process. Extreme floods events are expected to occur more frequently, and the integrated approach will help relevant departments implement corresponding measures.

4. Conclusions

A new integrated approach (L-WSCC) combined the various remote sensing factors and saturation-excess runoff was proposed at basin scale, which helped to clarify the difference between the distribution of water storage capacity in reality and the modelled runoff. Unlike the traditional WSCC method, the extended L-WSCC consisted of equations of the effect of different factors on runoff, which provided beneficial insights in understanding the behaviors of the rainfall-runoff process. L-WSCC enhanced the physical meaning of the original WSCC, and simulating flood in regions characterized by physical basis rather than calibration/empirical approach. L-WSCC was evaluated from the Misai basin to be very effective in flood simulation and the results showed that the physics-based approach is promising. The major findings of this paper are summarized as follows:

(1)

The application of L-WSCC in Misai basin demonstrated that the method proposed in this article performed well in hydrological simulating at the different timescales (mean NSE: 0.86 for daily scale, 0.82 for hourly scale; mean BE: 4.34% for daily scale, 11.02% for hourly scale). The results have proven to be robust and reliable, which were powerful evidence for the reasonability of the construction process of L-WSCC.

(2)

The contrasting results confirmed that L-WSCC performed favorably compared to the empirical WSCC method and identified novel L-WSCC with high confidence in flood simulation (mean NSE increased from 0.78 to 0.82, mean PE decreased from 21.66% to 12.74%, mean RE decreased from 12.18% to 11.02%). These results indicated that the L-WSCC can derive a set of feasible values for the runoff generation.

(3)

The correction of hydrologic responses resulting from L-WSCC is more obvious for relatively small flood event (increased by 15.96% of small floods in terms of mean NSE). The proposed method can more accurately capture the runoff process of small flood events according to various remote sensing factors, thereby improving the ability of WSCC to simulate small floods. Thus, the newly proposed model may better represent runoff generation that previously have been beyond the scope of the traditional WSCC method (which is considered to work well for large flood events simulation).

(4)

The results of sensitivity analysis showed that land-use and land cover (especially the change of vegetation) have the greatest impact on yield (mean ΔMAE: 131.38%), which demonstrated that land-use have a predominant control on runoff yield in the short term.

The benefit of the spatially distributed method (L-WSCC) indeed mainly lies in the broader range of applications such as spatial scenario simulations and the ability to runoff simulation in similar basins. Our suggested approach can serve as an alternative to regionalization techniques, when one is interested in runoff simulation at ungauged internal stations of a particular scale. Moreover, this new spatial description of runoff provided by L-WSCC may allow for a more realistic extension of the WSCC method to runoff related processes such as water quality model and sediment migration model.