Temporal covariance of spatial soil moisture variations: A mechanistic error modeling approach

When estimating field‐scale average soil moisture from sensors measuring at fixed positions, spatial variability in soil moisture leads to “measurement errors” of the spatial mean, which may persist over time due to persistent soil moisture patterns resulting in autocorrelated measurement errors. The uncertainty of parameters that are derived from such measurements may be underestimated when they are assumed to be independent. Temporal autocorrelation models assume stationary random errors, but such error models are not necessarily applicable to soil moisture measurements. As an alternative, we propose a mechanistic error model that is based on the spatial variability of the water retention curve and assumes a uniform water potential. We tested whether spatial soil moisture variability and its temporal covariance could be predicted based on (1) mean soil moisture, (2) water retention variability, and (3) (co)variances of the van Genuchten parameters using a first‐order expansion of the retention curve. The proposed models were tested in a numerical and a field experiment. For the field experiment, in situ sensor measurements and water retention curves were obtained in a field plot. Both experiments showed that water retention variability under a uniform water potential is a good predictor for spatial soil moisture variability, and that soil moisture errors are strongly correlated in time and neglecting them would be an incorrect assumption. The temporal error covariance could be predicted as a function of the mean moisture contents at two observation times. Further research is required to assess the impact of these temporal correlations on soil moisture predictions.


INTRODUCTION
With the development of cheap and autonomous soil monitoring sensors, real-time soil moisture data can play an important role in irrigation scheduling.The use of this information, however, remains a challenge mainly because the local sensor measurement may differ from the average soil moisture due to soil heterogeneity across the field.Even though the local measurement may correctly represent the local soil water content, in the context of irrigation scheduling that is based on averaged or mean soil moisture contents, its deviation from the average or mean water content can be interpreted as a "measurement error."In addition, when measuring soil moisture with sensors at fixed positions, a deviation from the mean soil moisture in the root zone may persist over time resulting in autocorrelated measurement errors.The uncertainty of variables or parameters that are derived from such time series may be underestimated if these measurement errors are interpreted as independent errors.Soil moisture measurements are used to calibrate models, for example, soil water balance models, that simulate soil moisture (Greenwood et al., 2010), which is a key variable for plant growth but also for microbial activity and biogeochemical cycles.These models typically simulate vertical soil moisture profiles or average moisture contents in soil layers.The simulated moisture contents therefore represent fieldscale average moisture contents at different depths or in soil layers.Bayesian methods can be used to parameterize soil water balance models based on soil moisture measurements, their uncertainty, and prior information about the model parameter distributions (Scharnagl et al., 2011(Scharnagl et al., , 2015)).The uncertainty of the measurements in this context refers to the deviation of the measurements from the field-scale average soil moisture.The deviation between model predictions and measurements is used to calculate likelihood functions where the weighting by the measurement uncertainty is crucial.In theory, the uncertainty of each measurement point and the correlation between all measurement errors need to be accounted for.In practice, however, error models are used to represent this uncertainty and correlation.The simplest model assumes that measurement errors are uncorrelated and do not vary over time, and that the variance of the measurement errors corresponds with the variance of deviations between the model predictions and measurements.This model neglects the measurement error correlation and the fact that the soil moisture variability may vary with the mean soil moisture content.More accurate models may consider autocorrelated errors and a relation between the error and mean soil moisture content.The parameters of these models, that is, hyperparameters, may be derived from the deviation between model simulations and measurements, or from measurements directly.In this paper, we present another approach where we use a simple

Core Ideas
• Soil moisture is spatially variable, and spatial variations are correlated in time.• Due to spatial variability, the estimated average soil moisture from a limited number of measurements is uncertain.• The error of the spatial mean estimated from a limited number of sensors at fixed positions is correlated in time.• A mechanistic error modeling approach is presented based on the spatial variability of the water retention curve.• Such a model can predict measurement variability and error covariance in time as a function of mean soil moisture.
mechanistic model to link this error and correlation to the spatial variability of the soil water retention curve.Soil properties including bulk density, saturated hydraulic conductivity, and van Genuchten (vG) water retention curve parameters (van Genuchten, 1980) are subjected to both spatial and temporal variations, which may lead to spatial and temporal variability of soil moisture in a field.Hence, accounting for the spatial distribution of these soil properties is crucial for irrigation management purposes (Feki et al., 2018).In addition to the heterogeneity of soil physical and hydraulic properties, field-scale soil moisture variability can result from a number of other reasons including variations in crop growth (canopy cover and root water uptake), variations in soil management (tillage, nutrients) and irrigation (nonuniform drip or overhead application), nonuniform precipitation, lateral water redistribution due to topography, and variations in microclimate or evaporative demand due to variations in aspect and slope (Wilson et al., 2004).Although small-scale heterogeneity was found to have a minor impact on root water uptake and long-term field-scale soil water balance computations (Schlüter et al., 2013), soil heterogeneity might still impact interpretation of soil monitoring such as soil moisture sensor measurements, as well as the short-term and local soil water balance.
Soil moisture variability has been widely studied in relation to mean soil moisture (Bell et al., 1980;Famiglietti et al., 2008;Irmak et al., 2022;Rosenbaum et al., 2012), and in relation to the spatially variable soil properties, vegetation cover, and topography as the immediate causes of this soil moisture variability (Albertson & Montaldo, 2003;Lawrence & Hornberger, 2007;Manns et al., 2014;Pan & Peters-Lidard, 2008;Schlüter et al., 2013;Teuling & Troch, 2005;Vereecken Vadose Zone Journal et al., 2007).Teuling and Troch (2005) successfully modeled soil moisture variability in a heterogeneous soil at field scale up to small catchment scale by accounting for variability in the saturated hydraulic conductivity, maximum land cover, porosity, and pore size distribution.This model was further developed by Lawrence and Hornberger (2007), who found that variability at low, mid-range, and high mean soil moisture levels is mainly determined by the wilting point, (unsaturated) hydraulic conductivity, and the porosity, respectively.
The uncertainty of repeated soil moisture measurements in a field is determined not only by soil moisture variability itself, but also by the persistency of local soil moisture measurement "errors" over time, which manifest themselves as a consistent under-or overprediction of the spatial mean or average soil moisture.The statistical measure of this error persistency is the autocorrelation.The default assumption, however, is that measurement errors are not correlated with each other, which generally results in an underestimation of the uncertainty of predictions that are based on repeated measurements with autocorrelated errors over time.Correlation removal strategies such as autoregressive models (AR, AR(I)MA) assume that measurement errors are temporally stationary random processes, that is, the covariance of errors depends only on the time lag between the measurements, and could be used when spatial soil moisture variations are a result of a spatially variable but temporally stationary random process such as rainfall.However, temporal correlations of local soil moisture deviations that are due to spatial variations in soil properties, vegetation cover, and topography are not stationary in time but process based (Doherty & Welter, 2010).This is because spatial variation in the properties and processes that generate this variability is persistent over time so that deviations in soil moisture at a certain location and time from the overall mean will have an imprint that does not vanish over time.Other ways to describe autocorrelation are required for soil moisture time series that are influenced by such temporally invariant but spatially variable properties and processes.
In this paper, we focus on how to model soil moisture variability and covariance in fields with spatially variable soil properties.Since we focus on irrigated and intensively managed crops, we assume that the vegetation cover is quite uniform across the field and that the considered fields are flat so that we may neglect the effect of topography and variable vegetation on soil moisture variability.Even though water retention curve variability has been described and discussed by several researchers (Cameron, 1978;Canone et al., 2008;Dohnal et al., 2006;Montzka et al., 2017;Shouse et al., 1995), this knowledge is rarely used to predict soil moisture variability, let alone soil moisture error correlation.A widely applied technique to represent retention curve variability introduced by Miller and Miller (1956) is the scaling technique that scales local retention curves to a reference curve, and it is useful in gaining insight into soil moisture variability in heterogeneous soils (Montzka et al., 2017;Sadeghi et al., 2016).Qu et al. (2015) provide another perspective on soil moisture variability and its relation with water retention variability, as they derived a first-order approximation of the vG model to predict soil moisture variability as a function of mean soil moisture based on the vG model parameters and their variability.
We propose a mechanistic error model that is based on the variability of the soil water retention curve and assumes a uniform soil water potential in the field (i.e., uniform in the horizontal direction, but variable in depth and time).We tested whether soil moisture variability and autocorrelation could be predicted as a function of the mean soil moisture based on (1) the variability of the soil water retention curves or (2) information on the variability and covariance of the vG parameters using a first-order expansion of the retention curve (Qu et al., 2015).The approaches and the assumptions on which these methods are based are described in the first section of this paper (Theoretical framework).Then, these methods were applied and validated both in a numerical and a field experiment.In the numerical experiment, the assumptions about a uniform soil water potential and spatial variability in soil moisture that is generated only by spatial variability in the soil hydraulic properties were tested using a virtual dataset of soil moisture measurements derived from water flow simulations in a heterogeneous field with a variability in soil hydraulic properties and root water uptake.The advantage of this approach is that the soil variability, based on a Miller-Miller scaling (Miller & Miller 1956), is known so that we can verify these assumptions.Second, the model was evaluated under real field conditions.For this evaluation, in situ soil moisture sensor measurements (TEROS 10; METER Group Inc.) as well as soil water retention curves were obtained at different locations in a field plot.

THEORETICAL FRAMEWORK: SOIL MOISTURE VARIABILITY AND "MEASUREMENT ERROR" AUTOCOVARIANCE
Sensor measurements and their errors, that is, their deviations from the (unknown) field-averaged soil moisture content at a certain depth and time  (  =   − θ ), are autocorrelated in time as each sensor measures soil moisture at a fixed position.In this study, the covariance between two "measurement errors" at the same location at two different moments is examined.Two measurements at two different locations are assumed to be independent as long as they are sufficiently far from each other (i.e., we do not consider spatial autocovariance and assume that the distance is larger than the horizontal correlation length [  ] of the soil properties).Different approaches are considered to describe and  1.To simplify the notation, the covariates will not be shown in the following but only the variables that are used to predict the covariance.First, standard deviations and autocovariances are derived directly from observations.For the autocovariances, these include the nonstationary temporal autocovariance (T-NS; Equation 1) and the stationary temporal autocovariance (T-S; Equation 2).The latter corresponds to the classic temporal covariance computation, that is, the covariance only depends on the time interval between the measurements.The nonstationary temporal autocovariance (T-NS) will serve as a reference for the stationary temporal autocovariance model (T-S), both of which are expressed as a function of time.
Additionally, a soil moisture-based autocovariance (SM model; Equation 3) is computed directly from observations.Instead of calculating the autocovariance of measurements at two different time points, the autocovariance is calculated for pairs of measurements that are observed at a certain pair of mean water contents.In this approach, it is assumed that the autocovariance between a pair of measurements is the same when the pair of mean water contents is the same, irrespective of the times at which the measurements are carried out.When longer time series of moisture contents are available and the time course of the mean soil moisture content crosses several times the same moisture content level, the degrees of freedom of a covariance model that uses mean soil water contents as independent variables is considerably reduced compared with a nonstationary temporal covariance model.Standard deviations are also calculated directly based on the observations as a function of mean soil moisture using the SM model.These standard deviations and autocovariances based on mean soil moisture will serve as the reference for the modeled standard deviations and autocovariances that are derived from the spatial variation of the soil moisture retention curve (RC, RC-FO, RC-FO-uncorr).
Next, standard deviations and autocovariances are estimated from fitted vG water retention curves.It is assumed that the water pressure heads are spatially uniform at a certain depth in the field.Based on this assumption, standard deviations and autocovariances of soil moisture are derived from the spatial variability of the vG water retention curves (RC; Equation 4), or from the variance and covariance of the spatially variable vG water retention parameters using a first-order expansion of the vG water retention function.Two first-order error models are considered: one that considers correlations between vG parameters (RC-FO; Equation 5) and one that assumes uncorrelated parameters (RC-FO-uncorr; Equation 6).

Temporal autocovariance
The nonstationary temporal covariance (T-NS) between measurement errors at time   and   is defined as follows: where  is the number of sensors,  m (  , ) is the measured soil water content at location  and time   , and θm (  ) represents the mean soil water content of all sensor measurements at time   .

Vadose Zone Journal
The classic temporal covariance of measurement errors made with regular time steps Δ can be computed for time lag Δ assuming a stationary random process in time using Equation (2).This covariance computation will be referred to as stationary temporal autocovariance, or T-S.
where  is the number of sensors and  is the number of measurements made over time. m (  , ) and  m (  − Δ, ) are the soil water contents measured at location  and at times   and   − Δ, respectively, and θm (  ) represents the mean soil water content of all  sensors at time   .Notice that the assumption of stationarity implies that the covariance is only a function of the time lag between two measurements but not of the measurement times themselves.

Soil moisture-based autocovariance
As soil moisture measurement errors are expected to depend on the soil moisture level itself, the covariance between two measurement times   and   can also be calculated as a function of the measured mean soil moisture θ and θ (SM; Equation 3).This covariance computation corresponds to T-NS but expresses covariance as a function of mean soil moisture instead of time and will serve as the reference for the modeled autocovariances as a function of mean soil moisture (RC, RC-FO, RC-FO-uncorr).
(3) where  is the number of sensors and  m ( θ , ) is the soil moisture measured by sensor  at time   when the mean measured soil moisture is θ .Notice that here, the covariance is a function of the mean soil moisture at two measurement times.
As measurement data are very discontinuous, a Gaussian filter with standard deviation  = 10, corresponding to 10 soil moisture steps, is applied to smooth the covariance surface and to interpolate over a wide range of soil moisture, with a larger standard deviation resulting in a wider and smoother filter.Depending on the size of the bins with a maximum width of 0.005 m 3 •m −3 , this may correspond to a maximum standard deviation of 0.05 m 3 •m −3 .When this covariance function is smoothed, pairs of measurements with very similar pairs of mean soil water contents will have similar covariance, irrespective of the time differences between the different measurements.
F I G U R E 1 Reference retention curve in black and an overview of main van Genuchten (vG) parameter (, ,  r ,  s ) effects.

Predicting soil moisture variability and autocovariance based on soil water retention variability
Field heterogeneity can be described by the variability of the soil water retention curve throughout a field (Figure 1).In practice, the retention curve variability can be obtained through repeated sampling of undisturbed soil cores.The reference retention curve corresponds to the reference vG parameter set fitted to the average of the soil water contents measured in all soil cores at a given pressure head.
Soil moisture variability can be predicted as a function of mean soil moisture based on the variability of the soil water retention curve.Assuming a uniform soil water potential, the mean soil moisture and its standard deviation can be calculated for each soil water potential.In addition to soil moisture variability, soil moisture covariance can also be estimated based on this water retention curve variability.The (co)variance between two soil water potentials, or two mean soil water contents, can be calculated using Equation (4).This approach will be referred to as RC.
where  rc is the number of water retention curves,  wrc (ℎ, (rc)) is the soil water content of water retention curve rc at water potential ℎ, with (rc) being the parameter set of the water retention curve, and with mean soil water content θwrc (ℎ) of all water retention curves at water potential ℎ.

Predicting soil moisture variability and autocovariance derived from a first-order expansion of the water retention curves
Soil moisture variability and covariance can be derived analogous to the derivation of Qu et al. (2015) using a first-order expansion of the vG function assuming a horizontally uniform water potential (Supporting Information Appendix B).
The expression is based on the variability of vG parameters , , and  s and accounts for vG parameter correlations, while variability in and correlations with  r are neglected.This method will be referred to as RC-FO (Equation 5).Another model is derived from this, now assuming the vG parameters are uncorrelated, and will be referred to as RC-FO-uncorr (Equation 6). (5) where   parameters are a function of mean vG parameters and soil water potential ℎ  (see Supporting Information Appendix A, Equations A1-A5) with  0 =  0 (ℎ  ), and   ,   , and   s represent the standard deviation of the respective vG parameters (Qu et al., 2015).

Theoretical approach: Virtual soil
The data used in the first part of this study included daily soil water contents (m 3 •m −3 ) at 0.1 m depth in a two-dimensional virtual soil: a vertical plane of 5 × 5 m, covered with grass with uncompensated root water uptake (see Supporting Information Appendix D, Figure D1).The data were made available from a study on soil heterogeneity (Schlüter et al., 2013;Schlüter, Vogel, et al., 2012) where soil water flow was simulated using 10-year weather data comprising daily values for precipitation and potential evapotranspiration of Magdeburg (eastern Germany), with alternating wet and dry periods.The lateral resolution of the data at 0.1 m depth was 0.02 m (250 data points over a distance of 5 m).The profile consisted of a heterogeneous silty cover and a sandy subsoil with embedded loam lenses (USDA classification).The silty cover layer consisted of a plow horizon down to 0.3 m and an untilled layer below.Plowing in the topsoil layer resulted in isotropic subscale heterogeneity, with horizontal and vertical correlation lengths  , = 0.04 m.As the distance between sample points is required to be larger than the horizontal correlation lengths to be considered independent, a sampling distance of 0.08 m was adopted.
The heterogeneity in the silty cover was characterized by a Miller-Miller scaling of the water retention and conductivity function (Miller & Miller, 1956), defined by a reference water retention function ℎ * m () and a reference conductivity function  * ().These hydraulic reference functions corresponded to a characteristic pore size  * with scaling factor (⃗ ) = (⃗ )∕ * , where ⃗  = (, ) is the position vector, and parameters of the hydraulic reference functions are shown in Table 2.A 0.08m-deep, loose seed bed with additional macroporosity was described by a bimodal retention function for which the fractions are denoted in the table.The plow pan was represented by a thin compacted soil layer perforated by macropores.
The roots were modeled either as a deterministic, depthdependent (Z) or variable, structure-dependent (S) root distribution, both of which were static.The deterministic root distribution (Z) was a function of depth but was independent of soil structure, resulting in spatially uniform roots and water uptake at a certain depth.In contrast, the variable root distribution (S) was a function of both depth and soil structure.Moreover, the root distribution function  r was also defined as a function of the scaling factor  for the structure-dependent root model (S) (Equation 8), as fine-textured areas indicated by a low ln() value exhibit a higher resistance against root penetration (Schlüter et al., 2013).
To reproduce the soil heterogeneity, the ln() normal distribution was sampled (n = 250), resulting in 250 different water retention curves for the 0.1 m depth (silty cover).Based on the sampled retention curves, the mean soil water content and its standard deviation were calculated under the assumption of a constant soil water potential.
Since the soil heterogeneity was here solely due to variability in , the variability in the other vG parameters, as well as all vG parameter covariances, was zero.As a result, there was no difference between two of the predictive models, RC-FO and the uncorrelated equivalent RC-FO-uncorr.Both prediction models Equations ( 5) and ( 6) could thus be reduced to Equation (9) (RC-FO).
T A B L E 2 van Genuchten (vG) parameters of the virtual soil (Schlüter et al., 2013).

Case study approach: Validation data
In two repeated field trials (2020 and 2021), field variability was assessed during two growing seasons of leek (July to December) based on undisturbed soil cores, while high-frequency dielectric capacitance soil moisture sensors (TEROS 10) were installed at different field locations at 15 cm depth to capture soil moisture variability, along with sensor measurement error autocorrelation.By applying the manufacturer's calibration for mineral soils to convert the raw sensor output in millivolts to volumetric soil moisture content (m 3 •m −3 ) (Equation 10), an absolute accuracy of 0.03 m 3 •m −3 could be obtained (METER Group, 2018).
The sensors were linked to a datalogger with a communication module for real-time online data access.The raw sensor measurements were aggregated to daily observations.
The experimental design consisted of four locations (25 m apart), each consisting of two or three embedded sublocations (2 m apart) in 2020 and 2021, respectively, with one sensor installed at each sublocation in the root zone at 15 cm depth.The field was a loamy sand to sandy loam soil and was located in Sint-Katelijne-Waver, Belgium.The duration of the growing season and the amount of rainfall during this growing season were 175 days and 330.4 mm in 2020, and 149 days and 334.7 mm in 2021.The field was irrigated with sprinklers with a total amount of 50 mm (four irrigation events) in 2020 and 30 mm (two irrigation events) in 2021.

Validation metrics
The predictive models of soil moisture standard deviations and covariances were compared to observed soil moisture standard deviations and covariances.The root-mean-square error (RMSE) (Equation 11) and relative root-mean-square error (rRMSE; Equation 12) are used as a quantitative metric to evaluate prediction skill.
where   is the ith observed value, ŷ is the ith predicted value, and  is the total number of data points.The rRMSE is computed as the RMSE divided by the mean observed value ȳ.
The variable  may represent the standard deviation or the covariance.
While the RMSE is the standard deviation of the residuals and is stated in the same units as the corresponding variable, rRMSE is a relative measure of variation in accuracy, which is easier to interpret.Model accuracy is considered excellent when rRMSE < 10%, good if 10% ≤ rRMSE < 20%, fair if 20% ≤ rRMSE < 30%, and poor if rRMSE ≥ 30% (Despotovic et al., 2016).

RESULTS AND DISCUSSION
4.1 Theoretical approach: Virtual soil

Soil heterogeneity
The sampled water retention curves and the mean soil moisture curve are shown in Figure 2. The Miller-Miller scaling of soil water potential ℎ resulted in a vertical shift of the water retention curve, which corresponds to a scaling of the vG shape parameter .Due to the lognormal character of the scaling, the distribution of  was positively skewed causing the arithmetic mean value for  (0.017 cm −1 ) to be different from the reference value (0.016 cm −1 ).The  parameter, the residual soil moisture  r , and the saturation point  s were assumed constant; hence, the retention curve was not compressed, that is, the slope of the curve remained unchanged and curves did not intersect.Based on this retention curve variability, soil moisture variability could be computed assuming a constant soil water potential (Figure 2, right axis).
The aggregated and sampled soil moisture time series illustrate the soil moisture variability and the autocorrelation of such time series at fixed locations (Figure 3).

Soil moisture variability
Soil moisture variability was predicted based on retention curve variability (RC) as illustrated in Figure 2 and based on the first-order vG model approximation (RC-FO) (Figure 4).The two predictive models diverge exponentially with increasing soil moisture.The steep increase in RC-FO at high soil water content is due to the  0 parameter becoming infinitely large for a small soil water potential ℎ (see Supporting Information Appendix A, Equation A1).The predictive models were compared with observed soil moisture variability (SM) for the depth-dependent (Z) and structure-dependent (S) root distribution scenario (Figure 4).The observed soil moisture ranged from 0.07 to 0.46 m 3 •m −3 over the 5 m distance in the 10-year period.
Figure 4 shows that all observed data points are on or above the predicted variability curve, as this curve is the minimal expected variability that results from soil heterogeneity.In this virtual soil, additional soil moisture variability can be due to plant variability for the case of structure-dependent roots (S), subsoil variability, and hydraulic nonequilibrium, that is, the water potentials are not necessarily uniform at a certain depth (Schlüter, Vanderborght, et al., 2012;Vogel et al., 2010).The RC was the better variability prediction compared to RC-FO, with an rRMSE of 20.71% (fair accuracy) and 40.28% (very poor accuracy), respectively.The rRMSE of RC was 11.79% and 25.84% for the depth-dependent (Z) and structuredependent (S) root distribution, respectively.The larger error for the latter (S) was expected due to the additional variability in root water uptake as a result of the variable root distribution and root water uptake.As RC-FO and RC are similar at low to mid-range soil moisture contents, RC-FO might still be a good predictor for irrigation purposes for which the wet range is of less interest.All (r)RMSE values are summarized in Supporting Information Appendix C (Table C1).

Soil moisture autocovariance
Next to soil moisture variability, soil moisture autocovariance was also calculated based on observed data (T-S, T-NS, SM) and estimated based on two predictive models (RC, RC-FO).
First, the stationary temporal autocorrelation of the time series (T-S) was computed using Equation (2).The T-S model was expected to be a poor estimator, as local soil moisture deviations are not a stationary random process.On Figure 5, the temporal autocovariance matrix shows a diagonal pattern, as each lag time difference has a fixed correlation, starting from the maximum covariance ( = 1) for a lag Δ = 0 to a covariance equal to 0 for a lag Δ = ∕2, which corresponds to 5 years.The T-S of the depth-dependent (Z) and structure-dependent (S) root distribution scenario is nearly identical, except for the maximum autocovariance, which is higher for the latter (S) (not shown).This model prediction can be compared to the nonstationary temporal autocovariance (T-NS) based on the observed time series (Figure 5C), where correlation does not depend on the time lag, but the measurement times themselves.It is clear that T-S (Figure 5B) does not agree with the reference observed autocovariance (T-NS; Figure 5C), which is confirmed by a high rRMSE of 87.09%, indicating that a certain time lag may have a different correlation depending on the measurement times.The covariance surface plotted versus measurement times (Figure 5C) is a very irregular surface, which does not correspond with a smooth surface versus time.This suggests that it is difficult to describe the covariance as a simple function of the times of the two measurements.Hence, assuming stationarity is proven inappropriate to estimate autocorrelation for a process-driven time series such as soil moisture.
Second, the autocovariances between pairs of deviations from mean soil water contents (SM) were calculated using Equation (3) with a maximum bin width of 0.005 m 3 •m −3 .This matrix will be our reference observed autocovariance as a function of mean soil moisture and is shown for both rooting scenarios in Figure 6.In contrast to the autocovariance surface plotted versus measurement times (Figure 5C), the autocovariance surface plotted versus the average soil moisture contents is a more smooth surface that could be represented by a bivariate function.This indicates that the autocovariance is rather related to the average soil moisture status of the field at two different times than to the times of the measurements or time difference.These covariances largely correspond to correlations close to 1, while minimum correlations of 0.65 and 0.4 were obtained for the depthdependent (Z) and structure-dependent (S) root distribution, respectively, for a mean soil moisture content near saturation point (0.46 m 3 •m −3 ) in combination with lower soil moisture contents (see Supporting Information Appendix D, Figure D2).This means that soil moisture variability at high mean soil moisture contents is less dependent on soil moisture variability at mean soil moisture contents lower than 0.40 m 3 •m −3 , that is, a large positive soil moisture deviation at one location in a dry to medium wet soil will not necessarily go hand in hand with such a large positive soil moisture deviation at that same location near saturation.These correlations confirm the expected high autocorrelation in soil moisture time series at fixed locations.
Third, soil moisture covariance was estimated based on soil moisture variability derived from water retention variability (RC).The covariance between two soil water potentials, or two mean soil water contents, was calculated using Equation (4), with a bin width of 0.01-0.0005m 3 •m −3 .As the retention curves were the same for both the depth-dependent (Z) and structure-dependent (S) root distribution, the same RC model applies to both scenarios.Figure 7A shows the covariance matrix that corresponds to an overall high correlation with a minimum correlation of 0.85 for a mean soil moisture content near saturation point (0.46 m 3 •m −3 ) in combination with lower soil moisture contents (see Supporting Information Appendix D, Figure D3).When comparing the RC prediction to the observed autocovariances (SM), the RC model is found to be a fair predictor (rRMSE(Z) = 29.63%),but tends to underestimate autocovariance (Figure 8, left), which is mainly because hydraulic nonequilibrium was not taken into account, since the RC model assumes uniform water potentials at a given depth.Interestingly, while high soil moisture levels showed good performance, the performance was lower for covariances between similar low soil moisture levels.This underestimation is even larger for the observed autocovariances of the structure-dependent (S) rooting scenario (Figure 8, right), where root density variability, resulting in variability in root water uptake, has a significant impact.Again, the performance was mostly low for covariances between similar low soil moisture levels, which is the range where the variable rooting has the largest impact.In practice, when applying this method to irrigated fields, soil moisture is not expected to reach such low soil moisture levels, and hence, good performance is assumed.In contrast, the correlation is mostly overestimated in both scenarios because of a general underestimation of soil moisture variability, as is shown in Figure 4.
Finally, soil moisture covariance was estimated based on the first-order approximation of the vG function (RC-FO) (Figure 7B).From Equation ( 9), the correlations could be F I G U R E 8 Scatterplots of predicted (RC) versus observed (SM) autocovariance by comparing the covariance of observed pairs of deviations from the mean soil moisture to the predicted covariance derived from retention curve variability for the depth-dependent (Z) and structure-dependent (S) root distribution.
expected to equal 1.In analogy with the predicted variance, soil moisture covariance becomes infinitely large at high soil moisture due to the  0 parameter.This results in a poor prediction (rRMSE >> 30%) with an overprediction of covariances at high soil moisture contents, while covariance is underpredicted at low soil moisture similar to the RC prediction.All (r)RMSE values are summarized in Supporting Information Appendix C (Table C2).

4.2
Validation: Case study approach

Soil heterogeneity
The water retention measurements and vG fitted curves are shown in Figure 9.For each soil water content, the standard deviation of the soil water content between the different water retention curves at a constant water potential was calculated and is shown by the green curve.In the literature, the shape parameter  is often assumed to be constant, or in other words, the soil is assumed to be macrosimilar (Sadeghi et al., 2016).In reality, the curves intersect (Figure 9), due to variability in the saturated soil moisture content in combination with variability in the  and/or  parameter.One should take note that no significant variability was observed in the residual water content  r since this parameter was always near zero.The vG model parameters seem strongly correlated (Figure 10).In literature, there is no consensus on correlation between vG parameters.Some report no significant correlation (Smith & Diekkrüger, 1996;Zhu & Mohanty, 2003), while others proclaim substantial correlations among vG parameters (Guber et al., 2004;Vereecken et al., 2010).Guber et al. (2004) stated that vG parameters  and  relate to the large and small aggregate contents in the soil, respectively.This supports the observed negative correlation between these two parameters in this study.On top of that, Guber et al. (2004) found that a decrease in saturated water content  s in silt loam and loam soils was accompanied by an increase of small aggregate contents.This is a plausible argument for the positive correlation between  s and  and the negative correlation between  s and  that was observed in 2021.However, the latter correlation was found to be positive in 2020.This could be explained by the difference in bulk density between the two years, indicating soil compaction in 2021, resulting in a decreased saturated water content and increased  with a higher degree of compaction.

Soil moisture variability
The standard deviation was predicted as a function of the mean soil water content and was validated with sensor measurements.Soil moisture variability was calculated for 2020 and 2021 based on the variability of the soil water retention curve (RC), and both the correlated (RC-FO; Equation 5) and the uncorrelated first-order model (RC-FO-uncorr; Equation 6).The three predictive models are shown together with the validation sensor data for both years (Figure 11).The range of the measurements was limited, having a wet growing season in 2021, and limited sensor measurements in 2020, mostly in the dry range.It should also be noted that the two growing seasons on the same field could not be described by the same model predictions since the soil heterogeneity was found to be temporally variable.Others also reported temporal variability in soil hydraulic properties in a field, both seasonal and interannual, such as bulk density, soil porosity, saturated hydraulic conductivity, and vG parameters (Alletto & Coquet, 2009;Jirků et al., 2013;Kargas et al., 2016).The predictive models showed reasonable results for the case study, both visually (Figure 11) and numerically, where RC showed the best predictions with an overall rRMSE of  16.09%, followed by RC-FO (rRMSE = 22.23%).Meanwhile, the uncorrelated first-order model-based variability (RC-FO-uncorr) showed poor predictions (rRMSE = 57.88%)as the variability was severely overpredicted.The difference between the first-order approximation models demonstrates the importance of vG model parameter covariances.The prediction skill of the three predictive models is shown for both growing seasons individually in Table 3.Even though the two growing seasons together resulted in a fair range of soil moisture levels, saturation point was never captured in the daily sensor measurements.Hence, the predictive models could not be validated for wet soil moisture observations and the results must be interpreted with caution.However, this should not be an issue as high soil moisture levels are of less importance for irrigation scheduling purposes.
A note of caution is due here since the observed soil moisture variability could also be caused by random sensor measurement noise, although one would expect higher observed standard deviations if this would be the case.It is also notable that the measured soil moisture variability sometimes even fails to reach the predicted variability, even though this prediction is the minimal expected variability.A possible explanation for this might be that the spread and number of sensors could not capture all soil heterogeneity in the field.Another possible explanation for this is an inherent measurement error of the sensor devices, which would result in a deviation of the mean soil moisture level.

Soil moisture autocovariance
In this section, the autocovariance model outputs are discussed and compared with each other.The soil moisture-based covariances derived from observations cannot be displayed as heatmaps since insufficient observation pairs covering the entire range of soil moisture combinations were available but are used as reference at the end of this section in a one-on-one comparison where each observation of a soil moisture pair is compared to the corresponding model estimation of that pair.First, the covariance matrices were calculated based on the water retention curves (RC) for both 2020 and 2021 and are shown as heatmaps in function of volumetric soil moisture (Figure 12A,B).Both case studies show similar patterns of positive and negative covariances.The covariances are generally larger for the 2021 case study.In 2020, negative covariances exist between high soil water contents close to saturation and moderate to low soil water contents.In contrast, positive covariances exist for all soil moisture contents when the moisture contents in the pairs does not differ more than 0.1 m 3 •m −3 or when both moisture content values are lower than 0.3 m 3 •m −3 .In Figure 9A, an intersection of the curves is visible around pF 1.7, equivalent to 0.34 m 3 •m −3 , and is characterized by a local minimum of the standard deviation.This intersection is also visible on the heatmap as a transition from positive to negative covariances.In 2021, also negative covariances exist between high soil water contents close to saturation and moderately high to low soil water contents, but the range of high soil water contents that are negatively correlated with lower ones is much narrower than that in 2020.
Another difference with 2020 is the larger covariance for intermediate water contents.In Figure 9B, an intersection of the curves is visible around pF 1.43, equivalent to 0.38 m 3 •m −3 , and is characterized by a local minimum of the standard deviation.This intersection is visible on the heatmap as well, as a transition from positive to negative covariances.In the dry range (pF > 4, or θ < 0.1 m 3 •m −3 ), covariances are negligibly small.Next, the covariance matrices were calculated based on the vG first-order approximation with and without considering vG parameter correlations (RC-FO and RC-FO-uncorr) for both 2020 and 2021 and are shown as heatmaps in function of volumetric soil moisture (Figure 12C-F).Again, both case studies show similar patterns of positive and negative covariances.It is already shown in Figure 11 that the standard deviation becomes much larger for the FO models for soil moisture contents larger than 0.4 m 3 •m −3 in 2020 and 0.3 m 3 •m −3 in 2021.These overestimations are reflected in the autocovariance matrices as well and can be identified as the dark area in the upper right corners.
The temporal autocorrelation of the time series (T-S and T-NS) was not computed for the case study, as T-S was already proven inappropriate in the theoretical approach.
Predicted and observed soil moisture covariances are compared visually and numerically.The performance of RC and RC-FO is shown in Figure 13.The model predicts a smaller range of variation in covariance than the observations, but the bias of the model-predicted covariance was small.The scatterplot reveals that the largest underestimation was observed for covariance predictions of low soil moisture levels, as was also observed in the virtual case (Figure 8), while covariance predictions at similar high soil moisture levels show good performance, and the model overpredicted the covariance between the high and low moisture contents.RC-FO showed the best predictions with a fair overall rRMSE of 28.04%, followed by RC (rRMSE = 34.38%).As expected, RC-FO-uncorr severely overpredicts soil moisture covariance (rRMSE > 100%), analogous to its soil moisture variability prediction.The prediction skill of the three predictive models is summarized for both growing seasons individually in Table 4.
As the ranges of observed moisture contents during each of the two growing seasons were small, the predictive models could only be validated for a limited set of soil moisture pairs; hence, the results must be interpreted with caution.

Soil moisture variability
The results of the retention curve-based models to estimate soil moisture standard deviation as a function of mean soil moisture content show that soil moisture variability can be predicted based on soil heterogeneity or, more specifically, soil properties and their variability, which can be quantified via undisturbed soil samples.This also indicates that soil heterogeneity is the main cause of soil moisture variability measured with sensors at multiple fixed locations across a field.However, as other factors might influence this variability, this prediction is only the minimum expected soil moisture variability and serves as a baseline.On top of soil heterogeneity, also crop and root growth variability can play a significant role, as was shown in this study by the virtual soil case with the structure-dependent (S) root distribution, where soil moisture variability at low mean soil moisture was around two times higher than for the deterministic (Z) root distribution scenario.
Next to soil and crop heterogeneity, other field factors such as hydraulic nonequilibrium, topography, soil management, nonuniform precipitation, and irrigation (Wilson et al., 2004)   can also contribute to soil moisture variability.Additionally, when measuring soil moisture with soil moisture sensors, an additional inherent measurement error must be kept in mind as well, as this is no part of the actual soil moisture variability in the field.Soil moisture variability can be derived from the variability of the soil water retention curves at a constant soil water potential (RC).RC is a good estimator (10% < rRMSE < 20%) for soil moisture variability assuming a horizontally uniform root distribution, but caution must be taken when predicting soil moisture variability in a field with large plant variability.The first-order approximation of the vG function taking vG parameter covariances into account (RC-FO) was found to be a good predictor for soil moisture variability at low to mid-range soil moisture levels, similar to RC, but it overestimates variability at high soil moisture levels.When vG parameter correlations were not considered (RC-FO-uncorr), soil moisture variability was overestimated.Hence, the difference between the two first-order approximation models emphasizes the importance of the correlation between vG model parameters.
The predicted variability curve based on the water retention curves (RC; Figure 4) is similar to other models and measurement data that are found in literature (Albertson & Montaldo, 2003;Famiglietti et al., 2008;Lawrence & Hornberger, 2007;Manns et al., 2014;Pan & Peters-Lidard, 2008;Rosenbaum et al., 2012;Schlüter et al., 2013;Teuling & Troch, 2005;Vereecken et al., 2007).For example, the following relation between standard deviation and mean soil moisture was derived by Famiglietti et al. (2008), as the coefficient of variation (CV =   ) in function of mean soil moisture can be assumed to follow an exponential curve: where  1 and  2 are fitting parameters determined by the CV model fit.This empirical soil moisture prediction curve has a similar shape compared to the predicted soil moisture based on the retention curve variability (RC), with zero variability at the minimum mean soil moisture, a maximal variability at mid-range soil moisture levels, and a decrease in variability toward high mean soil moisture levels.Such a model could be used as a hypermodel to estimate soil moisture variability.

Autocorrelation
Autocorrelated errors are expected in soil moisture time series measured by sensors at fixed locations.The measurement data in this study support this as the soil moisture measurement errors exhibited high autocorrelations.The virtual study even shows that if soil structure and its heterogeneity in the field remain the same, this autocorrelation could persist over years of time.
Soil moisture covariance can be predicted as a function of the spatial averaged soil moisture based on the variability of the soil water retention curves at a constant soil water potential (RC).RC was found to be a poor predictor (rRMSE ≥ 30%), as autocovariance might be underpredicted at extreme soil moisture levels because not all causes of heterogeneity were included in the model.The vG first-order approximation model (RC-FO) resembles the RC model at low mean soil water contents, but overestimates autocovariances at high soil water contents, resulting in misleading improved prediction results (20% < rRMSE < 30%).When vG parameter correlations are not considered (RC-FO-uncorr), soil moisture covariance was severely overpredicted (rRMSE > 100%).
Despite the low performance of the covariance prediction, which is mostly due to underprediction at low soil moisture levels, both the simulation data and field measurements, along with the model, show significant correlations and clearly indicate that neglecting autocorrelation would be an incorrect assumption.Although neglecting autocorrelation is the classic assumption in hydrological inverse modeling, data assimilation, and validation, the results of our study demonstrate that soil moisture covariance should be acknowledged.Ignoring or underestimating this covariance could result in an underestimation of the uncertainty of model predictions by a model that is calibrated to a time series of soil moisture measurements.Furthermore, the classic approach describing temporal autocorrelation of a soil moisture time series assumes stationarity and was proven inappropriate.An approach that describes the autocovariance as a function of the average soil moisture contents at the two different measurement times is more appropriate.In contrast to the covariance, the RC model overestimates the correlation between soil moisture measurements.This implies that the underestimation of the covariance by the model is due to an additional source of soil moisture variation that is not correlated to the estimated variation based on the soil water retention curve.

Assumptions and limitations
When predicting soil moisture variability and covariance using the proposed error models, some important assumptions and limitations need consideration.First of all, the error models predict soil moisture variability based on soil heterogeneity only and do not take into account the variability in root water uptake, even though it is known that both are linked closely (Clark et al., 2003;Wang & Smith, 2004).This can result in an underestimation of the actual (co)variance of soil moisture in a field, as was shown in this study for the virtual soil, when comparing soil moisture variability for the scenario of a depth-dependent (Z) and structure-dependent (S) root distribution (Figure 4).Teuling and Troch (2005) did consider land cover variability, in contrast with our simulations and assumptions.When plant variability is large and needs to be taken into account, that is, in a field with large plant variability or on a larger scale, the error model should be extended by including spatial variability in canopy cover (Lawrence & Hornberger, 2007).However, in irrigated fields, it is assumed that the canopy is uniform and closed.Additionally, variations in slope and aspect generate variations in microclimate and consequently variations in transpiration demand.Hence, this approach is valid for a uniform microclimate and flat topography, which is an acceptable assumption in our study area.
Second, hydraulic equilibrium is assumed, that is, soil water potential is assumed horizontally uniform in a field.This means that the error models do not consider soil moisture variability resulting from preferential flow and nonequilibrium (Jarvis, 2007;Vogel et al., 2010).The prevailing scenario where nonequilibrium conditions may occur is infiltration of water into relatively dry, heterogeneous soils due to an irregular water front and low hydraulic conductivities (Vogel et al., 2010).
Assuming hydraulic equilibrium implies that variability in hydraulic conductivity is not considered, nor is the correlation between the saturated hydraulic conductivity and vG model parameters.It has, however, been shown that considerable variation might exist in the unsaturated hydraulic conductivity function (Vogel et al., 2023), and moreover, the saturated hydraulic conductivity and shape parameter  are correlated (Zhu & Mohanty, 2003).Qu et al. (2015) also found that ln( s ) is the second most sensitive parameter to affect soil moisture variability.Therefore, hydraulic equilibrium at a certain depth is an important assumption, which could result in an underestimation of soil moisture variability.However, Schlüter et al. (2013) already showed that the horizontal variability in soil water potential was generally low for the virtual soil at 0.1 m depth, as used in this study (CV < 0.5), which ratifies this assumption.Additionally, in this study, the predicted soil moisture variability of the virtual soil without root water uptake variability (Z), based on the retention curve variability (RC), was only slightly underestimating the observed soil moisture variability, indicating only a small impact of soil water potential variability.Qu et al. (2015) also assumed hydraulic equilibrium and showed that this is an acceptable assumption for variability predictions.
Third, variability in the residual water content was not included in this study as it was not observed in our data.Qu et al. (2015) also assumed a constant residual soil water content.However, if the residual water content varies, its variability dominates soil moisture variability at low water contents and overrules violations of the hydraulic equilibrium assumption.
Finally, the error models assume soil heterogeneity to be static, even though other field studies showed seasonal changes in soil structure (Alletto & Coquet, 2009;Jirků et al., 2013;Kargas et al., 2016;Schlüter et al., 2011).This could mean that the error models become less accurate during the growing season and would be more difficult to parameterize.This temporal variability exists on two timescales, namely, (1) in between growing seasons (months-years) and (2) within a growing season (days-weeks) (Jirků et al., 2013;Kargas et al., 2016).Differences in soil properties in between growing seasons may be due to soil management, such as tillage and soil nutrients, wet conditions at the start of the growing season, erosion, and land use (crop rotation, cover crop).Jirků et al. (2013) found that climatic conditions, mostly during the winter and spring, played the leading role.This temporal variability also occurred in this study, as a significant change was observed between water retention curves in 2020 and 2021, even though they were measured on the same field, in the same period (July).The largest variation was found in the saturated soil water content and hence the bulk density, which could be expected to increase due to compaction as a result of the wet conditions in 2021.These findings are in line with literature, as soil properties, in particular bulk density, saturated hydraulic conductivity, and vG model parameters, were found to be subjected to temporal variation that often has a greater impact compared to their spatial variation (Feki et al., 2018).
On the other hand, soil properties may also be dynamic during a growing season.Factors such as soil organisms, nutrients, root formation, erosion, and periods of drought or water excess may have an impact on soil aggregates, water retention and conductivity, bulk density, pore distribution, and preferential flow (Alletto & Coquet, 2009;Jirků et al., 2013;Kargas et al., 2016;Schlüter et al., 2011).Another important and well-known factor is hysteresis of the water retention curve.More specifically, assuming a closed-loop hysteresis, the shape parameter  is the only variable during wetting and drying cycles ( drying ≤  wetting ) (Dohnal et al., 2006).However, Dohnal et al. (2006) found that the impact of hysteresis of the water retention curve on soil water dynamics is not significant and ignoring hysteresis did not lead to a substantial deviation of the soil water model predictions.Qu et al. (2015) implied that soil moisture variability is most sensitive to the pore distribution represented by shape parameter  in the vG model, which is expected to remain constant during wetting and drying cycles (Dohnal et al., 2006).However, shape parameters  and  relate to large and small aggregate contents in the soil, respectively (Guber et al., 2004), and were found to correlate negatively, which implies that both may have a significant impact on soil moisture variability.

Application in practice
The errors of which we calculated the covariance refer to the errors of the "true" field-scale averaged water contents at a certain depth that are estimated from measurements by a single sensor that stays at a fixed position.This error covariance can be reduced by installing more sensors.If the sensors are sufficiently far apart, the spatial covariance between the measurements is close to zero and the measurements by the different sensors can be assumed to be independent.The covariance of the errors of the true spatial mean that are estimated from the sample mean can then simply be obtained from dividing the single measurement error covariance by the number of sensors or the sample size n.
In practice, only a limited number of sensors are installed in the field, mainly due to the cost and management of the sensors, resulting in estimates of (co)variance that might not always be accurate.Based on the desired marginal error, an equation was proposed by Gilbert (1987) and Ott and Longnecker (1988) to determine the number of point measurements (N) that is required to accurately describe the mean value and spatial variability (Brocca et al., 2010): where Z is taken from the standard Z table for  = 0.05 (95% confidence interval),  is the standard deviation, and ME is the marginal error.For the case study, low soil moisture variabilities were predicted, indicating only two to four sensors would be required to obtain an ME = 10% of the mean value.However, these numbers vary significantly depending on the soil heterogeneity in a field.
To apply the proposed RC model in practice, accurate and detailed soil data are required, that is, a set of water retention curves, preferably resulting from multiple undisturbed soil samples throughout the field.As an alternative, retention curves could be approximated by using soil map information and/or scaling methods.The first-order approaches (RC-FO) do not necessarily require a set of water retention curves but might be used when parameters could be derived from other soil information.Finally, for both models, historical field data, for example, from previous years or from databases, will probably not suffice.We observed that soil water retention curves varied from year to year due to variations in soil compaction after tillage.Therefore, real-time data are preferred.
In this study, the cost of extra yearly soil core sampling was lower than extra sensor modules, considering a 3-year sensor usage period and repeated sampling.Additionally, installing more sensors in a field is not desired in practice by farmers, as it complicates field operations.However, the applicability of the method should be assessed on a case-by-case basis, depending on the goal, location, costs, budget, and so forth.
The RC and RC-FO methods also have the potential to be upscaled to a remote sensing resolution.When using in situ soil moisture data for validating or assimilating remote sensing data, the (co)variability within a pixel is unknown and measurement error autocorrelation is neglected.With the proposed soil heterogeneity-based method, one could compute the (co)variability of the mean soil moisture estimates based on the soil types within this pixel, while retention curves can be derived from soil maps and pedotransfer functions, analogous to the method of Montzka et al. (2017).However, on this scale, the assumption of horizontally uniform water potential head may not be a valid assumption, and the plant variability will be much larger, along with other factors such as topography, groundwater depth, and spatial differences in precipitation and irrigation, possibly resulting in a severe underestimation of soil moisture (co)variability.Nonetheless, it is likely preferable to consider an underestimated temporal covariance of the error of the estimated mean soil moisture from in situ sensors than to completely disregard it.

CONCLUSIONS AND OUTLOOK
This study provides a new perspective on the relation between soil moisture variability and soil heterogeneity and sheds light on the covariance between soil moisture deviations at fixed locations in a field as obtained with soil moisture sensors.In this study, a mechanistic error model was proposed to predict soil moisture variability and covariance based on the variability of the vG soil water retention curve under the assumption of a uniform soil water potential in the field.The proposed models were validated both in a theoretical and a practical setting.With such an error model, sensor measurements and their uncertainty can be interpreted with more accuracy.Using the sensor measurements directly to estimate their autocovariance could be an alternative to using an error model, but only if sufficient sensors are placed in the field.In practice, only a limited number of sensors are installed in the field, which might not give us an accurate estimate of (co)variance.
The main conclusions are as follows: 1.The variability of soil water retention curves can be used to predict soil moisture variability assuming a uniform water potential.Using this variability directly (RC) resulted in good predictions (10% < rRMSE < 20%), while the vG first-order approximation model (RC-FO) predicts similar variabilities at low soil moisture but is unreliable at high soil moisture levels.2. The predicted soil moisture variability resulting from soil heterogeneity is the minimum expected variability in a field.On top of this, plant variability may have a significant impact on soil moisture variability, as well as other factors such as nonuniform irrigation, subsoil, nonequilibrium, and so forth.3. Measured soil moisture time series demonstrate autocorrelation of the deviation of a local sensor measurement from the field-scale average.The classic approach describing temporal autocorrelation of soil moisture deviations assuming stationarity was proven inappropriate (rRMSE = 87%).Instead, a mechanistic approach was presented, as the variability of soil water retention curves can also be used to predict soil moisture covariance.This approach predicts the autocovariance as a function of the average soil moisture contents at the two observation times and not as a function of the time lag between observations and could describe the simulated and measured covariances better than a stationary covariance function.However, using RC variability directly resulted in poor predictions (rRMSE ≥ 30%), as covariance was mostly underpredicted.The vG first-order approximation model (RC-FO) somewhat overpredicted covariances resulting in misleading fair predictions (20% < rRMSE < 30%) and is unreliable for high soil moisture levels, where covariance was strongly overpredicted.4.This study emphasizes the importance of the correlation between vG model parameters.When predicting variability with the vG first-order approximation model, it is clear that soil moisture variability as well as covariance is overpredicted when vG parameter covariances are not considered.
The results of this study contribute to more reliable utilization of in situ soil sensing data.A realistic representation of spatial soil moisture variability has the opportunity to improve predictions of hydrological land surface processes such as evapotranspiration and runoff.When soil moisture measurements are used to calibrate such hydrological models, it is necessary to take into account the uncertainty of the measurements as well as the correlation between all measurement errors.The proposed mechanistic error models are a foundation to estimate soil moisture (co)variability based on soil heterogeneity, more specifically on the variability in the vG water retention curve in a direct approach (RC; Equation 4) as well as based on a first-order approximation (RC-FO; Equation 5).The error models result in covariance matrices in function of mean soil moisture, from which one can easily derive the measurement covariance matrix of measurement time series in that field.To parametrize these models, water retention analyses of multiple undisturbed soil cores are required to know the variability in the vG retention curve.
In conclusion, with the proposed mechanistic error models, we can better understand and interpret soil moisture sensor data measured at fixed locations and soil moisture variability resulting from soil heterogeneity.However, the mechanistic error model tended to underestimate the variance and autocovariance.Further research is needed to evaluate the effect of this underestimation on prediction uncertainty of models that are calibrated using in situ soil moisture measurements.In order to derive the autocovariance with a mechanistic error model, water retention curves should be measured on a large number of soil cores, and this should be done every year since soil properties are dynamic.Further research is needed to develop approaches that would either derive autocovariances from direct soil moisture measurements or combine information about measured variability with information about water retention curves, which could also be inferred from other proxies that are easier to measure or map.However, when insufficient information about soil variability is available or when the accuracy of the estimated covariance functions is insufficient, then a covariance function that is a bivariate function of average soil moisture contents could be tested with hyperparameters that could be derived from deviations between predictions by process models and observations.The functional form of this covariance function could be based on covariances that are estimated from the variability of retention curves.

AU T H O R C O N T R I B U T I O N S
Marit G. A. Hendrickx: Conceptualization; data curation; formal analysis; funding acquisition; investigation; methodology; project administration; resources; software; supervision; validation; visualization; writing-original draft.Jan Diels: Conceptualization; data curation; funding acquisition; methodology; project administration; resources; supervision; writing-review and editing.Pieter Janssens:

F
I G U R E 2 Soil heterogeneity: pF (black, left axis) and soil moisture standard deviation (green, right axis) predicted based on the variability of the soil water retention curve (RC) as a function of mean soil moisture.F I G U R E 3 Ten-year time series of the mean soil moisture at 0.1 m depth in the virtual field with the depth-dependent root distribution (A) and soil moisture deviations from the mean soil moisture shown for eight independent sample locations (B).

F
I G U R E 4 Soil moisture standard deviation as a function of mean soil moisture: based on retention curve variability RC (solid line), based on the first-order van Genuchten (vG) model approximation RC-FO (dashed line), and observed variability (SM) for the depth-dependent (Z) and structure-dependent (S) root distribution.

F
I G U R E 5 T-S: Temporal correlation function (A) and autocovariance matrix (B) for the depth-dependent (Z) root distribution.Time lags Δ > ∕2 are undefined.In comparison, the nonstationary temporal autocovariance matrix (T-NS) is also shown for the depth-dependent (Z) root distribution (C).

F
I G U R E 6 SM: Covariance matrices in function of mean soil moisture contents based on virtual observed soil moisture variability for the depth-dependent (Z) and structure-dependent (S) root distribution, after applying a Gaussian filter with a standard deviation of 10 soil moisture steps.F I G U R E 7 RC: Covariance matrix in function of mean soil moisture contents based on retention curve variability (A).RC-FO: Covariance matrix in function of mean soil moisture contents based on the first-order approximation of the van Genuchten function (B).The soil moisture range of the observations of the depth-dependent (Z) root distribution is indicated (---).

F
Measured soil heterogeneity: pF (black, left axis) and soil moisture standard deviation (green, right axis) as a function of mean soil moisture (RC) for 2020 (A) and 2021 (B).F I G U R E 1 0 Correlation heatmap of the van Genuchten (vG) model parameters of 2020 (A) and 2021 (B).

1 2
The covariance matrices in function of mean soil moisture contents resulting from the predictive autocovariance models (RC, RC-FO, and RC-FO-uncorr) for 2020 (left) and 2021 (right).T A B L E 4 Prediction skill of soil moisture covariance models (ACOV: RC, RC-FO-uncorr, and RC-FO) for two growth periods (2020 and 2021): root-mean-square error (RMSE) and relative root-mean-square error (rRMSE) of covariances, and RMSE of correlations.

F
Scatterplots of predicted (RC [A] and RC-FO [B]) versus observed (SM) covariances for 2020 by comparing the covariance of observed pairs of deviations from the mean soil moisture to their corresponding predicted covariances.
Overview of error models to estimate the standard deviation (STDEV) and autocovariance (ACOV) of soil moisture measurements.
T A B L E 1
T A B L E 3