Climate-driven soil suction variation using a natural-order Fourier series approach

. In unsaturated soil practice, the near-ground-surface moisture flow is commonly evaluated using the 1-dimensional (1D) suction diffusion equation. This study presents the application of a Natural-Order Fourier Series (NOFS) approach for representing monthly climate-driven soil suction variations near the ground surface, which is a boundary condition that is generally difficult to adequately model. The NOFS incorporates an algorithmic selection criterion to optimize the order of the Fourier series to improve the capture of the seasonal shifts and extreme climate periods, while maintaining acceptable computation efficiency. The climate-soil suction interaction is expressed using published empirical relationships between monthly rainfall, temperature, and soil index properties. An example and validation study of the proposed NOFS selection approach is presented using measured soil properties and historical weather data at a study site in Denver, Colorado USA. Key findings from performance and stability studies of five locations in the United States associated with differing climate regions are discussed. Limitations and recommendations for implementation are also included. The proposed NOFS approach for capturing climate-driven changes in suction near the ground surface can be efficiently implemented in unsaturated soil numerical analyses that are governed by moisture-dependent mechanical soil behavior and can help improve the computation time, stability, and performance associated with stochastic simulations.


Introduction
The adequate quantification of the moisture flux at the ground surface and the transient state of the moisture flow through the vadose zone continues to be a challenging problem in unsaturated soil analyses. The near-ground-surface moisture flow is commonly evaluated using the 1-dimensional (1D) diffusion equation derived by Mitchell in 1979 [1], which is a closed-form solution of the unsaturated flow equation [2]. The 1D diffusion equation can be used to model soil suction as a function of depth and time, provided that the boundary conditions can be established and there are negligible impacts from groundwater. The long-term minimum and maximum suction through the suction profile create the suction envelope, which is commonly used in foundation and pavement design around the globe [3,4,5,6].
There have been many successful studies which sought to measure the parameters and boundary conditions involved with the 1D diffusion equation and long-term suction envelopes. Several of those research efforts were able to relate the magnitude and depth of equilibrium suction to the long-term average (typically over 20 to 30 years) Thornthwaite Moisture Index (TMI) [5,7], referred to by the authors as the equilibrium TMI (TMIeq). Perera et al. [8] and Zapata et al. [9] developed useful models to relate the equilibrium suction of pavement subgrade with soil index properties and TMIeq. Vann & Houston [10] updated and developed several * Corresponding author: austin.olaiz@asu.edu new models for uncovered sites by relating TMIeq to equilibrium suction, magnitude of equilibrium suction, and the maximum expected change in climate-driven surficial suction using a compilation of literature studies [11][12][13][14][15][16][17][18] and an extensive field/lab investigation of sites in the United States (US) spanning differing climatic regions.
The surface flux boundary condition for unsaturated soils exposed to the natural environment is referred to herein as the climate-driven variation in surficial suction. A large portion of unsaturated soil analyses for pavement and shallow foundation design incorporate only one steady-state change in surficial suction, which is assumed to conservatively represent the largest shift from the initial to final suction state. For time-varying analyses, periodic functions are commonly applied in practice and have been shown to provide relatively simple and effective representations of the climatedriven variation in surficial suction, provided that the seasonal frequency and amplitudes are relatively stable year to year [3]. Unfortunately, in many scenarios periodic fits do not capture the minimum and maximum values of time-series data, which for climatic modelling may result in the extreme weather events not being adequately represented. Such extreme weather events like prolonged periods of drought and high heat, or consecutive intense rainfall events, can significantly impact infrastructure built with or constructed on unsaturated soils, due to the highly nonlinear mechanical response of soil to changing moisture states.
In 1979, Mitchell also proposed that an improved representation of the irregular patterns of the surface flux boundary condition could be achieved using a Fourier analysis [1]. Follow up studies by Aubeny & Long [18] and Olaiz et al. [19] preliminary demonstrated how the Fourier series can improve the representation of the climate-driven variation of surficial suction compared to a periodic fit. This improvement can be observed in Fig. 1, adapted from [19], which presents the estimated monthly suction near the ground surface over a 20-year period, which resulted in a better representation of the shifts from dry (high suction) to wet (low suction) periods and a significant increase of the adjusted-R 2 from 0.33 with a periodic fit to 0.73 with an 8th-order Fourier fit, as well as a decrease in the sum of square errors from 13.3 to 5.4, respectively. The improved goodness-of-fit with increased Fourier series order (k) is promising, however there are considerable limitations when applying a fixed order fit to time series data. Depending on the duration of the series and the number of differing frequencies and/or amplitudes, a significant number of terms may be required to encompass both the seasonal shifts and the extreme values (maximum and minimum) of the individual monthly data.
For most deterministic analyses in geotechnical engineering, limitations due to computational efficiency are near obsolete thanks to modern computer performance; however, stochastic models which generally need more input parameters and require thousands of simulation iterations can still push computation times past a feasible level for industry acceptance.
This study aimed to develop a new practical and automated approach to better capture the sporadic seasonal shifts and the extreme climatic events using an optimized Fourier series fit for use in 1D unsaturated soil-moisture flow numerical models involving timevarying climate-driven surface flux boundary conditions [1]. The authors strived to develop an approach that sufficiently and conservatively balanced increases in goodness-of-fit and the computational effort needed to achieve such fit.

Methodology
For infrastructure built upon or with unsaturated soils, incorporating all seasonal fluctuations and extreme values in the surface flux boundary condition of a numerical model will generally provide a more reasonable and conservative representation of the natural climate-driven soil response. Generally, the goodness-of-fit of a Fourier series representing discrete time series data improves as the order (number of sine and cosine terms) increases. Although any order can be chosen for a Fourier series fit to time-series data, the computational complexity also tends to increase as the order increases.
This study presents the application of an algorithmic Fourier analysis with an automated order selection criterion, referred to as the Natural-Order Fourier Series (NOFS), for modelling monthly climate-driven variation in surficial soil suction, which considers both the key seasonal shifts and the extreme climate events. The NOFS selection criterion consists of three statistical evaluations that evaluate the goodness-of-fit of the overall function, the initial condition, and the extreme values (minimum and maximum) of the time-series data set. Soil data from the extensive Arizona State University and National Science Foundation (ASU-NSF) Expansive Soil Study (NSF #1462358) [10,20,21] was used to test the performance and sensitivity of the proposed NOFS approach.

Climate-Soil Suction Soil Relationships
The time-varying relationship between climate and soil suction adopted herein for pavement subgrades (covered) and undeveloped (uncovered) sites use empirical models proposed by Zapata et al. [9] and Olaiz [22], respectively. The variation in soil suction for pavement subgrades are expected to exhibit lower amplitudes, but possibly similar fluctuations in frequency, compared to that of an uncovered surface (ceteris paribus), assuming the uncovered soil surface is exposed and continuously impacted by environmental effects (temperature, rainfall, humidity, wind, etc.).
Zapata et al. [9] improved on the work by Perera et al. [8] who developed empirical relationships between in-situ moisture content, soil suction, equilibrium TMI, and index properties of pavement subgrades from extensive field/lab testing at 44 sites across the US. The study included samples from below the pavement and from the shoulder of the roadways; however, only the relationship between TMIeq and the suction of the samples from below the pavement provided statistically significant results. The TMI-P200/wPI model was developed to estimate the matric suction ( ) in the subbase and/or subgrade material below the pavement section using the long-term average TMI (TMIeq) and the weighted Plasticity Index (wPI), which is the product of the Plasticity Index (PI) and the percent finer than the #200 sieve ( 200 P ) in decimal form: where,  is the matric suction of the soil in kPa; and α, β, γ, and δ are regression constants. For cases with materials that exhibit wPI < 0.5 (generally coarsegrained subbases), the regression constants are expressed as: Although the TMI-P200/wPI model was originally developed using a long-term average TMIeq, the model has potential to characterize seasonal climate-soil interaction near the ground surface (or just below the pavement surface), by substituting the long-term TMIeq in Equation (1) with a shorter-term average TMI, referred to by the authors as the seasonal TMI (TMIsn). The seasonal TMI adopted herein represents a 1-year average TMI calculated on a monthly basis and has been previously demonstrated as a useful parameter for describing time-varying climate-suction relationships [5,19,22,23].
The application of the TMI-P200/wPI model to seasonal estimations of soil suction can be considered valid given the understanding that: 1) for most scenarios, soil suction near the ground surface naturally deviates from an equilibrium suction state as seasons change and during various weather events (captured by TMIsn) but generally trends back towards the equilibrium state overtime (represented by TMIeq,), and 2) the soil suction at the ground surface (or just below the pavement surface) will be most affected by shortterm climate events but the suction changes throughout the rest of active zone below will lag behind those at the surface due to the relatively long equilibration time inherent in unsaturated soil moisture flow. As such, it can be assumed that if TMIsn began to stabilize over a longer period, the updated monthly calculations of TMIeq, would begin to trend towards TMIsn.
To further explore the applicability of using shortterm climatic parameters, such as TMIsn, to quantify seasonal variations in soil suction, Olaiz [22] statistically evalauted the upper 0.5 meters of soil borings from the ASU-NSF Expansive Soil Study database, which included data mining of over 100 new and historical geotechnical field/lab studies at sites associated with expansive soils throughout Arizona, Colorado, New Mexico, Texas, and Oklahoma (USA where, ψ is in pF (log cm of water), i represents the current month, and α, β, and γ are the regression coefficients based on the wPI, presented in Table 1.

Proposed Natural-Order Fourier Series (NOFS) Algorithm
The proposed NOFS approach for optimizing the equation fit to the sporadic climate-suction time-series data uses a selection criterion of three statistical evaluations which check various aspects of the Fourier series fit to the time-varying suction data. Note that the criteria thresholds presented herein are based on evaluations of the sites used in this study and are preliminary recommendations only. The threshold values may require a site-specific sensitivity study prior to implementation. The first NOFS check to optimize the representation of the monthly variation of suction near the ground surface uses the Mean Absolute Deviation (MAD) of the Fourier fit equation to the monthly suction data. This threshold provides the goodness of fit of the Fourier series equation to the monthly suction data. It is expected that the rate of decrease in the MAD initially starts of high as the order of the Fourier series moves away from a near-periodic function (k<4) to more asymmetrical representations. Based on the iterative optimization process performed during this study, it appears that the MAD of the Fourier series fit tends to fall below 0.05 pF once the order of the series has increase past approximately 25% of the maximized fit (i.e., order is equivalent to number of data points in time series). Furthermore, the order of the Fourier fit associated with a 0.05 pF MAD threshold value did not govern the selection of the order of the fit for any of the locations explored as part of this study. As such, a MAD threshold value of 0.05 pF is adopted.
The second NOFS check to optimize the representation of the monthly variation of suction near the ground surface uses the Absolute Error of the Fourier fit to the suction data at the initial monthly time (t) step (i.e., t=0). Errors associated with the equation representing the climate-driven variation in suction are commonly caused by high discrepancies in the initial point (starting suction value). Others have sought to correct for this issue using a phase shift of the Fourier series equation [3]. This phase shift process can be relatively simple to implement on a case-by-case basis but can lead to difficulties in reaching threshold values of the other Fourier fit criteria presented herein. As such, an evaluation of the absolute error at t=0 (AEt=0) for each order of the Fourier fit is recommended. Similar to the MAD evaluation, AEt=0 generally tends to decrease with increasing terms; however, there is commonly higher volatility in this check compared to the MAD. For the sites explored as part of this study, AEt=0 generally falls below 0.05 pF around 33% of the maximum fit and below 0.1 pF around 20% of the maximized fit.
Based on the iterative optimization effort performed, along with the recommendations from Lytton [24], an acceptable variability when attempting to determine the equilibrium soil suction in the subsurface is 0.1 pF, although error ranges up to 0.4 pF have been used. As such, the absolute error of the Fourier fit to the suction data at the initial month (t = 0) less than 0.1 pF was adopted.
The third NOFS check to optimize the representation of the monthly variation of suction near the ground surface uses the residual error of the Fourier series to the estimated suction at the extreme values (REEV). The residual error of the maximum and minimum data points in the time series can be used to improve the representation of extreme climate events. The residual error at the extreme values were quantified in a conservative direction which forces the Fourier equation fit to capture or be greater in absolute magnitude than the extreme values (i.e., the equation fit at the month associated with the maximum suction must be nearly equal to or greater than that maximum value, and the fit at the minimum suction must be nearly equal to or less than that minimum suction value).
This approach reduces the possibility of missing the crucial effects of extreme weather scenarios in the computation model. Although the desired threshold for the residual error at the extreme values should be zero (i.e., the model either mimics or over predicts the extreme clime events), such threshold results in orders of the Fourier series closer to the maximized fit and were not reproducible with the data evaluated in this study. The best fit is not always the most useful fit. As such, a more practical threshold 0.05 pF was adopted for the residual error at the extreme values (calculated in the conservative direction).

NOFS Selection Example
An example of the automated algorithm for the NOFS approach to optimize the equation fit to the climatedriven variation in soil suction is explored using Denver ASU-NSF Expansive Soil Study site. An average of the two uncovered borings from the Denver, CO (USA) site (DEN-2-U-N and DEN-3-U-N) was used to generate this example. Refer to Vann [20] for the detailed location and soil information at this study site. Thirty years of historical climate data (i.e., 360 data points) prior 8/2020 including the monthly average temperate and monthly rainfall were gathered from the nearby National Oceanic Atmospheric Administration (NOAA) weather station (ID: USW00023066). TMIsn was calculated for each month and the climate-suction relationship for uncovered sites, Equation (8), was used to generate empirical estimates of the time-varying suction at the ground surface (Fig. 2).

Fig. 2.
TMIsn and variation in surficial suction [22] for the Denver study site (8/1990 -8/2020). Fig. 3 presents the three NOFS criteria checks for the Denver study site using thresholds of MAD ≤ 0.05 pF, AEt=0 ≤ 0.1 pF, and REEV ≤ 0.05 pF. The Mean Absolute Deviation decreased from over 0.17 pF at k=8 to less than 0.05 pF around k=60. The Absolute Error at Initial Month was initially sporadic but stabilized below the threshold of 0.1 pF around k=75. The Residual Error at Extreme Value with a threshold of 0.05 pF governed the order selection at k=81. As such, the optimized minimum number of terms which provided a natural representation of the monthly data, as defined by the selection thresholds, was 81 (i.e., 81 sine and cosine terms).
For the study site evaluated, the algorithmic NOFS process provides a useful approach for producing a representative equation for the climate-soil suction relationship that encompasses the key natural amplitudes, frequencies, and extreme values while maintaining feasible computational efficiency by minimizing number of sine and cosine terms. Fig. 4 illustrates three different Fourier series fits to the 360 months of estimated surficial suction data to provide a visual comparison between a low order fit with k=8, the NOFS fit of k=81, and the maximum order fit of k=180. Fig. 4. Comparison of differing order fits to the variation in surficial suction for the Denver study site (8/1990 -8/2020).
Although the lower order Fourier fit (k=8) in Fig. 4 provides a reasonable representation of the yearly variations, it lacks the inclusion of several seasonal cycles and extreme events such as: 1) the full extent of the multiple seasonal cycles from relatively high (dry) to relatively low (wet) suction periods between 8/1990 and 1/1997, 2) the full magnitude of the relatively large shift from a wet season in 08/2001 to the extreme dry period (maximum suction) in 08/2002, and 3) the full magnitude of the largest shift in the data set from an extreme dry period (near maximum suction) in 1/2013 to the extreme wet period (minimum suction) in 7/2014. The first missed climate pattern of the lower order fit could provide misleading analyses results as fit indicates a relatively stable climate for the first few years of the design life, which could result in a conclusion that the soil reaches an equilibrium state much sooner than actual.
The second and third missed climate pattern identified may be more crucial than the previous, especially for analyses of climate-driven volume change, shear strength and water balance. The severe shift from a prolonged period in the extreme dry condition to the extreme wet state in such a short period of time would significantly affect the mechanical behaviour of the soil. Such extreme shifts in climate patterns could cause a soil to surpass its limiting design state but are not generally captured by standard practices for engineering foundation, pavement, slope, and evapotranspirative (ET) cover design. For example, a structure/pavement founded on very fine-grained soils could exhibit wide and deep desiccation cracks due to the prolonged drought, which provide direct paths for the rainfall infiltration during the extreme shift towards the wet period, resulting in moisture-driven changes in soil volume at depths much greater than those anticipated during design leading to heave-induced structural distresses. Additionally, an ET cover could reach maximum transient storage capacity due the quick shift from extreme dry to wet, surpassing the designed water balance scenarios, and resulting in runoff flow greater than that anticipated in the design.
Furthermore, the NOFS and the maximized order Fourier series fit in Fig. 4 appear to be visually equivalent, although the NOFS misses a few high frequency, low amplitude periods which fortunately are not typically of concern in unsaturated soil analysis due to the relatively long equilibration times associated with unsaturated hydraulic conductivity.

NOFS Performance Testing
The performance of the proposed NOFS approach for representing climate-driven surficial suction was evaluated through a stability and sensitivity study using five locations (Table 2), strategically chosen to fall within differing climate regions, as defined by Smith [25], Fityus [26], and the Australian Standards for Residential Slabs and Footings (AS2870) [6]. The TMI-P200/wPI climate-suction model for pavement subgrades was used to estimate the monthly variation of surficial soil suction. The stability and sensitivity of the NOFS fit to the estimated monthly suction data associated with differing locations, initial conditions (starting season), and duration of analysis (5,15, and 30 years) were explored. In general, the performance evaluation of the proposed NOFS approach indicated the following:  Differing climate regions may have differing optimal NOFS selection criteria threshold values, but the optimal thresholds generally hold true for most scenarios within a given climate region.  The month/season that the model initiated at had no noticeable impact on the estimations.  The soil suction estimated using the NOFS approach for both the pavement subgrade and uncovered soil sites did not appear to be sensitive to the duration of analysis for 15 and 30 years (i.e. the time series data used in the Fourier fit begins at least 15 years prior); however, when the duration of analysis dropped to 5 years, the estimated suction for various months differed from those estimated with 15 and 30 year durations.  Improvements in the goodness-of-fit are recorded for 15 and 30 years as the number of terms in the Fourier series increased. Significant differences in the goodness-of-fit between a standard periodic fit and the NOFS are observed as the duration of analysis increases from 5 to 30 years. Furthermore, the following additional assumptions and limitations associated with models presented herein:  The term covered is used herein when referring to pavement structures although there is still a threedimensional lateral moisture flow movement which is not accounted for in the NOFS model. This lateral movement would be directly affected by the width of the roadway and location of interest (e.g., pavement edge, wheel paths, centerline).  The Zapata et al. [9] model which relates long-term TMI values to soil suction under pavements was produced under the assumption that the soil suction measured below the pavement at the time of the field sampling corresponded to the equilibrium state of moisture reached sometime after construction of the roadway.  The error of the regression fit of the empirical models were not included in this deterministic approach for modelling monthly changes in the soil suction profiles.

Conclusions
An algorithmic model for the Natural Order Fourier Series fit of time-varying surface suction data was developed as part of this study to improve the representation of the unsaturated soil-climate interaction while attempting to avoid overfitting of empirical estimated data (compounding errors). The NOFS provides an improved representation of the sporadic climate-driven changes in soil suction by capturing key seasonal shifts and the extreme climate events. The proposed NOFS approach can be effectively implemented in geotechnical numerical modelling of unsaturated soils such as the design of ET covers, foundation/subgrade design on problematic soils, slope stability analyses in dry temperate to arid regions, and geo-forensic evaluations involving moisture-dependent soils. The algorithmic NOFS selection approach can significantly improve the computation time, stability, and the overall performance associated with the repetitive analyses required for stochastic simulations.