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

doi:10.4236/as.2011.24051

Paper Menu >>

Journal Menu >>

Vol.2, No.4, 392-405 (2011)

doi:10.4236/as.2011.24051

Agricultural Scienc es

Spatial and temporal variability of soil freeze-thaw

cycling across Southern Alberta, Canada

Andrew J. Phillips1,2*, Nathaniel K. Newlands1

1Environmental Health, Agriculture and Agri-Food Canada, Lethbridge Research Centre, Lethbridge, Canada;

*Corresponding Author: ajphilli@ucalgary.ca

2Department of Geomatics Engineering, Schulich School of Engineering, University of Calgary, Calgary, Canada.

Received 9 September 2011; revised 21 October 2011; accepted 30 October 2011.

ABSTRACT

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 through-

out the winter, freeze-thaw cycles generally oc-

cur in one of two categories; seasonal or winter

cycles. Current soil vegetation atmosphere mod-

els (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 typi-

cally lack regional context. A nested statistical

analysis was applied to identify the optimal set

of environmental variables (via a stepwise re-

gression 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 tem-

perature identified Chinook-driven freeze-thaw

events, with such interaction varying signifi-

cantly across the region and by soil depth. Cli-

mate-soil interactions were most significant

predictors of soil temperature during winter

months. The seasonal freeze-thaw cycle is es-

timated to vary between 112 - 131 days, consist-

ing 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 in-

teraction terms were considered. Our findings

highlight the importance of soil-water assume-

ptions within complex ecosystem SVAT models

in resolving regional-driven climatic trends. Alon g -

side improved representation of regional trends

aimed at reducing model-based uncertainty,

such efforts are expected to, in tandem, help

advance the geostatistical design, and imple-

mentation of agroenvironmental monitoring sys-

tems that combine in-situ and satellite/remote-

sensing derived estimates of near-surface soil

moisture.

Keywords: Freeze-Thaw; Soil Temperature;

Agro-Ecosystem Modeling; Regional Climate;

Soil Science

1. INTRODUCTION

Soil freeze-thaw cycles affect plant growth, greenhouse

gas emissions, soil structure, stability, microbe popula-

tions and nutrients. Freeze-thaw cycles have been shown

to increase grassland above-ground productivity, while

causing root damage [1]. Significant increases in nitrous

oxide (N2O) emissions, a major greenhouse gas that con-

tributes to global warming from agricultural systems, has

been documented during thawing events, with the great-

est emissions occurring during the spring thaw [2,3].

Freeze-thaw variability regulates the cumulative release

of N2O by controlling how water, in its different states,

physically blocks emission from the soil surface. At

depth, soil temperature, that regulates biological micro-

bial 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 freeze-

thaw 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 influ-

ence is also expected to increase, given the longer ab-

sence of insulating snow cover (i.e., longer growing sea-

sons at polar latitudes) that accompanies a rise in global

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

393393

mean air temperature due to global warming [7]. A better

understanding of freeze-thaw cycling behavior in rela-

tion to climate change, and its influence on soil condi-

tions can improve adaptive, multi-objective decision-

making 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 tem-

perature and soil freezing and thawing. These models are

often driven by boundary conditions based on a few cli-

mate parameters and are not parameterized with a re-

gional-specific context using statistical relationships be-

tween observed data and relevant related variables.

Models that do not consider regional-specific aspects of

climate interactions may produce unrealistic model out-

put, 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 re-

sources may reduce measurement uncertainty.

Several existing SVAT models consider soil freeze-

thaw 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 Day-

cent is deterministic. These models are all based on the

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

thaw algorithms.

Ecosystem Model Type Key Parameters Assumptions Upper Boundary Conditions

Ecosystem-atmosphere

simulation scheme

(EASS)

Stochastic

Latent heat of fusion,

equilibrium total water

potential, ice potential,

absolute temperature

1-D soil layer model,

explicit thermal

separation of

vegetation & soil,

water & energy

balance coupled,







11, 0,1

lim

min,,0

,0 1

uns wspp

PE PIz

Iz l









 





11, 0,1

= soil moisture flow fluxes between

layers l − 1 and l; l and l + 1, re-

spectively

Pun = precipitation rate at ground surface

Es = evaporation from soil surface

ρw = density of water

Ps = percolate water supplied to underlying soil

layers

Ip = potential infiltration rate before ponding

Ilim = limiting value to infiltration rate after

ponding

zp = depth of surface pond

Daily Century-Model

(DAYCENT) Deterministic

Diurnal temperature extremes,

depth, thermal diffusivity,

time step,

1-D, soil thermal

diffusion constant

with depth, diffusion

with no advection,

sinusoidal pattern of

soil surface

temperature

0.5





















Rd = diurnal temperature range at depth

Rs = diurnal temperature range at soil surface

d = depth

α = thermal diffusivity

p = time step

HYDRUS 1D Stochastic

Apparent thermal conductivity of

soil, temperature, time, depth,

coefficient of apparent thermal

conductivity, apparent thermal

conductivity of water,

oxygen uptake rate

1-D, convection &

conduction heat

transport, coupled

energy, vapour

& water transport

2π7π

sin 12

TTA p



 





T0 = temperature at soil surface

T = average temperature at soil surface

A = amplitude of temperature sine wave

pt = time period of temperature cycle

Simultaneous

Heat And Water

(SHAW)

Stochastic

depth, soil thermal conductivity,

soil temperature, density of water

& ice, specific heat capacity of

water & ice, liquid vapour flux,

total solutes present per mass of

soil, volumetric heat capacity of

soil, time, density of ice, latent heat

of fusion, soil vapour density, latent

heat of vaporization, vapour flux

1-D layered system,

convective heat

transfer by liquid

and latent heat

transfer by vapour,

layered system

RHLEG 

H = sensible heat flux

Rn = net all-wave radiation

E = total evapotranspiration

Lv = latent heat flux

G = soil heat flux

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

394

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 lay-

ers. 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 bound-

ary condition driving variables. All of these models as-

sume that one or two climate variables are the main

drivers of freezing and thawing in the soil and other cli-

mate variable effects are negligible. The problems with

current soil freeze-thaw models is that they do not con-

sider 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 like-

wise, 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 por-

tion 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 de-

pending 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 tempera-

ture.

Winter freeze-thaw cycles are events driven by synop-

tic weather conditions that bring unseasonably warm tem-

peratures to a region during a normally cold winter. Per-

sistent wind contributes greatly to the intensity and fre-

quency 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 con-

denses 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 dissi-

pates the region quickly returns to normal winter tem-

peratures. Climate-land interaction considerations pro-

vide 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 Os-

cillation (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 ma-

jor 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 re-

gional soil temperatures are difficult to quantify without

knowledge about details on timing of the events associ-

ated with regional impacts.

In this paper, we investigate the influence of climate-

soil 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 en-

vironmental variables for improved tracking of freeze-

thaw variability in time and space. Our statistical predic-

tions 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 tem-

perature boundary conditions for soil water and solute

vertical transport models. Our findings illustrate the im-

portance of regional characterization of freeze-thaw

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

395395

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 cover-

crop productivity as well as agro-ecosystem resilience.

Many of the world’s most productive agricultural re-

gions are situated near mountainous areas. The magni-

tude 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 en-

vironmental variables, contributes to reducing obser-

vational and modeling uncertainty, enhancing the ability

of decision-makers to identify and gauge potential vari-

ability in freeze-thaw cycling. In turn, agricultural deci-

sion-makers may better assess both short- and long-term

risks to agricultural crops, in relation to freeze-thaw ac-

tivity and management practices for the use-efficiency of

irrigation, fertilizers and herbicides (i.e., amount and

timing of applications).

2. METHODS & DATA

2.1. Data

The region of interest for our study was Southern Al-

berta, 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 pre-

cipitation 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 sta-

tions 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 sta-

tions (AAFC stations). Day lengths for each of the sta-

tions was computed using United States Naval Observa-

tory (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 geo-

spatial models, a statistical analysis was first under-

taken to identify leading variables, while awaiting more

fine scale data from the study region.

Figure 1. Long term climate monitoring stations distributed

within study region of Southern Alberta, Canada (n = 17).

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

Climate Stations Data Interval Variables Data Source

Barnwell, Bodo, Brocket, Champion,

Del Bonita, Foremost, Hussar, Morrin,

Olds, Oyen, Schuler, Stettler, Wrentham

daily

Min, Max & Mean air temperature,

Soil Temperature & Moisture at

5, 20, 50, 100 cm, 10 m Wind speed

Agriculture Drought

Monitoring Station (AGDM)1

Lacombe, Lethbridge,

Onefour, Stavely hourly

Min, Max & Mean air temperature

Soil Temperature & Moisture at

5, 20, 50, 100 cm 10 m Wind speed

Agriculture Agri-Food

Canada (AAFC)2

All Climate stations daily Day Length United States Naval

Observatory (USNO)3

1Available at: www.agric.gov.ab.ca/app116/stationview.jsp; 2Available internally through AAFC; 3Available at: aa.usno.navy.mil/data/docs/RS_OneYear.php.

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

396

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, mini-

mum air temperature, soil temperature at 5 cm, soil moisture at 5 cm and day length.

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

397397

2.2. Methods

To identify the best or more representative statistical

model to capture freeze-thaw observed historical vari-

ability, we evaluated in a nested fashion, a wide array of

possible models comprising different explanatory, pre-

dictor 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 data-

sets. Generalized linear models (GLM) assuming a

Gaussian distribution for each random variable were

applied to the dataset to find the best explanatory (de-

pendent) variable and set of predictor variables (inde-

pendent). Soil temperature and moisture at 5 cm, accu-

mulated 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 inde-

pendent variables in the GLM, with no variable interac-

tions being considered. Based on the lowest Akaike In-

formation Criterion (AIC) and resultant residual plots,

soil temperature at 5 cm was selected as the best ex-

planatory 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 simpli-

fies a linear model by reducing the AIC to find the best

independent variables and coefficients. Mean air tem-

perature, soil moisture at 5 cm, accumulated precipita-

tion, wind speed, and day length (both with and without

interactions) were considered as predictors for soil sur-

face temperature (i.e., surface depth of 5 cm). Examina-

tion of the residual plots showed unsatisfactory results,

so a second stepwise model selection was performed

next using minimum and maximum air temperature in-

stead of the mean. The new selected independent vari-

ables 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 tempera-

ture, 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 tempera-

ture 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 cor-

relation values. Spline surfaces were created from the

potted correlation values using varied degrees of smooth-

ing. 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 (geo-

statistical) 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 dif-

ference 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.

Table 3. Summary of the criterion used to judge the best dependent variable in a GLM, all models had 1121 residual deviance de-

grees of freedom and 1126 null deviance degrees of freedom.

Selected dependent variable (y),

n = 4 AIC BIC

Residual

Deviance

Null

Deviance

y = Precipitation 4872 4907 4913 5860

y = Mean Air Temperature 6205 6240 16,038 137,995

y = Soil Moisture at 5 cm 9500 9535 298,494 440,386

y = Soil Temperature at 5 cm 5386 5421 7754 97,091

y = Wind Speed at 10 m 6709 6745 25,093 30,854

y = Day Length 3914 3949 2101.1 9537

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

398

Figure 3. The soil temperatures (black) at 5 cm depth for the 2006/2007, 2007/2008 and 2008/2009

winters (November 1 to April 1). Each winter is labeled and separated by a vertical gridline.

In order to understand the temporal variability be-

tween 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 accumu-

lated precipitation were also computed to better charac-

terize inter-annual variation in freeze-thaw cycling be-

havior. 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 de-

noising with two different threshold frequencies (0.1 and

0.05) to remove the high-frequency (i.e., noise) variabil-

ity. The low frequency component of the filtered sea-

sonal 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 tem-

perature time series was created by averaging the 5 cm

soil temperature for each day between November 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.

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

399399

3. RESULTS AND DIS CUSSION

Mean daily air temperature alone was unable to ade-

quately describe observed trends in near surface-soil tem-

perature, as Chinook events disrupt climate-soil interac-

tion 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 best-

fitting seasonal and winter models. The best fitting mod-

els for the seasonal and winter datasets both used maxi-

mum, 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 inde-

pendent variables were considered (Table 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.

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.

Model-data partition Interactions AIC F-stat Terms DF Adjusted R²

seasonal no 5272 2404 6 1120 0.928

yes 4308 651 57 1169 0.971

winter no 1560 346 6 450 0.819

yes 1446 55 58 398 0.873

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

400

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 quan-

tile-quantile plot in Figure 4 shows that both the sea-

sonal and winter models need improved prediction ac-

curacy 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 tempera-

ture in Eq.1 had a p < 0.001. Eq.1 can be used to esti-

mate the soil temperature at 5 cm in the study region.

Eq.2 estimates winter soil temperatures between No-

vember 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 interac-

tions between climate variables that drive soil tempera-

ture throughout the year and during the winter season.

max min

minmax min

5.4 0.23430.0474

0.053220.007781

0.00787750.58.

TD MDPTT

PT MT T M

PMW

 





(1)

max minmaxmin

maxmin

0.1237 0.1877

0.009223.

TTT PT

TTM



 (2)

Ts = soil temperature at 5 cm

Tmax = daily max temperature

Tmin =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 individ-

ual 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 tem-

perature 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 precipi-

tation (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 pre-

viously shown to drive seasonal climate variation. Southern

Alberta receives very little winter precipitation; the av-

erage 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 pre-

cipitation 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 mois-

ture 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 sig-

nificant predictive ability for soil temperature when in-

teracting 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 interac-

tion term in the winter model was maximum air tem-

perature, 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 Chi-

nook 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 infil-

trate into the thawed upper soil layers. Maximum air

temperature and wind speed were not as strongly corre-

lated as wind speed and soil moisture during the winter,

0.66 at Onefour to 0.89 at Oyen with an average correla-

tion of 0.71 (Table 5). Figure 5 shows significant spatial

variability in the correlation between maximum air tem-

perature and wind speed. The various fitted spline sur-

faces 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

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

401401

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.

Station

name

Correlation of wind speed and

max air temperature

Predicted correlation at

withheld station

Absolute error of predicted

correlation

Error rank of predicted

correlation

Barnwell 0.70 0.700 0.000 17

Bodo 0.72 0.924 0.204 1

Brocket 0.76 0.778 0.018 14

Champion 0.70 0.744 0.044 6

Del Bonita 0.76 0.718 0.042 9

Foremost 0.66 0.602 0.058 5

Hussar 0.66 0.704 0.044 7

Lacombe 0.68 0.681 0.001 16

Lethbridge 0.74 0.734 0.007 15

Morrin 0.71 0.667 0.043 8

Olds 0.69 0.708 0.018 13

Onefour 0.60 0.636 0.036 10

Oyen 0.89 0.707 0.183 2

Schuler 0.71 0.771 0.067 3

Stavely 0.76 0.727 0.033 12

Stettler 0.68 0.713 0.033 11

Wrentham 0.64 0.699 0.057 4

Mean error: 0.052

Figure 5. Predicted strength of spatial occurrence of freeze-

thaw, 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.

Chinook event is a complex interaction between two or

more climate variables. Locations with high correlation

values between wind speed and maximum air tempera-

ture were expected to have an increased occurrence of

Chinook events and freeze-thaw cycles. Tab le 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 re-

corded, the absolute difference between observed corre-

lation value and the predicted values were calculated and

ranked. Barnwell, Lethbridge and Lacombe had the low-

est absolute difference between observed and predicted

correlation values, this can be explained by these sta-

tions’ relative close proximity to surrounding stations

with similar wind speed and maximum temperature cor-

relation values (Tabl e 5 ). The results of the leave-one-

out validation show that Oyen, Bodo and Schuler had the

largest errors (Ta b l e 5) and therefore were the most in-

fluential 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 sta-

tions 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-eco-

systems, additional monitoring stations likely need to be

placed near such locations where higher observed corre-

lation 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

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

402

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).

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 sen-

sitive 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 freez-

ing 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 tempera-

tures 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 (Tab le 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 maxi-

mum air temperature; this would explain the larger num-

ber 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 per-

sistence 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 mini-

mum air temperatures were the lowest of the three win-

ters. The 2008/2009 winter had a warmer mean soil tem-

perature, mean and maximum magnitude of soil freeze at

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).

Winters

Variables related to freeze-thaw 2006/2007 2007/2008 2008/2009

mean max air temperature (˚C) 1.302 0.533 –0.716

mean min air temperature (˚C) –9.923 –11.557 –12.172

mean soil temperature at 5 cm depth (˚C) –2.013 –2.894 –2.867

soil freezing days 114 112 127

Freeze-thaw cycles at 5 cm depth 4 1 5

estimated freeze-thaw cycles at 1 cm depth 18 20 12

average magnitude of soil freeze at 5 cm depth (˚C) –3.253 –4.469 –3.916

coldest soil temperature at 5 cm depth (˚C) –7.848 –10.325 –9.546

day of coldest soil temperature at 5 cm depth 12-Jan-07 29-Jan-08 26-Jan-09

first freeze 02-Nov-06 19-Nov-07 20-Nov-08

final thaw 10-Mar-07 10-Mar-08 31-Mar-09

length of winter freeze-thaw cycle 128 112 131

average lag time between freeze-thaw cycles 32 112 26

accumulated precipitation (mm) 52.415 36.943 69.118

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

403403

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 condi-

tion calibrated with a regional-specific context. All fitted soil temperature fit within the observed soil

temperature at 5 cm depth 95th percentiles.

5 cm than the previous winter (Tab le 6). More accumu-

lated 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 regional-

specific 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 con-

dition 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-spe-

cific climate data instead of generic boundary conditions

that require only one or two climate variable inputs.



13.247sin 0.0030.244

20.945sin3.496,

,Apr. 1Nov. 1

Apr. 1Nov. 1,

Nov. 1Apr. 1

iin

Tt t

























 













(3)

T = temperature [˚C]

t = time [s]

λ = decay rate

n = degree of B-spline

Pi = knot (control point)

bi = B-spline coefficient

4. CONCLUSIONS

Our findings illustrate the importance of regional

context in modeling freeze-thaw cycles. Our results in-

dicate that the interaction between air temperature and

wind drives variability in soil temperatures and if in-

cluded when specifying boundary conditions or cali-

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

404

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 consid-

ering 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 time-

periods.

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 none-

theless 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 sig-

nalling 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 ex-

ample 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 geostatis-

tical design, and implementation of agro-environmental

monitoring systems that combine in-situ and satellite/

remote-sensing derived estimates of near-surface soil

moisture.

5. ACKNOWLEDGEMENTS

The authors would to acknowledge R. Scott Erickson (Agriculture

and Agri-Food Canada) and Prof. Steven Liang (University of Calgary)

for helpful editorial comments that improved an earlier version of this

manuscript. This work was made possible by funding from the Sus-

tainable Agriculture Environmental Systems (SAGES) Growing For-

ward Program, Agriculture and Agri-Food Canada (Federal Govern-

ment of Canada).

REFERENCES

[1] Kreyling, J., Beierkuhnlein, C., Pritsch, K., Schloter, M.

and Jentsch, A. (2008) Recurrent soil freeze-thaw cycles

enhance grassland productivity. New Phytologist, 177,

938-945. doi:10.1111/j.1469-8137.2007.02309.x

[2] Chang, C. and Hao, X. (2001) Source of N2O emission

from a soil during freezing and thawing. Phyton, 41,

49-60.

[3] Matzner, E. and Borken, W. (2008) Do freeze-thaw

events enhance C and N losses from soils of different

ecosystems? A review. European Journal of Soil Science,

59, 274-284. doi:10.1111/j.1365-2389.2007.00992.x

[4] Bullock, M.S., Larney, F.J., Izaurralde, R.C.S. and Feng,

Y. (2001) Overwinter changes in wind erodibility of clay

loam soils in Southern Alberta. Soil Science Society of

America, 65, 423-430. doi:10.2136/sssaj2001.652423x

[5] Walker, V.K., Palmer, G.R. and Voordouw, G. (2006)

Freeze-thaw tolerance and clues to the winter survival of

a soil community. Applied and Environmental Microbiology,

73, 1784-1792. doi:10.1128/AEM.72.3.1784-1792.2006

[6] Henry, H.A.L. (2008) Climate change and soil freezing

dynamics: Historical trends and projected changes. Cli-

matic Change, 87, 421-434.

doi:10.1007/s10584-007-9322-8

[7] IPCC (2007) Climate change 2007: Synthesis report. In:

Pachauri, R.K and Reisinger, A., Eds., 4th Assessment

Report of the Intergovernmental Panel on Climate

Change IPCC, Geneva, 104.

[8] Eitzinger, J., Parton, W.J. and Hartman, M. (2000) Im-

provment and validation of a daily soil temperature sub-

model for freezing/thawing. Soil Science, 165, 525-534.

doi: 10.1007/s10584-007-9322-8

[9] Chen, B., Chen, J.M. and Ju, W. (2007) Remote sensing-

based ecosystem-atmosphere simulation scheme (EASS)—

Model formulation and test with multiple-year data.

A. J. Phillips et al. / Agricultural Sciences 2 (2011) 392-405

405405

Ecological Modelling, 209, 277-300.

doi: 10.1016/j.ecolmodel.2007.06.032

[10] Flerchinger, G.N. and Saxton, K.E. (1989) Simultaneous

heat and water model of a freezing snow-residue-soil

system. I. Theory and development. Transactions of the

American Society of Agricultural Engineers, 32, 565-

571.

[11] Flerchinger, G.N. (2000) The simultaneous heat and wa-

ter (SHAW) Model: Technical documentation. Northwest

Watershed Research Center, USDA Agricultural Research

Service, Boise.

http://afrsweb.usda.gov/SP2UserFiles/Place/53620000/S

hawDocumentation.pdf

[12] Šimůnek, J., Sejna, M., Saito, H., Sakai, M. and Genuch-

ten, M.T.V. (2009) The HYDRUS-1D software package

for simulating the one-dimensional movement of water,

heat, and multiple solutes in variably-saturated media.

PC Progress, Prague.

http://www.pc-progress.com/en/Default.aspx?h1d-downl

oads

[13] Nkemdirim, L.C. (1996) Canada’s chinook belt. Inter-

national Journal of Climatology, 16, 441-462.

doi:10.1002/(SICI)1097-0088(199604)16:4<441::AID-J

OC21>3.0.CO;2-T

[14] Nkemdirim, L.C. (1991) Chinooks and winter evapora-

tion. Theoretical and Applied Climatology, 43, 129-136.

doi: 10.1007/BF00867470

[15] Hurrell, J.W. (1996) Influence of variations in extra-

tropical wintertime teleconnections on northern hemis-

phere temperature. Geophysical Research Letters, 23,

665-668. doi: 10.1029/96GL00459

[16] Bonsal, B.R., Prowse, T.D., Duguay, C.R. and Lacroix,

M. P. (2006) Impacts of large-scale teleconnections on

freshwater-ice break/freeze-up dates over Canada. Jour-

nal of Hydrology, 330, 340-353.

doi:10.1016/j.jhydrol.2006.03.022

[17] Smith, N.V., Saatchi, S.S. and Randerson, J.T. (2004)

Trends in high northern latitude soil freeze thaw cycles

from 1988 to 2002. Journal of Geophysical Research,

109, 1-14. doi:10.1029/2003/2003JD004472

[18] Baker, D.G. and Ruschy, D.L. (1995). Calculated and

measured air and soil freeze thaw frequencies. Journal of

Applied Meteorology, 34, 2197-2205.

doi:10.1175/1520-0450(1995)034<2197:CAMAAS>2.0.

CO;2

[19] Henry, H.A.L. (2007) Soil freeze-thaw cycle experiments:

Trends, methodological weaknesses and suggested im-

provements. Soil Biology and Biochemistry, 39, 977-986.

doi: 10.1016/j.soilbio.2006.11.017

[20] Groffman, P.M., Driscoll, C.T., Fahey, T.J., Hardy, J.P.,

Fitzhugh, R.D. and Tierney, G.L. (2001) Colder soils in a

warmer world: A snow manipulation study in a northern

hardwood forest ecosystem. Biogeochemistry, 56, 135-

150. doi:10.1023/A:1013039830323