Spatial and temporal variability of soil freeze-thaw cycling across Southern Alberta, Canada

Soil freeze-thaw cycles play an important role in all aspects of agro-ecosystems, such as crop productivity, the evolution of the soil matrix, including trace-gas emissions. In regions that experience synoptic weather conditions throughout the winter, freeze-thaw cycles generally occur in one of two categories; seasonal or winter cycles. Current soil vegetation atmosphere models (SVAT’s) often include a heat-transport soil freeze-thaw algorithm, but lack detail on complex interactions between the main driving variables. Boundary conditions for these models are often based only on a few climate variables and typically lack regional context. A nested statistical analysis was applied to identify the optimal set of environmental variables (via a stepwise regression selection procedure) to track soil freeze-thaw dynamics. Historical data collected between the years 2006-2009, for 17 long-term climate stations distributed across southern Alberta Canada was utilized. Cross-correlation between wind speed and maximum air temperature identified Chinook-driven freeze-thaw events, with such interaction varying significantly across the region and by soil depth. Climate-soil interactions were most significant predictors of soil temperature during winter months. The seasonal freeze-thaw cycle is estimated to vary between 112 131 days, consisting of 12 20 winter cycles (1 cm depth), and 1-5 winter cycles (5 cm depth) with average lag time of 26 112 days. Freeze-thaw prediction was greatly improved when higher-order climate interaction terms were considered. Our findings highlight the importance of soil-water assumeptions within complex ecosystem SVAT models in resolving regional-driven climatic trends. Alongside improved representation of regional trends aimed at reducing model-based uncertainty, such efforts are expected to, in tandem, help advance the geostatistical design, and implementation of agroenvironmental monitoring systems that combine in-situ and satellite/remotesensing derived estimates of near-surface soil moisture.


INTRODUCTION
Soil freeze-thaw cycles affect plant growth, greenhouse gas emissions, soil structure, stability, microbe populations and nutrients.Freeze-thaw cycles have been shown to increase grassland above-ground productivity, while causing root damage [1].Significant increases in nitrous oxide (N 2 O) emissions, a major greenhouse gas that contributes to global warming from agricultural systems, has been documented during thawing events, with the greatest emissions occurring during the spring thaw [2,3].Freeze-thaw variability regulates the cumulative release of N 2 O by controlling how water, in its different states, physically blocks emission from the soil surface.At depth, soil temperature, that regulates biological microbial activity, is also affected by freeze-thaw variability.Freeze-thaw cycles reduce soil aggregate size and make soils more susceptible to erosion [4].There is evidence that soil microorganism populations adapt in order to survive in regions that experience multiple winter freezethaw cycles [5].Recent findings indicate that climate change may increase the frequency of freeze-thaw cycles [6].Alongside anticipated changes in the frequency of freeze-thaw events, the amplitude or strength of influence is also expected to increase, given the longer absence of insulating snow cover (i.e., longer growing seasons at polar latitudes) that accompanies a rise in global mean air temperature due to global warming [7].A better understanding of freeze-thaw cycling behavior in relation to climate change, and its influence on soil conditions can improve adaptive, multi-objective decisionmaking that includes crop selection, use of herbicides, irrigation and nutrient application timing, as well as leaching and erosion management.There are a number of soil-vegetation-atmosphere models (SVAT's) that have been developed to simulate soil water flow, soil temperature and soil freezing and thawing.These models are often driven by boundary conditions based on a few climate parameters and are not parameterized with a regional-specific context using statistical relationships between observed data and relevant related variables.Models that do not consider regional-specific aspects of climate interactions may produce unrealistic model output, underestimating the impact of freeze-thaw dynamics on agricultural land and its productivity.Along with modeling accuracy improvements, related advances in rapid automated monitoring of soil and climate at fine spatial scales and the integration of multiple sensor resources may reduce measurement uncertainty.
Several existing SVAT models consider soil freezethaw cycling, such as the Daily Century Model (Daycent) [8], Ecosystem-Atmosphere Simulation Scheme (EASS) [9], Simultaneous Heat And Water (SHAW) [10,11], and HYDRUS [12] (Table 1).All these models aim to model heat transport through a one-dimensional soil profile.SHAW, HYDRUS and EASS are stochastic, while Daycent is deterministic.These models are all based on the physical properties of heat and water transport through porous media.Cycling is represented through upper boundary conditions that use climate information to drive soil freeze-thaw cycling of the underlying soil layers.HYDRUS uses mean temperature and Daycent uses diurnal temperature ranges in the boundary conditions as the main climate drivers of soil freeze-thaw dynamics.SHAW uses net-solar radiation balance and EASS uses precipitation and infiltration rates as main climate boundary condition driving variables.All of these models assume that one or two climate variables are the main drivers of freezing and thawing in the soil and other climate variable effects are negligible.The problems with current soil freeze-thaw models is that they do not consider interactions between climate variables and the boundary conditions have little regional context to drive the models.
Freeze-thaw cycles can be divided into two categories: seasonal cycles and winter cycles.Seasonal cycles are characterized by initial freezing of the soil water in the late fall or early winter when maximum and minimum air temperatures stay below 0˚C consistently, and likewise, by thawing when soil temperature stays above freezing.The soil profile freezes from the top down, the depth of freeze and rate are affected by the duration and intensity of the cold weather as well as the insulation capabilities of the soil and snow cover.The thawing portion of a cycle begins when the air temperatures increase to above 0˚C for long enough to thaw the soil water.Thawing also occurs from the soil surface down; and if the warm weather persists long enough, the entire frozen portion will thaw, otherwise a return of cool weather conditions will re-freeze the upper portion of the soil profile.Such warm weather allows for partial thawing of the frozen soil profile or even a complete thawing depending of the extent of the initial freeze, duration and intensity of the warm weather.While the physical state of soil water is determined by soil temperature, other variables such as air temperature, soil moisture content, snow cover, precipitation and wind alter soil temperature.
Winter freeze-thaw cycles are events driven by synoptic weather conditions that bring unseasonably warm temperatures to a region during a normally cold winter.Persistent wind contributes greatly to the intensity and frequency of winter freeze-thaw cycles.In Canada, such wind-driven freeze-thaw cycles are most pronounced as Chinooks that occur in Canada almost exclusively within southern and central Alberta.Chinook winds belong to a family of mountain winds that blow in regions where long mountain chains are more or less perpendicular to the prevailing wind.Chinooks are established when moist Pacific air is blown westward by prevailing winds [13].As the air mass cools, it then rises over the Rocky Mountains of British Columbia.The moist air mass condenses and precipitates over the mountain ranges.The drier Pacific air mass reaches the leeward side of the mountains and warms at the dry adiabatic lapse rate as it descends onto the prairies with high wind speeds [13].The warm moisture-deprived air has the ability to quickly melt and evaporate any snow cover, while thawing the upper layers of the soil.Temperatures during a Chinook event can increase air temperatures by more than 20˚C in a few hours [14].A Chinook can last only a matter of hours to more than 24 hours.After the Chinook dissipates the region quickly returns to normal winter temperatures.Climate-land interaction considerations provide an explanation of the occurrence probability of dominant winter freeze-thaw events.
Freeze-thaw cycling is also influenced by long-term climate teleconnections, as large-scale atmospheric and oceanic circulations, such as the El Nino/Southern Oscillation (ENSO), Pacific Decadal Oscillation (PDO), Pacific North America (PNA), and North Pacific (NP) pattern.These large scale circulations have been shown to account for a large portion of cold season temperature trends in North America [15].Teleconnections have an influence on air temperature and fresh water freeze and break-up dates in Canada [16].Air temperature is a major driving force of northern latitude freeze and thaw cycles [17].Teleconnections influencing seasonal air temperature trends in North America will, in turn, also have a driving effect on seasonal soil freeze-thaw cycles.Nonetheless, the influence of teleconnections on regional soil temperatures are difficult to quantify without knowledge about details on timing of the events associated with regional impacts.
In this paper, we investigate the influence of climatesoil interaction/s on the intensity and frequency of soil freeze-thaw dynamics in time and space.For our study region (southern Alberta, Western Canada) we apply a nested, multivariate statistical modeling approach to identify the best set of variables to predict freeze-thaw cycling based on available data.In this way, our aim was to devise statistical models that capture the influence of interactions between climate, soil and other relevant environmental variables for improved tracking of freezethaw variability in time and space.Our statistical predictions of the timing, intensity, and inter-annual variability of soil freeze-thaw cycling are important for improving the reliability of SVAT models, the identification of most significant driving variables and interaction terms.Also, we obtain best-fit models that provide surface soil temperature boundary conditions for soil water and solute vertical transport models.Our findings illustrate the importance of regional characterization of freeze-thaw variability.We obtain correlation maps for identifying where additional climate stations might be located to improve the monitoring of freeze-thaw activity and its potential impact on growing-season and winter covercrop productivity as well as agro-ecosystem resilience.Many of the world's most productive agricultural regions are situated near mountainous areas.The magnitude and frequency of freeze-thaw cycling is anticipated to increase as a result of future global projected climate changes-yet, its potential regional variability is not well understood [7].For this reason, better tracking of freeze-thaw activity in time and space, in relation to environmental variables, contributes to reducing observational and modeling uncertainty, enhancing the ability of decision-makers to identify and gauge potential variability in freeze-thaw cycling.In turn, agricultural decision-makers may better assess both short-and long-term risks to agricultural crops, in relation to freeze-thaw activity and management practices for the use-efficiency of irrigation, fertilizers and herbicides (i.e., amount and timing of applications).

Data
The region of interest for our study was Southern Alberta, Canada, which experiences winter freeze-thaw cycles due to warming effects of the Chinook winds, and longer term seasonal freeze-thaw cycles that are strongly influenced by teleconnections.Historical maximum, minimum, mean air temperatures, wind speed and precipitation for 17 climate stations distributed across the study region were obtained (Figure 1).These stations recorded soil moisture and soil temperature at various soil depths.A seasonal time series was obtained for each station from 1 May 2006 to 31 May 2009 inclusive over three winters (Figure 2).Station elevations ranged from 767 m to 1310 m above mean sea level.The climate stations are operated by Agriculture Drought Monitoring (AGDM) and Agriculture and Agri-Food Canada (AAFC).The datasets were inspected for errors, involving the conversion of hourly to daily measurements for 4 stations (AAFC stations).Day lengths for each of the stations was computed using United States Naval Observatory (USNO) day length calculator and the latitude and longitude of each station.
Table 2 provides a brief summary of the data obtained for each station.Due the coarse spatial resolution of the climate stations, before considering more advanced geospatial models, a statistical analysis was first undertaken to identify leading variables, while awaiting more fine scale data from the study region.

Methods
To identify the best or more representative statistical model to capture freeze-thaw observed historical variability, we evaluated in a nested fashion, a wide array of possible models comprising different explanatory, predictor variables and assumed interactions between them (i.e., model and variable selection step).The daily mean of each variable was calculated from the 17 station datasets.Generalized linear models (GLM) assuming a Gaussian distribution for each random variable were applied to the dataset to find the best explanatory (dependent) variable and set of predictor variables (independent).Soil temperature and moisture at 5 cm, accumulated precipitation, mean air temperature, wind speed, and day length were each tested as a dependent variable in a GLM.The remaining variables were tested as independent variables in the GLM, with no variable interactions being considered.Based on the lowest Akaike Information Criterion (AIC) and resultant residual plots, soil temperature at 5 cm was selected as the best explanatory variable (Table 3).
A linear regression model was constructed with soil temperature at 5 cm as the explanatory variable.A nested or step-wise model selection was employed that simplifies a linear model by reducing the AIC to find the best independent variables and coefficients.Mean air temperature, soil moisture at 5 cm, accumulated precipitation, wind speed, and day length (both with and without interactions) were considered as predictors for soil surface temperature (i.e., surface depth of 5 cm).Examination of the residual plots showed unsatisfactory results, so a second stepwise model selection was performed next using minimum and maximum air temperature instead of the mean.The new selected independent variables were minimum, maximum air temperature, soil moisture at 5 cm, accumulated precipitation, wind speed, and day length, with and without interactions between predictors.
The values for each climate variable from November 1 to April 1 of each winter were extracted from the mean seasonal dataset.Figure 3 shows the winter dataset for three years (i.e., 2006-2007, 2007-2008, 2008-2009).A stepwise linear model selection strategy was also employed to select the best set of predictor and interaction terms during the winter season for soil temperature at 5 cm.The predictors tested were maximum, minimum air temperature, soil moisture at 5 cm, accumulated precipitation, wind speed, and day length as well as their interaction.
To examine climate-soil interaction variability across the study region, the statistical models were fitted to data collected from all stations.Geostatistical spline-based interpolation was then performed to interpolate between stations and to generate a representative "correlation map" of freeze-thaw variability across the study region.The correlation between wind speed and maximum temperature was most representative of freeze-thaw activity and was interpolated from best-fit predictions obtained at each station location.High correlations between wind speed and maximum air temperature were linked to Chinook-driven freeze-thaw events, such that locations with high correlation values between wind speed and maximum air temperature would also experience more freeze-thaw cycles relative to locations with lower correlation values.Spline surfaces were created from the potted correlation values using varied degrees of smoothing.Using a spline surface smoothing tension of 0.7, a leave-out-one validation was performed.Each station was individually removed from the spline-surface (geostatistical) calculation, and the predicted correlation at the location of each dropped or withheld climate station was then recorded.Model cross-validation error of the withheld correlation was calculated as the absolute difference between the actual correlation value and the withheld correlation value.The leave-out-one validation errors were then ranked in descending order to highlight the stations with the greatest and least influence on the correlation surface.In order to understand the temporal variability between the three winters, the average soil freeze-thaw characteristics across southern Alberta were com-puted for each winter from November 1 to April 1.From this data, a set of best-fit statistics was calculated for each winter, this included the mean maximum, minimum air temperature, mean soil temperature at 5 cm, total soil freezing days, freeze-thaw cycles at 5 cm, estimated freeze-thaw cycles at 1cm, the mean magnitude of soil freeze, the maximum magnitude and day of soil freeze at 5 cm.
The number of soil freezing days, day of first freeze and final thaw, length of seasonal freeze-thaw cycle, mean lag between 5 cm freeze-thaw cycles and accumulated precipitation were also computed to better characterize inter-annual variation in freeze-thaw cycling behavior.The length of each winter's seasonal freeze-thaw cycle was computed by subtracting the time between the initial freeze and final thaw of each winter.The mean magnitude of soil freeze was computed as the average of the negative 5 cm soil temperatures.The number of soil freezing days each winter was the total of days each winter where the average daily soil temperature at 5 cm was below 0˚C.Soil freeze-thaw cycles at 5 cm depth were calculated by counting the number of times the soil temperature at 5 cm dropped below 0˚C then rose above freezing again.The soil temperature at 5 cm was the shallowest soil measurement available and therefore shallow freeze-thaw cycles caused by daily temperature extremes were not captured at this measurement depth.Freeze-thaw cycles at the 1 cm soil depth were estimated from daily maximum and minimum air temperatures.One freeze-thaw cycle was defined by whenever the daily maximum air temperature is greater than or equal to 6.3˚C and the daily minimum is less than or equal to -3.5˚C [18].
The seasonal time series for soil temperature at 5 cm was filtered using Fourier domain smoothing and denoising with two different threshold frequencies (0.1 and 0.05) to remove the high-frequency (i.e., noise) variability.The low frequency component of the filtered seasonal time series represents the signal as the seasonal trend in soil temperatures at 5 cm.Each filtered time series of soil temperature was then fitted with a two damped sinusoidal model.An average winter soil temperature time series was created by averaging the 5 cm soil temperature for each day between November

RESULTS AND DISCUSSION
Mean daily air temperature alone was unable to adequately describe observed trends in near surface-soil temperature, as Chinook events disrupt climate-soil interaction via rapid increases in diurnal temperature change.For this reason, minimum and maximum air temperature better tracked temperature range and freeze-thaw cycling.
Figure 4 shows the distribution of residuals for the bestfitting seasonal and winter models.The best fitting models for the seasonal and winter datasets both used maximum, minimum air temperature, soil moisture at 5 cm, accumulated precipitation and wind speed at 10 m as predictors for soil temperature at 5 cm.Best-fitting models were obtained when interactions between independent variables were considered (Table 4).T Nested/stepwise model selection for the winter dataset did not improve the prediction skill of the winter model with interactions as much as the seasonal model with interactions.While more interactions were considered in the best fit winter model, the greater number of predictor terms in the best-fitting winter model is likely due to stronger climate fluctuations during winter.The quantile-quantile plot in Figure 4 shows that both the seasonal and winter models need improved prediction accuracy at soil temperature extremes.Eq.1 represents the reduced form of the best fitting statistical model for the seasonal time series.All predictor terms of soil temperature in Eq.1 had a p < 0.001.Eq.1 can be used to estimate the soil temperature at 5 cm in the study region.Eq.2 estimates winter soil temperatures between November 1 and April 1 in the study region and represents the reduced form of the best fitting statistical model for the winter time series.All the terms used for predicting soil temperature in Eq.2 had a p < 0.01.Eqs.1 and 2 are regional-specific and therefore capture complex interactions between climate variables that drive soil temperature throughout the year and during the winter season.T s = soil temperature at 5 cm T max = daily max temperature T min =daily min temperature M = soil moisture at 5 cm W =wind speed at 10 m P =precipitation D =day length Maximum and minimum air temperature, as individual predictors of soil temperature in both the winter and seasonal models, were significant (p < 0.01), with their interaction also significant (winter p < 0.001, seasonal p < 0.05).The higher significance level of climate-soil interactions during winter indicates that diurnal temperature range is critical for predicting soil temperature during the winter season.In the seasonal model the only significant predictors of soil temperature were: day length (p < 0.001), soil moisture (p < 0.01) and precipitation (p < 0.05).The high significance level of day length in the seasonal model was expected because in Northern mid-latitudes the length of day has been previously shown to drive seasonal climate variation.Southern Alberta receives very little winter precipitation; the average accumulated precipitation for the region was 138 mm in total over the three winters.Whereas the average accumulated precipitation between May 1, 2006 and May 31, 2009 for the study region was 1078 mm.In the winter model, precipitation did not significantly (p > 0.1) contribute to the prediction of soil temperature, but precipitation was significant (p < 0.05) in predicting soil temperature in the seasonal model.Soil moisture was only a significant (p < 0.01) predictor of soil temperature in the seasonal model.Likely precipitation and soil moisture individually play a more significant role in driving soil temperature during the spring and summer when precipitation events are more frequent and soil water is unfrozen.Despite the low amount of precipitation during the winter, precipitation and soil moisture still had significant predictive ability for soil temperature when interacting with other climate variables in the winter model.Interacting terms were significant predictors of soil temperature in the winter model, indicating that winter soil temperatures are driven by complex interactions between climate variables.A notable significant interaction term in the winter model was maximum air temperature, minimum air temperature and wind speed (p < 0.05), this interaction had a positive influence on soil temperature.The maximum, minimum air temperature and wind speed interaction may be capturing the Chinook events; these events can increase daily maximum, minimum air temperatures associated with high wind speed [13].The specific interaction terms that may be related to long term synoptic climate events, such as teleconnections nonetheless require more investigation as the complex nature of how teleconnections influence climate and freeze-thaw events over large areas.
Between November 1 and April 1, wind speed and soil moisture were strongly correlated across all stations, with correlation ranging between 0.81 at Barnwell and Stavely to 0.96 at Oyen, with a mean of 0.88.The strong winter co-occurrence of soil moisture and wind speed in southern Alberta is best explained as driven by Chinook events.Chinooks bring warm-dry-high-speed winds that have the ability to melt snow and thaw the upper layers of soil [14], snow melt that is not evaporated can infiltrate into the thawed upper soil layers.Maximum air temperature and wind speed were not as strongly correlated as wind speed and soil moisture during the winter, 0.66 at Onefour to 0.89 at Oyen with an average correlation of 0.71 (Table 5).Figure 5 shows significant spatial variability in the correlation between maximum air temperature and wind speed.The various fitted spline surfaces of this correlation pattern of climate-soil variability reveal a similar spatial pattern to the observed regional pattern of Alberta's Chinook signal from 1951-1990 [13].The similar spatial pattern between Alberta's average Chinook signal and the maximum air temperature and wind speed correlation surface supports the idea that a Openly accessible at http://www.scirp.org/journal/AS/401 401 Table 5. Best-fit model cross-validation (leave-one-out) statistical results: Maximum air temperature and wind speed correlations at each station before and after each station is withheld from the spatial interpolation, with the absolute error of predicted correlation value error and error rank of each station.The maximum air temperature and wind speed correlation values were used as an indicator potential freeze-thaw cycling due the relationship between Chinook events, unseasonably warm temperatures and wind speed.Chinook event is a complex interaction between two or more climate variables.Locations with high correlation values between wind speed and maximum air temperature were expected to have an increased occurrence of Chinook events and freeze-thaw cycles.Table 5 shows the results of the leave-one-out cross-validation, the pre-dicted wind speed and maximum temperature correlation value at the location of each withheld station was recorded, the absolute difference between observed correlation value and the predicted values were calculated and ranked.Barnwell, Lethbridge and Lacombe had the lowest absolute difference between observed and predicted correlation values, this can be explained by these stations' relative close proximity to surrounding stations with similar wind speed and maximum temperature correlation values (Table 5).The results of the leave-oneout validation show that Oyen, Bodo and Schuler had the largest errors (Table 5) and therefore were the most influential stations in calculation of the spline surfaces (Figure 6).Oyen, Bodo and Schuler are all in the most eastern part of the study region with a large spatial gap between these monitoring stations and the nearest stations to the west.For example, there is a larger variation in correlation pattern at Hussar (0.66) and Oyen (0.89).
In order to better track climatic influences on agro-ecosystems, additional monitoring stations likely need to be placed near such locations where higher observed correlation between climate and soil interaction is detected.
The area between Hussar and to the east Oyen, Bodo and Schuler is a prime candidate for an additional climate monitoring station that would improve the tracking of Chinook related climate interactions in eastern Alberta.The winter of 2006/2007 was the mildest of the three winters; it had the warmest mean soil temperature, mean maximum and minimum air temperatures (Table 6).The average and the maximum magnitude of soil freeze were warmer than the following winters.Consistently warmer air temperatures would have led to repetitive thawing and refreezing of the upper 1 cm of the soil leading to an estimated 18 freeze-thaw cycles at 1 cm soil depth.One explanation for there being only 4 freeze-thaw cycles at 5 cm depth is soil temperature at 1 cm depth will be sensitive to daily fluctuations in air temperature, whereas soil at 5 cm depth will need prolonged exposure to warm temperatures to thaw [19].The difference between the length of the seasonal freeze-thaw cycle and soil freezing days within the cycle show that when the soil thawed at 5 cm depth it remained unfrozen for multiple days.The 2006/2007 winter received the second most winter precipitation that may have provided some insulation to the soil [20] preventing the soil from freezing to the 5 cm depth.The 2006/2007 winter had the earliest date of initial soil freeze and final thaw at 5 cm, cold temperatures may have come early with warmer temperatures coming earlier as well.The date of maximum magnitude of soil freeze being in mid-January instead of later in the month provides evidence of an earlier winter.The winter of 2007/2008 had the lowest amount of accumulated winter precipitation; this may have contributed to the low average soil temperature, the coldest minimum soil temperature and the low number of freeze-thaw cycles at 5 cm depth (Table 6).The lack of winter precipitation may have led to minimal snow cover to insulate the soil.Without much snow cover the soil surface would be more susceptible to air temperature extremes [20].The 2007/2008 winter experienced the largest range between the mean minimum air temperature and the mean maximum air temperature; this would explain the larger number of 1 cm depth freeze-thaw cycles when compared to the 2008/2009 winter.A contributing factor to the lack of freeze-thaw cycles at 5 cm depth during the 2007/2008 winter may be that the positive temperatures lacked persistence and intensity to thaw the soil much below 1 cm depth, leading to an estimated 20 cycles near the surface and yet only 1 at the 5 cm depth.An interesting aspect of the 2007/2008 winter is that the soil freezing days and length of the seasonal cycle are the same, supporting the idea that the warm weather was never able to thaw the soil to 5 cm.The date of coldest soil temperature, initial soil freeze and final thaw indicate that the 2007/2008 winter was shorter, with a late start and earlier end than the preceding and following winters.The coldest of the winters was 2008/2009; the mean maximum and minimum air temperatures were the lowest of the three winters.The 2008/2009 winter had a warmer mean soil temperature, mean and maximum magnitude of soil freeze at  5 cm than the previous winter (Table 6).More accumulated winter precipitation is likely the main reason for warmer soils, snow cover may have insulated the soil from colder air temperatures [20].The 2008/2009 winter had the least number of soil freeze-thaw cycles at 1 cm (12) and the most at 5 cm (5); this is likely caused by a few intense warming events with a long enough duration thaw deeper soil layers instead of many short-duration low-intensity warming events.The 2008/2009 winter had the latest day of initial soil freeze and final thaw at 5 cm but was still the longest winter as indicated by the soil freezing days (Table 6).Figure 7 shows the fitted B-spline curve fitted to the average winter soil temperature time series for southern Alberta.The B-spline fitted curve had a R² = 0.99 when compared to the actual average winter soil temperature at 5 cm depth and was used as the winter season regionalspecific portion of Eq.3.The B-spline cure fitted tends to overestimate the peaks and under estimate troughs observed in the average winter soil temperature time series.The regional-specific soil surface boundary condition is represented by Eq.3 and can be used to drive freeze-thaw algorithms using regional-specific climate data.Fitting sinusoidal curves to a temporally filtered seasonal signal provides the regional-specific seasonal boundary condition for seasonal soil temperature trends in Eq.3.Fitting the B-spline curve the average winter time series provides a regional-specific winter boundary condition that can estimate winter freeze-thaw cycles.The proposed model was developed using regional-specific climate data instead of generic boundary conditions that require only one or two climate variable inputs.

CONCLUSIONS
Our findings illustrate the importance of regional context in modeling freeze-thaw cycles.Our results indicate that the interaction between air temperature and wind drives variability in soil temperatures and if included when specifying boundary conditions or cali-brating to historical climate, may significantly improve the prediction skills of models framed at many different spatial and temporal scales, e.g.soil water and solute vertical transport or ecosystem scale models.By considering a nested statistical modeling approach to capture regional-specific freeze-thaw variability, we utilize a set of best-fit measures (leave-one-out cross validation), and geostatistical spine interpolation, to better identify the influence of measurement and model-based uncertainty on our model predictions of freeze-thaw variability over our study region and both annual-scale and winter timeperiods.
The spatial analysis of the correlation of wind speed and maximum air temperature revealed areas of southern Alberta where the correlation between these two variables varies most strongly.This kind of spatial analysis that reveals the most influential climate monitoring stations on model predictions and can be used to infer where additional climate monitoring stations could be located to better track climate signals and freeze-thaw dynamics within the study region (southern Alberta).The interwinter variability revealed that freeze-thaw cycling behaviour varies differently in the winter months compared to the full annual seasonal cycle.From this result, we infer that, to more reliably predict and track freeze-thaw activity, its spatial and temporal variability needs to be better represented.Often historical data is limited, but nonetheless statistical modeling can be performed and used to guide how new monitoring data can, in the future, be collected to improve SVAT model prediction skill.Specifically, our findings suggest that freeze-thaw signalling within southern Alberta is strongest during winter months at several key monitoring stations.Where climate-soil correlation is most statistically significant are guiding regions where monitoring activity could be sampled at higher frequency and at greater resolution using networks of wireless soil moisture sensors.
Our findings also provide a statistical soil surface boundary condition that captures both the seasonal and winter-scale variability of soil temperature.This boundary condition can be incorporated in SVAT and soil water/ solute vertical transport models.Our work illustrates how models can be given regional context, by calibrating them to a representative statistical model that captures main climatic trends, regional variation, and climate-soil interaction terms.In such a statistical representative model, we show here that freeze-thaw cycling is one example of a climate-soil interaction term.Using this statistical approach, regional calibration of models that require tracking of soil water movement can be improved, especially in regions where freeze-thaw activity is strong.Use of a single, simplified statistical model also potentially offers a way to simplify climate-related inputs to complex SVAT and agro-ecosystem models.
Future modeling work will explore how the spatial pattern and inter-annual variation of soil freeze-thaw cycling dynamics is regulated by climate teleconnections and synoptic climate events.Related work will look at the integration of fine scale climate and soil data with multiple sensor resources to help advance the geostatistical design, and implementation of agro-environmental monitoring systems that combine in-situ and satellite/ remote-sensing derived estimates of near-surface soil moisture.

=
soil moisture flow fluxes between layers l − 1 and l; l and l + 1, respectively P un = precipitation rate at ground surface E s = evaporation from soil surface ρ w = density of water P s = percolate water supplied to underlying soil layers I p = potential infiltration rate before ponding I lim = limiting value to infiltration rate after ponding z p = depth of surface pond Daily Century

Figure 1 .
Figure 1.Long term climate monitoring stations distributed within study region of Southern Alberta, Canada (n = 17).

Figure 2 .
Figure 2. Average daily climate and soil conditions from May 1, 2006 to May 31, 2009, data was collected from long term climate stations (n = 17) in southern Alberta.The 95% confidence interval is shown for maximum, minimum air temperature, soil temperature at 5 cm, soil moisture at 5 cm and day length.
1 and April 1 for the 2006/2007, 2007/2008 and 2008/2009 winters.The average winter 5 cm soil temperature time series was fitted with a B-spline curve tracking winter freeze-thaw variability.

Figure 4 .
Figure 4. Statistical distribution quantiles (5%, 95%) and pattern of residuals for the seasonal model (left) and the winter model (right).The standardized residuals lie within two standard deviations of the observations on the scale-location plot of both models.Both the seasonal and winter models loose accuracy at extremes as shown in the quantile-quantile plots.

Figure 5 .
Figure 5. Predicted strength of spatial occurrence of freezethaw, based on station monitoring data (n = 17) for maximum air temperature and wind speed (10 m) using thin-plate spline interpolation under varying degrees of smoothing.

Figure 6 .
Figure 6.The results of the leave-one-validation for the three stations with the largest absolute error of the predicted strength of freeze-thaw spatial occurrence (spline tension = 0.7).

Figure 7 .
Figure 7. B-spline curve fitted to average soil temperature (5 cm depth, R² = 0.99) in southern Alberta between November 1 and April 1.Fitted soil temperature B-spline curve represents the boundary condition calibrated with a regional-specific context.All fitted soil temperature fit within the observed soil temperature at 5 cm depth 95th percentiles.

Table 1 .
A summary of leading parameters, assumptions and corresponding upper boundary conditions for current soil water freezethaw algorithms.

Table 2 .
Summary of the data collected for each climate station in southern Alberta (n = 17).

Table 3 .
Summary of the criterion used to judge the best dependent variable in a GLM, all models had 1121 residual deviance degrees of freedom and 1126 null deviance degrees of freedom.

Table 4 .
Soil temperature predicted over a seasonal and winter time series, with and without interaction term considered, where p < 0.001 for all models.DF = Degrees of freedom.

Table 6 .
The average climate variability of the southern Albertan winters (November 1 to April 1) for each winter.Each reported variable is an average value computed using all the southern Alberta climates stations (n = 17).