Evaluation of Soil Water Index of distributed Tank Model in a small basin with field data

The objective of the present study was to determine the spatial behaviour of the Soil Water Index (SWI) by applying a distributed version of the Tank Model (D-Tank Model) to the Araponga river basin (5.26 ha) in southern Brazil and to verify its reliability through the comparison to soil moisture estimated with the measured water-tension values and the water retention curve. The study area has a monitoring system for rainfall, discharge (5min interval), and soil-water tension (10-min interval). The simulation results showed that the D-Tank Model has a reliable performance. The correlation between SWI and HAND 5 was reasonable (r=0.6) meanwhile that between SWI and the Topographic Wetness Index was high (r=0.88). The comparison between the spatially distributed values of the SWI and soil moisture confirmed the high potential of the SWI derived from the D-Tank Model to be applied for predictions related to hydrological and environmental sciences.


10
Soil moisture is very important in the energy and water exchanges between pedosphere and atmosphere, which exert direct influence on the processes of infiltration, drainage, evapotranspiration, surface runoff, among others (Entin et al., 2000). Therefore, soil moisture is widely used as a variable in a lot of environmental, hydrological, meteorological and agricultural studies (Walker et al., 2004). It is well known that the spatio-temporal monitoring of soil moisture supports better management of water resources, prediction of floods and droughts, and so on. In spite of its importance, representing the spatial variation patterns 15 of soil moisture is challenging (Mälicke et al., 2019), and this variable is not regularly monitored due to the high cost (Brocca et al., 2017).
Similarly to the TWI, Rennó et al. (2008) found that the elevation difference to the nearest stream was correlated to soil water content distribution and proposed a new model, Heigh Above the Nearest Drainage (HAND). And both TWI and HAND are considered indicators of flood-prone areas (Nardi et al.,2006;Manfreda et al., 2015;Zheng et al.,2018;Tavares da Costa et al.,2019), and for landscape classification (Gharari et al.,2011;Loritz et al.,2019).
Physically-based distributed models have shown suitability for studying soil moisture and runoff response (Castillo et al.,2003;30 Loritz et al.,2019). Furthermore, some distributed hydrological models, such as WBMGB proposed by Saldanha et al. (2012) which divides the analyzed area into square cells of 10-km resolution, spatially calculate the water balance in the soil and consequently deal with soil moisture conditions. Sheikh et al. (2009) developed the BEACH model, a simple two-layer soil water balance model that provides spatially distributed initial soil moisture content. They emphasized the importance of soil moisture in hydrological modeling. Rahmati et al. (2018) also concluded that it is important to evaluate the basin physical 35 characteristics like soil moisture, in order to properly understand the hydrogeomorphic processes therein.
As soil moisture condition results from a water cycle, it is quite reasonable to treat it in the context of hydrological models simulating the water balance. One of the classic hydrological models is the Tank Model proposed by Sugawara (1961Sugawara ( , 1995, which is originally lumped and deterministic. This model was well evaluated and recommended by World Metereological Organization (WMO) (1975,1992) and Franchini and Pacciani (1991). Since it is computationally simple and generates good 40 results of hydrograph estimation, this model has been applied to various issues, for example, landslide prediction (Kobaishi and Suzuki,1987;Shuin et al.,2013;Nie et al.,2017); debris flow investigation (Takahashi and Nakagawa, 1991); flood forecasting (Tingsanchali and Gautam, 2000); Paddy field management (Chen et al., 2003); water quality (Kato, 2005); and sediment yield estimate (Lee and Singh, 2005), and also modified to its distributed version Huang et al.,2007).
The Tank Model was initially considered a storage model, composed of a series of vertical reservoirs (tanks), which schemat-45 ically represents the stratification of the soil layers from the surface to the base. The hypothesis based on this model structure is that the sum of the water storage calculated in these tanks would be a representative measure of the real condition of soil moisture in the basin. Then, the Soil Water Index and Tank Moisture Index derived from the Tank Model were proposed by Okada et al. (2001) and by Lindner and Kobiyama (2009), respectively. Thus, the use of the Tank Model as well as these indexes can be very useful in hydrological studies and, consequently, in natural disasters, water resources and basin managements (Saito 50 and Matsuyama,2012; Chen et al.,2013;Oku et al.,2014;Mukhlisin et al.,2015;Saito and Matsuyama,2015;Jun,2016;Chen et al.,2017;Matlan et al.,2018). However, the relationship between the water storage of each tank and the real soil moisture condition measured in a real basin has never been evaluated.
The objective of the present study was, therefore, to propose a distributed version of the Tank Model, to take into account the model and parameter uncertainties to compare the Soil Water Index (SWI) of Okada et al. (2001) with soil moisture indicators 55 (TWI and the HAND topology), and to compare the results with the soil moisture obtained in a small experimental basin.

Study Area
The study area is Araponga river basin (5.26 ha) which is a small experimental basin located in the rural area of Rio Negrinho municipality of the Santa Catarina state, southern Brazil. In this basin, the mean slope of the main channel and drainage density 60 are 0.30 m.m −1 and 9.62 ·10 3 · m −1 respectively . According to Alvares et al. (2013), the climate of the region is classified as Cfb -Humid subtropical, oceanic climate, without dry season, with temperate summer.
The basin is covered with the Mixed Ombrophilous Forest (Araucaria Forest or Brazilian Pine Forest) forming Montana (altitudes between 400 and 1000 m). The predominant soil in the region is the Cambisol, moderate and prominent with medium and clayey texture. (Santa Catarina, 1986). From field survey, Mota (2017) identified two soil layers and denominated the top one as A and the subsurface layer as B. The layer A is 0.2-m thick and characterized as a dark brown soil with 78% sand, 18% silt, and 4% clay. The layer B is 1.2-m thick and characterized as a dark yellow soil with 19% sand, 28% silt, and 53% clay).
Figure 1(a) shows the hypsometric map of the basin and, also the location of the rain gauge, the discharge weir and 9 tensiometer-batteries. The rainfall and discharge data were automatically measured in 5-min interval and, soil water tension data every 10 minutes.

70
According to Mota et al. (2017), the tensiometer sensors were installed by considering the soil depth and profile along the three hillslopes selected for the monitoring. In each hillslope, there were three points of tensiometer-batteries: near the river (higher slope with a mean soil depth of 0.5 m), at the middle hillslope (smaller slope with soil depths from 1.0 to 1.9 m), and near the basin divisor (smaller slope with soil depths from 1.0 to 1.9 m). Therefore, considering the variability of soil depth along each hillslope, there were adopted two configurations for the depths' distribution of the soil water tension sensors.

75
The type II configuration is characterized with the points very near the river (tensiometers A12, B12 and Z12), and the type I configuration with the other points (tensiometers A3, A4, B3, B4, Z3 and Z4). For the type II, 2 sets with 3 sensors were

Used Data
The topographic information of the basin was extracted from the digital elevation model (DEM) of the Araponga basin at 1:5000 scale of 1-m resolution that was generated through the topographic survey carried out by Mota et al. (2017).
The precipitation and discharge data used in the present study were obtained during the period from March 2011 to December 2015 as showed in Figure 2. The evapotranspiration was calculated from the measurements at the Feio meteorological station, 85 located 3 km distant from the Araponga basin, using the same procedure of Lindner and Kobiyama (2009) which utilized the Penman modified method of Doorenbos and Pruitt (1977).
By using the equation of Van Genuchten (1980), Mota (2017) estimated the parameters of the water retention curves for the two soil-layers of the study area from 12 undisturbed soil samples (6 at 10-cm depth and 6 at 50-cm depth) collected near  the Richards pressure plate until the drainage stabilized and then their water contents were determined. The RETC program (Van Genuchten et al., 1991) was used to fit the retention curves from pressure head and soil water content (Table 1). simulations to a certain basin area A (km 2 ), (Sugawara, 1995) suggested the following equation: Since A = 0.0524 km 2 in the present study, T is 0.034 h, approximately 2.1 min. As the field measurement was conducted every 5 minutes, the value of 5 minutes was adopted for T.
The hydrological and the water tension data series presented many periods of no-data caused by technical problems with the sensors. For this reason, it was not possible to use the entire series of data obtained by the monitoring. Thus, there were 100 selected 5 rainfall events for calibration and 2 for validation which all lasted 3 days. Table 2 shows the period, the mean, cumulative, measured and calculated values of precipitation and flow and the mean and cumulative values of evapotranspiration for the calibration events: I (occurred between January 14th to 16th, 2012), II (occurred between April 28th to 30th, 2012), III (occurred between June 4th to 6th, 2012), IV (occurred between June 20th to 22nd, 2013), V (occurred between March 8th to

Tank Model
According to Sugawara (1995), the number of tanks used in the modeling depends on the basin size, the soil type and use, and the time interval used in the simulation. In the case of large basins, and long-term simulations with daily data, there is a strong influence of the base flow, so a 4-tanks structure is used to represent the hydrological processes. In the case of small basins,

110
where short rainfall events will be evaluated, the greatest contribution to the flow will come from the surface runoff because the time response of the basin is too fast, and a few tanks are needed. For this reason, a 2-tanks structure was chosen for simulating the rainfall-discharge processes and estimating of soil moisture (Figure 3 (a)).
Then, the water balance equations were considered. Two versions of the Tank Model were applied in the study area. At first, the lumped Tank Model was created in order to obtain the parameters to be used in the distributed model. Secondly, apart from 115 the parameters encountered in the first step, the distributed Tank (D-Tank) Model was established.
The application of the lumped model was carried out by considering the structure composed of two tanks (reservoirs), and the equation for the flow calculation was obtained as follows: where S is the height of water stored in each tank (mm); qs2, qs1 and qs3 represent the surface, sub-surface and groundwater flow, respectively (mm/5min); qb1 is the percolation from the upper to the lower tank (mm/5min); a1, a2, and b1 are the coefficients of discharge; and a0 the coefficient of percolation; and HA1 and HA2 are the heights of the upper tank holes.
The present study calls the flow and percolation coefficients of all the holes, and the heights of the lateral and vertical holes of each tank as model parameters. These parameters´values are usually obtained through calibration. According to Ishihara 130 and Kobatake (1979), they are related to the soil type and use as well as to geological characteristics of the basin. For the application of D-Tank Model, the study area was divided into square cells, with 2-m resolution. In each of them, the lumped Tank Model with the equations (2) to (5) was applied.
The flow traveled through the cells depending on slope, where the upper tank of one cell passed to the upper tank of the next cell (Figure 3 (b)). Thus, each of the tanks cell received flow generated by neighboring cells, which calculation is as follows: where Q is the total flow received by the adjacent tank; q is the flow calculated by lateral orifice; i represents the tank (1: upper tank, 2 lower tank); j is the lateral hole; and c is the number of cells that send outflow directly to the considered cell. The water storage in each tank was equal to that used by Kato et al. (2005), and determined by the following equations: where P is the precipitation; Ct is the total number of covered cells; Q1 and Q2 are the total flow received by the upper and lower tank, respectively; and t is the time step.
The present work adopted the D-infinity method proposed by (Tarboton, 1997) to determine flow direction between the 145 neighboring cells. This method calculates the flow direction (represented by an angle between 0 and 2π rad), which is determined by the steepest direction of the eight triangular facets formed on a grid of 3x3 cells, centered on the pixel of interest.

Soil Water Index
As aforementioned, there are two indexes derived from the Tank Model which indicate the soil moisture condition in a basin: Soil Water Index -SWI (Okada et al., 2001), and the Tank Moisture Index -TMI (Lindner and Kobiyama, 2009). Both 150 were analyzed for their applicability to prediction and understanding of different natural disasters, the first one focusing on landslides (Okada et al., 2001), and the second one on flood and drought occurrences (Lindner and Kobiyama, 2009). Because of its simplicity, the SWI was adopted in the present study.
According to Okada et al. (2001), the SWI is defined as the sum of the storage heights (S1, S2) of the Tank Model, which indicates the soil moisture condition of the basin in the lumped version and of each cell in the distributed version, i.e.,

155
SW I = S1 + S2 2.5 Topographic Wetness Index and HAND Beven and Kirkby (1979) proposed the Topographic Wetness Index (TWI) to quantify topographic control on hydrological processes. This index is calculated by where α is the flow accumulation in a catchment area; and tanβ is the slope. In the present study, the flow accumulation was calculated by using D-infinity algorithm. The slope was calculated in degrees, and after transformed into radians.
The Heigh Above the Nearest Drainage (HAND) proposed by Rennó et al. (2008) was calculated using the software Ter-raViewHidro. The generated HAND topology was reclassified aiming to better represent the basin characteristics.

165
In order to perform an automatic calibration and uncertainty analysis of model parameters, the present work applied the automatic calibration algorithm Differential Evolution Adaptive Metropolis (DREAM(ZS)) proposed Laloy and Vrugt (2012) and Vrugt (2016). DREAM (ZS) is used for predict uncertainty in hydrological models (Cunha David et al., 2019), and is used for calibration of soil moisture distributed models (Linde and Vrugt, 2013). In the present study a generalized likelihood function (GL) proposed by Schoups and Vrugt (2010) was used for the inference of the hydrological model parameters. The residual model used a homoscedastic Gaussian likelihood distribution.
The last 7500 sets of parameters sampled with the DREAM (ZS) algorithm were used to represent the uncertainty associated with the parameter values and to create the probabilistic streamflow simulations. The present study evaluated the performance 175 of the models through two functions: the Nash-Sutcliffe efficiency coefficient (NASH) proposed by Nash and Sutcliffe (1970) and the root mean square error (RMSE): where Qobs(t) and Qsim(t) are the observed and simulated flows at the time t, respectively; and (Qobs) is the mean value 180 of the observed flow along the horizon t = 1 to T. The NASH values vary from −∞ to 1. The RMSE is always positive, and RMSE = 0 means the perfect adjustment. These two functions were used also for evaluating the model validation.

Application Procedure
The first step of the work was the MATLAB implementation of the Tank Model and D-Tank Model which is similar to that proposed by Kato et al. (2005). At this step, the algorithm was also implemented to calculate the flow direction in the study 185 area.
After, when the distributed version was implemented, it was necessary to divide the study area into square cells with 2-m side, because the available DEM had 1-m resolution, and in the present study it was adopted 2-m resolution. This procedure was At the following step, 5 rainfall-events were selected for the calibration of the Tank Model through the DREAM(ZS) algorithm. For each event, a set of parameters was chosen, and this set was also applied to each cell of the D-Tank Model.
This transfer of parameters was validated as follows: for each parameter found in the five events of the previous step, its average was calculated to obtain a single set of parameters to be used for the model validation with two rainfall events.
As aforementioned, the parameters used in the validation step were those corresponding to the mean values of the parameters 195 obtained in the calibration. Two periods were selected for validation of the distributed model, both with 3-days duration.
For the events of the validation period, the SWI values were compared through the linear correlation coefficient (r) with the calculated soil moisture values. Furthermore, the SWI values were also compared with TWI and HAND values.

200
The DREAM (ZS) was performed within a range pre-established by the user with the minimum and maximum limits of the decision variables which are the Tank parameters. Those limits are listed in Table 3. Figure 4 presents the predictive uncertainty for each event selected for calibration, which permits to understand how the uncertainty in the model parameters translates into Tank Model uncertainty. The model seems to be able to match the peak of the hydrographs, but the representation of the smaller flow values could be better. The events I, IV and V possess less 205 uncertainty than events II and III. The set of parameters selected after the automatic calibration for each of the 5 analyzed events is listed in Table 4. These values were used later by the D-Tank Model.  the study area has uniform land-use and that precipitation and evapotranspiration were also considered uniform over the whole basin area.
In general, the lumped model presented a reliable performance in relation to the observed hydrograph, both in magnitude and in time, especially at the peaks. There was an underestimation event occurred during the period June 4th to 6th, 2012 (Event III), in which a peak flow was observed from 0.43 mm/5min, and other occurred between April 28th and 30th, 2013 (Event II), 215  which had an observed peak flow of 0.095 mm/5min. In the other events, the maximum flows were close to 0.2 mm/5min, and the calibration was considered efficient in adjusting the hydrograph peaks.
In relation to the lower flows, the lumped model corresponded adequately to three events (I, II and III), and in the others, it underestimated the minimum flows.
The calculated flow by the distributed version also showed good correspondence, maintaining the good peaks adjustment.

220
There was a tendency to overestimate the flows of intermediate values ( Figure 5).

225
The Tank Model simulation resulted in a satisfactory water balance for the two events analyzed, the RMSE close to zero and NASH very close to 1 (Table 6).    In the similar way, the D-Tank Model showed a good performance for the event VI, with RMSE also close to zero, and NASH  These results show that the distributed model was adequate for the application of the following steps (SWI generation and point comparison with the soil moisture values estimated from the tensiometer data).  The spatial variation of the SWI values at three different moments (the beginning of the event, the peak, and the last moment) of the two validation events VI and VII were elaborated in Figure 9. In both events the SWI behavior was very similar, especially at the beginning. It resulted from the fact that the initial storage heights of the two reservoirs were the same in both events. In 245 the SWI maps related to the peak of the hydrograph, there is also a similarity between the values in the two events, which is consistent, since both had the same maximum flow rate in 0.15 mm/5min. In general, the points of the highest SWI values were located along the drainage network. However, high SWI values were also observed at certain points of hillslopes. Furthermore,

Soil Water Index
the SWI values were larger in the area closer to the basin exit, which is reasonable.
higher, reaching 251 mm, while most of the basin was showing its values of up to 193 mm. In the event VII, the maximum value was slightly smaller (214 mm), however the remaining area of the basin was wetter, represented by a SWI of 205 mm.
Between the two events, the spatial behavior of the SWI was similar, where wetter areas at any time kept their wetness under other circumstances (Figure 8). In relation to the upstream section, the SWI values in the section closer to the outlet were slightly larger, which is coherent.
255 Vachaud et al. (1985) observed this behavior of soil moisture and proposed the concept of temporal stability of soil moisture, i.e., there is a high probability that a moist condition for a moment will remain itself at other times. By using geostatistical tools, Gonçalves et al. (1999) observed the persistence in time of moisture distributions. This is an interesting result because the confirmation of the temporal stability of the soil moisture distribution allows reducing the sampling number for monitoring or estimating the soil water storage.

260
The spatial variation of the SWI can be observed also in Figure 10, where two cross sections were selected for SWI evaluation at the peak of the event VI. In the cross-section 1, there is a 60-m elevation difference between the highest point and the river, although these points are about 200-m distant.
The SWI value in this profile was almost constant, with a slight increase of 20 mm in the point located in the river. In the cross section 2, the SWI behavior in a profile closer to the outlet is presented. With a difference of 30 m of elevation between 265 the highest point and the river, but distant from one another only by 80 m, the SWI variation throughout the section was slightly larger, and, as expected, the highest SWI value occurred at the river.

Linear Correlation between SWI and Soil Moisture
For linear correlation analysis between SWI and the only measured values of soil moisture, the present study used only the tensiometers' data which were considered minimally coherent in relation to the precipitation of the event. Therefore, the 270 tensiometer of cell A12 was excluded from all the comparisons.
In the event VI, a strong correlation between the calculated SWI and the soil moisture in all the analyzed cells and depths.
They are all below 0.7 (Table 7). Figure 11 shows the linear relation between SWI and Soil Moisture calculated by the tension of to the sensor A3, tensiometer location at 10 cm of depth during the Event VI.
Significant correlations between the SWI and the soil moisture were also obtained for the event VII (Table 8). However, 275 in this event, those correlations were more significant in the uppermost layers of the soil. It can be said that all the analyses demonstrated that a SWI behavior is similar to that of the soil moisture.
TWI and HAND were generated ( Figure 12) to be used as soil moisture indicators. Both maps were compared to SWI results ( Figure 9) in order to verify a correlation between them. A visual analysis permits to verify the similarity of spatial behavior between TWI and SWI maps. The HAND topology also shows the similarity with TWI (the lower values are located in the 280 areas more close to the drainage, that are usually more moist). The results of the correlation analysis show the high correlation of SWI with the TWI values (in both events and in all 3 different hydrograph moments) (Table 9). HAND presented smaller values of r, but still presented significant correlation with SWI.
The correlation between SWI and TWI could be explained for the structure of the D-Tank Model, where similarly as the 285 TWI, the accumulated flow has a determinant influence on the generated results. Thus, it can be said that SWI has a major relation with flow accumulation and slope than the vertical distance to the stream.

Linear Correlation of Soil-Moisture Spatial Distribution with SWI, TWI and HAND
For spatial analysis of linear correlation between the Soil Moisture and SWI, TWI and HAND, the present study used only the soil moisture data estimated with the tensiometers located at 30-cm depth, because of the two reasons: (i) these were the 290 sensors that generated more coherent results at all the analyzed events; and (ii) Ehret (2014) reported that the mean moisture of the first 30-cm soil layer was a robust indicator of the slopes initial conditions prior to a rainfall event.
Several studies used geostatistics for estimate soil moisture patterns from point observations (Bardossy and Lehmann,1998;Western et al.,1998;Western and Bloschl,1999;Perry and Niemann,2008;Yang et al.,2017). The Inverse Distance Weighting (IDW) method was used for interpolating the soil moisture values calculated by tensiometer information, 295 and the distributed soil moisture maps were generated. However, in this case, only the covered area by the sensors was treated because of the limitation of this method ( Figure 13).
The different patterns of spatial distribution of soil moisture can be observed between two analyzed events (Event VI and Event VII). Though the total values of rainfall of the Event VI and VII were similar (104.59 mm and 139.12 mm, respectively), the values of the initial soil water content were slightly different (0.55 cm 3 /cm 3 in Event VI and 0.41 cm 3 /cm 3 in Event VII), which might result in the differences in the peak flow between these events (0.15mm in Event VI and 0.09mm in Event VII). Zehe et al. (2007) and Uber et al. (2018) found significant effects of initial soil moisture condition on discharge values.  Furthermore, Zehe and Sivapalan (2009) demonstrated how the initial soil moisture condition influenced on the infiltration process and consequent runoff generation.
The distribution features of the soil moisture had high correlations with SWI, TWI, and HAND (Table 10).   According to , soil moisture patterns depend on different topographic properties at different times.The authors found that TWI can explain soil moisture variability at catchment scale. Minet et al. (2011) analyzed soil moisture patterns in terms of runoff response, and concluded that the use of TWI for evaluating the soil moisture pattern is a suitable technique. Observing Table 10 more in detail, it can be said that the SWI generated by D-Tank Model performed slightly better than TWI. Hence, the use of SWI could be more recommended for the soil moisture representation. Loritz et al. (2019) dis-310 cussed hydrologic similarity in different catchments comparing TWI and HAND performance, and concluded that despite the similarities between these indices, they represent different hydrological aspects, where the modified-HAND version performed better than TWI and HAND approaches. Contrary to their study, the present study showed that HAND was the worst index to represent soil moisture patters. Downer and Ogden (2003) also used a distributed conceptual model to estimate soil moisture, concluding that physically-315 based hydrologic models can be used to make predictions of peak discharge and soil moisture. The present study can infer that a distributed conceptual hydrological model can be used to predict soil moisture patterns in catchment scale.

Conclusions
By using the D-infinity method of (Tarboton, 1997), the present study constructed a distributed version of Tank Model and called it D-Tank Model. Adopting the SWI proposed by Okada et al. (2001), the D-Tank Model was applied in order to 320 evaluate the soil moisture distribution in the Araponga river basin in the rural area of Rio Negrinho municipality of the Santa Catarina state, southern Brazil. Then the performance of the D-Tank Model and the spatially-distributed SWI values were verified with rainfall, discharge and soil moisture data.
Based on the obtained results, it can be concluded that the D-Tank Model had satisfactory performance for the studied period. Hydrographs between the observed and calculated flow showed a good fitting. The spatially distributed values of 325 the SWI, generated from the D-Tank Model, performed also satisfactorily, representing well the spatial variation of the water storage within the basin, which were qualitatively confirmed with the soil moisture data. Even compared with TWI and HAND, the SWI values presented better fit with the values of linear correlation between the analyzed maps.
The spatially-distribution analysis of the soil moisture was used to validate SWI, which shows good correlation. Thus, the use of SWI could be recommended for the soil moisture representation. However, the D-Tank Model and its derived SWI 330 values should be checked with other basins which are different from the Araponga basin in terms of the basin size, forest type, land-use, among others.
Code and data availability. The full analysis scripts are published on Github (https://github.com/mvsofia/dtankmodel). The data is available upon request.
Author contributions. The methodology was developed by SM, and supervised by MK. The data was provided by AM. All code was devel-335 oped by SM. The manuscript was written by SM, with contributions of MK in the Introduction and Discussions, and AM in the Introduction.