Time Stability of Soil Water Content

1. Introduction

Soil water is a key variable controlling water and energy fluxes in soils (Vereecken et al., 2007). It is necessary for plant and vegetation growth and development. Research has indicated that soil water content (SWC) varies both in space and time. Variations in both space and time present a substantial challenge for applications such as precision agriculture and soil water management.

Since the contribution of Vachaud et al. (Vachaud et al., 1985), a large body of research has indicated the presence of time stability of SWC (Biswas and Si, 2011c; Comegna and Basile, 1994; Grayson and Western, 1998; Hu et al., 2009; Martínez-Fernández and Ceballos; 2005; Mohanty and Skaggs, 2001), which means that the spatial pattern of SWC does not change with time at a certain probability. According to this concept, if a field is repeatedly surveyed for SWC, there is a high probability that a location with certain wetness characteristics (i. e., wet, dry, intermediate) will maintain those characteristics on subsequent occasions. Time stability has also been extended to describe the characteristics of SWC at point scales. A location will be regarded as time stable provided it can estimate the average SWC of an area.

Time stability of SWC has been observed at a large variety of scales ranging from plot (Pachepsky et al., 2005) to region (Martínez-Fernández and Ceballos, 2003) and related studies cover a range of investigated areas, sampling schemes, sampling depths, investigation periods, and land uses (Biswas and Si; 2011c; Brocca et al., 2009; Cosh et al., 2008; Hu et al., 2010a;Tallon and Si, 2004; Vachaud et al., 1985). As a result, a variety of methods have been developed to evaluate time stability of SWC, each with its own advantages and disadvantages. Time stability is usually used to characterize time persistence of the spatial pattern of SWC between measurement occasions, either at the measurement scale or at different scales (Biswas and Si, 2011c; Cosh et al., 2006; Kachanoski and De Jong, 1988; Vachaud et al., 1985). On the other hand, SWC at time stable locations can be used to estimate the spatial average SWC of a given area (Grayson and Western, 1998). Therefore, quick and accurate estimation of soil water content at a field or large areas may be possible with only one representative location. Areal estimation of SWC from point source data has the potential to substantially reduce both the capital and labor costs involved in estimating average SWC, making the method appealing to a wide range of disciplines.

Time stability of SWC is controlled by various factors including soil, topography, vegetation, and climate (Brocca et al., 2009; Gómez-Plaza et al., 2001; Grayson and Western, 1998; Hu et al, 2010a; Tallon and Si, 2004). Information on the controls of time stability provides an essential insight into the mechanisms of soil water movement and storage. In addition, it is also important to identify the time stable location for average SWC estimation a priori.

The goal of this chapter is to provide a comprehensive review on the study of time stability of SWC. In doing so, we have introduced the concept, methodology, application, and the controlling factors step by step.

2. Concept

Soil water content varies with time and space. However, if a field is repeatedly surveyed for SWC, there is a high probability that a location with a certain relative wetness condition (i. e., wet, dry, and intermediate) at a given time will remain at the same relative wetness condition at other times. This phenomenon was first explained as time stability by Vachaud et al. (1985). This can be defined as the “time invariant association between spatial location and classical statistical parametric values” of SWC, most often the mean.

Time stability of SWC can be divided into two types, one is to describe the overall similarity of the spatial pattern between measurements, and the other is to describe the time-invariance of the relative SWC of a given location.

Repeated survey indicated that some locations are always wetter or drier than the average SWC of an area, resulting in the preservation of their ranks. Hence, the spatial pattern of SWC measured at one time will be similar to those measured on subsequent occasions. Statistically, this kind of time stability indicates that the orders, or ranks, of SWC at different locations do not change over time at some probability. If all locations were to maintain their rank on subsequent measurement occasions, then the spatial pattern of SWC will be identical. From this point, Chen (2006) augured that the term of “rank stability” or “order stability’’ should be better than “time stability” because “the stability is the order or rank of a soil property at different spatial points that does not change at some probability”.

The most widely used method to describe time stability of SWC spatial patterns is Spearman’s rank correlation analysis (Vachaud et al., 1985), while Pearson correlation analysis has also been widely used (Cosh et al., 2004). This is because Pearson correlation measures strictly linear relationships, while Spearman Rank correlation also measures nonlinear correlation between SWC values measured at two different occasions. Time stability of the spatial pattern of SWC is scale dependent because of the interaction between measurement locations (Kachanoski and De Jong, 1988). Scale dependence of the time stability was also investigated by different methods including spatial coherency analysis (Kachanoski and De Jong, 1988), wavelet coherency analysis (Biswas and Si, 2011c), and multivariate empirical mode decomposition (Rehman and Mandic, 2010). A detailed introduction of these methods is given in Section 3.

For time stability at point scales, SWC at each location is usually scaled in terms of the field average (Vachaud et al., 1985). This indicates that at a given location, the change in scaled SWC exhibits little change over time. Many indices are available to examine the degree of time stability at point scales. These indices include standard deviation of relative difference (Mohanty and Skaggs, 2001; Schneider et al., 2008; Vachaud et al., 1985), root mean square error (Bosch et al., 2006; Jacobs et al., 2004), width of the 90% empirical tolerance interval of relative water content (Guber et al., 2008), chi-squared statistic (Guber et al., 2008), root-mean-squared differences (Guber et al., 2008), temporal coefficient of variability (Starr, 2005), and mean absolute bias error (Hu et al., 2010a, 2010b). One of the applications of time stability at point scales is to identify catchment average soil moisture monitoring (CASMM) site (Grayson and Western, 1998) or benchmark location (Tallon and Si, 2004) for average SWC evaluation. Detailed introduction of these indices and applications is given in Section 3.

Generally, these two types of time stability are correlated. If the spatial pattern is time stable, there is a larger possibility of time stability at point scales. On the other hand, if more points are time stable, there is a larger possibility of time stability of spatial pattern. However, these two types of time stability are also distinguished for the following three aspects: (1) time stability of spatial pattern is used to evaluate the time stability of an area, while time stability of relative SWC is used to evaluate the time stability of a point. Therefore, time instability of spatial pattern does not mean the absence of time stable points (Grayson and Western, 1998; Schneider, 2008). On the other hand, if no points are time stable, it does not mean that spatial pattern is time unstable for any two measurement occasions; (2) they have different evaluation criteria as mentioned above, and (3) they have different applications. Time stability of spatial pattern is used to qualitatively describe the similarity of spatial distribution of SWC between different measurement occasions for the better understanding of soil water related process and influencing factors (Lin; 2006), while time stable locations can be selected to estimate the average SWC for an area (upscaling) by directly using the soil water content of the time stable location (Brocca et al., 2009; Cosh et al., 2008; Grayson and Western, 1998; Martínez-Fernández and Ceballos, 2005) or by considering the offset between the mean value and that of the time stable location (Grayson and Western, 1998; Hu et al., 2010b; Starks et al., 2008).

3. Methodology

There are different methods for analyzing time stability. We have documented the methods in two sections: (1) time stability of spatial patterns and (2) time stability of points.

3.1. Time stability of spatial patterns

3.1.1. Pearson correlation analysis

Pearson correlation coefficient is sensitive only to a linear relationship between two variables. A Pearson correlation coefficient, r_j,j’, between two spatial series of soil water content measured at time j and j’ can be defined by

$$r_{j,j^{\text{'}}}=\frac{{\displaystyle \sum _i(SW{C}_j\left(i\right)-⟨SW{C}_j⟩)(SW{C}_{j^{\text{'}}}\left(i\right)-⟨SW{C}_{j\text{'}}⟩)}}{\sqrt{{\displaystyle \sum _i{(SW{C}_j\left(i\right)-⟨SW{C}_j⟩)}^2}}\sqrt{{\displaystyle \sum _i{(SW{C}_{j^{\text{'}}}\left(i\right)-⟨SW{C}_{j\text{'}}⟩)}^2}}}$$E1

where ${ SWC}_j\left(i\right)$ and $⟨{SWC}_j⟩$ are soil water content at location i and spatial average at a given time j, respectively. ${ SWC}_{j^{\text{'}}}\left(i\right)$ and $⟨{SWC}_{j\text{'}}⟩$ are soil water content at location i and spatial average at another time $j^{\text{'}}$, respectively. The resulting coefficients refer to the correlation of the spatial patterns of SWC from one time to another. It is expected that closely correlated patterns have a $r_{j,j\text{'}}$ near one, while uncorrelated patterns are indicated by r_j,j’ values near zero (Cosh et al., 2004). Student’s t is generally used to test the significance of Pearson correlation coefficient. It can be implemented by some statistical software, such as Excel, SPSS, MathCad, Matlab, and SAS.

3.1.2. Spearman’s rank correlation analysis

Spearman's rank correlation coefficient is a non-parametric measure of statistical dependence between two variables. It assesses how well the relationship between two variables can be described using a monotonic function. It is the Pearson correlation between the ranks of one series and the ranks of another series. Because ranking linearize some of the nonlinear relationships, it is sensitive to nonlinear relationships. Let R_ij be the rank of soil water content ${SWC}_j \left(i\right)$ and R_i^’_j the rank of ${SWC}_{j\text{'}} \left(i\right).$ The Spearman rank correlation coefficient, r_s, is calculated by

$$r_s=1-\frac{6{{\displaystyle \sum _{i=1}^n\left(R_{ij}-R_{i{j}^{\text{'}}}\right)}}^2}{n(n^2-1)}$$E2

where n is the number of observations. A value of r_s=1 corresponds to identity of rank for any sites, or perfect time stability between time j and j’. The closer r_s is to 1, the more stable the spatial pattern will be. Student’s t can be used to test the significance of r_s (Zar, 1972). Software such as Excel, SPSS, MathCad, Matlab, and SAS can implement this analysis.

3.1.3. Average spatial coefficient of determination and temporal coefficient of variability

If perfectly time stable soil water content pattern exists, the ratio of ${SWC}_j \left(i\right)$ to $⟨{SWC}_j⟩$, can be a time independent scaling factor, $r \left(i\right),$ which is expressed as (Starr, 2005)

$$\frac{SW{C}_j\left(i\right)}{⟨SW{C}_j⟩}=r(i)$$E3

Complete time stability is usually disturbed by a series of processes and therefore Eq. (3) needs an extra term to fit a more general situation (Starr, 2005)

$$\frac{SW{C}_j\left(i\right)}{⟨SW{C}_j⟩}=r(i)+e_j(i)=s_j(i)$$E4

where $s_j \left(i\right)$ is a scaled factor depending on both time and location, and $e_j \left(i\right)$ is an additional term that accounts for random measurement errors, random sampling variability, and any true deviations from time stability. Assuming that $e_j \left(i\right)$ and $r \left(i\right)$ are independent, then the variance of $s_j \left(i\right)$ over space at any time, ${\gamma }_i^2\left[s_j \left(i\right) \right]$, can be written as

$${\gamma }_i{}^2\left[s_j(i)\right]={\gamma }_i{}^2\left[r(i)\right]+{\gamma }_i{}^2\left[e_j(i)\right]$$E5

where ${\gamma }_i^2\left[r \left(i\right) \right]$ and ${\gamma }_i^2\left[e_j \left(i\right) \right]$ are the variances of $r \left(i\right)$ and $e_j \left(i\right),$ respectively over space. The ${\gamma }_i^2\left[s_j \left(i\right) \right]$ can be calculated by averaging the calculated variances for each observation. Then, the average spatial coefficient of determination, $R_s^2$, which describes the proportion of the total variance explained by the time stability model of Eq. (3) can be written as

$$R_s^2=\frac{{\gamma }_i{}^2\left[r(i)\right]}{⟨{\gamma }_i{}^2\left[s_j(i)\right]⟩}$$E6

where $⟨{\gamma }_i^2\left[s_j \left(i\right) \right]⟩$ is the mean value of ${\gamma }_i^2\left[s_j \left(i\right) \right]$ over time. Given $e_j \left(i\right)$ with an average value of zero over time, $r\left(i\right)$ can be approximated by time average of $s_j \left(i\right),⟨s_j \left(i\right) ⟩$, i. e.,

$$r(i)=⟨s_j(i)⟩$$E7

The value of $R_s^2$ approaches one for a complete spatial pattern described by Eq. (3) and zero for a situation where no variances can be explained by the time stability model. Therefore, greater $R_s^2$ means stronger time stability of spatial pattern. However, Starr (2005) pointed out it is not suitable for the situation where no spatial variability of SWC and completely time stability exists.

To deal with this issue, temporal coefficient of variability, ${CV}_t$, was developed to describe the degree of time stability (Starr, 2005) and can be written as

$$C{V}_t=100\frac{1}{n}{\displaystyle \sum _{i=1}^n\frac{{\sigma }_t\left[s_j(i)\right]}{⟨s_j(i)⟩}}$$E8

where ${\sigma }_t\left[s_j \left(i\right) \right]$ is the standard deviation of $s_j \left(i\right)$ over time, $⟨s_j \left(i\right) ⟩$ is the time average of $s_j \left(i\right).$ The smaller ${CV}_t$ value indicates stronger time stability. The ${CV}_t$ approaching zero indicates a perfect time stable pattern.

3.1.4. Empirical orthogonal function method

Empirical orthogonal function (EOF) analysis is a decomposition of a data set in terms of orthogonal basis functions which are determined from the data. It is the same as performing a principal components analysis on the data, except that the EOF method finds both time series and spatial patterns. The basic principle of EOF method is to partition a series into time-invariant spatial patterns (EOFs) of SWC and coefficients (ECs) which vary temporally but is constant spatially (Perry and Niemann, 2007). The original spatial series of SWC can be obtained by taking the sum of product of EOFs by ECs. A limited number of EOFs that present significant spatial variation of SWC can be selected and used to identify the dominant factors determining the spatial pattern.

First, the spatial anomalies of SWC are computed. The spatial anomaly at location $i$ and time $j$, ${ z}_j \left(i\right),$ can be calculated as

$$z_j\left(i\right)=SW{C}_j\left(i\right)-\frac{1}{n}{\displaystyle \sum _{i=1}^{n}SW{C}_j\left(i\right)}$$E9

Considering the spatial anomalies, the spatial covariance $v_{j,j\text{'}}$ at time $j$ and $j\text{'}$ can be calculated as

$$v_{j,j^{\text{'}}}=\frac{1}{n}{\displaystyle \sum _{i=1}^n{z}_j}\left(i\right)z_{j^{\text{'}}}\left(i\right)$$E10

where ${ z}_{j\text{'}} \left(i\right)$ is the spatial anomaly at the same location i but time $j\text{'}$.

To consider all the measurement times, the matrix of spatial anomalies, Z, and its spatial covariance, V, can be written as

$$Z=\left(\begin{array}{ccc}{z}_{11}& \dots & z_{1m}\\ ⋮& \ddots & ⋮\\ z_{n1}& \cdots & z_{nm}\end{array}\right)$$E11

and

$$V=\frac{1}{n}{Z}^{T}Z$$E12

respectively, where m is the number of sampling times.

The next step is to diagonalize the spatial covariance matrix V by finding its eigenvectors, E, and eigenvalues, L. Mathematically, they should satisfy the equation

$$VE=LE$$E13

where E is an m×m matrix that contains the eigenvectors as columns

$$E=\left(\begin{array}{ccc}{e}_{11}& \cdots & e_{1m}\\ ⋮& \ddots & ⋮\\ e_{m1}& \cdots & e_{mm}\end{array}\right)$$E14

and L is an m×m matrix that contains the associated eigenvalues along the diagonal and zeros at off-diagonals

$$L=\left(\begin{array}{ccc}{l}_{11}& \cdots & 0\\ ⋮& \ddots & ⋮\\ 0& \cdots & l_{mm}\end{array}\right)$$E15

The eigenvectors in E represent the weights applied to each component in V to diagonalized V. This transformation is a rotation of the original axes in multi-dimensional space, where each dimension corresponds to a sampling time. Axes are orthogonal to each other and explain different amount of covariance in the spatial anomaly dataset. The eigenvalues contained in L represent the variance that occurs in the direction of each newly identified axis. After diagonalization of V, E and L are arranged accordingly to keep the eigenvalues in L sorted in a descending order. Therefore, the first axis explains the most covariance in the spatial anomaly dataset, and the second axis explains the second most covariance and so on. Each axis represents a direction in the multidimensional space. The total variance of the spatial anomaly data is the sum of the diagonal values in L and is equal to the total variance of the original spatial anomalies. Therefore, the portion of the variance, $P_j$, that the j^th axis can explain is

$$P_j=\frac{l_{jj}}{{\displaystyle \sum _{k=1}^m{l}_{kk}}}$$E16

The relative importance of each axis in explaining the variability of spatial anomalies can be defined by the relative magnitude of the associated eigenvalue in L.

Each spatial anomaly can be described in terms of the new variable axis. A matrix F containing the coordinates of the spatial anomalies on the new axis can be obtained by projecting the spatial anomalies onto the rotated axis. Mathematically, this operation is

$$F=ZE$$E17

Each column in F is EOF, and corresponding column in E is expansion coefficient (EC). From the aspect of time stability application, usually a limited number of EOFs that explains significant amount of the spatial variability of the original SWC series is selected to present the underlying time-invariable spatial pattern. North et al. (1982) indicated a set of selection criterion. An EOF is considered significant provided the lower confidence limit of its eigenvalue is larger than the upper confidence limit of the eigenvalue of the next most important EOF (North et al., 1982).

3.1.5. Spectral coherency analysis

Spectral coherency is used to measure the similarity between two spatial series in a frequency domain, which can be converted to a spatial scale. This method assumes the spatial series to be stationary and linear. Therefore, it can be used to find the scale specific information on time stability of spatial pattern between different measurements.

Spectral coherency analysis involves calculation of the power spectrum of a variable $V_t^2\left(f_K\right)$ which estimates the observed variance as a function of spatial scale.

$$V_t^2\left(f_K\right)={\left(2h+1\right)}^{-1}{N}^{-1}{\displaystyle \sum _{l=-h}^h\left[{\left({\displaystyle \sum _{j=1}^{N}SW{C}_t\left(j\right)\cos \left(2\pi f_{K+1}j\right)}\right)}^2+{\left({\displaystyle \sum _{j=1}^{N}SW{C}_t\left(j\right)\sin \left(2\pi f_{K+1}j\right)}\right)}^2\right]}$$E18

where $f_K=K/N$, K = 0, 1, 2, …, N/2 cycles h^-1, and h is the smoothing coefficient determining the degree of averaging of adjacent independent frequencies and the degrees of freedom for the individual spectral variance estimates. The value of h cannot be too large because power spectrum of a variable may incorporate more bias from smoothing of the spectra, although larger value of h will result in smaller variance of the estimate (Gómez-Plaza et al., 2000).

The covariance between ${ SWC}_j \left(i\right)$ and ${SWC}_{j\text{'}} \left(i\right)$ can be estimated as a function of spatial scale by calculating sample cross spectrum $V_{j,j\text{'}}\left(f_K\right)$

$$V_{j,j^{\mathrm{\text{'}}}}\mathrm{ }\left(f_{{}_K}\right)\mathrm{ }={\mathrm{ }(2h+1)\mathrm{ }}^{-1}{N}^{-1}\sum _{l=-h}^h{F}_j{\mathrm{ }\left(f_{K+1}\right)\mathrm{ }}^{\mathrm{*}}{F}_{j^{\mathrm{\text{'}}}}\mathrm{ }\left(f_{K+1}\right)\mathrm{ }$$E19

where the asterisk (*) is a complex conjugate of

$$F_j\left(f_{{}_K}\right)={\displaystyle \sum _{i=0}^{N-1}SW{C}_j}\left(i\right)\left[\cos \left(2\pi f_{K}i\right)-i\sin \left(2\pi f_{K}i\right)\right]$$E20

$$F_{j^{\text{'}}}\left(f_{{}_K}\right)={\displaystyle \sum _{i=0}^{N-1}SW{C}_{j^{\text{'}}}}\left(i\right)\left[\cos \left(2\pi f_{K}i\right)-i\sin \left(2\pi f_{K}i\right)\right]$$E21

where $f_K=K/N$, K = 0, 1, 2, …, N/2 cycles h^-1, and i = (-1)^0.5.

Spectral coherency function $R_{j,j\text{'}}^2\left(f_K\right)$ can then be calculated in a similar manner to calculate coefficient of determination (Kachanoski and De Jong, 1988)

$$R_{j,j^{\text{'}}}{}^2\left(f_{{}_K}\right)=\frac{{\left|V_{j,j^{\text{'}}}\left(f_K\right)\right|}^2}{V_j{}^2\left(f_K\right)V_{j^{\text{'}}}{}^2\left(f_K\right)}$$E22

The coherency function $R_{j,j\text{'}}\left(f_K\right)$ estimates the proportion of the spatial variance of ${SWC}_j \left(i\right)$ which can be explained by the spatial variance of ${SWC}_{j\text{'}} \left(i\right),$ as a function of spatial scale. Thus, spectral coherency is useful to measure time stability of soil water content as a function of spatial scale.

The significance test method for spectral coherency consist of parametric methods based on an assumed theoretical distribution and nonparametric method such as reshuffling, bootstrapping, and bagging method. We encourage readers to see Si (2008) for a more detailed account on significance test.

3.1.6. Wavelet coherency analysis

Wavelet analysis, differing with spectral analysis, can identify localized features of soil processes (Si, 2008). It is suitable to reveal the scale and location specific time-persistence of spatial pattern of SWC between sampling occasions assuming the soil water content system is linear (Biswas and Si, 2011c). Wavelet coherency of two spatial series can describe the time stability of spatial pattern. It requires the calculation of wavelet coefficients for each of the two data series and associated cross-wavelet spectrum. Many publications on the introduction of wavelet analysis (Farge, 1992; Kumar and Foufoula-Georgiou, 1993, 1997) and wavelet coherency (Biswas and Si, 2011c; Grinsted et al., 2004; Si, 2008; Si and Zeleke, 2005) can be found. Here, we will present the basic procedure to calculate wavelet coherency.

First, wavelet coefficient, $W_i^Y\left(s\right)$, is calculated with the continuous wavelet transform (CWT) for a SWC series of length n (Y_i, i = 1, 2, …, n) with equal incremental distance $\delta x$. This can be defined as the convolution of Y_i with the scaled and normalized wavelet using fast Fourier transform (Si and Zeleke, 2005; Torrence and Compo, 1998):

$$W_{{}_i}^Y\left(s\right)=\sqrt{\frac{\delta x}{s}}{\displaystyle \sum _{j=1}^n{Y}_j}\psi \left[\left(j-i\right)\frac{\delta x}{s}\right]$$E23

where $\psi \left[\right]$ is the mother wavelet function and s is the scale. The mother wavelet functions include Morlet, Mexican hat, Harr, and others (Si, 2008). Depending on the purpose, different mother wavelet functions can be selected. The Morlet wavelets allow us to detect both location dependent amplitude and phase for different frequencies in the spatial series (Torrence and Compo, 1998), which can be written as

$$\psi \left[\eta \right]={\pi }^{-1/4}{e}^{i\overline{\omega }\eta -0.5{\eta }^2}$$E24

where $\omega$ is dimensionless frequency and $\eta$ is dimensionlessspace. The wavelet is stretched in space (x) by varying its scale ($s$), so that $\eta =s/x$. The Morlet wavelet (with $\omega$ =6) is a good choice for feature extraction purpose like identifying the scales and locations, since it provides a good balance between space and frequency localization.

The wavelet coefficients $W_i^Y\left(s\right)$ can be expressed as a + ib where a and b are the real and imaginary components of $W_i^Y\left(s\right)$. For the polar form of complex numbers, $W_i^Y\left(s\right)= \left|W_i^Y\left(s\right)\right|\left(cos\theta +isin\theta \right)$, where $\theta =arctan\frac{b}{a}$ is called the phase or argument of $W_i^Y\left(s\right)$. The wavelet power spectrum is defined as ${\left|W_i^Y\left(s\right)\right|}^2$ and the local phase is defined as the complex argument of $W_i^Y\left(s\right)$.

After calculating the wavelet spectra $W_i^Y\left(s\right)$ and $W_i^Z\left(s\right)$ corresponding to two SWC spatial series Y and Z, respectively for two different times, cross wavelet power spectrum $W_i^{YZ}\left(s\right)$ at scale s and location i can be defined as

$$W_i^{YZ}\left(s\right)=\left[W_i^Y\left(s\right)\right]\overline{\left[W_i^Z\left(s\right)\right]}$$E25

and the wavelet coherency of two spatial series can be written as

$${\rho }_i^2\left(s\right)=\frac{{\left|{\overleftrightarrow{W}}_i^{YZ}\left(s\right)\right|}^2}{{\left|{\overleftrightarrow{W}}_i^Y\left(s\right)\right|}^2{\left|{\overleftrightarrow{W}}_i^Z\left(s\right)\right|}^2}$$E26

where $$\overleftrightarrow{(\cdot )}$$ is a smoothing operator, which can be written as

$$\overleftrightarrow{W}=S{M}_{scale}\left[S{M}_{space}\left(W\right)\right]$$E27

where ${SM}_{space}$ denotes smoothing along the wavelet scale axis and ${SM}_{scale}$ smoothing in spatial distance. The following smoothing function is the normalized real Morlet wavelet and has a similar footprint as the Morlet wavelet

$$\frac{1}{s\sqrt{2\pi }}\exp \left(-\frac{{\tau }^2}{2s^2}\right)$$E28

where $$\tau$$ denotes location. Therefore, the smoothing along locations can be written as

$$S{M}_{space}\left(W,s,\tau \right)={\displaystyle \sum _{k=1}^m\left\{W\left(s,\tau \right)\frac{1}{s\sqrt{2\pi }}\exp \left[-\frac{{\left(\tau -x_k\right)}^2}{2s^2}\right]\right\}}|{}_s$$E29

The Fourier transform of Equation (28) is $\mathrm{ }\mathrm{e}\mathrm{x}\mathrm{p}\left(-2s^2{\omega }^2\right)$, where $\omega$ is the frequency. Eq. (29) can be implemented using Fast Fourier Transform (FFT) and Inverse FFT (IFFT) based on convolution theorem and is written as

$$S{M}_{space}\left(W,s,x\right)=IFFT\left\{FFT\left[W\left(s,\tau \right)\right]\left[\exp \left(-2s^2{\omega }^2\right)\right]\right\}$$E30

The smoothing along scales can be written as:

$$S{M}_{scale}\left(W,s_k,x\right)=\frac{1}{2p+1}{\displaystyle \sum _{j=k-p}^{k+p}\left[S{M}_{space}\left(W,s_j,x\right)\prod \left(0.6s_j\right)\right]}|{}_x$$E31

where $P$ is the number of terms on each symmetrical half of the window, and П is the rectangle function. The factor of 0.6 is the empirically determined scale decorrelation length for the Morlet wavelet (Si and Zeleke, 2005; Torrence and Compo, 1998).

There are many methods available to test the statistical significance of wavelet coherency. The Monte Carlo simulation or reshuffling method is among the suggested ones. A detailed description of these methods can be found in Si (2008). Matlab codes for calculating wavelet coherency are available at URL: http://www.pol.ac.uk/home/research/waveletcoherence/ (Grinsted et al., 2004).

3.1.7. Multivariate empirical mode decomposition

Empirical mode decomposition (EMD), extracts oscillations from the soil water content series into a finite and often small number of intrinsic mode functions (IMFs) according to the energy associated with different space scales (Huang et al., 1998; Huang and Wu, 2008). Unlike spectral and wavelet methods, EMD does not call for any assumption of the data and works directly in the spatial domain with the basis completely derived from the data. Therefore, it is intuitive, direct, a posteriori and adaptive (Huang et al., 1998; Huang and Wu, 2008). The locality and adaptivity of EMD can deal with different types of spatial series including non-stationary and nonlinear (Huang et al., 1998). Therefore, EMD has a great potential to find the underlying scale of spatial series of soil moisture without imposing any mathematical assumption on the measured data (Biswas and Si, 2011b).

Multivariate empirical mode decomposition (MEMD) is the multivariate extensions of standard EMD (Huang et al., 1998). An important step for MEMD is the computation of the local mean because local extrema are not well defined directly for multivariate spatial data. Moreover, the notion of “oscillatory model” defining an IMF is rather confusing for multivariate spatial data. To deal with these problems, Rehman and Mandic (2010) produced a multiple m-dimensional envelopes by taking projections of multiple inputs along different directions in an m-dimensional space.

Assuming ${\left\{\mathrm{v}\left(s\right)\right\}}_{s=1}^{\mathrm{S}}=\left\{{\mathrm{\upsilon }}_1\left(s\right),{\mathrm{\upsilon }}_2\left(s\right),\dots ,{\mathrm{\upsilon }}_{\mathrm{m}}\left(\mathrm{s}\right)\right\}$ being the m spatial data sets of soil water content as a function of space (s), and ${\mathrm{x}}^{{\theta }_k}=\left\{{\mathrm{x}}_1^k,{\mathrm{x}}_2^k,\dots ,{\mathrm{x}}_{\mathrm{m}}^k\right\}$ denoting a set of direction vectors along the directions given by angles ${\theta }^k=\left\{{\theta }_1^k,{\theta }_2^k,\dots ,{\theta }_{\mathrm{m}-1}^k\right\}$ on a m-1-dimensional sphere ($k$ is the number of direction used to calculate the projections and envelope curves). Then, IMFs of the m spatial data sets can be obtained by MEMD as follows:

1. Choose a suitable point set for sampling on an m-1-dimensional sphere. This can be done by sampling unit hyperspheres (m-spheres) based on both uniform angular sampling methods and quasi-Monte Carlo-based low-discrepancy sequences.

2. Calculate a projection, denoted by ${\left.p^{{\theta }_k}\left(t\right)\right\}}_{s=1}^{\mathrm{S}}$, of the spatial data sets ${\left\{\mathrm{v}\left(s\right)\right\}}_{s=1}^{\mathrm{S}}$ (S is the number of sampling point for each time) along the direction vector${\mathrm{ }\mathrm{x}}^{{\theta }_k}$, for all $k ($ the whole set of direction vectors), giving ${\left.p^{{\theta }_k}\left(s\right)\right\}}_{k=1}^K$ as the set of projections.

3. Find the spatial instants $s_i^{{\theta }_k}$ corresponding to the maxima of the set of projected data sets ${\left.p^{{\theta }_k}\left(s\right)\right\}}_{k=1}^K$.

4. Interpolate $\left[s_i^{{\theta }_k},\mathrm{v}\left(s_i^{{\theta }_k}\right)\right]$ to obtain multivariate envelope curves ${\left.{\mathrm{e}}^{{\theta }_k}\left(s\right)\right\}}_{k=1}^K$ for all the data sets considered.

5. For a set of $K$ direction vectors, the mean $m\left(s\right)$ of the envelope curves is calculated as:

$$m\left(s\right)=\frac{1}{K}\sum _{k=1}^K{\mathrm{e}}^{{\theta }_k}\left(s\right)$$E32

1. Extract the “detail” $d\left(s\right)$ using $d\left(s\right)=\mathrm{v}\left(s\right)-m\left(s\right)$. If the “detail” $d\left(s\right)$ fulfills the stoppage criterion for a multivariate IMF, apply the above procedure to v $\left(s\right)-d\left(\mathrm{s}\right)$, otherwise apply it to $d\left(s\right)$. The stoppage criterion for multivariate IMFs is similar to that proposed by reference (Huang et al., 2003).

MEMD has the ability to align “common scales” present within multivariate data. Each common scale is manifested in the common oscillatory modes in all the variates within an m-variate IMF. After MEMD analysis, scale and location specific time stability of SWC can be easily identified by comparing and calculating spearman’s rank correlation coefficient between the IMFs with the same numerical numbers for different measurements using Eq. (2). The exact scale for each IMF can be obtained from the instantaneous frequencies by Hilbert transformation with IMF (Huang et al., 1998). The instantaneous frequencies of soil water content can be converted to period (1/frequency), which was further converted to the spatial scale after multiplying the period with the sampling interval. MEMD can be completed using Matlab program such as that written by Rehman and Mandic (2009) (http://www.commsp.ee.ic.ac uk/~mandic/research/emd.htm). This method requires that soil water content measurement should be obtained in intervals with equidistance.

3.2. Time stability at points

Time stability at point scales are usually obtained in terms of average SWC estimation quality. Direct and indirect methods can be used to estimate average SWC by the measurement at time stable locations. The direct method estimates average SWC directly by measuring SWC at a time-stable location (Brocca et al., 2009; Cosh et al., 2008; Grayson and Western, 1998; Pachepsky et al., 2005). The indirect method estimates average SWC by considering the offset between the average SWC and the measurement value at a time-stable location (Grayson and Western, 1998; Starks et al., 2006). All the indices are listed in Table 1. Each index has its own advantages and disadvantages, and these are discussed in the application (Section 4.2). For all indices, smaller value of time stability index at a given location indicates stronger time stability and also higher quality of average SWC estimation by that location.

Index	Formula and explanation	Average SWC estimation Method^†	Reference
Standard deviation of relative difference, $SDRD\left(i\right)$	where $$SDRD\left(i\right)=\sqrt{\frac{1}{m-1}{\displaystyle \sum _{j=1}^m{\left({\delta }_j\left(i\right)-⟨{\delta }_i⟩\right)}^2}}$$ and $${\delta }_j\left(i\right)=\frac{SW{C}_j\left(i\right)-⟨SW{C}_j⟩}{⟨SW{C}_j⟩}$$	Direct and indirect	Vachaud et al., 1985
Root mean square error, $RMSE\left(i\right)$		Direct	Jacobs et al., 2004
Standard deviation of relative SWC, $\sigma \left({\beta }_i\right)$ ^‡	$$⟨{\delta }_i⟩=\frac{1}{m}{\displaystyle \sum _{j=1}^m{\delta }_j}\left(i\right)$$; where $$RMSE(i)={\left({⟨{\delta }_i⟩}^2+SDRD{\left(i\right)}^2\right)}^{1/2}$$ and $$\sigma \left({\beta }_i\right)=\sqrt{\frac{1}{m-1}{\displaystyle \sum _{j=1}^m{\left({\beta }_j\left(i\right)-⟨{\beta }_i⟩\right)}^2}}$$	Indirect	Pachepsky et al., 2005
Width of the 90% empirical tolerance interval of relative water content, $T\left(i\right)$	where $\mathrm{ }{\beta \left(i\right)}_{P=0.95}$ and ${\beta \left(i\right)}_{P=0.05}$ are the relative SWC values at cumulative probability of 95% and 5%, respectively.	Indirect	Guber et al., 2008
Chi-squared statistic, ${\chi }^2\left(i\right)$	where $\mathrm{ }{\sigma }_j$ is the standard deviation of SWC at observation time j.	Indirect	Guber et al., 2008
Root-mean-squared differences, $D\left(i\right)$		Indirect	Guber et al., 2008
Mean absolute bias error, $MABE\left(i\right)$		Indirect	Hu et al., 2010a, 2010b

Table 1.

Time stability indices at point scales and their formulas

†Direct refers to estimating average SWC directly by measuring SWC at a time-stable location and indirect refers to estimating average SWC by considering the offset between the mean value and the measurement value at a time-stable location.

‡ $\sigma \left({\beta }_i\right)=SDRD\mathrm{ }\left(i\right).\mathrm{ }$

4. Application

4.1. Time stabilityof spatial pattern

Information on the spatial pattern of soil water content with time is important to understand the changes in hydrological processes and in making predictions with physical and statistical models. In this section, we have documented the use of each method introduced above to characterize the similarity of the spatial pattern of SWC.

Spearman’s rank correlation analysis is the most widely used method to examine the similarity of the spatial patterns of SWC. Numerous studies have used this method and have confirmed the presence of time stability of the spatial pattern (Biswas and Si, 2011c; Hu et al., 2009, 2010b; Martínez-Fernández and Ceballos, 2003; Vachaud et al., 1985; Vivoni et al., 2008). In the first of this kind, Vachaud et al. (1985) measured soil water storage of 0 to 1.0 m depth at 17 neutron access tubes located at 10 m intervals for 24 occasions over a time span of two and a half years. According to the cumulative probability function for the two extreme conditions, i. e., the wettest on 7 September 1982, and the driest on 17 August 1983, many locations maintained their ranks (Figure 1). For example, location #1 always had the maximum soil water storage, location #4 had the minimum soil water storage irrespective of the environmental conditions (e. g., wet or dry). Spearman’s rank correlation coefficients of the soil water storage between measurements were also calculated to examine the degree of time stability of the spatial pattern. As an example, the correlation between seven measurements representing the range of the water storages recorded over the entire observation period is presented in Table 2. All rank correlation coefficients were highly significant at the 0.1% two-tailed test, indicating strong time stability of the spatial pattern of soil water storage. On the other hand, time instability of spatial pattern between measurements was also reported (Comegna and Basile, 1994; Mohanty and Skaggs, 2001). Although it is very popular among researchers, Spearman’s rank correlation analysis should be viewed only as a statistical tool measuring the degree of concordance between two rankings (Vachaud et al., 1985). It may be questionable if differences between measured values are smaller than experimental uncertainties. This may be the case in situations where either the probability density function is very uniform or the experimental determinations are very crude.

Figure 1.

Cumulative probability function of soil water storage at 0-1.0 m measured in the most dry, and the most wet conditions. Numbers refer to measuring locations. (Reproduced from Vachaud et al. (1985)).

Dates	7 Sep 1982	28 Oct; 1981	2 Dec 1982	29 Jul 1981	18 Aug 1981	16 Jul 1982	25 Aug 1981
Average storage (mm)	428.1	425.77	424.77	382.06	354.72	347.23	344.46
(1)	1
(2)	0.953	1
(3)	0.941	0.953	1
(4)	0.988	0.961	0.946	1
(5)	0.953	0.922	0.882	0.968	1
(6)	0.863	0.789	0.843	0.836	0.882	1
(7)	0.860	0.824	0.824	0.794	0.863	0.939	1

Table 2.

Matrix of Spearman’s rank correlation coefficient between series of soil water storage measurements obtained on the 17 locations on seven dates. (Reproduced from Vachaud et al. (1985))

As an alternative to Spearman’s rank correlation analysis, Pearson correlation coefficients were also used to describe time stability of the spatial pattern (Cosh et al., 2004). It was observed that the Pearson correlation analysis revealed a pattern similar to that of the Spearman’s rank correlation analysis (Cosh et al., 2004). Note that the tendency of soil moisture sensors to occasionally report erroneous measurements may negatively affect the value of Spearman’s rank correlation coefficient, while the Pearson correlation coefficient is less sensitive to this problem. This is because the total number of measurements will not affect Pearson correlation coefficient as it affects the rank correlation coefficient (Cosh et al., 2006, 2008). However, both Spearman’s rank correlation and Pearson correlation analysis are only suitable to describe how much spatial pattern at one measurement time can be preserved at another time (Rolston et al., 1991). If we want to get the information on the similarity of spatial pattern over multiple measurement occasions, mean value of correlation coefficients can be used (Hu et al., 2009).

With the aim to characterize the degree of time stability of spatial pattern among multiple data series, Starr (2005) suggested two indices, i. e., average spatial coefficient of determination (Rs²) and temporal coefficient of variability (CV_t). The authors reported that 26 to 76% (with mean of 47%) of the observed variability was explained by the time stable pattern. The reported CV_t value ranged from 11% to 18% (with mean of 16%) for different sampling transects. On the other hand, the Rs² may not be a good measure of how time stable a SWC pattern is at a field that is uniform and completely time stable. This is because the Rs² will be zero in this case (Starr, 2005).

Perry and Niemann (2007) used the EOF method to characterize the time stable spatial pattern of SWC. They identified two dominant spatial patterns of soil water variability at Tarrawarra catchment (Figure 2). The first spatial pattern was more important for medium wet periods and was associated with lateral redistribution of soil water. The second spatial pattern served to modify the first spatial pattern at some days which was associated to evapotranspiration. The importance of each spatial pattern on each measurement occasion can be reflected by the EC values (Figures 2c and 2d). There was a large decrease in EC1 weight during the dry season, which agrees well with the finding of Western et al. (1999) that soil moisture patterns on moderately wet occasions show a strong topographic influence. The EOF method can also be used to predict the SWC distribution for unobserved times if the average SWC and the temporal varying coefficient (ECs), which reflect the relative importance of different spatial patterns can be obtained. Average SWC can be estimated by many methods such as linear regression (Western et al., 1999), a dynamic multiple linear regression of soil moisture against topographic attributes (Wilson et al., 2005), and by measurement from representative location (Grayson and Western, 1998). Examples of EOF method can also be found in other study (Ibrahim and Huggins, 2011; Jawson and Niemann, 2007; Joshi and Mohanty, 2010; Korres et al., 2010; Yoo and Kim, 2004). Different numbers, usually one to four, of significant EOFs were obtained to explain most of SWC variations (usually more than 70%). The time-invariant spatial patterns were usually related to topography-related factors (Yoo and Kim, 2004), sand content (Jawson and Niemann, 2007), mixed influences of rainfall, topography, and soil texture (Joshi and Mohanty, 2010), elevation and wetness index under wet conditions and soil properties (ECa and bulk density) in dry conditions (Jawson and Niemann, 2007).

Figure 2.

a) EOF1 and (b) EOF2 computed from the spatial anomalies of soil moisture. The bars in (c) and (d) are EC1 and EC2, respectively. The spatial average and variance of the soil moisture on each sampling date are also shown in (c) and (d) for comparison. (Reproduced from Perry and Niemann (2007))

Time stability of soil water content is usually scale dependent because of the interactions of soil moisture between measurement locations (Kachanoski and De Jong, 1988). The first effort to reveal the scale specific time stability was made by Kachanoski and De Jong (1988) who used spectral coherency to explore the time stability of soil water storage in a rolling landscape (Figure 3). They observed that for both wet and drying periods, the spectral coherency was high at large scales (>40 m). For smaller scales (<40 m), however, the spectral coherency was high for the drying period but greatly weakened during the alternating wetting and drying period. Therefore, their results clearly showed the loss of time stability at small scales (<40 m). Gómez-Plaza et al. (2000) also used the spectral coherency to explore the scale dependence of time stability in burned and unburned areas in a semi-arid environment. They clearly observed that soil moisture patterns were maintained in all cases at large scales (transect scale), while soil moisture tended to become time unstable, depending on the range of scale and section of transect considered at medium scales.

Figure 3.

Coherency spectrum for soil water storage for the drying (June 28-July 12) and recharge (July 12-October 1) period. (Reproduced from Kachanoski and De Jong (1988))

Spectral coherency was verified to be robust to look for the scale dependency of time stability. Its multi-variate extension is also feasible (Si, 2008), although no publications have been found on this point. However, spatial coherency analysis assumes the stationarity in the measured data series over space, hence may not capture localized features of time stability (Biswas and Si, 2011a).

Figure 4.

The inter-season wavelet coherencies between the soil water storage measured during (a) spring (2 May 2008) and summer (23 August 2008), (b) spring (2 May 2008) and fall (22 October 2008), and (c) summer (23 August 2008) and fall (22 October 2008) for the whole soil profile. Cross sectional view of the transect with landform elements is shown at the top. The X-axis indicates distance along the transect; the Y-axis indicates the scale in meter; the solid black line indicates 5% significance level; the color bar indicates strength of correlation, and the direction of arrow indicates the phase relationship or type of correlation. (Reproduced from Biswas and Si (2011c)).

Wavelet coherency can deal with the non-stationary series, it is therefore suitable to identify both scale and location specific time stability (Biswas and Si, 2011a). An excellent application of wavelet coherency in time stability analysis was made by Biswas and Si (2011a). They identified scales and locations of time stability between the spatial patterns of soil water storage in a hummocky landscape in central Saskatchewan, Canada. As Figure 4 shows, strong time stability at all scales and locations existed between summer and fall, when the vegetation was present in field. The spatial pattern in spring was different from that of summer or fall. The large-scale (> 72 m) time stability was present over the whole transect. The medium-scale spatial patterns were time stable in large depressions that had fewer landform elements. The small-scale time stability was mostly random.

Wavelet coherency can also be used to identify the correlation between different variables at different scales and positions (Lakshmi et al., 2004; Si and Zeleke, 2005). One limitation of wavelet coherency analysis is that it is based on the assumption of linearity of systems. Moreover, the application of wavelet coherency on time stability analysis was only confined to two spatial series. Extension of this method to multiple spatial series is expected to understand the scale and location specific time stability when soil water content from different seasons are considered together.

MEMD is a data-driven method that can deal with both non-stationary and nonlinear series. Recently, Hu et al. (unpublished data) (2012a) applied MEMD to time stability analysis using soil water storage measurements from 20 occasions in a hummocky landscape in central Saskatchewan, Canada. At first, multiple spatial series of soil water storage were separated into different IMFs. The IMFs for the 0-20 cm layer on a selection of measurement occasions is given in Figure 5. Each IMF corresponds to a specific scale common to all the spatial series considered. Spearman’s rank correlation analysis was then conducted for each scale. According to their study, the dominant scale of variations of soil water storage was 104-128 m (IMF5), for 0-20 cm, 20-40 cm, 40-60 cm, and 60-80 cm, and 315-412 m (IMF7) for 80-100 cm, 100-120 cm, and 120-140 cm. Time stability generally increased with scale and it was the strongest at the dominant scale of variation. At different scales, the time stability generally ranked in the order of intra-season>inter-annual>inter-season. Therefore, the study of Hu et al. (2012a) verified that scale specific time stability of soil moisture can be identified using MEMD method.

According to the study reviewed above, time stability analysis is developing from two spatial series to multiple series, from one scale to multiple scales. With the advance of methodology, a deeper understanding of spatial patterns of soil moisture can be expected. Future work should make full use of various methods, especially those most recently developed, such as EOF, wavelet coherency, and MEMD to gain greater insight into the spatial pattern dynamics at different scales.

4.2. Identify time stable locations for average soil water estimation

One of the most important applications of time stability is to identify a representative location for average soil water content estimation for a given area. Identification of representative locations cover different land uses, different climate zones, and places with various topographical attributes (Gao et al., 2011).

Figure 5.

Intrinsic mode functions (IMFs) and residues based on multivariate empirical mode decomposition using soil water storage series (top) of the surface layer (0-20 cm). The horizontal axis is the distance between a sampling location and the origin of the transect. The vertical solid bar in Y-axis shows the scale of IMFs. (Reproduced from Hu et al. (2012a)).

Various time stability indices for use at the point scale have been used to identify representative locations. Most studies have been aimed at finding locations where soil water content can directly represent the mean of an area (Brocca et al., 2009; Cosh et al., 2008; Grayson and Western, 1998; Martínez-Fernández and Ceballos; 2005). In this way, SWC at the time stable location should approximate the average SWC of the study area. By extension, point scale time stability representing the areal average also means that the mean relative difference of soil water content, $⟨{\mathrm{\delta }}_{\mathrm{i}}⟩$, at the stable location i should be close to zero at all times. Practically, it is difficult to find a location where the mean relative difference is the close to zero at all times. In this situation, the time stability index SDRD is widely used to judge the persistence of the stable location. Therefore, the mean relative difference and associated SDRD should be considered together. Successful application of this method can be found in many publications (Brocca et al., 2009; Cosh et al., 2006, 2008; Gómez-Plaza et al., 2000; Grayson and Western, 1998; Tallon and Si, 2004). In the study of Grayson and Western (1998), although the overall spatial pattern of soil moisture was time unstable, some locations were time stable and can represent the average SWC of a landscape. In the Tarrawarra watershed (Figure 6), the mean relative differences for locations 1, 4, 12, and 13 were close to zero, and their SDRD were also relatively small, so these locations were time stable and could be used as representative location for average SWC estimation. In the study of Tallon and Si (2004) which was conducted on a chernozemic soil on the Canadian prairies (75 km northeast of Saskatoon, SK, Canada) time stability of spatial patterns of soil moisture for different depths (30, 60, 90, 120 and 160 cm) were observed. With the exception of one location that occurred at two separate depths, the locations representative of average SWC differed with soil depths.

Figure 6.

Ranked mean relative differences of soil water content for 0-60 cm measured by neutron probe in the Tarrawarra watershed. Also shown is one standard deviation error and site numbers (Reproduced from Grayson and Western (1998)).

Although mean relative difference and SDRD have been used successfully in the past in identifying representative locations for average SWC, the question becomes a matter of which index is more important when the location with mean relative difference closest to zero disagrees with that with the minimum SDRD. Most researchers believed that the SDRD was more important (Brocca et al., 2009; Cosh et al., 2006, 2008; Gómez-Plaza et al., 2000; Grayson and Western, 1998), while mean relative difference was also emphasized (Martínez-Fernández and Ceballos, 2003).

Instead of considering two indices simultaneously, one index, RMSE, combines both indices into one, thus making the identification of time stable location more objective. By this index, drier locations were usually found to have lower RMSE values than the wetter points (Cosh et al., 2004, 2006). This suggests that improved consistency will result from the selection of slightly drier sampling locations (Cosh et al., 2004). Operationally, it facilitates the selection of representative locations. However, it may confuse the concept of time stability, because for locations with large absolute values of relative difference (i. e. extremely dry or wet conditions) but very small SDRD, they may have large value of RMSE although they are completely time stable.

Another method for estimating average soil water content with time stable location involves a consideration of the constant offset (mean relative difference) for time stable location (Grayson and Western, 1998). Then $⟨{SWC}_j⟩$ can be estimated by

$${\beta }_j\left(i\right)=\frac{SW{C}_j\left(i\right)}{⟨SW{C}_j⟩}$$E33

In this situation, time stable locations are usually defined by the value of SDRD (Grayson and Western, 1998). However, there are two limitations in using SDRD as a criterion. First, the value of the SDRD may be affected to some extent by the magnitude of the mean relative difference; second, the value of the SDRD cannot be directly related to the estimation error of the mean value (Hu et al., 2010b).

Indices that can consider the estimation error are needed for identifying representative time stable locations. One example was conducted by Guber et al. (2008). They applied ${\chi }^2\left(i\right)$ and $D\left(i\right)$ to rank the utility of sampling locations to estimate average SWC by considering the estimation error. Another alternative is the mean absolute bias error (MABE), which was developed directly to relate to the relative error of mean estimation (Hu et al., 2010a, 2010b). Although MABE was positively correlated to SDRD in some cases (Gao et al., 2011; Hu et al., 2010a, 2010b), the lowest value of SDRD did not coincide with the lowest value of MABE (Figure 7). Therefore, MABE outperformed SDRD in terms of representative location identification for indirect estimation of average SWC (Hu et al., 2010b).

The degree of representativeness of time stable locations may also change with season or year due to the change of the relative importance of different factors affecting time stability (Guber et al., 2008; Martínez-Fernández and Ceballos, 2003). Schneider et al. (2008) showed that the selected representative locations were appropriate to predict average soil water content of the sites for multiple years, although the time stability characteristics of some locations varied between years. An interesting study conducted by Guber et al. (2008) indicated that the utility of a location to estimate the estimation accuracy of field average SWC was not influenced by measurement frequency.

Figure 7.

Rank ordered mean absolute bias error (MABE) and standard deviation of relative differences (SDRD) for the soil water storage of 0-1.0 m measured with neutron probe from October, 2004 to September, 2006. Also shown is location number (Reproduced from Hu et al. (2010b))

Many indices have been available to identify the representative location for average SWC estimation. However, few works have been conducted to compare their performances with a few exceptions (Gao et al., 2011; Hu et al., 2010a, 2010b, 2012b). MABE was reported to be a better index than SDRD in terms of average soil water storage estimation by the indirect method (Hu et al., 2010a, 2010b). More recently, with the soil water storage datasets from a transect in a Canadian prairie area and a watershed on the Chinese Loess Plateau, Hu et al. (2012b) evaluated the performance of different time stability indices in terms of the average soil water storage estimation quality judged by different goodness-of-fit indices. Their results showed that MABE, χ², D, and CV_t outperformed SDRD and T. If root mean squared deviations (RMSD) and Nash-Sutcliffe coefficient of efficiency (NSCE) were used as goodness-of-fit indices, D was the best time stability index. If absolute mean difference (BIAS) and absolute bias relative to mean (RBIAS) were adopted, MABE was the best time stability index. Their results also showed that average soil water storage estimation by the indirect method was more accurate than that by the direct method. It is necessary to compare their performances under other different hydrological backgrounds to obtain a more suitable index for identifying representative location. Furthermore, all these methods still need an intensive measurement of soil water content over time before a representative location can be identified. According to Martínez-Fernández and Ceballos (2005) about one year (a complete seasonal cycle) of measurements is needed to determine the representative location. Therefore, there is an urgent need to look for the possible definable features of the terrain and soils that could be used to determine the representative locations a priori (Martínez-FernándezandCeballos, 2005).

5. Controlling factors

Time stability of soil water content is controlled by factors that consistently influence soil moisture distribution in a similar way at all times. Most studies focus on factors influencing the underlying spatial pattern of soil moisture, and relatively few studies are associated with factors influencing the degree of time stability at point scales. We will review the controlling factors of the two types of time stability below.

5.1. Controlling factors of time stability of spatial pattern

Controlling factors of underlying spatial pattern are usually defined by comparing spatial patterns of time stability and various factors (Hu et al., 2010a; Perry and Niemann, 2007; Schneideret al., 2008). Related studies referred to different climate conditions, from semi-arid (Hu et al., 2010a) to humid environments (Perry and Niemann, 2007) and different land uses, such as rangeland (Mohanty and Skaggs, 2001), grass (Jacobs et al., 2004), forest (Lin, 2006), and agricultural land (Starr et al., 2005).

Soil properties such as soil texture are usually found to influence the spatial pattern of time stability. For example, Vachaud et al. (1985) found that the scaled factor of soil water content were significantly correlated to the silt plus clay contents on a field in St. Martin d’Heres, France. Mohanty and Skaggs (2001) observed that the sandy loam soil produced the best time-stable soil moisture pattern and silty loam soil produced an intermediate level of time stability in the Little Washita watershed of the southern part of the Great Plains of the US. Jacobs et al. (2004) pointed out that lower sand content may produce more time stable patterns in the Iowa. Grant et al. (2004) found a significant correlation between clay content, coarse fragment and the rank of mean relative difference of soil water storage in the Reynolds Mountain East catchment in the Owyhee Mountains. Starr (2005) found that coarser particle size classes were generally drier than adjacent finer particle size class soil on a potato farm located near Houlton, Maine. Hu et al. (2010a) recognized that sand and silt content can explain about 41.4% - 65.5 % of the spatial pattern of mean relative difference of soil water content. Except for soil texture, other soil properties such as organic matter content (Hu et al., 2009), bulk density (Jacobs et al., 2004), and soil thickness (Zhu and Lin, 2011) have also been recognized as the controlling factors for time stability of spatial pattern. However, soil properties may also have no influences on the spatial pattern of time stability. For example, Tallon and Si (2004) reported that clay content showed the least amount of control of spatial patterns on a chernozemic soil on the Canadian prairies (75 km northeast of Saskatoon, SK, Canada). On ungrazed sites, Schneider et al. (2008) found that neither bulk density, organic carbon nor sand and clay content could explain the time-stable characteristics of the sampling points. Comegna and Basile (1994) and Starr (2005) both observed the absence of the controlling influence of soil texture and they attribute it to the homogeneous distribution of soil texture in their study areas.

Topographic characteristics have also been recognized to influence the time stability of spatial patterns. Gómez-Plaza et al. (2000) showed that when the factors affecting soil moisture were limited to topographical position or local topography at the transect scale, spatial patterns of soil moisture presented time stability. Joshi et al. (2011) observed that fields with flat topography (LW21) showed the worst time stability features compared to the fields having gently rolling topography (LW03 and LW13) in Oklahoma. The influence of topography on the time stability of spatial pattern may be scale dependent. For example, although the Pearson correlation coefficient between soil water storage and elevation was very small, Biswas and Si (2011a) identified that strong scale-specific correlations existed between soil water storage and elevation at different scales using wavelet coherency, which contributed to the time stability of soil water storage at different scales. Kachanoski and De Jong (1988) attributed to the loss of time stability of soil water storage at scale smaller than 40 m to the role of spatial pattern of surface curvature.

Usually, soil properties and topographical properties jointly control the time stability of soil moisture spatial pattern. According to Joshi et al. (2011), soil properties (i. e., percent silt, percent sand, and soil texture) and topography (elevation and slope) jointly affected the spatiotemporal evolution and time stability of soil moisture at both point and footprint scales in the Little Washita watershed, Oklahoma and in the Walnut Creek watershed, Iowa. Zhu and Lin (2011) also observed that both soil and terrain controlled soil moisture variation at different seasons and soil depths at farm scale.

The relative role of soil and topographic properties in controlling the spatial variability and time stable pattern of soil moisture is usually related to the dominant hydrological process taking place for different soil water conditions between wet and dry periods. The SWC patterns existing between wet and dry periods are referred to as preferred states (Graysonet al., 1997; Western et al., 1999). When dry conditions dominate, evapotranspiration is greater than precipitation, and local controls dominate mainly by vertical movement of soil water with no connection between adjacent points. In this situation, differences in vegetation and soil properties such as texture would be responsible for the spatial pattern of soil moisture (Graysonet al., 1997). During wet states, evapotranspiration is less than precipitation and non-local controls govern soil moisture distribution. These controls are related to upslope topography and include catchment area, aspect, depth and soil profile curvature (Gómez-Plaza et al., 2001). Western at al. (1999) attributed these differences in local and non-local controls to the high degree of organization during wet periods that consist of connected bands of high soil water content in drainage lines. Soil water in these drainage lines is laterally redistributed by both overland and subsurface flow. Using EOF analysis, Perry and Niemann (2007) found that underlying spatial patterns of soil moisture were controlled by local soil properties in wet and dry conditions and topographic characteristics during intermediate conditions. Their results agreed with the findings of Western et al. (1999) that soil moisture patterns on moderately wet dates show a strong topographic influence, one that seems to diminish on both dry and very wet days. Note that soil porosity, rather than topography, will dominate the spatial distribution of soil moisture under very wet conditions (Western et al., 1999). Non-local control can also be found in the study of Biswas and Si (2012). In their study, the coefficient of determination between soil water storage and elevation increased more than eightfold after shifting the spatial series of soil water storage by a length equal to that of the existing slope. Although local and non-local control was widely recognized, it is still difficult to make a very clear limit between them, most likely due to the possible influences of topography on soil properties

Because of the different controls in different soil water conditions, time stability of soil water content tends to be stronger under the similar soil water conditions than those with quite different soil water conditions. For example, Martinez-Fernández and Ceballos (2003) observed that the periods with the lowest time stability coincided with situations involving the transition from dry to wet. Gao et al. (2011) revealed that the time stability of root zone soil water (0–60 cm) was higher in either dry or wet season than that including both, and soil water exhibited very low time stability during the transition period from dry to wet.

In addition to soil and topography, other factors such as vegetation (Pachepsky et al., 2005; Schneideret al., 2008), soil depth (Huet al., 2009) can also influence the time stability of soil moisture. For example, Gómez- Plaza et al. (2000) found that the existence of vegetation could destroy the time stability of spatial patterns of soil moisture. Mohanty and Skaggs (2001) found varying degrees of time stability depending on vegetation cover as well as topography. Generally, time stability of spatial pattern is strongest in deep soil layers than in shallow layers (Cassel et al., 2000). However, stronger time stability of spaital pattern at shallower depth (20 cm) than deeper soil depths was also observed (Hu et al., 2009, 2010a).

Although many factors have been revealed to influence the time stability of spatial pattern, no single, or even combined, factor can explain all time stable patterns. This suggests that some influencing factors may still remain undiscovered. In addition, random measurement errors also contribute to some variability of the time stable pattern (Starr, 2005). According to the study of Starr (2005), random measurement error can contribute 20% of the variability of underlying spatial pattern. Therefore, to better understand the spatio-temporal variability and time stability of spatial pattern of soil moisture, controlling factors should be made clearer when more accurate measurement of soil water content should also be expected.

5.2. Controlling factors of time stability at the point scale

Revealing the controlling factors of time stability at the point scale is crucial in identifying time stable locations for average SWC estimation. Soil and topography are usually found to influence the time stability of soil water content at point scales. Brocca et al. (2009) observed that the time stable positions were linked to topographic characteristics, primarily the upslope drainage area but also elevation and slope. Cosh et al. (2008) found that almost 60% of the SDRD variability can be credited to the soil type parameters. Soil parameters as major controlling factors is a finding that is in agreement with those of Hu et al. (2010a) who found that soil texture can significantly (P<0.05) affect the stability of soil water content. Gao et al. (2011) observed that time-stable locations corresponded to relatively high clay contents, mild slopes and planar surfaces, which agreed with the previous studies, where both soil type and topography influenced time stability at point scales (Jacobs et al., 2004; Joshi et al., 2011; Mohanty and Skaggs, 2011). However, poor relationships between time stable locations and soil and topographic properties were also found (Tallon and Si, 2004). Land use was observed to have no influences on time stability of soil water content at point scales (Hu et al., 2010a; Joshi et al., 2011).

Time stability at point scales was usually associated with the soil water conditions of that location. For example, Martinez-Fernández and Ceballos (2005) found that time stability as characterized by SDRD was usually stronger in dry locations than wet locations for all depths, and they attribute it to the predominance of the sandy fraction with the resultant weak ability to retain water in dry locations. Similar results were obtained in other studies (Bosch et al., 2006; Jacobs et al., 2004; Martínez-Fernández and Ceballos, 2003). However, the relationship of time stability and soil water condition is likely influenced by the time stability index used. In the LYMQ watershed on the Loess Plateau in China, Hu et al. (2010a) also found stronger time stability in drier locations using SDRD, while they observed stronger time stability in wetter locations when using MABE.

Time stability at point scales depends on different periods and soil depths. Generally, stronger time stability was found in wetter periods (Guber et al., 2008; Zhao et al., 2010; Zhou et al., 2007). Zhou et al. (2007) believed that the lateral redistribution of soil moisture in the wet/dry transition period and uneven evapotranspiration and root water uptake in drying process was the reasons for the relatively higher time stability during wetting periods. Guber et al. (2008) related the weaker time stability in the May–June period to the active vegetative growth of corn. Zhao et al. (2010) attributed the high temporal stability under wet conditions to an enhanced capillary movement of water from the subsoil to the topsoil. Time stability at point scale was usually found to be stronger at deeper depths observed by different time stability indices (Guberet al., 2008; Hu et al., 2010a, 2010b; Starks et al., 2006). This was generally in agreement with the relationship of time stability of spatial pattern with depth (Cassel et al., 2000).

Due to the complex influences of various factors on the time stability at the point scale, future work should be devoted to quantitatively define the contribution of various factors on time stability of soil water content. Models of the relationship of time stability and environmental factors should also be developed for the purpose of identification of time stable positions with known soil and topographic information. In addition, because of the different relationships between factors and different time stability indices (Hu et al., 2010a), attention should also be paid to the selection of time stability index for average SWC estimation in a field.

6. Conclusion

We reviewed the studies on time stability of soil water content, including its concept, methodology, application, and controlling factors. The following conclusions can be drawn: (1) Two types of time stability of soil water content can be classified; one is to describe the overall similarity of soil water content spatial patterns between different occasions, and the other is to describe the time invariance of the relative soil water content with time at point scales; (2) Time stability can be used to analyze the time persistence of spatial pattern of soil water content. Spearman’s rank correlation, Pearson correlation analysis, average spatial coefficient of determination, and temporal coefficient of variability can be used to characterize the time stability of spatial pattern for two spatial series. Empirical orthogonal function can be used to extract the same spatial patterns of soil moisture for multiple spatial series. Spectral coherency, wavelet coherency, and multivariate empirical mode decomposition can be used to identify the scale specific time stability of a spatial pattern. However, the spectral coherency method assumes spatial series are stationary and linear, wavelet methods assume the spatial series to be linear, and both of them are only applicable to the case of two spatial series. For multivariate empirical mode decomposition, it can deal with non-stationary and nonlinear systems with multiple spatial series; (3) Time stability can be used to identify the most time-stable locations for average SWC evaluation for a given field. Standard deviation of relative difference, root mean square error, standard deviation of soil water content, width of the 90% empirical tolerance interval of relative water content, and mean absolute bias error can be used to identify the most time-stable location. However, the performance of these indices should be compared in different environmental backgrounds to look for the most suitable index; (4) Time stability of soil water content can be influenced by many factors, including soil, topography, vegetation, and climate. It also shows non-local control and local control, depending on the soil water conditions. Knowledge of the controlling factors for time stability of spatial pattern and that at point scale is important for understanding the soil water processes, soil water management, and identification of time-stable locations for average SWC evaluation. Future work should quantify the influences of different factors on time stability, and make deeper understanding of nonlocal and local control on time stability of soil water content; (5) Time stability of soil water content is also depth-wise and scale specific. Usually, time stability of soil moisture increased with depth and scale. Time stability analysis can not only facilitate our understanding to soil water related processes but also greatly reduce the sampling work for soil water content in fields.

Acknowledgements

The project was funded by Natural Science and Engineering council (NSERC) of Canada.