Articles | Volume 25, issue 8
Hydrol. Earth Syst. Sci., 25, 4357–4372, 2021
Hydrol. Earth Syst. Sci., 25, 4357–4372, 2021

Research article 10 Aug 2021

Research article | 10 Aug 2021

Daily soil temperature modeling improved by integrating observed snow cover and estimated soil moisture in the USA Great Plains

Daily soil temperature modeling improved by integrating observed snow cover and estimated soil moisture in the USA Great Plains
Haidong Zhao1, Gretchen F. Sassenrath1,2, Mary Beth Kirkham1, Nenghan Wan1, and Xiaomao Lin1 Haidong Zhao et al.
  • 1Department of Agronomy, Kansas State University, Manhattan, KS, USA
  • 2Southeast Research and Extension Center, Kansas State University, Parsons, KS, USA

Correspondence: Xiaomao Lin (


Soil temperature (Ts) plays a critical role in land–surface hydrological processes and agricultural ecosystems. However, soil temperature data are limited in both temporal and spatial scales due to the configuration of early weather station networks in the USA Great Plains. Here, we examined an empirical model (EM02) for predicting daily soil temperature (Ts) at the 10 cm depth across Nebraska, Kansas, Oklahoma, and parts of Texas that comprise the USA winter wheat belt. An improved empirical model (iEM02) was developed and calibrated using available historical climate data prior to 2015 from 87 weather stations. The calibrated models were then evaluated independently, using the latest 5-year observations from 2015 to 2019. Our results suggested that the iEM02 had, on average, an improved root mean square error (RMSE) of 0.6 C for 87 stations when compared to the original EM02 model. Specifically, after incorporating the changes in soil moisture and daily snow depth, the improved model was 50 % more accurate, as demonstrated by the decrease in RMSE. We conclude that, in the USA Great Plains, the iEM02 model can better estimate soil temperature at the surface soil layer where most hydrological and biological processes occur. Both seasonal and spatial improvements made in the improved model suggest that it can provide a daily soil temperature modeling tool that overcomes the deficiencies of soil temperature data used in assessments of climatic changes, hydrological modeling, and winter wheat production in the USA Great Plains.

1 Introduction

A reliable estimate of soil temperature (Ts) is useful for understanding agricultural ecological systems, hydrological processes, and land–atmosphere interactions (Lembrechts et al., 2020; Qi et al., 2016; Zhang et al., 2018) due to the fact that Ts governs physical, chemical, and biological processes of the soil and interactions between the atmosphere and land surface (Smith, 2000; Soong et al., 2020). In particular, Ts has been widely used for a better understanding of changes in soil moisture (Lakshmi et al., 2003), the ecosystem carbon balance (Goulden et al., 1998), and the nitrogen mineralization process (Persson and Wirén, 1995), although a larger prevalence of air temperature observations are available as a soil temperature proxy. From a practical perspective, Ts is critical for agricultural system models, such as the crop environmental resource synthesis (CERES) models, to assess the impacts of extreme climate on crop production and stress tolerance, thereby allowing producers to better prepare for proactive and reactive field management (Bergjord et al., 2008; Persson et al., 2017; Williams et al., 1989). Frequent extreme climate events, such as spring freezes and summer heat stress, can impact winter wheat (Triticum aestivum L.) growth and development, reducing grain yields by more than 7 % in the USA winter wheat belt (Tack et al., 2015; Paulsen and Heyne, 1983). These effects are also modulated through land–surface interaction processes (Hillel, 1998; Araghi et al., 2017).

To improve the accuracy of crop management modeling, a bare soil temperature (Ts) at the 10 cm depth, a standard soil temperature variable, has commonly been considered as being a more direct and useful variable than air temperature (Ta) measured at 1.5 or 2 m height in crop phenology (Onwuka and Mang, 2018), plant photosynthesis and soil respiration (Meyer et al., 2018; Wu and Jansson, 2013), plant nutrient uptake (Yan et al., 2012), and estimate of crop production (Araghi et al., 2017; Hillel, 1998). There are many Ts modeling techniques, mostly based on the land–surface interaction process (Qi et al., 2019; Yener et al., 2017). Most Ts models are rooted in theories of soil heat exchange and surface energy balance (Rankinen et al., 2004; Nobel and Geller, 1987; Chalhoub et al., 2017). The theory-based simulation for surface energy balance usually includes solar radiation (incoming and outgoing), infrared radiation (absorbed and reflected), turbulent flux energy (latent heat and sensible heat), and net ground heat flux through the ground surface into soil layers thermodynamically (Mihalakakou et al., 1997; Chalhoub et al., 2017). Obviously, the energy-balance-based model usually requires more detailed near-surface and soil variables, such as turbulent flux quantities (sensible heat flux and latent heat flux), to make the model reliable and accurate; however, determining quality turbulent flux quantities is not a trivial task (Dhungel et al., 2021; Kutikoff et al., 2021). In addition, seasonal variations in soil thermal conductivity and underestimates of actual evapotranspiration usually lead to overestimated surface soil temperatures (Bittelli et al., 2008). Therefore, simpler empirical models with fewer dynamic processes for Ts prediction have been explored (Zheng et al., 1993; Plauborg, 2002; Liang et al., 2014; Badache et al., 2016; Kang et al., 2000). However, these empirical models might result in relatively large estimated errors of over 2 C due to the lack of details about physical processes, such as uncertainties of the soil volumetric heat capacity and thermal conductivity (Badía et al., 2017). For example, the volumetric heat capacity was higher for a clay soil (1.48–3.54 MJ m−3C−1) than for a sand soil (1.09–3.04 MJ m−3C−1) when the soil moisture content was between 0 to 0.25 kg kg−1 (Abu-Hamdeh, 2003). Currently, the USA Department of Agriculture (USDA) provides a high-resolution Gridded Soil Survey Geographic (gSSURGO) database product (, last access: 23 July 2021) that includes static soil physical property data at 10 km resolution. The gSSURGO data facilitate Ts modeling, especially for better performance in large-scale Ts modeling due to its spatial variations in soil properties and soil moisture. These data sets have been widely used in the estimation of root zone soil water content (Miller et al., 2018) and subsurface hydrologic properties (Dirmeyer and Norton, 2018). The empirical model proposed by Plauborg (2002) performed better than energy-balance-based models when applied in the USA Great Plains for the last 5 years. Due to the lack of information about static soil properties on a large scale 1 or 2 decades ago, either over- or underestimates of Ts occurred, which, in turn, leads to large deviations in the assessment of crop stress and crop production (Gupta et al., 1990; Stone et al., 1999).

Recent studies have shown that estimated soil temperature usually deviates from observed soil temperature in the winter due to snow cover, frozen soil, and wide spatial and temporal heterogeneity in frozen soil properties (Nagare et al., 2012; Zhang et al., 2008; Rankinen et al., 2004). The impact of snow cover on soil temperature has been investigated (Rankinen et al., 2004) and is partially accounted for by incorporating correcting factors in land surface modeling and ecosystem models (Zhang et al., 2008) and soil and water assessment tools (SWAT; Qi et al., 2019). For both empirically and physically based soil temperature modules embedded in SWAT, the predictions of soil temperature in regions with thick snow cover seldom agree with field measurements in winter (Qi et al., 2019).

In the USA Great Plains, there has been increasing interest in improving hydrological process modeling of surface water and groundwater due to the Ogallala Aquifer's depletion in recent decades (Haacker et al., 2019). However, the automated weather station networks that observed soil temperature were not commissioned in this region until the late 1980s and early 1990s (Brock and Crawford, 1995). Not only were there few continuous observations for Ts earlier than the 1990s, these automated weather station networks also had limited stations in each state of the USA Great Plains. Such a lack of reliable soil temperature data both spatially and temporally makes the long-term assessment of water resources, crop phenology, and crop production modeling difficult.

The objectives of this study include the following: (1) to develop a robust Ts model using limited surface climate variables by integrating soil moisture and snow depth observations, (2) to demonstrate the error contributions in soil temperature modeling, and (3) to evaluate the performance of an improved model to predict Ts compared to current models. The data sets and methods are described in Sect. 2. Section 3 provides modeling results, and conclusions are presented in Sect. 4.

Figure 1Specific soil textures (a) and soil bulk density (b) at 87 weather stations in the USA winter wheat belt, including the states of Nebraska (NE), Kansas (KS), Oklahoma (OK), and part of Texas (TX) in the USA Great Plains.

2 Data sets and methods

2.1 Weather stations and data sets

The spatial domain of this study covers the winter wheat belt in the USA Great Plains, comprising the states of Nebraska (NE), Kansas (KS), Oklahoma (OK), and part of Texas (TX), where soil texture and bulk density vary (Fig. 1). In this study, three surface climate data sets were obtained from (1) the Automated Weather Data Network (AWDN) (, last access: 23 July 2021), commissioned in the 1980s for Nebraska and Kansas; (2) the Oklahoma Mesonet (OK Mesonet), a daily climate data source for Oklahoma, which started in the 1990s (, last access: 23 July 2021); and (3) the Soil and Climate Analysis Network (SCAN), which gives daily climate observations (, last access: 23 July 2021) that we selected for Texas due to limited quality data available in its automated weather station network. The number of selected stations was 26 in NE, 8 in KS, 44 in OK, and 9 in TX. The selection of these 87 stations was based on the completeness of climate data and data length (at least longer than periods of continuous 15 years). In addition to the weather station data sets, soil data sets providing soil attributes and characteristics were obtained from the standard USDA-NRCS Soil Survey Geographic (gSSURGO) database product (, last access: 23 July 2021) from which soil bulk density (ρb; grams per cubic centimeter), soil organic matter (fOM; percent), sand (fsa; percent), clay (fcl; percent), silt (fsl; percent) contents, soil porosity (Ø; percent), and soil surface albedo (α; –) were used for all weather stations. Note that all symbols and corresponding descriptions for variables used in this study are listed in the Table A1 (see the Appendix). The snow depth data were taken from the daily Global Historical Climatology Network (GHCN; Menne et al., 2009; Lin et al., 2017). Detailed data set sources and data variables used in each data set are shown in Table A2.

2.2 Soil temperature models

2.2.1 Empirical model

There are two common soil temperature models, i.e., empirical and process based. After examining both types of models for our study region, the current empirical model was selected because it was more accurate than the process-based model in this area. Plauborg (2002) developed a statistical soil temperature (Ts; degrees Celsius) model based on the current and previous 2 d air temperatures (Ta; degrees Celsius), annual and semi-annual cycles in the soil temperature fluctuations, and a daily soil temperature offset at a specific site, as shown in Eq. (1) (called EM02, thereafter) as follows:

(1) T s , j = γ + α 0 T a , j + α 1 T a , j - 1 + α 2 T a , j - 2 + β 1 sin ( ω j ) + δ 1 cos ( ω j ) + β 2 sin ( 2 ω j ) + δ 2 cos ( 2 ω j ) ,

where γ is an offset constant (degrees Celsius), and coefficients α0, α1, and α2 are dimensionless. The units of the coefficients β1, β2, δ1, and δ2 are in degrees Celsius. The j and ω denote day of the year and annual frequency (2π/365 d or 2π/366 d in leap years) in an annual soil temperature signal.

Figure 2Effects of soil moisture on (a) thermal conductivity (λ) and (b) soil thermal diffusivity (k), as obtained by Eqs. (7)–(11).


2.2.2 Improved empirical model

The improved model, based on the EM02, was developed through the following three steps: (1) prolonging the time window of Ta to include 1 extra prior day Ta; (2) constructing a new fictive environmental temperature (Tenv; degrees Celsius), defined as a function of air temperature and surface skin temperature (Tsfc; degrees Celsius; Williams et al., 1984) by utilizing Tenv to replace the original Ta; and, most importantly (3), incorporating site-specific daily soil thermal diffusivity and snow depth. This improved empirical model (iEM02) can be described by Eqs. (2)–(6) as follows:

(2) T s , j = γ + α 0 T env , j + α 1 T env , j - 1 + α 2 T env , j - 2 + α 3 T env , j - 3 + β 1 sin ( ω j ) + δ 1 cos ( ω j ) + β 2 sin ( 2 ω j ) + δ 2 cos ( 2 ω j ) × f D S , j × DR eff , j


where a fictive environmental temperature (Tenv) is assumed to be the weighted mean of air temperature (Ta) at 2 m and surface temperature (Tsfc). The β refers to the weighting coefficient, which defines the relative weight of the air temperature. This weighted fictive temperature will help weigh surface cooling and heating due to radiative and convective process (Dolschak et al., 2015). The Tsfc in Eq. (2) was estimated iteratively from the 3 d running average of daily air temperature (Ta), daily maximum temperature (Tmax; degrees Celsius), and daily solar radiation (Rs; megajoule per square meter per day, hereafter MJ m−2 d−1). The α denotes soil surface albedo (–) and initial Tsfc,j-1 was set as annual mean Ta in Eq. (3). The constant of 33.5 is an empirical constant (MJ m−2 d−1; Williams et al., 1984). The function of snow cover on the jth day is given as f(DS,j) and was introduced based on the work of Rankinen et al. (2004). The fS and DS are empirical soil heat damping parameters (meters) and snow depth (meters). The damping ratio of soil at the soil depth of h (h=0.1 m in this study) is DReff,j (Rosenberg et al., 1983). The weighting coefficient for the damping ratio (–) is k0. The p represents the period (365 or 366 d in leap years) in an annual cycle. The thermal diffusivity ks,j (square meters per second) is equivalent to thermal conductivity (λ; watts per meter per kelvin) divided by volumetric heat capacity (C; joule per cubic meter per kelvin) and reflects both the ability of soil to transfer heat and to change temperature when the heat is supplied or dissipated (Fig. 2). The estimate of thermal conductivity (λ) and volumetric heat capacity (C) can be described by Eqs. (7)–(11) (Lu et al., 2014) as follows:


where λdry is oven-dried soil thermal conductivity derived from a linear function of soil porosity (Ø; percent). Both b1 and b2 are the shape factors of the λ curve that are estimated by soil texture components. Soil water content is defined as θj on the jth day (cubic centimeters per cubic centimeter) and was calculated by the soil water balance model (Chalhoub et al., 2017). Briefly, the iEM02 operates on a daily time step as daily soil moisture is a function of soil moisture storage capacity (θ*; millimeters), 24 h precipitation (P; millimeters), and Penman–Monteith reference evapotranspiration (ET0; millimeters) and is estimated by Eqs. (12)–(15) as follows:


where θr and θs define residual and saturated volumetric soil water contents (cubic centimeters per cubic centimeter). θs is assumed to be equal to soil porosity while βd,j is a weighting coefficient for the difference between ET0 (Allen et al., 1998) and P on the jth day (–). The initial soil water content (θj−1) is assumed to be equal half of soil porosity.

Figure 3Percentage increments of soil temperature modeling improvement in iEM02, as determined by RMSE changes -100RMSEimproved-RMSEoriginalRMSEoriginal (a) after introducing an air temperature of Ta,j-3, (b) after substituting the air temperature Ta with a fictive environmental temperature (Tenv), and (c) after integrating the impacts of soil thermal diffusivity and snow cover. The color bar was coded by the improved percentage of iEM02 against the EM02 model.

Figure 4Spatial variations in the improved empirical model (iEM) coefficients are as follows: (a–d) for α0, α1, α2, and α3; and (e–h) for β1, δ1, β2, and δ2. Panel (i) shows the snow damping ratio (fs) and (j) the soil damping ratio coefficients (k0). The color bar defines the values of the model's coefficients.

Figure 5The 1:1 plots of absolute mean errors between the complete model and reduced model, where one independent variable term was removed in the improved empirical model (iEM02). Panels (a–d) compare results with vs. without α0 in Nebraska (NE), Kansas (KS), Oklahoma (OK), and Texas (TX), respectively; panels (e–h) compare results with vs. without β1 in NE, KS, OK, and TX, respectively; and panels (i–l) compare results with vs. without δ1 in NE, KS, OK, and TX, respectively. The ratio of the root mean square error (RMSE) shown includes both RMSE values calculated from complete and reduced models, respectively. The color bar indicates the number of observed data points.


Climate observation data prior to the year 2015 were selected to calibrate the iEM02 for each station. For NE, KS, and OK, daily soil temperature observations at each station had at least 10 years of daily time series for calibrations. Data sets from TX had at least 4 years available for calibrations. Climate variables used for calibration included air temperature, precipitation, snow depth, solar radiation daily observations, and the site's static soil property. The optimal parameter values for each weather station were estimated when a minimum root mean square error (RMSE) between estimated and observed soil temperature was achieved. These parameters, for all 87 stations, are listed in Table A3.

2.3 iEM02 evaluation

In the data sets selected, all 87 station observations were longer than 15 years, except for the stations located in Texas. The last 5-year observations (2015 to 2019) were used to independently conduct model validation for all 87 stations. The metrics used to evaluate model performance were RMSE and mean absolute error (MAE). Soil temperature modeling improvement was evaluated by relative RMSE changes -100RMSEimproved-RMSEoriginalRMSEoriginal and by intercomparison between the fully complete model and the reduced model.

Figure 6Spatial distribution of the mean absolute error (MAE) (a, c) and RMSE (b, d) for an empirical model (EM02; a, b) and improved modelEM02 (iEM02; c, d). The color bar defines values of MAE (in degrees Celsius) and RMSE (in degrees Celsius).

3 Results and discussion

3.1 Improved empirical model (iEM02)

The iEM02 was evaluated from 2015 to 2019 for 87 weather stations. Soil temperature modeling, using different soil textures, was improved in different ways in the iEM02 model (Fig. 3). The improvement in soil temperature modeling was indicated by relative RMSE changes that were different across sites. The weather stations located in NE and KS, as well as TX, showed less improvement after introducing the air temperature of Ta,j-3 compared to OK (Fig. 3a). The soil types in OK are more clay and silt compared to NE and KS (Fig. 1). However, the improvement, when using the fictive environmental temperature, was significant in northern areas of NE and KS (sandy soil) but not in the southern area of OK and part of TX (clay and silt soil; Fig. 3b). Overall, latitude-dominated air temperature should play a role in improving estimated soil temperature. Most of the 87 stations achieved a 15 % to 40 % improvement in simulated soil temperature by introducing air temperature Ta,j-3 and replacing Ta with Tenv. This improvement was in agreement with a previous study (Dolschak et al., 2015). By incorporating changes in soil moisture and daily snow depth, additional improvements in soil temperature simulation of up to 50 % could be achieved (Fig. 3c) compared to the original model (EM02). It should be noted that there were fewer stations available in KS and TX compared to NE and OK. Overall, integrating snow cover and soil moisture data in iEM02 improved the simulated soil temperature (Fig. 3). The daily soil temperature modeling could be further improved if high-resolution (e.g., 30 m and daily), satellite-based soil moisture and/or snow cover products become available, for example, products based on the Soil Moisture Active Passive (SMAP) or Sentinel satellites (Das et al., 2019).

Figure 7Seasonal comparison between the estimated and observed soil temperatures. (a–d) The empirical model (EM). (e–h) The improved empirical model (iEM02). RMSE was calculated as the root mean square error between estimated and observed soil temperature. N refers to the sample size, and the gray line represents the 1:1 line. The color bar describes the number of data points.


Figure 8Daily soil temperature comparison between the observed (gray line), original model (EM02; green line), and improved model (iEM02; blue line) in western (>100 W), central (between 97 and 100 W), and eastern (<97 W) Nebraska (a–c) and Kansas (d–f) during the winter wheat growing seasons from 2015 to 2019. RMSE is the root mean square error (degrees Celsius). The values in brackets refer to both the RMSE values of the original model (green numbers) and the improved model (blue numbers) during the periods of October–November, December–February, and March–June, respectively. Shaded areas highlight the winter season (December–February).


3.2 iEM02's parameters

The parameters described in iEM02 for each weather station are indicative of soil temperature sensitivities for each independent variable in Eq. (1), although, strictly speaking, they are not mathematical sensitivities (Fig. 4 and Table A2). For Tenv, the current day Tenv was the most weighted, as expected (Fig. 4a). The parameters of Tenv for the prior day 1 to day 3 were relatively weak in terms of absolute magnitudes due to autoregression properties in the soil temperature (Fig. 4b–d). Interestingly, in the iEM02 model, the prior day 2 was negatively associated with soil temperature (Fig. 4c), which cannot be interpreted by soil physical processes but rather in a more autoregressive sense in which the soil temperature signals are superimposed. The periodic property embedded in iEM02 was two low-frequency components (i.e., semi-annual and annual signals). Obviously, the annual signal strength indicated by β1 and δ1 was 1 order of magnitude stronger than the semi-annual signal strengths in soil temperature (Fig. 4e–h). The result also suggested that the strong β1 and δ1 spatial contexts of the northern region (e.g., in Nebraska and Kansas) were differently weighted from those in the southern region (e.g., in Oklahoma and Texas). For the snow damping factor, the snow cover had a larger impact on soil temperature in the northern region when compared to the southern region (Fig. 4i). However, the soil damping ratio factor was relatively evenly distributed (Fig. 4j).

Figure 9The same as Fig. 8 but for western, central, and eastern Oklahoma (a–c) and Texas (d–f).


RMSE performance is shown in Fig. 5, when the iEM02 was a complete model vs. the reduced iEM02 model, where one independent variable term from the complete model was removed. When removing any one independent variable, the modeled soil temperature RMSE increased from 110 % to 130 % (Fig. 5), indicating a 20 % rise in RMSE. Specifically, the iEM02 model performance decreased (i.e., RMSE increased from 0.1 to 0.4 C) when the α0 term was removed (Fig. 5a–d). Unlike α0, removing the β1 term was not as sensitive and gave an increase of 0.1–0.2 C RMSE, on average, for all states in the region (Fig. 5e–h). However, it is clear that the iEM02 model was the most sensitive to δ1. With the removal of δ1 from the complete iEM02 model, the RMSE increased 0.3–0.4 C for all four states (Fig. 5i–l). Due to the location dependency of the above coefficients, further spatial interpolation of the iEM02 model would be beneficial to predict soil temperature for irrigated agricultural areas without weather stations in the USA Great Plains and to improve water and crop management modeling.

3.3 Spatial and temporal modeling performance

A graphical summary of how closely the modeled soil temperature agreed with the observed soil temperature for each weather station is shown in Fig. 6. Daily Ts estimated in the iEM02 model outperformed that in the original EM02 model for all 87 weather stations. For example, both MAE and RMSE were decreased, on average, by 0.6 C when the iEM02 model was used to estimate Ts. Individually, the improved model showed a less than 1.6 C RMSE for any individual station, but 16 % of the stations had a larger than 2 C RMSE in the original EM02. In addition, we compared the performance of iEM02 against a recent energy balance model (Chalhoub et al., 2017). Our prediction of Ts was improved by 1.2 C RMSE compared to the energy balance model (not shown).

Spatial distributions of RMSE showed that the majority of weather stations had better performance in Oklahoma, with a mean RMSE of 1.9 and 1.1 C for EM02 and iEM02, respectively, whereas Nebraska had a RMSE of 2.1 and 1.3 C for EM02 and iEM02, respectively. The different modeling performance was associated with the soil heat transport process and how frequent snowfall could be observed in Nebraska and Oklahoma. Similar results were presented in a recent study by Huang et al. (2017). On the other hand, the high quality of weather data from the Oklahoma Mesonet is considered to be the gold standard for the statewide weather network (Lin et al., 2016), thus ensuring the quality of both model calibrations and observed soil temperature in Oklahoma.

Seasonal Ts indicated that iEM02 modeling was mostly improved in the spring season, from 2 to 1.3 C RMSE (Fig. 7a), but the original model (EM02) showed that the uncertainty was in good agreement with the performance achieved in Plauborg (2002). All other seasons were improved in similar ways, from 1.8 to 1.2 or 1.3 C RMSE. The improvement for all seasons could be attributed to introducing soil diffusivity, which changed with daily soil moisture and snow cover, and this affected the soil thermal conductivity (Rankinen et al., 2004; Zhang, 2005). Moreover, although modeling wintertime soil temperature improved from 1.8 to 1.3 C RMSE, which was the same as in the summer (Fig. 7), the soil temperature located in more frequent snow-covered states (e.g., Nebraska and Kansas), was better improved when Tenv and snow depth were introduced into the model. Our findings confirmed those reported by Rankinen et al. (2004) and Dutta et al. (2018).

Since precipitation gradients exist in the USA Great Plains from western to eastern regions (Evett et al., 2020), three subregions were classified for each state as western (100 W; towards the west), central (between 97 and 100 W), and eastern (97 W; towards the east). Figure 8 displays the time series of EM02 modeled, iEM02 modeled, and observed soil temperatures only covering winter wheat growing seasons (1 October to 30 June) for four growing seasons from 2015 to 2019 (validation periods) in Nebraska and Kansas. All subregions in Nebraska and Kansas showed improvement when using the iEM02 model (Fig. 8). Similarly, the iEM02 improved the RMSE during four growing seasons in Oklahoma and Texas (Fig. 9). The EM02 model had the best performance in Oklahoma, with a mean RMSE of 1.0 C, while the mean RMSE in Kansas was 1.4 C in EM02. Soil temperatures estimated by iEM02 had approximately a 0.3 to 1.4 C RMSE (Figs. 8 and 9). In addition, larger improvements by iEM02 were observed in most subregions during wintertime, which would be beneficial for modeling winter wheat yields and potential yields (Persson et al., 2017).

4 Conclusion

The primary intention of this work was to develop an improved soil temperature model for the USA Great Plains that can predict soil temperature by using common weather station variables as inputs. The improved empirical model (iEM02) integrated soil thermal diffusivity and snow cover factors, and these significantly improved the estimate of soil temperature for the 87 weather stations in the USA Great Plains that were studied. Specifically, after incorporating changes in soil moisture and daily snow depth, the improved model showed a near 50 % gain in performance in terms of RMSE decrease when compared to the original model. The value of RMSE across 87 stations was 0.6 C lower, on average, than the original model from 2015 to 2019. We concluded that the iEM02 model can better estimate soil temperature at the surface soil layer where most hydrological and biological processes occur. Both seasonal and spatial improvements made in the improved model demonstrated the robustness of the iEM02 model, suggesting that this improved model can provide a reliable simulation of soil temperature to use in modeling hydrological processes and crop production in the USA Great Plains.

Appendix A

Table A1Table of symbols and corresponding descriptions used in this paper.

Download Print Version | Download XLSX

Table A2List of data sets used in this study, including the data source (Networks), state names (Coverage states), and specific data variables (Variables). Data sources include the Gridded Soil Survey Geographic (gSSURGO), the Automated Weather Data Network (AWDN) – High Plains Regional Climate Center (HPRCC), the Oklahoma Mesonet (OK Mesonet), the Soil Climate Analysis Network (SCAN), and the daily Global History Climatology Network (dGHCN). Weather stations from four states were located in the USA Great Plains, including Nebraska (NE), Kansas (KS), Oklahoma (OK), and Texas (TX). Climate data report the daily maximum (Tmax; degrees Celsius) and minimum air temperature (Tmin; degrees Celsius) at 2 m height, relative humidity (RH; percent), rainfall (prcp; millimeters), solar radiation (Rs; MJ m−2 d−1), wind speed at 2 m (WS; meters per second), and snow depth (Ds; millimeters). Soil data consist of the daily bare soil temperature at 10 cm depth (Ts; degrees Celsius), albedo of soil surface (α; –), organic matter content (fOM; percent), bulk density (ρb; grams per cubic meter), porosity (Ø; percent), sand (fsa), silt (fsl), and clay (fcl) content (percent).

Download Print Version | Download XLSX

Table A3List of model parameters for each weather station in the USA Great Plains. The location consists of latitude (Lat) and longitude (Long). There are 12 parameters in the improved EM model, including parameters of air temperature (β; –); parameters for current day to the previous 3 d of Tenv, including α0 (–), α1 (–), α2 (–), α3 (–), and constant offset γ (degrees Celsius); annual and semi-annual waves of sine and cosine functions parameters are β1, β2, δ1, and δ2 (degrees Celsius); and parameters for the snow depth damping factor (fS; meters) and the soil damping factor (k0; –). The bold font indicates that estimated coefficients are not statistically significant at 95 % confidence intervals.

Download XLSX

Code availability

MATLAB code is available upon request.

Data availability

Data used in this study are available from the links given in Sect. 2 (AWDN,, last access: 23 July 2021, HPRCC, 2021; OK Mesonet,, last access: 23 July 2021, Mesonet, 2021; and SCAN,, last access: 23 July 2021, NRCS, 2021).

Author contributions

HZ and XL designed the experiments, conducted simulations, analyzed the data, and wrote the paper. NW helped with the data analysis, result interpretation, and discussion. MBK and GFS provided suggestions, held a discussion, and wrote and revised the paper.

Competing interests

The authors declare that they have no conflict of interest.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


This study was supported, in part, by the USA Department of Agriculture, National Institute of Food and Agriculture (grant nos. 2016-68007-25066 and 2016-68007-25066), and the Kansas Crop Improvement Association, the USA Department of Agriculture National Institute of Food and Agriculture (hatch project no. 1018005; contribution no. 20-252-J). We appreciated Gerard Kluitenberg and Jesse Tack at Kansas State University, for providing helpful suggestions to improve the quality of paper. We thank Dallas Staley, for her outstanding contribution in editing and finalizing the paper. Her work continues to be at the highest professional level.

Financial support

This research has been supported by the USA Department of Agriculture (USDA) (grant no. 58-3090-5-009), USDA (NIFA grant no. 2016-68007-25066), USDA (ARS grant no. 58-3070-6-023), and Hatch fund (grant no. 1018005).

Review statement

This paper was edited by Lixin Wang and reviewed by Ray Anderson and one anonymous referee.


Abu-Hamdeh, N. H.: Thermal Properties of Soils as affected by Density and Water Content, Biosyst. Eng., 86, 97–102,, 2003. 

Allen, R. G., Pereira, L. S., Raes, D., and Smith, M.: Crop evapotranspiration – Guidelines for computing crop water requirements-FAO Irrigation and drainage paper 56, Fao, Rome, available at: (last access: 23 July 2021), 1998. 

Araghi, A., Mousavi-Baygi, M., Adamowski, J., Martinez, C., and van der Ploeg, M.: Forecasting soil temperature based on surface air temperature using a wavelet artificial neural network, Meteorol. Appl., 24, 603–611,, 2017. 

Badache, M., Eslami-Nejad, P., Ouzzane, M., Aidoun, Z., and Lamarche, L.: A new modeling approach for improved ground temperature profile determination, Renew. Energy, 85, 436–444,, 2016. 

Badía, D., López-García, S., Martí, C., Ortíz-Perpiñá, O., Girona-García, A., and Casanova-Gascón, J.: Burn effects on soil properties associated to heat transfer under contrasting moisture content, Sci. Total Environ., 601, 1119–1128, 2017. 

Bergjord, A. K., Bonesmo, H., and Skjelvåg, A. O.: Modelling the course of frost tolerance in winter wheat, Eur. J. Agron., 28, 321–330,, 2008. 

Bittelli, M., Ventura, F., Campbell, G. S., Snyder, R. L., Gallegati, F., and Pisa, P. R.: Coupling of heat, water vapor, and liquid water fluxes to compute evaporation in bare soils, J. Hydrol., 362, 191–205,, 2008. 

Brock, F. V. and Crawford, K. C.: The Oklahoma Mesonet_A Technical Overview, J. Atmos. Ocean. Tech., 12, 5–19,<0005:TOMATO>2.0.CO;2, 1995. 

Chalhoub, M., Bernier, M., Coquet, Y., and Philippe, M.: A simple heat and moisture transfer model to predict ground temperature for shallow ground heat exchangers, Renew. Energy, 103, 295–307,, 2017. 

Das, N. N., Entekhabi, D., Dunbar, R. S., Chaubell, M. J., Colliander, A., Yueh, S., Jagdhuber, T., Chen, F., Crow, W., and O'Neill, P. E.: The SMAP and Copernicus Sentinel 1A/B microwave active-passive high resolution surface soil moisture product, Remote Sens. Environ., 233, 111380,, 2019. 

Dhungel, R., Aiken, R., Evett, S. R., Colaizzi, P. D., Marek, G., Moorhead, J. E., Baumhardt, R. L., Brauer, D., Kutikoff, S., and Lin, X.: Energy Imbalance and Evapotranspiration Hysteresis under an Advective Environment: Evidence from Lysimeter, Eddy Covariance, and Energy Balance Modelling, Geophys. Res. Lett., 48, e2020GL091203,, 2021. 

Dirmeyer, P. A. and Norton, H. E.: Indications of surface and sub-surface hydrologic properties from SMAP soil moisture retrievals, Hydrology, 53, 36,, 2018. 

Dolschak, K., Gartner, K., and Berger, T. W.: A new approach to predict soil temperature under vegetated surfaces, Model. Earth Syst. Environ., 1, 32,, 2015. 

Dutta, B., Grant, B. B., Congreves, K. A., Smith, W. N., Wagner-Riddle, C., VanderZaag, A. C., Tenuta, M., and Desjardins, R. L.: Characterising effects of management practices, snow cover, and soil texture on soil temperature: Model development in DNDC, Biosyst. Eng., 168, 54–72,, 2018. 

Evett, S. R., Colaizzi, P. D., Lamm, F. R., O'Shaughnessy, S. A., Heeren, D. M., Trout, T. J., Kranz, W. L., and Lin, X.: Past, present, and future of irrigation on the US Great Plains, T. ASABE, 63, 703–729, 2020. 

Goulden, M., Wofsy, S., Harden, J., Trumbore, S. E., Crill, P., Gower, S., Fries, T., Daube, B., Fan, S.-M., and Sutton, D.: Sensitivity of boreal forest carbon balance to soil thaw, Science, 279, 214–217, 1998. 

Gupta, S. C., Radke, J. K., Swan, J. B., and Moncrief, J. F.: Predicting soil temperature under a ridge-furrow system in the U.S. Corm Belt, Soil Till. Res., 18, 145–165, 1990. 

Haacker, E. M., Cotterman, K. A., Smidt, S. J., Kendall, A. D., and Hyndman, D. W.: Effects of management areas, drought, and commodity prices on groundwater decline patterns across the High Plains Aquifer, Agr. Water Manage., 218, 259–273, 2019. 

Hillel, D.: Environmental soil physics: Fundamentals, applications, and environmental considerations, Academic Press, San Diego, CA, USA, 1998. 

HPRCC: AWDN, available at:, last access: 23 July 2021. 

Huang, Y., Jiang, J., Ma, S., Ricciuto, D., Hanson, P. J., and Luo, Y.: Soil thermal dynamics, snow cover, and frozen depth under five temperature treatments in an ombrotrophic bog: Constrained forecast with data assimilation, J. Geophys. Res.-Biogeo., 122, 2046–2063,, 2017. 

Kang, S., Kim, S., Oh, S., and Lee, D.: Predicting spatial and temporal patterns of soil temperature based on topography, surface cover and air temperature, Forest Ecol. Manage., 136, 173–184, 2000. 

Kutikoff, S., Lin, X., Evett, S. R., Gowda, P., Brauer, D., Moorhead, J., Marek, G., Colaizzi, P., Aiken, R., Xu, L., and Owensby, C.: Water vapor density and turbulent fluxes from three generations of infrared gas analyzers, Atmos. Meas. Tech., 14, 1253–1266,, 2021. 

Lakshmi, V., Jackson, T. J., and Zehrfuhs, D.: Soil moisture–temperature relationships: results from two field experiments, Hydrol. Process., 17, 3041–3057, 2003. 

Lembrechts, J. J., Aalto, J., Ashcroft, M. B., De Frenne, P., Kopecky, M., Lenoir, J., Luoto, M., Maclean, I. M. D., Roupsard, O., Fuentes-Lillo, E., Garcia, R. A., Pellissier, L., Pitteloud, C., Alatalo, J. M., Smith, S. W., Bjork, R. G., Muffler, L., Ratier Backes, A., Cesarz, S., Gottschall, F., Okello, J., Urban, J., Plichta, R., Svatek, M., Phartyal, S. S., Wipf, S., Eisenhauer, N., Puscas, M., Turtureanu, P. D., Varlagin, A., Dimarco, R. D., Jump, A. S., Randall, K., Dorrepaal, E., Larson, K., Walz, J., Vitale, L., Svoboda, M., Finger Higgens, R., Halbritter, A. H., Curasi, S. R., Klupar, I., Koontz, A., Pearse, W. D., Simpson, E., Stemkovski, M., Jessen Graae, B., Vedel Sorensen, M., Hoye, T. T., Fernandez Calzado, M. R., Lorite, J., Carbognani, M., Tomaselli, M., Forte, T. G. W., Petraglia, A., Haesen, S., Somers, B., Van Meerbeek, K., Bjorkman, M. P., Hylander, K., Merinero, S., Gharun, M., Buchmann, N., Dolezal, J., Matula, R., Thomas, A. D., Bailey, J. J., Ghosn, D., Kazakis, G., de Pablo, M. A., Kemppinen, J., Niittynen, P., Rew, L., Seipel, T., Larson, C., Speed, J. D. M., Ardo, J., Cannone, N., Guglielmin, M., Malfasi, F., Bader, M. Y., Canessa, R., Stanisci, A., Kreyling, J., Schmeddes, J., Teuber, L., Aschero, V., Ciliak, M., Malis, F., De Smedt, P., Govaert, S., Meeussen, C., Vangansbeke, P., Gigauri, K., Lamprecht, A., Pauli, H., Steinbauer, K., Winkler, M., Ueyama, M., Nunez, M. A., Ursu, T. M., Haider, S., Wedegartner, R. E. M., Smiljanic, M., Trouillier, M., Wilmking, M., Altman, J., Bruna, J., Hederova, L., Macek, M., Man, M., Wild, J., Vittoz, P., Partel, M., Barancok, P., Kanka, R., Kollar, J., Palaj, A., Barros, A., Mazzolari, A. C., Bauters, M., Boeckx, P., Benito Alonso, J. L., Zong, S., Di Cecco, V., Sitkova, Z., Tielborger, K., van den Brink, L., Weigel, R., Homeier, J., Dahlberg, C. J., Medinets, S., Medinets, V., De Boeck, H. J., Portillo-Estrada, M., Verryckt, L. T., Milbau, A., Daskalova, G. N., Thomas, H. J. D., Myers-Smith, I. H., Blonder, B., Stephan, J. G., Descombes, P., Zellweger, F., Frei, E. R., Heinesch, B., Andrews, C., Dick, J., Siebicke, L., Rocha, A., Senior, R. A., Rixen, C., Jimenez, J. J., Boike, J., Pauchard, A., Scholten, T., Scheffers, B., Klinges, D., Basham, E. W., Zhang, J., Zhang, Z., Geron, C., Fazlioglu, F., Candan, O., Sallo Bravo, J., Hrbacek, F., Laska, K., Cremonese, E., Haase, P., Moyano, F. E., Rossi, C., and Nijs, I.: SoilTemp: A global database of near-surface temperature, Global Change Biol., 26, 6616–6629,, 2020. 

Liang, L. L., Riveros-Iregui, D. A., Emanuel, R. E., and McGlynn, B. L.: A simple framework to estimate distributed soil temperature from discrete air temperature measurements in data-scarce regions, J. Geophys. Res.-Atmos., 119, 407–417,, 2014. 

Lin, X., Pielke Sr, R. A., Mahmood, R., Fiebrich, C. A., and Aiken, R.: Observational evidence of temperature trends at two levels in the surface layer, Atmos. Chem. Phys., 16, 827–841,, 2016. 

Lin, X., Harrington, J., Ciampitti, I., Gowda, P., Brown, D., and Kisekka, I.: Kansas trends and changes in temperature, precipitation, drought, and frost-free days from the 1890s to 2015, J. Contemp. Water Res. Educ., 162, 18–30, 2017. 

Lu, Y., Lu, S., Horton, R., and Ren, T.: An Empirical Model for Estimating Soil Thermal Conductivity from Texture, Water Content, and Bulk Density, Soil Sci. Soc. Am. J., 78, 1859–1868,, 2014. 

Menne, M. J., Williams Jr., C. N., and Vose, R. S.: The US Historical Climatology Network monthly temperature data, version 2, B. Am. Meteorol. Soc., 90, 993–1008, 2009. 

Mesonet: OK Mesonet, available at:, last access: 23 July 2021. 

Meyer, N., Welp, G., and Amelung, W.: The temperature sensitivity (Q10) of soil respiration: controlling factors and spatial prediction at regional scale based on environmental soil classes, Global Biogeochem. Cy., 32, 306–323, 2018. 

Mihalakakou, G., Santamouris, M., Lewis, J., and Asimakopoulos, D.: On the application of the energy balance equation to predict ground temperature profiles, Solar Energy, 60, 181–190, 1997. 

Miller, K., Luck, J., Heeren, D. M., Lo, T., Martin, D., and Barker, J.: A geospatial variable rate irrigation control scenario evaluation methodology based on mining root zone available water capacity, Precis. Agricult., 19, 666–683, 2018. 

Nagare, R. M., Schincariol, R. A., Quinton, W. L., and Hayashi, M.: Effects of freezing on soil temperature, freezing front propagation and moisture redistribution in peat: laboratory investigations, Hydrol. Earth Syst. Sci., 16, 501–515,, 2012. 

Nobel, P. S. and Geller, G. N.: Temperature modelling of wet and dry desert soils, J. Ecol., 75, 247–258, 1987. 

NRCS: Soil Climate Analysis Network (SCAN) Data & Products, available at:, last access: 23 July 2021. 

Onwuka, B. and Mang, B.: Effects of soil temperature on some soil properties and plant growth, Adv. Plants Agric. Res., 8, 34–37, 2018. 

Paulsen, G. M. and Heyne, E. G.: Grain production of winter wheat after spring freeze injury, Agron. J., 75, 705–707,, 1983. 

Persson, T. and Wirén, A.: Nitrogen mineralization and potential nitrification at different depths in acid forest soils, in: Nutrient uptake and cycling in forest ecosystems, Springer, 55–65, 1995. 

Persson, T., Bergjord Olsen, A. K., Nkurunziza, L., Sindhöj, E., and Eckersten, H.: Estimation of Crown Temperature of Winter Wheat and the Effect on Simulation of Frost Tolerance, J. Agron. Crop Sci., 203, 161–176,, 2017. 

Plauborg, F.: Simple model for 10 cm soil temperature in different soilswith short grass, Eur. J. Agron., 17, 173–179, 2002. 

Qi, J., Li, S., Li, Q., Xing, Z., Bourque, C. P.-A., and Meng, F.-R.: A new soil-temperature module for SWAT application in regions with seasonal snow cover, J. Hydrol., 538, 863–877, 2016. 

Qi, J., Zhang, X., and Cosh, M. H.: Modeling soil temperature in a temperate region: A comparison between empirical and physically based methods in SWAT, Ecol. Eng., 129, 134–143, 2019. 

Rankinen, K., Karvonen, T., and Butterfield, D.: A simple model for predicting soil temperature in snow-covered and seasonally frozen soil: model description and testing, Hydrol. Earth Syst. Sci., 8, 706–716,, 2004. 

Rosenberg, N. J., Blad, B. L., and Verma, S. B.: Microclimate: the biological environment, John Wiley & Sons, New York, NY, USA, 1983. 

Smith, K. A.: Soil and environmental analysis: physical methods, revised, and expanded, Marcel Dekker, New York, 2000. 

Soong, J. L., Phillips, C. L., Ledna, C., Koven, C. D., and Torn, M. S.: CMIP5 models predict rapid and deep soil warming over the 21st century, J. Geophys. Res.-Biogeo., 125, e2019JG005266,, 2020.  

Stone, P., Sorensen, I., and Jamieson, P.: Effect of soil temperature on phenology, canopy development, biomass and yield of maize in a cool-temperate climate, Field Crops Res., 63, 169–178, 1999. 

Tack, J., Barkley, A., and Nalley, L. L.: Effect of warming temperatures on US wheat yields, P. Natl. Acad. Sci. USA, 112, 6931–6936,, 2015. 

Williams, J., Jones, C., and Dyke, P. T.: A modeling approach to determining the relationship between erosion and soil productivity, T. ASAE, 27, 129–144, 1984. 

Williams, J. R., Jones, C. A., Kiniry, J. R., and Spanel, D. A.: The EPIC Crop Growth Model, T. Am. Soc. Agricul. Eng., 32, 497–511, 1989. 

Wu, S. H. and Jansson, P.-E.: Modelling soil temperature and moisture and corresponding seasonality of photosynthesis and transpiration in a boreal spruce ecosystem, Hydrol. Earth Syst. Sci., 17, 735–749,, 2013. 

Yan, Q., Duan, Z., Mao, J., Li, X., and Dong, F.: Effects of root-zone temperature and N, P, and K supplies on nutrient uptake of cucumber (Cucumis sativus L.) seedlings in hydroponics, Soil Sci. Plant Nutr., 58, 707–717,, 2012. 

Yener, D., Ozgener, O., and Ozgener, L.: Prediction of soil temperatures for shallow geothermal applications in Turkey, Renew. Sustain. Energ. Rev., 70, 71–77,, 2017. 

Zhang, T.: Influence of the seasonal snow cover on the ground thermal regime: An overview, Rev. Geophys., 43, RG4002,, 2005. 

Zhang, T., Shen, S., Cheng, C., Song, C., and Ye, S.: Long-Range Correlation Analysis of Soil Temperature and Moisture on A'rou Hillsides, Babao River Basin, J. Geophys. Res.-Atmos., 123, 12606–12620,, 2018. 

Zhang, Y., Wang, S., Barr, A. G., and Black, T.: Impact of snow cover on soil temperature and its simulation in a boreal aspen forest, Cold Reg. Sci. Technol., 52, 355–370, 2008. 

Zheng, D., Hunt Jr., E. R., and Running, S. W.: A daily soil temperature model based on air temperature and precipitation for continental applications, Clim. Res., 2, 183–191, 1993. 

Short summary
This study was done to develop an improved soil temperature model for the USA Great Plains by using common weather station variables as inputs. After incorporating knowledge of estimated soil moisture and observed daily snow depth, the improved model showed a near 50 % gain in performance compared to the original model. We conclude that our improved model can better estimate soil temperature at the surface soil layer where most hydrological and biological processes occur.