Applying transfer function-noise modelling to characterize soil moisture dynamics: a data-driven approach using remote sensing data

Abstract The increasing availability of remotely sensed soil moisture data offers new opportunities for data-driven modelling approaches as alternatives for process-based modelling. This study presents the applicability of transfer function-noise (TFN) modelling for predicting unsaturated zone conditions. The TFN models are calibrated using SMAP L3 Enhanced surface soil moisture data. We found that soil moisture conditions are accurately represented by TFN models when exponential functions are used to define impulse-response functions. A sensitivity analysis showed the importance of using a calibrated period which is representative of the hydrological conditions for which the TFN model will be applied. The IR function parameters provide valuable information on water system characteristics, such as the total response and the response times of soil moisture to precipitation and evapotranspiration. Finally, we encourage exploring the possibilities of TFN soil moisture modelling, as predicting soil moisture conditions is promising for operational settings.

In general, three methods exist for estimating soil moisture at various spatiotemporal scales: in situ (Dobriyal et al., 2012;Susha Lekshmi et al., 2014), remote sensing (Fang and Lakshmi, 2014;Petropoulos et al., 2015;Zhuo and Han, 2016) and hydrological modelling (Vischel et al., 2008;Zhuo and Han, 2016). In situ soil moisture sensors provide accurate information at local scales, since soil-specific calibration procedures can be performed. However, in situ sensors typically have limited spatial coverage (Susha Lekshmi et al., 2014). Remote sensing and hydrological modelling are alternative sources for providing spatially distributed soil moisture information at larger scales. Remotely sensed soil moisture information is often retrieved using active and passive microwave sensors (Petropoulos et al., 2015). The temporal coverage of remote sensing is limited in comparison with in situ methods, as satellite imagery is only available during satellite overpasses. In addition, only surface soil moisture can be retrieved by remote sensing due to sensor capabilities (Zhuo and Han, 2016). Furthermore, vegetation dynamics, surface roughness, and satellite sensor uncertainty significantly affect remote sensing retrievals (Petropoulos et al., 2015;Benninga et al., 2019).
Hydrological modelling provides a means to estimate soil moisture at various spatiotemporal scales (Vereecken et al., 2016;Brocca et al., 2017). The complexity of unsaturated zone models ranges from simple conceptual lumped models to complex integrated physically-based distributed models. Water managers often regard model accuracy as a limiting factor for the application of hydrological modelling for decision-making in operational water management (Pezij et al., 2019a). The accuracy of hydrological models is partly based on which dataset is used for calibration purposes. Data assimilation schemes can be applied to update model simulations with available observations, although such schemes often require significant computational power (Liu et al., 2012;Weerts et al., 2014;Pezij et al., 2019b). Furthermore, physically-based unsaturated zone models are often based on the Richards equation, which is highly non-linear and therefore poses challenges for numerical solutions ( � Simůnek et al., 2003;Vereecken et al., 2016). In addition, Richards-based models (e.g. SWAP or Hydrus) are generally developed for local applications (Van Dam et al., 2008;� Simůnek and van Genuchten, 2008). Significant computational power is required to apply such models at regional scales (Van Walsum and Groenendijk, 2008). Therefore, the application of Richards-type soil water flow models is not trivial in operational water management at catchment scales.
Data-driven modelling methods are suitable alternatives for processbased modelling (Todini, 2007;Solomatine and Ostfeld, 2008), especially when large amounts of data are available, such as in the Netherlands. Data-driven methods for the prediction of soil moisture conditions are promising (Kolassa et al., 2018;Cai et al., 2019). The availability of new high-resolution remotely sensed soil moisture data offers new opportunities for data-driven modelling methods (Petropoulos et al., 2015). However, data-driven methods, such as neural networks and random forests, are often considered black box models. It can be challenging to interpret and apply the results of data-driven methods (Padarian et al., 2020). Therefore, although suitable for predictions, such data-driven methods are limitedly suitable to increase our understanding of hydrological processes and characteristics.
In this study, we present a data-driven framework based on transfer function-noise (TFN) modelling for describing soil moisture dynamics and its controls. TFN modelling is a data-driven method to model an observed time series by applying a linear transformation of one or more deterministic input series known as stress series (Von Asmuth et al., 2002). The stress series are transformed using one or more impulse-response (IR) functions. The strengths of the framework are: (1) TFN modelling is a fast and easy-to-construct data-driven alternative for complex process-based models. (2) The IR functions are set up using only observational data. TFN modelling does not need prior assumptions on system characteristics, which is an interesting property since no model structure is expected to work best everywhere (Peterson and Western, 2014). (3) The IR functions contain information on the response of a water system to input stresses such as precipitation and therefore increase our understanding of hydrological characteristics and processes. (4) Due to the stochastic nature of the included noise model, TFN models can model system dynamics which are not well explained by physical laws (Von Asmuth et al., 2002).
The applicability of TFN modelling for groundwater studies has been shown extensively (Yi and Lee, 2004;Bakker et al., 2007;Manzione et al., 2010;Fabbri et al., 2011;Obergfell et al., 2013;Sutanudjaja et al., 2013;Peterson and Western, 2014;Zaadnoordijk et al., 2018;Bakker and Schaars, 2019). These studies show that TFN modelling can be used to describe groundwater dynamics. In addition, the IR functions contain valuable information on groundwater system characteristics, such as response times. A TFN model is suitable for short-term soil moisture predictions using short-term meteorological predictions. Also, depending on the time period covered by the meteorological input data, one could use a TFN model to construct historical surface soil moisture time series. Commonly applied TFN modelling tools for groundwater modelling in the Netherlands are Menyanthes (Von Asmuth et al., 2012) and Pastas (Collenteur et al., 2019a).
We are interested in the applicability of TFN modelling as an innovative data-driven method for explicitly calculating soil moisture dynamics, which has not been studied before. Our main focus period is the year 2018. The soil moisture drought caused by the 2018 European heat wave significantly impacted water management and agricultural practices (Vogel et al., 2019), which shows the importance of retrieving up-to-date soil moisture information. We want to show that TFN modelling is a fast and relatively simple method for soil moisture modelling. Also, as data-driven models are often black boxes, we want to show that TFN models can help in identifying controls on soil moisture dynamics. We address the following research question in this study: to what extent can TFN modelling describe and predict soil moisture dynamics using remotely sensed soil moisture information? This paper is organised as follows: Section 2 describes the research methodology. Section 3 presents the results, which are discussed in Section 4. Finally, conclusions are drawn in Section 5.

Transfer function-noise modelling
The time series modelling approach applied in this study originates from TFN Autoregressive Moving Average (ARMA) modelling. TFN models are often applied in hydrological applications, since they are fast and provide accurate predictions (Yi and Lee, 2004;Von Asmuth et al., 2008;Manzione et al., 2010;Fabbri et al., 2011;Peterson and Western, 2014). The main concept is to interpret the output of an observed system as a combination of analytically defined impulse responses of a linear system, which can be described by transfer functions (Kennaugh and Moffatt, 1965;Box and Jenkins, 1970). Among others, Marco (1993) and Remesan and Mathew (2015) discuss the applicability of transfer functions in hydrology, such as in the Unit Hydrograph approach applied for estimating river discharge rates (Sherman, 1932). The use of transfer functions for modelling subsurface dynamics has been studied by for example Besbes and De Marsily (1984), Gehrels et al. (1994), and Von Asmuth et al. (2002). Ramirez-Beltran et al. (2008) showed the applicability of transfer functions for soil moisture modelling using in situ soil moisture measurements as validation data.
A TFN model has a deterministic component for the linear transformation of observations to a change of the output variable and an exponential noise model comparable to an auto-regressive lag 1 (AR1) model (Peterson and Western, 2014;Remesan and Mathew, 2015). The linear transformation is performed using transfer functions. The transfer functions, or IR functions, are not known a priori and have to be estimated. The IR functions can be estimated using the Predefined Impulse Response Function In Continuous Time (PIRFICT) method (Von Asmuth et al., 2002;Von Asmuth and Bierkens, 2005). The PIRFICT method assumes that IR functions have known analytical expressions.
Von Asmuth et al. (2002) show that the PIRFICT method overcomes the following limitations of estimating IR functions in regular TFN ARMA models: (1) PIRFICT allows the use of data with an irregular time interval, and (2) the iterative Box-Jenkins style model identification procedure is generally applied to define the number of parameters in ARMA models (Box and Jenkins, 1970). The Box-Jenkins model identification procedure can be knowledge-and labour intensive. Since the analytical expression of the IR functions has to be defined a priori, the order of the IR functions does not have to be defined using a Box-Jenkins model identification procedure (Von Asmuth et al., 2012). The PIRFICT method allows selecting a typical IR function for a specific type of series . The PIRFICT method is applied in this study.
A general form of a continuous-time TFN model, in this case formulated for soil moisture dynamics, is (Von : where hðt j Þ is the observed soil moisture state at time step t j [m 3 m À 3 ], N stress is the number of stress series which influence the soil moisture state [-], h i ðt j Þ is the change in the soil moisture state due to stress series i at time t j [m 3 m À 3 ], d is the baseline soil moisture state [m 3 m À 3 ], and n res ðt j Þ is a residual time series [m 3 m À 3 ]. The subscript j indicates the day number. A TFN model can have multiple stress series as input. The contribution of stress series i to the soil moisture state (h i ðt j Þ) is determined by solving a convolution integral in continuous time using an IR function, which describes the variation of the soil moisture state due to an individual stress series. h i ðt j Þ is defined as: where R i is the value of a stress series i [mm] at the times up to time step t j and Θ i ðt j Þ is the IR function of the corresponding stress series i. To solve Equation (2), Θ i ðt j Þ should be known. The type and shape of the functions depend on the type of stress and water system characteristics. We assume that the main drivers of soil moisture dynamics are precipitation (Entekhabi and Rodriguez-Iturbe, 1994;Vereecken et al., 2016) and evapotranspiration (Syed et al., 2004). Therefore, we use time series of precipitation and evapotranspiration as stress series for the TFN model. Von Asmuth et al. (2012) state that, independently of system properties, analytical expressions such as the scaled gamma function fit the behaviour of many hydrogeological systems. The scaled gamma step response function is a commonly applied IR function for precipitation and evapotranspiration stress series in groundwater TFN modelling (Besbes and De Marsily, 1984;Von Asmuth et al., 2012). The scaled gamma IR function is defined as: where A gam corresponds to the unit response of the soil moisture state due to the stress at time, t j ¼ 0 [m 3 m À 3 ], a gam is a shape parameter [day], n is a shape parameter [-], and ΓðnÞ is the gamma function of the form ðn À 1Þ! [-]. Typical behaviour of the gamma block response function is visualized in Fig. 1.
In addition to the gamma function, an exponential function is used. Particularly, we studied whether an exponential function is a better representation of the unit response of soil moisture to precipitation and evapotranspiration, as we expect that soil moisture shows a fast response to precipitation and evapotranspiration. The exponential IR function is defined as: where A corresponds to the unit response of the soil moisture state due to the stress at time t j ¼ 0 [m 3 m À 3 ] and a is a shape parameter [day]. Fig. 1 also shows typical behaviour of the exponential block response function. The parameters of the IR functions are unique for every location which is analysed . The output of a deterministic model will never match observations perfectly. Residuals of TFN models often show large autocorrelation (Collenteur et al., 2019a). Modelling the residuals helps to provide more accurate predictions of the soil moisture state at unobserved time steps. Therefore, the residuals n res ðt j Þ, defined in equation (1), are modelled using a noise model to satisfy a white noise requirement (Von Asmuth and Bierkens, 2005). The residuals are assumed white noise if all important system stresses are considered. The noise model, with exponential decay of the residuals, is formulated as: where υðt j Þ is white noise resulting from a random process for time step t j [m 3 m À 3 ], α is a decay parameter [day], n res ðt jÀ 1 Þ is the residual value at the previous time step [m 3 m À 3 ], and Δt j is the length of the time step [day], which is one day in this study. In addition, the noise model allows the application of a least squares objective function (Von Asmuth and Bierkens, 2005;Peterson and Western, 2014). The optimal parameter sets for the IR functions used in equations (1) and (2) are determined by minimizing the objective function. The A gam , a gam , and n parameters are estimated when the gamma function is used as IR function. The A and a parameters are estimated when the exponential function is used as IR function. More information on the parameter estimation procedure and the noise model can be found in Von Asmuth and Bierkens (2005).

TFN modelling library: pastas
To set up the TFN models, we use the open source library Pastas (Collenteur et al., 2019a, b), which is a Python 3 implementation of the TFN PIRFICT modelling approach described in section 2.1. More information on the Pastas library can be found at https://pastas.readthedocs. io.

Study area and data
We assessed the applicability of TFN modelling for soil moisture predictions in the Twente region in the eastern part of the Netherlands, see Fig. 2. The Twente area has been part of soil moisture studies (e.g. Dente et al., 2012;Benninga et al., 2019;Chen et al., 2019;Pezij et al., 2019a; Van der Velde et al., 2019) We will compare the results of the TFN modelling approach with the SMAP (Soil Moisture Active Passive) soil moisture retrieval studies of Colliander et al. (2017), Chan et al. (2018), and Kolassa et al. (2018), who included the Twente area in their assessment. The study area is situated in a temperate marine climate zone (Hendriks et al., 2014), has an elevation ranging between 3 and 85 m.a.s.l. and has an extent of approximately 40 km by 50 km. The main soil types are sandy and loamy sandy (W€ osten et al., 2013). The primary land use is agriculture. Annual precipitation varies between 800 and 850 mm (Kaandorp et al., 2018). Table 1 provides an overview of the data products used in this study. We use two types of stress series: precipitation and Makkink reference crop evapotranspiration. Makkink reference crop evapotranspiration (in the following evapotranspiration) describes the potential evapotranspiration from a reference surface covered with grass (Makkink, 1957). The Makkink method requires less input variables than the Penman-Monteith method. De Bruin and Lablans (1998) show that the reference crop evapotranspiration estimates obtained using the Makkink method correspond well to the reference crop evapotranspiration estimates obtained using the Penman-Monteith method for the Netherlands. Open source precipitation and evapotranspiration data from the Royal Netherlands Meteorological Institute (KNMI) are used (KNMI, 2018a, b). The precipitation data are based on radar data which are corrected using KNMI station data by applying ordinary kriging. The precipitation data have a spatial resolution of 1 km by 1 km. The evapotranspiration data are based on extrapolating KNMI station data using Thin Plate Spline (TPS) interpolation. The station data are calculated by KNMI using incoming shortwave radiation and mean daily temperature measurements at the KNMI stations. Hiemstra and Sluiter (2011) suggest that at least 15 stations should be used to interpolate daily evapotranspiration at a national scale. In total 32 stations are used for the interpolation in the Netherlands. Of these 32 stations, three stations are located in or near the Twente study area (Twenthe, Heino, and Hupsel). The evapotranspiration data are gridded at a spatial resolution of 10 km by 10 km.
We use the SMAP L3 Enhanced radiometer-only daily gridded soil moisture product to calibrate the parameters of the TFN models (Entekhabi et al., 2010;Chan et al., 2018;O'Neill et al., 2018). Also, we use the SMAP L3 Enhanced product to validate the TFN models. The SMAP Enhanced products are developed by applying the Backus-Gilbert optimal interpolation technique to the SMAP brightness temperature data (T A ) . The latter is posted on a grid with a spatial resolution of 36 km by 36 km (Chan et al., 2016). After various correction and calibration procedures of the regridded brightness temperatures, the currently operational SMAP baseline soil moisture algorithm is used to produce the SMAP L3 Enhanced surface soil moisture estimates. The SMAP L3 Enhanced product is gridded with a spatial resolution of 9 km by 9 km.
In general, SMAP soil moisture products perform well in the Twente region (Colliander et al., 2017). Chan et al. (2018) assessed the SMAP L3 Enhanced product using in situ soil moisture measurements. They found an unbiased root mean square error (uRMSE) of 0.056 m 3 m À 3 (using the SCA-V retrieval algorithm) for the Twente region. The SMAP data products were designed to meet a soil moisture retrieval accuracy of 0.040 m 3 m À 3 (uRMSE) (Entekhabi et al., 2010). The SMAP observations are available for the study area approximately every 2-3 days. The footprint of the SMAP L3 Enhanced product in the Twente region is visualized in Fig. 2. We have analysed the morning SMAP retrievals for the period January 1, 2016-January 1, 2019. The SMAP soil moisture retrievals are affected by, for example, sensor noise and uncertainties in surface roughness and vegetation changes. By applying the noise model with exponential decay of the residuals (equation (5)), we partly include these uncertainties in the TFN modelling approach.
Additionally, we use in situ soil moisture measurements from a monitoring network to assess whether the TFN models can describe soil moisture field conditions. The monitoring network, operating since 2009, is maintained by the ITC faculty of the University of Twente (Dente et al., 2012; Van der Velde, 2018; Van der Velde et al., 2019). Both volumetric moisture content and soil temperature are measured at 20 locations in the Twente region. Stations 1,2,3,4,5,11,12,13,14,15,16,17,18, and 19 cover agricultural grass fields. Stations 6, 7, 8, 9, and 10 cover maize fields. Station 20 is installed in a forest area. The station locations are shown in Fig. 2. The monitoring network consists of Decagon 5TM probes at 5, 10, 20, 40, and 80 cm soil depths and provide a reading every 15 min. We use daily averaged measurements at 5 cm soil depth, since remotely sensed soil moisture data are limited to surface soil moisture (Petropoulos et al., 2015;Benninga et al., 2019). Fig. 2 provides an overview of the SMAP L3 Enhanced footprint relative to the in situ locations in the study area.

General workflow
We set up a TFN model including corresponding IR functions for each location of the in situ soil moisture monitoring network in the study area. Fig. 3 shows the general research workflow to set up a TFN model. The workflow focuses on three main parts, indicated by the yellow boxes in the figure. First, the TFN models are set up and calibrated using SMAP data (SMAP calibration). Next, these models are validated using SMAP data for a different period (SMAP validation). Finally, we assess the applicability of the TFN models for estimating soil moisture at field scales using in situ measurements (Field validation). The numbers 1 to 16 in the following sections refer to the steps shown in Fig. 3.

SMAP calibration
First, (1) the input datasets, which are described in section 2.3, are selected. The datasets are split in a (2) calibration set and a (3) validation set. The calibration set is used to determine the parameter sets of the IR functions for the implementation of the TFN models. The validation set is used to assess the TFN model results. We use a calibration period which covers at least the response time of the hydrological system that we observe. We refer to Section 3.2 for the determination of the response time. Especially, we are interested in the predictive capabilities of the TFN model for the dry summer period of 2018. Therefore, the calibration set covers the period January 1, 2016-January 1, 2018, while the validation set covers the period January 1, 2018-January 1, 2019. Since more SMAP observations are becoming available, the calibration period can be continuously extended in an operational setting.
We assessed the influence of the calibration period length by performing a sensitivity analysis in which both the length and period of calibration set are varied. The calibration periods used for the sensitivity analysis are a summer period (  First, (4) the SMAP series of the calibration set are used as the observational calibration dataset for the TFN model. Then, (5) the stress series are defined for the calibration period. Subsequently, (6) a IR function is defined for each stress series. Both the observations and the stress series are added to a (7) Pastas model object. Then, Pastas applies (8) a least squares optimization approach to find (9) optimal parameter sets for the precipitation and evapotranspiration IR functions for each in situ location by solving equation (1). These sets lead to the best fit of the SMAP soil moisture observations for the calibration period. Finally, the sets are used to estimate soil moisture dynamics in the validation period.
We assessed the applicability of two functions to define the IR function of each stress series: a gamma and an exponential function. Table 2 lists the combinations. The explained variance percentage (EVP) is calculated to assess the applicability of each combination. The EVP is defined as: where σ 2 h is the variance of the SMAP soil moisture observations [ðm 3 m À 3 Þ 2 ] and σ 2 n is the variance of the TFN model residuals as defined in equation (1) [ðm 3 m À 3 Þ 2 ] (Von Asmuth et al., 2002). An EVP of 100% indicates a perfect simulation of the observations, since no residuals exist in that case. As a rule of thumb, one generally accepts the results of a TFN model if the EVP � 70% (Van Engelenburg et al., 2020).
Additionally, the noise series should not be significantly autocorrelated. Autocorrelation would indicate that the white noise assumption does not hold (Von Asmuth et al., 2002). We use the Ljung-Box test to determine whether the noise series shows significant autocorrelation (Ljung and Box, 1978).

SMAP validation
We use SMAP observations from the year 2018 to validate the TFN models. We define (10) the precipitation and evapotranspiration stress series for the period January 1, 2018-January 1, 2019. These series are used to (11) set up a Pastas model for the validation period. The (9) optimized parameter sets from the calibration set are applied to define the IR functions. Again, Pastas solves equation (1) using the defined IR functions, the optimized parameter sets, and the stress series, resulting in (12) predictions of soil moisture for the validation period. We use (13) the SMAP validation set to (14) assess the TFN model results using the unbiased root mean square error (uRMSE), bias, and Pearson correlation coefficient, see Appendix A.

Field validation
Furthermore, we are interested in the applicability of the TFN models at field scales compared to the regional scales represented by the SMAP observations. Therefore, we use (15) in situ soil moisture measurements from the soil moisture monitoring network to (16) evaluate the TFN model results at field scales using the uRMSE, bias, and Pearson correlation coefficient error metrics. The evaluation is performed for all in situ location for which measurement data are available for the validation period. The following eleven stations provide data for the 2018 validation period: 2, 4, 7, 9, 10, 11, 13, 14, 15, 16, and 17.

Selection of IR functions
First, we evaluated which combination of IR functions leads to the best fit of the TFN models in terms of EVP. Fig. 4 shows the spatially averaged EVP for the four function combinations as defined in Table 2. The GG and GE combinations lead to TFN models which cannot sufficiently explain soil moisture dynamics. Both the GG and GE combinations show spatially averaged EVP values lower than 50%. The EE and EG combinations show the best model behaviour, even though the EG combination leads to a rejection of the TFN model for one location, i.e. station 16 (67%). The EVP of the EE combination is consistently larger than 70% for all stations, exceeding the model acceptance criterion.
Since the EE combination consistently shows good accuracy, we use the exponential function for both the precipitation and evapotranspiration IR functions in the remainder of the study.
Because both the GG and GE combinations are rejected, the TFN models will be rejected when a gamma function is used for the precipitation stress series. So, although Von Asmuth et al. (2002) show that the gamma function is suitable to model the response of groundwater head to recharge using precipitation stress series, the gamma function is not the best choice for precipitation stress series when modelling surface soil moisture dynamics. The difference is less distinct for the evapotranspiration stress series. Either a gamma or exponential function leads to approximately similar results in terms of EVP.

Calibrated IR functions
The IR functions contain valuable information on the soil moisture response to the stress series. Fig. 5 shows the calibrated IR functions of the precipitation and evapotranspiration stresses for location 2. As expected, the precipitation stress series has a positive impact on the soil moisture state, while the evapotranspiration stress series decreases the soil moisture state. Also, the time scale of the precipitation stress series is smaller than the time scale of the evapotranspiration stress series. Furthermore, the initial response of soil moisture at day one is much larger for the precipitation stress series than for the evapotranspiration stress series. Similar observations hold for all individual locations. These findings are physically reasonable, since precipitation causes immediate spikes in soil moisture, while drydown due to evapotranspiration takes place on longer time scales.
In addition, Fig. 5 shows that the length of the precipitation IR function for location 2 is approximately 75 days, while the evapotranspiration IR function has a length of approximately 150 days. The length of the IR functions can be interpreted as the system response to that specific stress series. The calibration period of two year covers this length multiple times. Therefore, we can conclude that the calibration period is of sufficient length to estimate soil moisture dynamics in the Twente region.

Assessment of soil moisture modelling
We will show the TFN model results for the first in situ location for which a full field validation dataset is available, which is location 2. On the other hand, the TFN model underestimates soil moisture during July-August 2018. Finally, the SMAP observations overestimate soil moisture in the transition from summer to fall in 2018. Van der Velde et al. (2019) evaluated SMAP surface soil moisture data in the Twente region and stated that ''in the summer/fall of 2018, dry spells were ended by a sequence of substantial rain events that exposed the disparity in sampling depth between SMAP and the in situ sensors''. Shellito et al. (2016) and Benninga et al. (2018) found similar results.
We quantify the accuracy of the TFN models by calculating the uRMSE, bias and correlation coefficient error metrics with respect to the SMAP observations for each in situ location in the 2018 validation  period. The figure also shows the results for the 2016 validation period, which are discussed in Section 4.2. The dots in Fig. 7 (SMAP validation  2018) visualize the error metrics of the TFN models results with respect to the SMAP observations for the year 2018. The uRMSE varies between 0.059 and 0.070 m 3 m À 3 , the bias varies between 0.0040 and 0.019 m 3 m À 3 , and the correlation coefficient varies between 0.79 and 0.82 ½ À � for all locations. Recognizing that the TFN models do not consider the over-and underestimation of SMAP in frozen and dry conditions, the TFN models perform well in predicting SMAP surface soil moisture. An implication is that the calibrated TFN model can be used to estimate surface soil moisture and extend SMAP data if precipitation and evapotranspiration data are available.
In addition, Fig. 7 shows the uRMSE, bias and correlation coefficient of the TFN model results with respect to the in situ soil moisture measurements for eleven locations in the 2018 validation period (Field validation 2018). Although a fundamental difference exists in spatial scales represented by the SMAP satellite footprint and the in situ measurements, the TFN models accurately predict field scale soil moisture for seven out of eleven locations in terms of uRMSE, for four out of eleven locations in terms of bias, and for all locations in terms of correlation coefficient.

Verification of TFN modelling approach
The results show that TFN modelling using the PIRFICT method can be applied to predict surface soil moisture conditions in the Twente region using SMAP surface soil moisture remote sensing data as calibration set. As part of the TFN model verification, we assessed whether the noise series show autocorrelation using the Ljung-Box statistical test (Ljung and Box, 1978). The autocorrelation is assessed considering a significance level of 0.05 [-]. The Ljung-Box test shows that no significant autocorrelation is observed for all stations. Therefore, the white noise assumption holds for all stations.
It is generally known that the validity of data-driven models is often limited to the hydrological conditions for which the model is calibrated (e.g. Abrahart et al., 2010;Kornelsen and Coulibaly, 2014;Mount et al., 2016). Therefore, because of the data-based nature of TFN models, there is the risk to extrapolate results to situations for which no references are available in the calibration set (Von Asmuth et al., 2012). For example, the TFN model of location 2 does not capture the extremely dry summer period of 2018 well, as seen in Fig. 6. According to the TFN model, the volumetric moisture content drops to almost 0 m 3 m À 3 in that period. However, both the SMAP observations and the in situ measurements show that soil moisture has a physical lower limit of approximately 0.1 m 3 m À 3 . The TFN models do not identify the lower limit. This limitation might be explained by evapotranspiration reduction, which is a non-linear mechanism which reduces actual evapotranspiration when only low amounts of moisture are available. The linear nature of TFN models indeed decreases their accuracy during long dry periods, according to Peterson and Western (2014). Moreover, the Makkink method describes the potential evapotranspiration, which is the maximum evapotranspiration occurring when water is not a limiting factor. Actual evapotranspiration observations could be considered to improve the TFN model results.
These findings stress the need for a representative calibration dataset. In addition, the length of the calibration period may affect the TFN model results. The influence of the calibration period will be assessed in the next section.

Sensitivity of calibration period
We performed a sensitivity analysis on both the calibration set period and the length of the period based on the RMSE error metric with respect to the in situ measurements in the 2018 validation period.  calibration set, the TFN models can simulate the dry period of 2018 correctly. Thus, the TFN models can represent the drought period of 2018 when a representative calibration period is selected. We want to stress that the representativeness of the calibration period is defined by the characteristics of the soil moisture dynamics occurring in that period rather than the calendar date. One would have to assess the mean, variance, and other statistics to assess whether the calibration set is representative for the period of interest.

Results for 2016 validation period
The dry period in the summer of 2018 is an extreme event. Probably, the sensitivity analysis results are case-specific and thus not generalizable. To test the applicability of the TFN models in more general situations, we switched the calibration and validation set periods. The TFN models are calibrated for the period January 1, 2017-January 1, 2019 and validated for the year 2016. The triangles in Fig. 7 (SMAP validation  2016) show the uRMSE, bias, and correlation coefficient of the TFN models with respect to the SMAP observations when the year 2016 is used as validation set. The 2016 results have consistently higher accuracies than the SMAP 2018 validation results. The uRMSE varies between 0.042 and 0.052 m 3 m À 3 , the bias varies between 0.0061 and 0.023 m 3 m À 3 for the locations, and the correlation coefficient varies between 0.82 and 0.88 ½ À � for the locations. Furthermore, the stars in Fig. 7 (Field validation 2016) show the error metrics of the TFN models with respect to the in situ measurements when the year 2016 is used as validation set. The following stations have a full dataset for the 2016 validation period: 1, 2, 3, 7, 9, 16, and 18. No apparent change in accuracy is found for the 2016 validation results at field scales, which are represented by the in situ measurements.

Accuracy of TFN modelling estimates
This section discusses the accuracy of the TFN modelling estimates as shown in Fig. 7. Table 3 shows the performance of the soil moisture TFN modelling approach in comparison with other SMAP validation studies in the Twente study area. Colliander et al. (2017) and Chan et al. (2018) validated the SMAP and SMAP Enhanced soil moisture products, respectively. Kolassa et al. (2018) used the SMAP microwave observations and a neural network to retrieve surface soil moisture estimates. Since these studies considered spatially averaged results for the Twente study area, the TFN model results are also spatially averaged. Especially the TFN model validated for the year 2016 performs well as shown by the uRMSE, bias, RMSE, and correlation coefficient error metrics. Less extreme dry and wet events occurred in 2016 in comparison with 2018. The underestimation of the dry 2018 period has a significant impact on the model accuracy, which explains the increased accuracy of the TFN models in the 2016 validation period. The comparison with the data-driven approach by Kolassa et al. (2018) shows that SMAP TFN modelling for the year 2016 has a comparable accuracy. These results indicate that a calibration period including the extreme dry summer of 2018 will lead to more accurate TFN models for drought predictions,  especially at spatial scales similar to the SMAP satellite footprint.
The SMAP satellite retrievals have a design accuracy of 0.04 m 3 m À 3 in terms of uRMSE (Entekhabi et al., 2010), which is represented by the red line in Fig. 7. We want to stress that the accuracy of the original SMAP L3 Enhanced retrievals do not meet this design accuracy in the Twente study area (0.056 m 3 m À 3 ) . The upper panel of Fig. 7 shows the accuracy of the TFN model estimates with respect to the design accuracy. Although not meeting the design accuracy, SMAP validation 2018 and SMAP validation 2016 approach the accuracies found by Chan et al. (2018) in terms of the uRMSE, bias, and correlation coefficient error metrics. Also, the field validation 2016 estimates correspond well with the field observations, although local scale issues cause discrepancies for some stations. Large discrepancies can be found for the field validation 2018, which are mainly caused by large biases as shown in panel B of Fig. 7.

Spatial and temporal resolution of TFN soil moisture estimates
The spatial resolution of the TFN soil moisture estimates depends on the spatial resolution of the input data. Only considering the spatial resolution, the reference crop evapotranspiration data form the limiting factor for the TFN soil moisture estimates. These data have a resolution of 10 km by 10 km, see Table 1. However, evapotranspiration does not vary as much as precipitation over short distances (Dalezios et al., 2002;Hess et al., 2016). Therefore, we assume that the TFN soil moisture estimates are bound by the spatial resolution of the SMAP input data, which is 9 km by 9 km. The results shown in Table 3 and Fig. 7 reflect this statement. Table 3 shows that, at a regional scale, the accuracy of the TFN soil moisture estimates approximates the accuracy of other soil moisture products. We define a regional scale as the extent of our study area, which approaches the spatial resolution of the SMAP Passive soil moisture product (36 km by 36 km).
At local scales, the accuracy of the TFN soil moisture estimates mainly depends on whether the SMAP input data reflect local conditions, as shown in Fig. 7. The in situ soil moisture monitoring stations represent these local scales. Although the correlation shows good corresponding between the TFN soil moisture estimates and in situ data, the uRMSE and bias error metrics show large discrepancies for approximately half of the locations. Therefore, the spatial resolution of the input data is an important aspect which should be considered in TFN modelling. In addition, the SMAP soil moisture retrievals are limited to shallow soil depths (Petropoulos et al., 2015). As the IR functions are calibrated using the SMAP data, the TFN soil moisture estimates are also limited to shallow soil depths. Other methods, such as data-assimilation using process-based models (Liu et al., 2012;Pezij et al., 2019b), are necessary for the translation to soil moisture in deeper layers.
In contrast, the temporal resolution of the TFN soil moisture estimates only partially depends on the temporal resolution of the input data. The IR functions enable prediction of soil moisture if precipitation and reference crop evapotranspiration data are available. As both the precipitation and reference crop evapotranspiration are available on a daily basis, the temporal resolution of the TFN soil moisture estimates is also considered to be on a daily basis. As the SMAP satellite does not provide daily soil moisture observations, the TFN models can be used to fill the data gaps occurring on days without a satellite overpass, as discussed in Section 3.3. The SMAP input data have an availability of 2-3 days. This study shows that this temporal resolution is sufficient to develop IR functions for accurately estimating soil moisture. However, the soil moisture response due to some individual precipitation events might not be included in the calibration procedure if the SMAP satellite does not overpass the area of interest. Therefore, the accuracy of the TFN soil moisture estimates might be affected when heavy precipitation events are missed.

Water system characteristics
Generally, IR functions can provide information on characteristics of the system which is observed. Among others, these functions describe the unit step response and response time of groundwater dynamics when a groundwater system is observed Zaadnoordijk et al., 2018). The exponential function applied in this work has two parameters: A and a, see equation (4). The parameter A is related to the total change in soil moisture due to a unit stress. A large A indicates a large total change. The shape parameter a is related to the time scale on which a unit stress affects soil moisture. A large a indicates a slow response. Fig. 10 shows the calibrated parameters for the precipitation and evapotranspiration IR functions for all locations based on the 2016-2017 calibration set. We have assessed the relationship of the parameters with the following variables: longitude, latitude, soil elevation, vegetation type, and soil type. The elevation data are obtained from Actueel Hoogtebestand Nederland (2019) The vegetation characteristics are obtained from Van der Velde et al. (2019). The soil type characteristics are obtained from the BOFEK2012 dataset (W€ osten et al., 2013). The next sections will elaborate on the relationships found for parameters A, related to the total change, and a, related to the response time.

Total change
The change in total soil moisture volume due to precipitation is strongly correlated to longitude. The change is larger for stations in the western part of the study area. Longitude is strongly correlated with distance to the coast in the Netherlands. Daniels et al. (2014) found that, on average, the precipitation amount is higher near the coast. Similarly, the change in total soil moisture due to evaporation is correlated to latitude. The change is larger for stations in the north of the study area. In addition, the total change in soil moisture might be related to soil physical characteristics and vegetation type. However, the elevation, soil type and vegetation type do not show clear patterns with respect to the total change in soil moisture. These characteristics are quite homogeneous for the locations in the study area, with mainly sandy soils and grass vegetation. A valuable addition can be to study whether different IR functions and corresponding parameters are found in areas with other soil characteristics, such as peaty or clayey areas, or areas with different vegetation types.

Response time
A linear trend is seen in the spatial distribution of parameter a in relation with longitude. In general, if a large soil moisture response time scale due to precipitation is found at a specific location, the corresponding soil moisture response time scale due to evapotranspiration is relatively small and vice versa. The time scale of the precipitation IR function is larger for locations in the western part of the study area. Also, the subsurface in the western part of the study area contains thick sand layers. Precipitation infiltrates relatively easy in sandy layers, which causes a fast response of shallow soil moisture. Moraines of clay are found in the eastern part of the study area. Clay layers have low infiltration rates, which results in a slow increase of soil moisture content. However, the soil type characteristics do not show a clear relation with respect to the a parameter. The latitude, elevation, and vegetation type characteristics do not show clear patterns with respect to the soil moisture response time. More research is needed to generalize these findings.
One should be careful in relating these parameters to physical processes as the selection of a representative IR function is an assumption (Von Asmuth et al., 2012). For example, Fig. 5 shows that the time scale of the precipitation IR function is approximately 75 days. To our knowledge, this time scale cannot be directly connected to physical phenomena. More research on the IR function parameters is needed to increase the understanding of their physical meaning in terms of soil moisture.

Conclusions
We studied the applicability of transfer function-noise modelling (TFN) for describing and predicting soil moisture dynamics. TFN modelling is a fast alternative for process-based models, taking only seconds to simulate a full year of daily soil moisture conditions. TFN modelling is based on the assumption that soil moisture dynamics can be explained by linearly transforming precipitation and evapotranspiration stress series using impulse-response (IR) functions. The SMAP L3 Enhanced surface soil moisture product is used to calibrate the TFN models. We found that exponential functions describe the IR functions of both the precipitation and evapotranspiration stress series better than gamma functions.
TFN models describe soil moisture conditions well when comparing the TFN model results with the SMAP observations. In addition, the TFN model results were compared with in situ soil moisture measurements to assess the field scale applicability of TFN modelling. The accuracy of the TFN models mainly depends on the representation of the SMAP satellite product for that specific spatial scale.
A practical application for operational water management is that the TFN modelling approach can be used to estimate soil moisture dynamics using predictions of precipitation and evapotranspiration. The application is promising if sufficient calibration data are available, although one should be careful when interpreting results in extreme situations, since the TFN models do not consider the physical lower and upper limits of soil moisture. However, a sensitivity analysis and variation of the calibration and validation periods showed that selecting a suitable calibration period can significantly increase the TFN model capabilities in both regular and extreme situations. In addition, the IR function parameters potentially provide valuable information on water system characteristics, such as the total response and the response times of soil moisture to precipitation and evapotranspiration. Thus, the TFN models have value in explaining hydrological processes and characteristics, in contrast to other data-driven tools such as neural networks, which are black boxes. However, more research on the physical meaning of these parameters is needed to understand their applicability. Concluding, we consider the applicability of TFN modelling for explaining soil moisture dynamics promising and propose to explore the possibilities of TFN modelling for predicting soil moisture in operational settings.

Declaration of competing interest
The authors declare that they have no know competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. thank the Water Resources department of the ITC faculty of the University of Twente, in particular Rogier van der Velde and Harm-Jan Benninga, for sharing the in situ soil moisture measurements from their monitoring network. The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix A. Error measures
The Root Mean Square Error (RMSE) is defined as: