Using NDII pattern for a semi-distributed rainfall-runoff model in tropical nested catchments

A parsimonious semi-distributed rainfall-runoff model has been developed for flow prediction. In distribution, attention is paid to both timing of runoff and heterogeneity of moisture storage capacities within sub-catchments. This model 10 is based on the lumped FLEXL model structure, which has proven its value in a wide range of catchments. To test the value of distribution, the gauged Upper Ping catchment in Thailand has been divided into 10 sub-catchments, which can be grouped into 5 gauged sub-catchments where internal performance is evaluated. To test the effect of timing, firstly excess rainfall was calculated for each sub-catchment, using the model structure of FLEXL. The excess rainfall was then routed to its outlet using the lag time from storm to peak flow (TlagF) and the lag time of recharge from the root zone to the 15 groundwater (TlagS), as a function of catchment size. Subsequently, the Muskingum equation was used to route subcatchment runoff to the downstream sub-catchment, before adding to runoff of the downstream sub-catchment, with the delay time parameter of the Muskingum equation being a function of channel length. Other model parameters of this semidistributed FLEX-SD model were kept the same as in the calibrated FLEXL model of the entire Upper Ping basin, controlled by station P.1 located at the centre of Chiang Mai Province. The outcome of FLEX-SD was compared to: 1) observations at 20 P.1; 2) the results of the calibrated FLEXL model; and 3) the semi-distributed URBS model another established semidistributed rainfall-runoff model. FLEX-SD showed better performance than URBS, but a bit lower than the calibrated FLEXL model with NSE of 0.74, 0.71, and 0.76, respectively. Subsequently, at the level of the gauged internal subcatchments, runoff estimates of FLEX-SD were compared to observations and calibrated FLEXL model results. The results demonstrate that FLEX-SD provides more accurate runoff estimates at P.1, P.67 and P.75 stations which are located along 25 the main Ping River, compared to those provided by the lumped calibrated FLEXL model. The results were less good at 2 tributary stations (P.20 and P.21), where calibrated FLEXL output performed better, while performance was similar at one tributary station (P.4A). Overall, FLEX-SD performed better than URBS at 5 out of 6 stations except at P.21. Subsequently, the effect of distributing moisture storage capacity was tested. Since the FLEX-SD uses the same Sumax value the maximum moisture holding capacity of the root zone for all sub-catchments, FLEX-SD-NDII was set-up making use 30 of the spatial distribution of the NDII (the normalized difference infrared index). The readily available NDII appears to be a good proxy for moisture stress in the root zone, particularly during dry periods. The maximum moisture holding capacity in https://doi.org/10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c © Author(s) 2020. CC BY 4.0 License.

the root zone assumed to be a function of the maximum seasonal range of NDII values. The spatial distribution of this range among sub-catchments was used to calibrate the semi-distributed FLEX-SD-NDII model. The additional constraint by the NDII improved the performance of the model and the realism of the distribution. To test how well the model represents root zone soil moisture, the performance of the FLEX-SD-NDII model was compared to time series of the soil wetness index (SWI). The correlation between the root zone storage and the daily SWI appeared to be very good, even better than the 5 correlation with the NDII, because NDII does not provide good estimates during wet periods. The SWI, which is partly model-based, was not used for calibration, but appeared to be an appropriate index for verification.

Introduction
Runoff is one of the most important components of the hydrological cycle and can be monitored by the installation of a gauging station. Unfortunately, there are only a limited number of gauging stations available due to many complicating 10 factors such as topography, financial and human resources. Various rainfall-runoff models have been developed in gauged and ungauged catchments for runoff estimation at any required locations. Most rainfall-runoff models are categorised as lumped models, which can provide runoff estimates only at the site of calibration. These models include FLEXL, FLEX-Topo (Euser et al., 2015;Gao et al., 2014), NAM (Bao et al., 2011;Tingsanchali and Gautam, 2000;Vaitiekuniene, 2005; Yew Gan et al., 1997), SCS (Hawkins, 1990;Lewis et al., 2000;Mishra et al., 2005;Suresh Babu and Mishra, 2011;Yahya 15 et al., 2010), and many others.
To alleviate the limitation of lumped-rainfall-runoff models, URBS was developed as a semi-distributed nonlinear rainfall runoff routing model, which can account for the spatial and temporal variation in rainfall by separating a catchment into a series of sub-catchments (Mapiam and Sriwongsitanon, 2009). Therefore, URBS claims to provide runoff estimates not only at a gauging station but also at any required locations upstream (Carroll, 2004;Malone, 1999). URBS has been applied 20 successfully for real time flood forecasting in a range of catchments from small to very large basins in Australia and in many countries worldwide (Malone, 2006;Malone et al., 2003;Mapiam and Sriwongsitanon, 2009;Mapiam et al., 2014;Rodriguez et al., 2005;Sriwongsitanon, 2010). However, this model only addresses the distribution of travel times and does not address the effect of distributed storage capacities that affect the partitioning of moisture and hence the water balance.
Among the wide range of existing lumped-rainfall-runoff models, FLEXL has proven to be an adequate model for runoff 25 estimation in a wide range of catchments Fenicia et al., 2008;Gao et al., 2014;Kavetski and Fenicia, 2011;Tekleab et al., 2015). This model was further developed by Gao et al. (2016) to account for the spatial variability of landscape characteristics (FLEX-TOPO), useful for prediction in ungauged basins. Moreover, Sriwongsitanon et al. (2016) demonstrated that catchment-scale soil moisture content in the rootzone of vegetation computed from FLEXL is correlated with the remotely sensed Normalized Difference Infrared Index (NDII), which is widely used to monitor the equivalent 30 water thickness (EWT) in the root zone, for a number of gauging stations in the Upper Ping River Basin, especially during periods of moisture stress. https://doi.org /10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c Author(s) 2020. CC BY 4.0 License.
This study aims to utilize the fundamental model structure of FLEXL, include distributed time lags and channel routing as used in URBS, and include distributed root zone soil moisture capacity per sub-catchment so as to create a new parsimonious semi-distributed FLEX model for flood and flow monitoring within the sub-catchments of the gauged Upper Ping River Basin. Distribution of time lags is expected to improve hydrograph shape, particularly the timing and shape of the peaks, which would improve best-fit parameters, but it does not affect the partitioning of the hydrological fluxes or the water 5 balance. Since the root zone storage is the main control on flux partitioning, the distribution of the root zone moisture storage capacity would potentially have a larger impact on model performance. Therefore, the spatial variation of the NDII, as an indicator of root zone moisture stress, has been used to distribute moisture storage capacities among sub-catchments, while the model based SWI, as an estimator of moisture storage, was used for validation. The basin is dominated by well-forested, steep mountains in a generally north-south alignment (Sriwongsitanon and Taesombat, 2011). The areal average annual rainfall and runoff of the basin from 2003-2013 are 1,184 mm/yr and 277 15 mm/yr, respectively. The land use for the UPRB in 2013 can be classified into 6 main classes comprising forest, irrigated agriculture, rainfed agriculture, bare land, water body, and others, which cover approximately 77.40%, 3.11%, 12.54%, 1.99%, 1.23%, and 3.73% of the catchment area, respectively (Land Development Department, LDD). The landform of the UPRB varies from an undulating to a rolling terrain with steep hills at elevations of 1500-2000 m, and valleys of 330-500m (Mapiam and Sriwongsitanon, 2009;Sriwongsitanon, 2010). Chiang Dao district, north of Chiang Mai is the origin of the 20 Ping River, which flows downstream to the south to become the inflow of the Bhumibol Dama large dam with an active storage capacity of about 9.7 billion m 3 (Sriwongsitanon, 2010). The climate of the basin is dominated by tropical monsoons.
The southwest monsoon causes a rainy season between May and October and the northeast monsoon brings dry weather and low temperatures between November and April. Only 6,142 km 2 of the total area controlled by the runoff station P.1 (situated at the centre of Chiang Mai) is selected for this study (Fig. 1). 25

Rainfall data
Daily rainfall data from 48 non-automatic rain-gauge stations located within the UPRB and its surroundings from [2003][2004][2005][2006][2007][2008][2009][2010][2011][2012][2013] were used in this study. These data are owned and operated by the Thai Meteorological Department and the Royal Irrigation Department. These data have been validated for their accuracy on monthly basis using double mass curve and some inaccurate data were removed from the time series before spatially averaging using an inverse distance square (IDS) to be 30 applied as the forcing data of URBS, FLEXL, and FLEX-SD. https://doi.org/10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c Author(s) 2020. CC BY 4.0 License.

Runoff data
The Royal Irrigation Department (RID) operates 7 daily runoff stations in the study area between 2003 and 2013 as shown in Fig. 1. Catchment P.56A was rejected from the study because it is located upstream of Mae Ngat reservoir. Outflow data from the reservoir were used as input data in model calibration. Runoff data at the remaining 6 stations were used for the study since they are not affected by large reservoirs. The data have been checked for their accuracy by comparing them with 5 average rainfall data covering their catchment areas at the same periods. Table 1 presents the catchment characteristics and hydrological data for these 6 gauging stations in the UPRB. In this study, the catchments of these 6 stations were divided into 10 sub-catchments (see Fig. 1) with the areas ranging from 57 to 1,256 km 2 . High variation of catchment size is due to the proximity between the locations of these runoff stations and the outlets of the tributaries.

NDII Data 10
The Normalized Difference Infrared Index (NDII) is a ratio of the near-infrared (NIR) and shortwave infrared (SWIR) bands, centred at 859 and 1,640 nm, respectively, as shown in Eq. (1). In this study, the NDII was calculated using the MODIS level 3 surface reflectance product (MOD09A1), which is available at 500 m resolution in an 8-day composite of the gridded level 2 surface reflectance products. Atmospheric correction has been carried out to improve the accuracy and can be downloaded (1)

SWI Data
The near real-time Soil Water Index (SWI) is derived from the reprocessed Surface Soil Moisture (SSM) data derived from 20 the ASCAT sensor (Brocca et al., 2011;Paulik et al., 2014), which is a C-Band Scatterometer measuring at a frequency of 5.255 GHz in VV-polarisation (Paulik et al., 2014). The product makes use of a two-layer water balance model to describe the time series relationship between surface and profile soil moisture. This dataset of moisture conditions is available on a daily basis for eight characteristic time windows 1, 5, 10,15,20,40,60

FLEXL model
FLEXL is a lumped hydrological model comprising five reservoirs: a snow reservoir (Sw), an interception reservoir (Si), an unsaturated soil reservoir (Su), a fast-response reservoir (Sf), and a slow-response reservoir (Ss) (Gao et al., 2014). Excess rainfall from a snow reservoir, an interception reservoir, and an unsaturated soil reservoir is divided and routed into a fast-5 response reservoir and a slow-response reservoir using two lag functions. It includes the lag time from storm to peak flow (TlagF) and the lag time of recharge from the root zone to the groundwater (TlagS). Each reservoir has process equations that connect the fluxes entering or leaving the storage compartment to the storage in the reservoirs (so-called constitutive functions) . The water balance equations and constitutive equations for each conceptual reservoir are summarised in Fig. 2 and Table 2. The total number of model parameters is 11. Forcing data include daily 10 average rainfall and potential evaporation derived by the Penman-Monteith equation.

Snow reservoir
The snow routine, not very relevant in Thailand, can play an important role in areas with snow. When there is snow cover and the temperature (Ti) is above Tt, the effective precipitation is equal to the sum of rainfall (Pi) and snowmelt (Mi). The snowmelt (Mi) is calculated by the melted water per day per degree Celsius above Tt (FDD) (Eq. 2). The snow reservoir uses 15 the water balance equation, Eq. (3), where Swi (mm) is the storage of the snow reservoir.

Interception reservoir
Interception is more important in summer and autumn. The interception evaporation Eii was calculated by potential evaporation (Epi) and the storage in the interception reservoir (Sii), with a daily maximum storage capacity (Imax) (Eqs. 4, 5). The interception reservoir uses the water balance equation, Eq. (6), presented in Table 2. 20

Root zone reservoir
The root zone routine, which is the core of the hydrological models, determines the amount of runoff generation. In this study, we applied the widely used beta function of the Xinanjiang model (Ren-Jun, 1992) to compute the runoff coefficient for each time step as a function of the relative soil moisture. In Eq. (7), Cri indicates the runoff coefficient, Sui is the storage in the root zone reservoir, Sumax is the maximum moisture holding capacity in the root zone and β is the parameter 25 describing the spatial process heterogeneity of the runoff threshold in the catchment. In Eq. (8), Pei indicates the effective rainfall and snowmelt into the root zone routine; Rui represents the generated flow during rainfall events. In Eq. (9) Sui, Sumax and potential evaporation (Epi) were used to determine actual evaporation from the root zone Eai; Ce indicates the fraction of Sumax above which the actual evaporation is equal to potential evaporation, here set to 0.5 as previously https://doi.org/10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c Author(s) 2020. CC BY 4.0 License.
suggested by Savenije (1997) otherwise Eai is constrained by the water available in Sui. The unsaturated soil reservoir uses the water balance equation, Eq. (10), presented in Table 2.

Fast response reservoir
In Eq. (11), Rfi indicates the flow into the fast-response routine; D is a splitter to separate recharge from preferential flow.
Equations (12) and (13) were used to describe the lag time between storm and peak flow. Rft-i+1 is the generated fast runoff in 5 the unsaturated zone at time t -i + 1, TlagF is a parameter which represents the time lag between storm and fast runoff generation, clagF(i) is the weight of the flow in i -1 days before and Rfli is the discharge into the fast-response reservoir after convolution.
A linear-response reservoir, representing a linear relationship between storage and release, was applied to conceptualize the discharge from the surface runoff reservoir, fast response reservoirs and slow-response reservoirs. In Eq. (14), Qffi is the 10 surface runoff, with timescale Kff, active when the storage of the fast-response reservoir exceeds the threshold Sfmax. In Eq.
(15), Qfi represents the fast runoff; Sfi represents the storage state of the fast response reservoirs; Kf is the timescales of the fast runoff. The fast response reservoir uses the water balance equation, Eq. (16), presented in Table 2.

Slow response reservoir
In Eq. (17), Rsi indicates the recharge of the groundwater reservoir. Equations (18) and (19) were used to describe the lag 15 time of recharge from the root zone to the groundwater. Rst-i+1 is the generated slow runoff in the groundwater zone at time t -i + 1, TlagS is a parameter which represents the lag time of recharge from the root zone to the groundwater, clagS(i) is the weight of the flow in i -1 days before and Rsli is the discharge into the slow-response reservoir after convolution. In Eq. (20) , Qsi represents the slow runoff; Ssi represents the storage state of the groundwater reservoir; Ks is the timescales of the slow runoff. The slow response reservoir uses the water balance equation, Eq. (21), presented in Table 2. 20

URBS model
URBS was developed by Queensland Department of Natural Resources and Mines in 1990 based on the structures of RORB (Laurenson and Mein, 1990) and WBNM (Boyd et al., 1987). URBS is categorized as a semi-distributed rainfall-runoff model that can provide runoff estimates not only at the calibrated station but also at any required location upstream. The calibrated catchment area needs to be divided into sub-catchments to obtain different areal rainfall and different catchment 25 and channel travelling time. Table 3 presents 5 main processes used in URBS comprising the calculation of the initial loss, proportional loss, excess rainfall, catchment routing and channel routing. Excess rainfall is calculated separately between pervious and impervious areas. For the pervious area, URBS assumes that there is the maximum initial loss rate (ILmax) to be reached before any https://doi.org/10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c Author(s) 2020. CC BY 4.0 License. rainfall becoming the effective rainfall ( ). The initial loss (ILi) can be recovered when the rainfall rate (Ri) is less than the recovering loss rate (rlr) per time interval ( t) (see Eq. (22)).
Excess rainfall for each time step is calculated using Eq. (23) by weighting the excess rainfall between pervious and impervious area using a ratio of the cumulative infiltration (Fi) and the maximum infiltration capacity (Fmax). The recovering rate is included by simply reducing the amount infiltrated after every time step using the reduction coefficient (k t) as shown 5 in Eq. (24), and the pervious excess rainfall ( ) is calculated using the Eq. (25), where pr is the proportional runoff coefficient. The remaining water (1-pr) will infiltrate to the root zone storage (dFi) (see Eq. (26)). Excess rainfall is then routed to the centroid of any sub-catchment using a nonlinear reservoir relationship (Si = K ). The parameter m is the catchment non-linearity and K is the catchment travel time, which can be calculated for different sub-catchment using the multiplication between the catchment lag time coefficient (β) and square root of each sub-catchment area (A) (see Eq. (27)). 10 Thereafter, the outflow at the centroid of each sub-catchment is routed along a reach downstream of each sub-catchment using the Muskingum equation ( ℎ = Kch(XIi+(1-X)Qi)). The parameter X is the Muskingum coefficient and Kch is the channel travel time, which can be calculated for different sub-catchment using the multiplication between the channel lag coefficient (α) and the reach length (L) between the closest location in the channel to the centroid and the outlet of each subcatchment (see Eq. (28)). 15

Development of the semi-distributed FLEX model
The first step in distribution is to account for the timing of floods and the rooting of flood waves as a function of topographical factors. The resulting semi-distributed FLEX-SD model therefore is expected to better represent the shape of hydrographs, although it would not affect the partitioning of fluxes or the water balance. The root zone storage capacity is a 20 strong control on partitioning, affecting both runoff generation and evaporation. Therefore, distribution of this parameter would potentially affect overall model performance more strongly than merely the timing of the peaks. Therefore, in a second step, the NDII, as a proxy for moisture storage, is used to assess the distribution of moisture storage among subcatchments.

Accounting for distributed timing and channel-routing 25
FLEX-SD is set-up by applying lumped models for each sub-catchment, adding up to a semi-distributed model for a downstream calibration site. Therefore, the catchment area of any gauging station needs to be divided into sub-catchments.
Runoff estimates at each sub-catchment can be simulated using the structure of the original FLEXL by calculating different excess rainfall for each sub-catchment. The excess rainfall of each sub-catchment is routed to its outlet using the lag time from storm to peak flow (TlagF) and the lag time of recharge from the root zone to the groundwater (TlagS). In this study, 30 https://doi.org/10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c Author(s) 2020. CC BY 4.0 License.
TlagF and TlagS are calculated in hours instead of days to increase model performance. The lag time is distributed among sub-catchments using the following equations.
TlagS sub =TlagS√A sub A ⁄ where, Tlag is a lag time parameter for the entire catchment of a calibrated gauging station. The lag time of each sub-5 catchment (Tlagsub) is scaled by the square root of each sub-catchment area divided by the overall catchment area (A).
Runoff estimates from an upstream sub-catchment is later routed from its outlet to the outlet of a downstream sub-catchment using the Muskingum method (Eq. (31)) before adding to the runoff estimates of the downstream sub-catchment.
K sub = L sub (32)  10 where, α and X are the delay time parameter and the channel routing parameter for the entire catchment, respectively. The delay time parameter of each sub-catchment (Ksub) can be calculated by the multiplication between α and the main channel length of each sub-catchment as shown in Equation (32).

Accounting for distributed root zone storage at sub-catchment scale using NDII
The Normalized Difference Infrared Index (NDII) was used to estimate root zone storage capacity for each sub-catchment. 15 The NDII values, which are available at 8 day intervals, were found to correlate well with the 8-day average root zone moisture content (Su) simulated by FLEXL during the dry period in eight sub-catchments in the UPRB . The relation between NDII and Su can be described with a power function of the type: aNDII b +c, with c close to zero. The maximum value that Su can achieve is Sumax, the storage capacity of the root zone. The hypothesis is that the ecosystem creates sufficient storage to overcome a critical period of drought (Gao et al., 2014;Savenije and Hrachowitz, 20 2017). Every year has a maximum range of storage variation. If a sufficiently long NDII record is available, then the maximum of the annual ranges of the NDII should provide an estimate of the root zone storage capacity Sumax. Based on this assumption, the maximum difference between the annual maximum and minimum values of NDII b were used to calculate Sumaxi for each sub-catchment i (see Eq. (33) https://doi.org/10.5194/hess-2020-82 Preprint. Discussion started: 10 March 2020 c Author(s) 2020. CC BY 4.0 License.

Applications of URBS, FLEXL, FLEX-SD and FLEX-SD-NDII
URBS, FLEX-SD and FLEX-SD-NDII were calibrated at P.1 station located in the city of Chiang Mai. Since these models are semi-distributed rainfall-runoff models, they can provide runoff estimates in any required locations upstream of P.1 station, resulting in runoff estimates for P.4A, P.20, P.21, P.75 and P.67. As benchmarks for analysis, the calibrated FLEXL model was also used to estimate runoff at these 5 stations. The model parameters of the calibrated models were determined 5 using MOSCEM, which is a multi-objectives optimization algorithm. The evaluation criteria include the Nash-Sutcliffe Efficiency (NSE) and the Kling-Gupta Efficiencies for high flows, low flows, and the flow duration (KGEE, KGEL and KGEF), respectively.

Accuracy of runoff estimates simulated by URBS, FLEXL, FLEX-SD and FLEX-SD-NDII 10
The performance of the runoff estimates simulated by URBS, FLEXL, FLEX-SD and FLEX-SD-NDII is presented in Table   4. Figures 3, 4 and 5 present the output of the four models for a selection of stations (P.1, P.20, and P.67) compared to observations, as accumulated flows, hydrographs on logarithmic scale, and duration curves, respectively. In the appendix A the hydrographs of all 5 stations are presented in Figures A1, A2 and A3. Table 4 surprisingly shows that FLEX-SD and FLEX-SD-NDII provide better or equally good runoff estimates at P.1, P.67 15 and P.75 stations located along the main Ping River, compared to those provided by the calibrated FLEXL model. However, FLEXL shows better results in 2 out of 3 tributaries (P.20 and P.21) and similar results at P.4A. For a fair comparison, it should be realized that FLEXL was calibrated individually for each sub-catchment, while the other models were used in a predictive mode. In addition, overall, the FLEX-based models appear to perform better than URBS, especially with regard to low flow estimates where URBS showed very poor performance indicated by a very low KGEL compared to other models. 20 Looking more closely at three stations: P.1 on the main outfall, P.67 on the main stream, and P.20 as a major tributary, the following observations are made. affected by reservoir operation makes the calibration and validation difficult. However, close observation of Figure 4 indicates that FLEXL-SD-NDII follows the pattern of the observations more closely. In P.1 there is almost complete correspondence, especially in the dry season, where the other models are not capable of doing this. The strong constraint of the NDII on the root zone storage is the main reason for this.
The model parameters used in FLEXL, FLEX-SD and FLEX-SD-NDII are summarized in Table 5. The SD model provides 5 different values for TlagF (the time lag between storm and fast runoff generation), and TlagS (the lag time of recharge from the root one to the groundwater); the other parameters are kept the same as the calibrated values for P.1. Since TlagF and TlagS were designed to be related to the catchment area, the parameter values for each station are more reasonable compared to the values given by FLEXL. In the FLEX-SD-NDII model, in addition, the maximum moisture holding capacity in the root zone (Sumax) is varied between sub-catchments, making use of the spatial distribution of the NDII ranges. In general, 10 the FLEX-SD-NDII model provides lower Sumax estimates, constraining evaporation in the dry season (which provides more realistic recessions in Figure 4), but compensates for this reduction by a smaller  value, so as to limit excessive flood generation. Since these parameters jointly control Eq. (7), they can compensate for each other, leading to equifinality. If one of the parameters is constrained by additional information, as is the case here using the NDII, then this is no longer possible.
The performance with respect to best fit parameters may reduce in the process, but the model has gained realism and hence 15 predictive power.
In general, we see that FLEXL-SD-NDII shows the highest realism (illustrated clearly in Figure 4 and Figure A2 in the appendix A) but not a very good performance in the sub-catchment P.20, although still better than the other SD models. P.20 remains a difficult sub-catchment to predict due to its flow regulation and water consumption. Also we see that adding constraints to model calibration does not always improve best-fit performance, as compared to free calibration, but that 20 realism can be improved. To further test the realism of the models, in the following section the output of the models are compared to observations of NDII and the global scale SWI dataset for verification. Sriwongsitanon et al. (2016) suggested that NDII can be used as a proxy for soil moisture storage in hydrology. Therefore, the 8 day average NDII values were compared to the 8 day average root zone moisture storage (Sui) as calculated by FLEXL, 25 FLEX-SD and FLEX-SD-NDII. Table 6 shows the coefficients of the power relationships and the coefficients of determination (R 2 ) for the wet season, and the dry season for all sub-catchments. The table shows that the time series of NDII values correlate well with Su values during the dry season by giving R 2 value (average of all sub-catchments) of 0.71, 0.68, 0.77 for the 3 models respectively, with FLEX-SD-NDII performing best. During the wet season these correlations are much worse, resulting in average R 2 values of 0.36, 0.33, and 0.42 respectively, with FLEX-SD-NDII performing best. Fig.  30 moisture can only be observed by this remote sensing product when the vegetation is experiencing moisture stress. Hence correlations between root zone soil moisture and NDII are poor during the wet season.

The relationship between the average root zone soil moisture storage (Sui) and the average NDII and SWI
Because the FLEX-SD-NDII was constrained by the spatial variability of NDII ranges, the good correlation between Su and NDII during the dry season may not be surprising. Therefore, an additional test was done, testing the modelled Su values at daily time step with the daily SWI values, for all models. Figure 7 and Table 7 surprisingly show that the daily SWI values 5 correlate even better with the daily Sui values than the 8 day average NDII values for all sub-basins during the wet season with an R 2 of 0.75, although hysteresis effects are present. This result is much better than the wet season correlation on the basis of NDII in Table 6. One should realise, however, that the SWI is partly model based and that this may affect the good correspondence during the wet season. It can therefore be concluded that the NDII is a suitable parameter to constrain hydrological models during moisture recession, but that it works less well under wet conditions. The SWI, being partly 10 model-based, is less attractive as a model constraint, but does not suffer from a drawback during wet conditions, and hence serves well as an assessment criterion, particularly during wet conditions. As a result, the NDII appears to be useful to constrain hydrological models during dry conditions and both SWI and NDII appear to be useful to test model performance and to assess moisture states of river basins.

Conclusion 15
Most lumped rainfall-runoff models are controlled by a gauging station at the outfall on which it is calibrated. Runoff estimation at any location upstream requires indirect approaches such as model parameter transfer from gauged stations to ungauged locations, or applying relationships between model parameters and catchment characteristics to the ungauged locations. By using any of these approaches, uncertainty in runoff estimation for ungauged catchments is unavoidable. A semi-distributed hydrological model could offer a better alternative. FLEX-SD was introduced in this study based on the 20 model structure of FLEXL. The lag time from storm to peak flow and the lag time of recharge from the root zone to the groundwater have been distributed among sub-catchments using their catchment areas and reach lengths. Runoff routing between sub-catchment has been calculated using the Muskingum equation. FLEX-SD proved to be an effective model for runoff estimates not only at the calibrated stations but also at any location upstream, as compared to the original FLEXL, which required model calibration at each of the selected locations. FLEX-SD also proved to yield a better performance 25 compared to URBS, especially for low flow estimates.
However, FLEX-SD still used the same maximum moisture holding capacity in the root zone (Sumax) for all subcatchments. Therefore, FLEX-SD was further constrained by using NDII as a proxy for root zone moisture and to distribute the Sumax between sub-catchments, using the spatial distribution of annual NDII variation. We concluded that the maximum zone moisture storage (Sui) of all sub-basins controlled by runoff stations. The model parameters provided by the semidistributed FLEX models are more realistic compared to the original FLEXL since they are distributed according to catchment characteristics comprising catchment area, reach length, and remote sensing indices (NDII and SWI). A next step in the analysis is to account for diversity in landscape composition and related model structures among sub-catchments , which would allow for a distinction between the main rainfall-runoff mechanisms belonging to different 5 landscape types.   Snow Proportional Loss and Excess Rainfall