Abstract

In many practical applications snow depth is known, but snow water equivalent (SWE) is needed as well. Measuring SWE takes ∼20 times as long as measuring depth, which in part is why depth measurements outnumber SWE measurements worldwide. Here a method of estimating snow bulk density is presented and then used to convert snow depth to SWE. The method is grounded in the fact that depth varies over a range that is many times greater than that of bulk density. Consequently, estimates derived from measured depths and modeled densities generally fall close to measured values of SWE. Knowledge of snow climate classes is used to improve the accuracy of the estimation procedure. A statistical model based on a Bayesian analysis of a set of 25 688 depth–density–SWE data collected in the United States, Canada, and Switzerland takes snow depth, day of the year, and the climate class of snow at a selected location from which it produces a local bulk density estimate. When converted to SWE and tested against two continental-scale datasets, 90% of the computed SWE values fell within ±8 cm of the measured values, with most estimates falling much closer.

1. Introduction

As global temperatures rise, the world’s snow resources are predicted to change in significant ways (Hosaka et al. 2005; Christensen et al. 2007; Räisänen 2008; Deser et al. 2010). Long-term changes in global, regional, and local snow depth (hs), snow water equivalent (SWE), and extent will ultimately have major ramifications for ecosystem function, human utilization of snow resources, and the climate itself through feedback mechanisms like snow albedo (Barry 1996). Unfortunately, of the three snow metrics listed above, only extent [i.e., snow cover area (SCA)] is easily monitored using satellites. This monitoring, under way for several decades (Robinson 1993; Frei and Gong 2005), has shown that global SCA has been decreasing for the past 30 years (Lemke et al. 2007). SCA, however, is only an indirect measure of the world’s snow water resources (Brown 2000; Brown et al. 2003). To fully understand global snow water trends, the most fundamental metric to monitor is SWE, with depth a close second. For monitoring these metrics, two potential methods are available: passive microwave remote sensing and estimations based on direct measurements. Many issues make estimating snow mass from remote sensing problematic (König et al. 2001), leaving field and station measurements currently the primary means of inferring critical trends in snow resources.

Of the two fundamental measurements, depth is quicker and easier to measure in the field than SWE. No detailed estimates of the total number of depth and SWE measurements made worldwide is available, but what is available suggests that considerably more depths are collected than SWE measurements. For example, Environment Canada operates 1556 snow depth sites, but only 27 sites where SWE is measured. In the United States, over 700 snow pillows (Beaumont 1965; Johnson and Schaefer 2002), chiefly operated by the National Resource Conservation Service (NRCS), are used to make continuous SWE measurements, but conservative estimates (based on sales of sonic sounders) suggest that thousands of depth-monitoring stations are in operation. The U.S. Weather Service as well as the Swiss Service measure more depths than water equivalents.

The two measurements are related but not interchangeable. The decision to make one type of measurement versus the other (or a mix of both) is often determined as much by tradition as an objective assessment of accuracy, need, and available resources. To our knowledge, no formal cost–benefit analyses of measuring depth versus SWE have been conducted, but a rough estimate can be made of the relative time effort involved in making each type of measurement. During field campaigns we try to split our available time equally between the two. A census of records for the past decade indicates an approximate ratio of depth to SWE measurements of about 30 to 1. In a more objective comparison, our standard field protocol calls for 201 depth measurements and 10 SWE measurements at each station. These are done in the same amount of time giving a 20:1 ratio. Part of the reason for this disparity is that we use an efficient self-recording depth probe (U.S. Patent No. 5864059). In the U.S. alone approximately 6000 man hours are expended monthly during winter and spring in the collection of ∼20 000 SWE data, but based on the 20:1 ratio, the same work effort could produce nearly half a million depth measurements.

What if there was a reliable method of converting depth measurements into SWE? Such a method would allow for a significant expansion in the number of locations where SWE could be estimated without incurring much additional expense. This would improve our ability to assess worldwide snow water resources. Here we describe a method that uses historical station data and knowledge of snow climate classes (Sturm et al. 1995) to estimate the bulk density of a snowpack. The bulk density is then used to convert depth into SWE. Specifically, using a large (n = 25 688) training set of depth, density, and SWE measurements, we fit the data, including the climate class of each datum (Sturm et al. 1995), with a nonlinear analysis of covariance model (ANCOVA) using Bayesian methods. From the analysis an equation is derived that allows estimation of bulk density based on input of location, snow depth, and the day of the year. We examine the residuals to assess model fit and test the model against a large independent dataset to provide the reader with an understanding of model accuracy.

The method opens the possibility of converting thousands of depth observations that have been, or currently are being, collected in to SWE values. For example, it could be used to convert Environmental Technical Applications Center monthly observed global snow depth climatologies (Foster and Davy 1988) into SWE climatologies. It is also likely to become even more useful in the future if one promising new method of measuring snow depth becomes operational: airborne LiDAR (Hopkinson et al. 2004; Deems and Painter 2006; Deems et al. 2006; Minoru and Hiroshi 2006; Schaffhauser et al. 2008). Airborne LiDAR snow surveys are conducted much like gamma ray surveys (Peck et al. 1980) with a set of measurements made prior to the start of the snow season subtracted from a second set made after the snow has accumulated. Unlike gamma ray surveys, however, LiDAR surveys produce snow depth swath maps along the flight path that can contain literally millions of depth data. Accuracy depends on flight altitude, terrain, and GPS navigation, but is typically about ±15 cm in absolute terms, with considerably higher accuracy possible when frequent bench marks are included in the swath. Coupled with the depth-to-SWE conversion method described here, airborne LiDAR could greatly improve our ability to estimate local-to-hemispheric snow resources.

2. Background

At a given point, snow depth (hs) is related to SWE by the local bulk density (ρb):

 
formula

where depth is measured in centimeters, density in grams per centimeters cubed, ρw is the density of water (1 g cm−3), and SWE is measured in centimeters of water. It has long been observed (e.g., Dickinson and Whitely 1972; Steppuhn 1976) that the natural range of hs is many times greater than the range of ρb. For a northern Alaskan dataset (n = 5323), 95% of all bulk density values fell between 0.12 and 0.42 g cm−3 while 95% of all depths fell between 8 and 100 cm, a dynamic range 4 times greater than that of density. For the training dataset used in this paper the ratio is about 10:1. Based on these ratios, estimating the more conservative parameter (ρb) while directly measuring the more dynamic (and easier to measure) parameter (hs) is the most practical and potentially accurate method of estimating SWE.

The form of Eq. (1), along with the large dynamic range of hs, ensures that there will be a strong correlation between SWE and depth, but it is mute about the relationship between depth and bulk density. These two parameters are functionally related, albeit in a complex way. Steppuhn (1976) noted that Eq. (1) requires a covariance term between depth and density (C) when used to compute mean areal values (indicted by overbars):

 
formula

The covariance between depth and bulk density, however, tends to be negligible for a snow cover less than 80 cm deep (Pomeroy and Gray 1995), and at greater depths is only about 2.5% of SWE (data from this study). The error in neglecting C is therefore typically smaller than other potential error sources in Eqs. (1) or (2), so for simplicity, C is ignored.

Despite the fact that estimating ρb using Eq. (1) makes it possible to determine SWE from depth, there have been relatively few formal attempts to do so. A handbook for the NRCS field surveyors (Davis et al. 1970) provides qualitative rules for adjusting the local mean snow density upward or downward depending on wind, exposure, thaws, and the day of the year. Wilks and McKay (1996) developed a power-law version of Eq. (1) with the bulk density term (their “pseudodensity”) varying with climate. More recently, estimates of ρb have been produced using snow models that solve the surface energy balance and evolve a snowpack layer by layer (Brun et al. 1989; Liston and Elder 2006; Lehning et al. 2006), but these models require extensive meteorological input as well as the calculation of settlement of individual snow layers. They tend to be difficult to apply to large basins and regions because of their heavy data and computational demands.

3. Data and modeling

a. Data

The bulk density model was developed using a training set of 25 668 records of snow depth, density, and SWE from three countries and two continents (Table 1). These data span the climate classes of seasonal snow defined by Sturm et al. (1995): alpine, maritime, prairie, tundra, taiga, and ephemeral. The classes are defined by the general physical attributes of a snow cover (depth, density, type of snow layers, etc.), attributes that tend to remain relatively consistent within climate zones. The set combines research surveys done in Canada, Switzerland, and Alaska by the authors with data from monthly agency surveys done by the NRCS on snow courses in the western United States (Fig. 1a). Because few systematic measurements have been made in ephemeral snow it had to be excluded from the analysis. The data in the training set tend to be of high quality because a high percentage were taken during research campaigns.

Table 1.

Sources of SWE, depth, and bulk density data for the training dataset.

Sources of SWE, depth, and bulk density data for the training dataset.
Sources of SWE, depth, and bulk density data for the training dataset.
Fig. 1.

(a) The location of the sites in the training dataset (Western Hemisphere only) used to develop the model. Swiss data from a basin located at 47°04′N, 8°43′E are included in this set, and (b) the location of the sites in the test dataset (Meteorological Service of Canada 2000).

Fig. 1.

(a) The location of the sites in the training dataset (Western Hemisphere only) used to develop the model. Swiss data from a basin located at 47°04′N, 8°43′E are included in this set, and (b) the location of the sites in the test dataset (Meteorological Service of Canada 2000).

At NRCS stations, data are collected on the first of the month starting in November and continuing through May or June each year. The Canadian and Alaskan data were taken by the authors during over-snow traverses conducted in March and April when conditions were still well below freezing but the snowpack was at near-maximum depth (Table 1). The Swiss data were taken at Alptal, Switzerland (Staehli and Gustafsson 2006) on a weekly and biweekly schedule between November and April. All data used in the analysis are available online (see http://cdp.ucar.edu/cadis).

After developing the model, we located a large depth–density–SWE dataset from Canada (i.e., the test dataset; Fig. 1b; n = 226 009) against which the model was tested. This set, compiled in 2000 (Meteorological Service of Canada 2000) and updated through 2004 by R. Brown, comes from snow courses measured between 1935 and 2004 and is wholly different than the Canadian data used in model development. The practices and equipment used in collection of the large set varied widely and are poorly documented. Considering this, and the age of the data, they are probably not as accurate as the training set. Additionally, the set is heavily biased toward sites clustered along the U.S.–Canada border, many of which are classified as maritime sites despite having markedly lower depths than in the training set.

Both the training and test sets undoubtedly include erroneous depth and SWE data (cf. Carroll 1995). SWE measurements, in particular, are difficult to make accurately. Using a corer that is either plunged or twisted into the snow until it hits the ground, a snow core is removed from the snowpack and weighed in the corer using a spring balance calibrated to read out directly in SWE units, or bagged and weighed later using a digital balance. In either case, the measurement is gravimetric and returns the bulk density of the snow; though in one case the measurement is reported in SWE units. While the corer is touching the ground, the depth is read off graduated marks on the side of the corer with a nominal precision of ±0.2 cm, but an accuracy that is closer to ±1.0 cm.

At least four models of snow corer are in use in North America [the Federal or Mt. Rose, the Adirondack, the Eastern Snow Conference (ESC), and the Snow-Hydro] and all suffer from similar problems, the most common being a loss of snow. If the corer fails to retain a plug of soil, vegetation, and/or ice, snow will fall out the bottom of the corer when it is extracted from the snowpack (Turcan and Loijens 1975). The corer can also hit an ice layer, which will push snow out of the corer’s path, again resulting in a light sample. We have tried to quantify the magnitude of this undersampling by routinely comparing core-based density values to those determined from a density profile done in nearby snow pit. For measuring the latter, a 3-cm high steel box cutter with a volume of 100 cm3 is used with a digital balance. By integrating the individual layer densities measured in the pit, a bulk density can be computed and compared to the adjacent values obtained by coring. Our experience suggests that the pit measurements are more reliable and accurate, albeit more labor-intensive, than the core measurements. The results (spanning 8 yr, 5 field programs, and 3 types of corers: Table 2) suggest core-based bulk densities average 7.1% lower than layer-integrated values.

Table 2.

Difference (%) between bulk density determined by coring and obtained by integrating individual snow-layer densities. The average value is weighted by the number of snowpits (n).

Difference (%) between bulk density determined by coring and obtained by integrating individual snow-layer densities. The average value is weighted by the number of snowpits (n).
Difference (%) between bulk density determined by coring and obtained by integrating individual snow-layer densities. The average value is weighted by the number of snowpits (n).

Our results contradict prior studies (Work et al. 1965; Peterson and Brown 1975; Goodison 1978; Farnes et al. 1982) that found that coring results in a high SWE bias. This happens when excess snow is forced into a core tube, a situation that is said to happen more frequently with smaller diameter cutters (like the Federal Sampler with a cross-sectional area of 11 vs 30 cm2 for the ESC and Snow-Hydro samplers), deeper, denser snow (Peterson and Brown 1975), and poor cutter design (bad taper; fewer and/or dull teeth). The type of snow (nature and amount of depth hoar, number and thickness of ice layers) also matters. SWE errors also arise when a spring balance is used for weighing the cores (Bray 1973). Given our findings, and those of the previous studies, we conclude that within the training set (i) the SWE errors are small and (ii) they are a function of the equipment used, the care taken in collecting the measurement, and the nature of the snowpack (none of which is routinely recorded). They are effectively random in nature, and therefore cannot be corrected easily, so no correction is applied.

b. Modeling

While it would be possible to develop a model relating SWE directly to depth, it is more appropriate to model bulk density (as detailed above) and use that relationship to convert measured or assumed depths to SWE. Bulk density is a complex function of snow depth (hs), snow temperature (θ), snow deposition history (t), and the initial density of the individual snow layers (ρ0):

 
formula

thus Eq. (1) can be written as

 
formula

Without knowing explicitly the form of Eq. (3), it is clear that SWE is a complex and nonlinear function of hs.

To model the bulk density we used Bayesian statistical methods to develop a nonlinear ANCOVA model. Bayesian methods offer versatility when modeling complex systems (Congdon 2003; Gelman et al. 2004; see also http://www.bayesian.org/ and the journal Bayesian Analysis). When ρb is plotted as a function of depth, it exhibits nonlinearity in the mean response, a skewed error structure, and a nonconstant error variance (heteroscedasticity). These complexities can be handled using Bayesian methods, but are more problematic when using nonlinear least squares. The Bayesian methods also produce useful posterior probability distributions that define the uncertainties associated with the model predictions, as well as accounting for temporal and spatial dependencies. Three easily obtained predictor variables: snow depth, DOY (day of the year), and snow class, were used as input to the model. The snow class, as discussed below, was intended to capture environmental variables like temperature and initial density, which are well known for impacting the bulk density.

Many different functional forms of Eq. (3) (i.e., linear models with and without quadratic and cubic functions of depth and DOY; nonlinear functions) and all combinations of predictor variables were tested. To choose which form was best, as well as whether all three predictor variables were needed, we used an objective method of evaluation called the deviance information criterion (DIC; Burnham and Anderson 2002; Spiegelhalter et al. 2002; Wheeler et al. 2010). DIC scores become lower (in our case increasingly negative) as model accuracy goes up, but if, and only if, the improvement in accuracy exceeds a penalty assessed for increasing model complexity. A decrease in DIC score of 5–10 points implies a significant improvement in model performance. The final model that was selected had a DIC score ∼1630 lower than its closest competitor, suggesting that it was significantly better.

The variance in the bulk density data (when plotted against snow depth, DOY, and subdivided into climate classes) was modeled using both beta and normal distributions, but in the end, a beta distribution was selected because it could accommodate the skew in the data. Incorporating a term to address the nonconstant variability (heteroscedasticity) present in the training dataset improved model performance (decreasing DIC by almost 2000 points) and was also incorporated into the final model.

The software package WinBUGS (Lunn et al. 2000) was used to run the various models and determine the optimal result. Multiple Markov chains (using a range of initial values) were run for more than 100 000 iterations to simulate posterior distributions, with the first 50 000 iterations eliminated to allow for burn-in. We verified convergence by examining run histories, autocorrelation functions for simulated values, the Gelman–Rubin convergence statistic (Gelman and Rubin 1992) as modified by Brooks and Gelman (1998), and by performing additional diagnostic tests using the CODA (Plummer et al. 2006) libraries in Program R (R Development Core Team 2008).

Two potential types of errors affect the modeling. The first relates to the use of snow climate classes as a proxy for the physical processes that lead to the densification of a snow cover. Snow bulk density is known to increase with time and depth as the weight of the overlying snow compacts underlying layers (Kojima 1966). It is also moderated by temperature, solar radiation, the nature of the initial snow deposit, infiltration of meltwater, and the rate, sequence, and rapidity with which new snow layers are added to the snow on the ground. Bilello (1969), McKay and Findlay (1971), and McKay and Gray (1981) first demonstrated that these differences in bulk snow densification could be captured empirically using a climate classification; Sturm and Holmgren (1998) confirmed this finding. Nonetheless, the approach is general. The alternative would be to explicitly model compaction processes, as has been done in several physically based snow models (cf. Anderson 1976; Koren et al. 1999; Liston et al. 2007; see Rutter et al. 2009 for an extensive list). The problem is that these physical models require high-quality daily or even hourly weather and snowfall data because they must track individual snow layer settlement through time. Simpler quasi-physical models have been proposed based on time since snow deposition (Martinec 1966; Elder et al. 1991; Sturm and Holmgren 1998), but even these models are too complex and too computationally intensive to be used in global and regional applications (cf. Anderton et al. 2004). So while our model is general and potentially less accurate than explicit locally applied physical models, it is far easier to apply globally or regionally.

The second potential error is misclassification of the snow cover. We used a 1-km by 1-km global snow class map (G. E. Liston and M. Sturm 2010, unpublished manuscript) for classifying both the training and test datasets. The map is a refinement of the 0.5° by 0.5° grid map that appears in Sturm et al. (1995). For the training set the classification was checked using aerial photos (vegetation can indicate class) and satellite images, and in many cases there was local knowledge of the snow. Classification of the test set was based strictly on location (latitude–longitude) and likely includes greater classification errors.

4. Results

The training dataset (Table 3 and Fig. 2) has a mean density (ρb) of 0.312 g cm−3 and exhibits a nearly symmetric distribution, while depth (mean = 108.4 cm) and SWE (mean = 38.8 cm) have asymmetric distributions. Consistent with the findings of Skaugen (1999) and Bocchiola and Rosso (2007), these distributions can be fit using normal (density) and lognormal (depth and SWE) functions, respectively. The similarity between the SWE and depth distribution curves is indicative of the fact that SWE is more closely linked to depth than it is to bulk density, as previously discussed. A regression line of SWE versus depth [forced through (0, 0), see Fig. 3a] has a slope of 0.394 and an r2 value of 0.95 (Table 3). The correlation value is high, but it should be given the nature of Eq. (1). The same equation suggests that the slope of the regression should be equal to the mean bulk density (0.312 g cm−3), but it is not. Instead the slope is 26% higher because the functional relationship between depth and SWE is nonlinear. The regression line is skewed upward by the greater snow densities associated with deeper snow. This problem illustrates one of the motivations for modeling ρb rather than estimating SWE directly from depth.

Table 3.

Statistics for the training and test datasets with model errors also listed.

Statistics for the training and test datasets with model errors also listed.
Statistics for the training and test datasets with model errors also listed.
Fig. 2.

Empirical PDFs for (a) SWE, (b) depth, and (c) bulk density for the training dataset (filled bars) and the test dataset (black lines).

Fig. 2.

Empirical PDFs for (a) SWE, (b) depth, and (c) bulk density for the training dataset (filled bars) and the test dataset (black lines).

Fig. 3.

(a) SWE vs depth, (b) bulk density vs depth, and (c) SWE vs bulk density for the training dataset (n = 25 688). In (a) and (c) a line has been fit to the data. In (b) an exponential curve [Eq. (4)] has been fit.

Fig. 3.

(a) SWE vs depth, (b) bulk density vs depth, and (c) SWE vs bulk density for the training dataset (n = 25 688). In (a) and (c) a line has been fit to the data. In (b) an exponential curve [Eq. (4)] has been fit.

A plot of bulk density (ρb) versus depth (hs) [Fig. 3b; see also Fig. 3 in Pomeroy and Gray (1995)] shows considerably more scatter. There is a distinct increase in the minimum value of ρb (the floor of the data cloud in Fig. 3b) with increasing depth, a result that could have been anticipated because with increasing depth there is greater compaction. There is almost no increase in the maximum bulk density (ρmax) with depth; this too is consistent with previous work that has found that the maximum density of seasonal snow (in the absence of water infiltration) approaches 0.6 g cm−3 (Paterson 1981). A linear function could be used to relate bulk density to depth, but a better fit is achieved using a nonlinear function asymptotic to the maximum seasonal density:

 
formula

where ρ0 is the initial density and k is a fitting parameter. Unlike the SWE versus depth regression, the bulk density fit should not be forced through (0, 0) because even a vanishingly thin snow layer will have a nonzero density.

Within the training dataset there are distinctly different density patterns for each class of snow (Fig. 4). Maritime snow, for example, has the largest range (in both depth and density) of all the classes, forming a data cloud that has an axis oriented primarily in the horizontal direction, reflecting the large range of snow depths this class exhibits (Sturm et al. 1995). Deep maritime snow, always relatively close to its melting point, compacts rapidly to produce high bulk densities by late winter. In contrast, the tundra snow data cloud has a vertical axis, indicating a large density range, but little range in depth. Tundra snow never gets very deep, but its two major components, depth hoar and wind slab, bracket the full range of density found in seasonal snow (Benson and Sturm 1993). Variations in the relative percentages of these components can alter the bulk density dramatically without much variation in depth.

Fig. 4.

Different density patterns by snow classes for the training dataset. Contours in the main figure are for SWE (cm) and satisfy Eq. (1). The patterns are similar for the test dataset.

Fig. 4.

Different density patterns by snow classes for the training dataset. Contours in the main figure are for SWE (cm) and satisfy Eq. (1). The patterns are similar for the test dataset.

Within each snow class shown in Fig. 4 there are further variations in depth and density controlled by landscape, vegetation, microtopography, and time (e.g., Adams and Roulet 1982; Elder et al. 1991; Lapen and Martz 1996). These within-class variations are well illustrated by the tundra snow class. In the tundra regions, locations exposed to higher winds (in Alaska the arctic coastal region) tend to exhibit higher bulk densities and lower snow depths than more protected places like the foothills just north of the Brooks Range (Fig. 5: see also Bilello 1957, 1984). Here, we ignore these local, landscape-controlled density differences, though they could be introduced into the methodology in manner similar to that of Anderton et al. (2004), Watson et al. (2006), and Jonas et al. (2009).

Fig. 5.

Data from the tundra snow class in Alaska (Kuparuk River basin, March 1995) showing variation in the bulk density by landscape (upland vs coastal).

Fig. 5.

Data from the tundra snow class in Alaska (Kuparuk River basin, March 1995) showing variation in the bulk density by landscape (upland vs coastal).

The scatter in Figs. 3a,b suggest that in estimating the bulk density of snow using a statistical model there is likely to be a fairly large uncertainty. To place the model results in context, it is useful to look at the range of bulk density values observers might encounter if they were to take several accurate measurements at a single site or within a very small basin. This natural range is the base against which model-based estimates need to be compared. Unfortunately, there has been relatively little systematic effort to quantify it. A sense of its magnitude, however, can be obtained by analyzing data collected at 124 sites (10 measurements per site) during 2 of the Alaskan traverses that are included in the training set. Each of these sites was flat and nearly featureless (cf. Liston and Sturm 2002; Sturm and Liston 2003) limiting local variations because of topography or vegetation. The bulk density varied between 0.25 to 0.33 g cm−3 in 2000, and between 0.23 to 0.28 g cm−3 in 2002 (Fig. 6), a range that when converted to SWE implies a 27% and 20% variation, respectively, about the mean. In a more comprehensive survey that included five different locations and several snow types, Jonas et al. (2009) concluded that the local natural variation in bulk density ranged from 7% to 23% of the mean bulk density value, corresponding with a variation in SWE of 13%–80%, an amount that translates into 2.5–6 cm of SWE.

Fig. 6.

Box and whisker plots for bulk density variability in Northwest Alaska. In 2000, 49 stations were measured. In 2002, 79 stations were measured, including the same 49 measured in 2000. At each station, 10 bulk density measurements were made from which the range was computed. The shaded boxes bracket 25% and 75% of the data. The whiskers are at 90% and 10%.

Fig. 6.

Box and whisker plots for bulk density variability in Northwest Alaska. In 2000, 49 stations were measured. In 2002, 79 stations were measured, including the same 49 measured in 2000. At each station, 10 bulk density measurements were made from which the range was computed. The shaded boxes bracket 25% and 75% of the data. The whiskers are at 90% and 10%.

Finally, in preparation for testing the model results, training and test datasets were compared using standard statistical measures (Table 3). The mean values of depth, density, and SWE were higher for the training set by 35%, 12%, and 42%, respectively. We note that these differences are of the same order as the local variances shown in Fig. 6. The test set was also more scattered than the training set, probably because it is considerably larger, older, and contains more problematic data points.

5. Modeling the bulk density

The simplest method of converting snow depth to SWE is to replace ρb in Eq. (1) with the mean density of the training set (Table 3; 0.312 g cm−3). A refinement of this approach would be to first classify the sites being modeled and then use the mean densities of the individual classes. The mean density of the whole training set differs from the whole test set by only 12%, a small error considering other uncertainties in calculations of global snow reserves. The differences in snow class densities range from 1% to 19% (in all cases the training set values are higher than the test set), still reasonably close for regional and global modeling.

A more accurate approach, however, is to directly include the effect of snow depth [hs in Eq. (5)], the effect of snow aging (older snow is generally deeper and denser than younger snow) through the use of DOY, and to indirectly include the effects of climate (temperature and wind) through the use of snow classes. Equation (5) becomes

 
formula

where k1 and k2 are densification parameters for depth and DOY, respectively; ρmax, ρ0, k1, and k2 vary with snow class (Table 4); and i indicates the ith observation. Because the winter season in the Northern Hemisphere spans two calendar years, DOY runs from −92 (1 October) to +181 (30 June), with no 0 value. In this system, 1 February is DOY = 32, while 15 November is DOY = −47.

Table 4.

Model parameters by snow class.

Model parameters by snow class.
Model parameters by snow class.

Before we settled on Eq. (6) other functional forms relating bulk density to depth and DOY were evaluated using the DIC criteria. Several linear models with and without quadratic and cubic functions of depth and DOY produced nonphysical results and were discarded. Incorporating nonconstant variability (i.e., the heteroscedasticity exhibited in Fig. 3b) improved model performance, decreasing DIC by almost 2000 points, but had virtually no effect on the fitted curve and resulting estimates of k1 and k2. The final model [Eq. (6)] included all three predictor variables (snow depth, DOY, and snow class) and had a DIC score that was 1630 lower than the closest nonlinear competitor (which included only snow class and depth as input; Fig. 7). The chosen model explains about 45% of the variability in the data. By way of comparison, the adjusted r2 for a linear model with depth as the only input explains less than 24% of the variability (Fig. 7).

Fig. 7.

DIC scores plotted against adjusted r2 values for various models tried during the development of Eq. (6) and Table 4. For nonlinear models the value plotted against the y axis was the ratio of the difference between the variation in density data and the variation in the residuals of the model fit to the density data divided by the variation in the density data. Linear models produced nonphysical results (negative density) even with low DIC scores (i.e., model 6) and were ultimately not used. The arrow between models 7 and 8 suggests the impact of switching from homo- to hetero-scedastic variance models. The arrow from models 4 to 5 suggests the impact of changing from single to multiple slope models. The overall decrease in DIC illustrates the impact of increasing the number of parameters utilized in the model from one to three.

Fig. 7.

DIC scores plotted against adjusted r2 values for various models tried during the development of Eq. (6) and Table 4. For nonlinear models the value plotted against the y axis was the ratio of the difference between the variation in density data and the variation in the residuals of the model fit to the density data divided by the variation in the density data. Linear models produced nonphysical results (negative density) even with low DIC scores (i.e., model 6) and were ultimately not used. The arrow between models 7 and 8 suggests the impact of switching from homo- to hetero-scedastic variance models. The arrow from models 4 to 5 suggests the impact of changing from single to multiple slope models. The overall decrease in DIC illustrates the impact of increasing the number of parameters utilized in the model from one to three.

To quantify the model fit, bulk density residuals (ρobsρmodel) for the training set were used to calculate the root mean squared error (RMSE) as a function of snow depth:

 
formula

where n is the number of observations falling in depth class h1 to h2, and obs and mod indicate measured and modeled values, respectively. As expected, the predictive power of the model for snow depths less than 30 cm is relatively poor (Fig. 8), but the error diminishes quite rapidly as snow depth increases, quickly dropping to levels well below those of natural local variability (Fig. 6). For example for snow deeper than 50 cm, modeled bulk density errors are less than 0.066 g cm−3.

Fig. 8.

Residual values [Eq. (7)] after applying the model to the training dataset.

Fig. 8.

Residual values [Eq. (7)] after applying the model to the training dataset.

Application of the model over an area, region, or the Northern Hemisphere is simple: each snow depth station is classified using the snow class map of Sturm et al. (1995) or the newer 1-by-1 km map developed by G. E. Liston and M. Sturm (2010, unpublished data). Equation (6) is applied to the station snow depth using the date of observation and the appropriate parameters from Table 4. The resulting value of ρb is multiplied by the depth to produce the SWE for the station.

6. Testing the model

For each site and date of observation in the test dataset Eq. (6) was applied using the appropriate values from Table 4 based on the site classification. The estimated bulk density was compared to the observed value and two errors were defined:

 
formula

and

 
formula

Probability distribution functions for these errors are shown in Figs. 9a,b. Ninety percent of all model-estimated bulk densities fell within ±0.13 g cm−3 of the measured density. For comparison, site-specific field surveys based on multiple SWE cores commonly find local variations in density of ±0.07 g cm−3 (Fig. 6; see also Goodison 1978 and Table 4 in Jonas et al. 2009). SWE values produced from the model densities fell between −7 and +9 cm of the observed SWE: locally observed variations in SWE are on the order of half this range.

Fig. 9.

PDFs for (a) the bulk density error and (b) the SWE error as defined in Eqs. (8a) and (8b). The test dataset has a low bias in comparison to the training dataset and the model results (see text) that are marked on distributions by vertical dashed lines.

Fig. 9.

PDFs for (a) the bulk density error and (b) the SWE error as defined in Eqs. (8a) and (8b). The test dataset has a low bias in comparison to the training dataset and the model results (see text) that are marked on distributions by vertical dashed lines.

The biases between the training and test datasets listed in Table 3 and shown in Figs. 2a–c reappear in Figs. 9a,b. The test data are on average 0.04 g cm−3 lower in density than the training data and 2.2 cm lower in SWE (shown by vertical dashed lines in Figs. 8a,b). There is no definitive way to tell which dataset is better. However, given the uncertain heritage of many of the datum in the test set, we favor the training set as being more universal, in part because it is more regional in scope, and in part because its heritage is better known and it is known to include data of higher quality.

The probability of error as a function of snow depth has been contoured in Figs. 10a–c, removing the previously noted bias in order to focus on how the error varies with depth. The bulk density error is nearly constant across the entire range of depth (Fig. 10a), while the probable error in SWE increases as depth increases (Fig. 10b). The relative error in SWE, on the other hand, decreases with increasing depth (Fig. 10c). This is because the depth increases faster than the error. We examined whether the errors showed any systematic variation with other factors besides depth (and SWE) but found none. The errors are also of approximately the same magnitude across all latitudes, longitudes, and dates in the test dataset.

Fig. 10.

Contours of probability density for (a) the error in bulk density as a function of snow depth, (b) the error in SWE as a function of depth, and (c) the relative error in SWE (error/SWE) as a function of depth.

Fig. 10.

Contours of probability density for (a) the error in bulk density as a function of snow depth, (b) the error in SWE as a function of depth, and (c) the relative error in SWE (error/SWE) as a function of depth.

In a second test, model results were compared against data collected from three midlatitude basins in the United States where considerable snow research has been done. These are: Tokopah basin in the Sierra’s of California (Maritime or Alpine snow: see Sickman et al. 2001), Green Lakes basin in the Front Range of Colorado (Alpine snow; data courtesy of M. Williams), and the Fraser and North Park basins measured during the Cold Lands Processes Experiment (CLPX) in Colorado in 2001 and 2002 (Alpine and Prairie snow). The comparison (Table 5) shows reasonable agreement for estimated SWE values for all areas with the exception of the very shallow prairie snow of North Park. With the exception of the latter, the SWE error as a percentage of the observed value ranged from −1.5% to 26%, and was in all cases less than the standard deviation of the observed values.

Table 5.

Comparison of model results with measurements from Colorado and California.

Comparison of model results with measurements from Colorado and California.
Comparison of model results with measurements from Colorado and California.

7. Discussion

The strength of the model developed here is its simplicity and ease of use. It seems ideally suited for use in sparsely populated snow-covered areas like the Arctic and mountainous regions where access can be difficult. It is also useful where the computational cost of physically based modeling to estimate bulk density may not be warranted, or where meteorological station density may be too low to support such modeling. Given its simplicity, the performance of the model is surprisingly good. The reason the model works is because bulk density is a conservative value. Fully 50% of all bulk density measurements for seasonal snow lie in a narrow range between 0.244 to 0.375 g cm−3 (Fig. 2c); 80% lie between 0.194 and 0.435. Reasonable estimates of SWE can be derived even using a fixed bulk density value (either 0.312 g cm−3 for all snow, or the values in Table 3 for each snow class).

This approach of deriving SWE from snow depth is not new: NRCS instructions from the 1970s for converting aerial marker snow depths into SWE (Davis et al. 1970) incorporate several of the basic ideas presented here, though in a rule-based form. What is new is our use of snow climate classes to refine bulk density estimates, though the work of McKay and Findlay (1971), Wilks and McKay (1996), and Sturm and Holmgren (1998) presaged this part of the approach. The addition of the climate classes broadens the area over which estimates can be derived and improves the estimates.

The methodology applied here (using historical data to develop statistically based methods for estimating bulk density) has already been applied to a single snow class (Alpine snow: see Jonas et al. 2009) where snow depth, altitude, and DOY were used as input to a within-class statistical model. The same principle, that bulk density is conservative, served as the basis for this local study and ensured that the approach worked.

The ultimate value of the method is that it could potentially improve local, regional, and global estimates of snow resources at a time when budgets for operating traditional measurement monitoring networks have become more difficult to obtain and sustain. Using the model increases efficiency and cost effectiveness in snow monitoring by leveraging depth measurements (which can be taken 20 times as fast as SWE) into SWE. Dickinson and Whitely (1972), Steppuhn (1976), and Rovansek et al. (1993) have all proposed schemes that in some way do this, though here we have made the application of the method easier. Steppuhn (1976), writing more than 30 years ago, stated that:

“…statistically valid water equivalents [could] be derived from sampling where depth measurements outnumber[ed] those for density four-fold or more.”

We could not agree more, and put our model forward as a potential way to do this.

8. Conclusions

Here we present a method of converting snow depths to snow water equivalents (SWE). Using as input data snow depth, day of year, and station location (from which can be inferred the climate class of snow following Sturm et al. 1995), estimates of the bulk density of the snow are derived using a statistical model. The model was trained on historical depth–density–SWE data from the United States, Canada, and Switzerland, then tested on an extensive set of independent data from Canada and a second set from the western United States. The method works well because there is a satisfactory statistical relationship between bulk density and snow depth, time of year, and climate class of snow. It also works because density is a conservative variable that in nature is constrained between fairly narrow limits, while snow depth varies over a much greater range. Consequently, estimates of bulk density derived from historical data generally fall close to measured values. In a test against extensive Canadian data, 90% of the model-derived SWE estimates fell within ±8 cm of measured values, with relative errors in SWE not much higher than the range of SWE that would be encountered at a single site due to local real variability. The model is easy to apply and requires simple input data, making it useful in local to global snow water resource applications.

Recent efforts using remote sensing to monitor snow water resources have met with limited success. For example, direct satellite remote sensing methods for determining SWE are still not reliable enough (Cline et al. 1998; Mote et al. 2003) for tracking local, regional, or global snow mass or volume. As an alternative, our ability to measure snow extent has led to hybrid procedures where SCA has been used with a melt model to compute water volume (Cline et al. 1998; Molotch and Margulis 2008; Molotch 2009). In light of these limitations in our ability to measure and monitor local to global snow water resources, there remains a critical need for more station measurements of SWE. The method described here allows reasonable estimates of SWE to be made from simple observations of snow depth.

Acknowledgments

We express our thanks to several generations of NRCS, Canadian, and Swiss snow researchers, unknown to us by name, who painstaking collected the data used here. In our research programs, we thank Jon Holmgren, Eric Pyne, April Cheuvront, Wylie Bogren, Chris Polashenski, Paul Heflinger, Arvid Silis, Peter Toose, Andrew Rees, and others for the depth and SWE measurements they made in the field. James Sickman and John Melack provided snow data from Tokopah Basin. Mark Williams provided data from Niwot Ridge in Colorado. CLPX data were provided by Don Cline and Kelly Elder. Noah Molotch suggested the inclusion of Table 5 as a second test of model results. Ron Barry helped develop the Bayesian model. Ross Brown generously provided the test dataset. This work was supported by the U.S. Army Cold Regions Research and Engineering Laboratory and by the National Science Foundation, Office of Polar Program, Grants 0629279, 0632398, and 0632131.

REFERENCES

REFERENCES
Adams
,
W. P.
, and
N. T.
Roulet
,
1982
:
Areal differentiation of land and lake snowcover in a small sub-Arctic drainage basin.
Nord. Hydrol.
,
13
,
139
156
.
Anderson
,
E. A.
,
1976
:
A point energy and mass balance model of a snow cover.
NOAA Tech. Rep. NWS-19, 150 pp
.
Anderton
,
S. P.
,
S. M.
White
, and
B.
Alvera
,
2004
:
Evaluation of spatial variability in snow water equivalent for a high mountain catchment.
Hydrol. Processes
,
18
,
435
453
.
Barry
,
R. G.
,
1996
:
The parameterization of surface albedo for sea ice and its snow cover.
Prog. Phys. Geogr.
,
20
,
63
79
.
Beaumont
,
R. T.
,
1965
:
Mt. Hood pressure snow gauge.
J. Appl. Meteor.
,
4
,
626
631
.
Benson
,
C. S.
, and
M.
Sturm
,
1993
:
Structure and wind transport of seasonal snow on the Arctic slope of Alaska.
Ann. Glaciol.
,
18
,
261
267
.
Bilello
,
M. A.
,
1957
:
A survey of Arctic snow-cover properties as related to climatic conditions.
CRREL Research Rep. RR 39, 9 pp
.
Bilello
,
M. A.
,
1969
:
Relationship between climate and regional variations in snow-cover density in North America.
CRREL Research Rep. RR 267, 11 pp
.
Bilello
,
M. A.
,
1984
:
Regional and seasonal variations in the snow-cover density in the U.S.S.R.
CRREL Rep. RR 84-22, 70 pp
.
Bocchiola
,
D.
, and
R.
Rosso
,
2007
:
The distribution of daily snow water equivalent in the central Italian Alps.
Adv. Water Resour.
,
30
,
135
147
.
Bray
,
D. I.
,
1973
:
A report on the variability of snow water equivalent measurements at a site.
Proc. 30th Eastern Snow Conf., Amherst, MA, Eastern Snow Conference, 20–31
.
Brooks
,
S. P.
, and
A.
Gelman
,
1998
:
Alternative methods for monitoring convergence of iterative simulations.
J. Comput. Graph. Stat.
,
7
,
434
455
.
Brown
,
R. D.
,
2000
:
Northern Hemisphere snow cover variability and change, 1915–97.
J. Climate
,
13
,
2339
2355
.
Brown
,
R. D.
,
A. E.
Walker
, and
B.
Goodison
,
2000
:
Seasonal snow cover monitoring in Canada: An assessment of Canadian contributions for global climate monitoring.
Proc. 57th Eastern Snow Conf., Syracuse, NY, Eastern Snow Conference, 1–11
.
Brown
,
R. D.
,
B.
Brasnett
, and
D.
Robinson
,
2003
:
Gridded North American monthly snow depth and snow water equivalent for GCM evaluation.
Atmos.–Ocean
,
41
,
1
14
.
Brun
,
E.
,
E.
Martin
,
V.
Simon
,
C.
Gendre
, and
C.
Coleou
,
1989
:
An energy and mass model of snow cover suitable for operational avalanche forecasting.
J. Glaciol.
,
35
,
333
342
.
Burnham
,
K. P.
, and
D. R.
Anderson
,
2002
:
Model Selection and Multi-Model Inference: A Practical Information-Theoretic Approach.
2nd ed. Springer-Verlag, 496 pp
.
Carroll
,
S. S.
,
1995
:
Modeling measurement errors when estimating snow water equivalent.
J. Hydrol.
,
172
,
247
260
.
Christensen
,
J. H.
, and
Coauthors
,
2007
:
Regional climate projections.
Climate Change 2007: The Physical Science Basis, J. T. Houghton et al., Eds., Cambridge University Press, 847–940
.
Cline
,
D.
,
K.
Elder
, and
R.
Bales
,
1998
:
Scale effects in a distributed snow water equivalence and snowmelt model for mountain basins.
Hydrol. Processes
,
12
,
1527
1536
.
Congdon
,
P.
,
2003
:
Applied Bayesian Modelling.
Wiley, 457 pp
.
Davis
,
R. T.
, and
Coauthors
,
1970
:
SCS national engineering handbook: Section 22.
Snow Survey and Water Supply Forecasting, 1–210
.
Deems
,
J. S.
, and
T. H.
Painter
,
2006
:
LiDAR measurement of snow depth: Accuracy and error sources.
Int. Snow Science Workshop, Telluride, CO, Internatinoal Snow Science Steering Committee, 330–338
.
Deems
,
J. S.
,
S. R.
Fassnacht
, and
K. J.
Elder
,
2006
:
Fractal distribution of snow depth from Lidar data.
J. Hydrometeor.
,
7
,
285
297
.
Deser
,
C.
,
R.
Thomas
,
M.
Alexander
, and
D.
Lawrence
,
2010
:
The seasonal atmospheric response to projected Arctic sea ice loss in the late twenty-first century.
J. Climate
,
23
,
333
351
.
Dickinson
,
W. T.
, and
H. R.
Whitely
,
1972
:
A sampling scheme for shallow snowpacks.
IASH Bull.
,
17
,
247
258
.
Elder
,
K.
,
J.
Dozier
, and
J.
Michaelsen
,
1991
:
Snow accumulation and distribution in an alpine watershed.
Water Resour. Res.
,
27
,
1541
1552
.
Farnes
,
P. E.
,
N. R.
Peterson
,
B. E.
Goodison
, and
R. P.
Richards
,
1982
:
Metrication of manual snow sampling equipment by Western Snow Conference Metrication Committee.
Proc. 50th Annual Western Snow Conf., Reno, NV, Western Snow Conference, 120–132
.
Foster
Jr.,
D. J.
, and
R. D.
Davy
,
1988
:
Global snow depth climatology.
USAFETAC/Tech Note 88/006, USAF Environmental Technical Applications Center, 48 pp
.
Frei
,
A.
, and
G.
Gong
,
2005
:
Decadal to century scale trends in North American snow extent in coupled atmosphere-ocean general circulation models.
Geophys. Res. Lett.
,
32
,
L18502
.
doi:10.1029/2005GL023394
.
Gelman
,
A.
, and
D. B.
Rubin
,
1992
:
Inference from iterative simulation using multiple sequences.
Stat. Sci.
,
7
,
457
472
.
Gelman
,
A.
,
J. B.
Carlin
,
H. S.
Stern
, and
D. B.
Rubin
,
2004
:
Bayesian Data Analysis.
2nd ed. CRC Press, 668 pp
.
Goodison
,
B. E.
,
1978
:
Accuracy of snow samplers for measuring shallow snowpacks: An update.
Proc. 35th Annual Eastern Snow Conf., Hanover, NH, Eastern Snow Conference, 36–49
.
Hopkinson
,
C.
,
M.
Sitar
,
L.
Chasmer
, and
P.
Treitz
,
2004
:
Mapping snowpack depth beneath forest canopies using airborne Lidar.
Photogramm. Eng. Remote Sens.
,
70
,
323
330
.
Hosaka
,
M.
,
D.
Nohara
, and
A.
Kitoh
,
2005
:
Changes in snow cover and snow water equivalent due to global warming simulated by a 20 km-mesh global atmospheric model.
SOLA
,
1
,
93
96
.
doi:10.2151/sola.2005-025
.
Johnson
,
J. B.
, and
G. L.
Schaefer
,
2002
:
The influence of thermal, hydrologic, and snow deformation mechanisms on snow water equivalent pressure sensor accuracy.
Hydrol. Processes
,
16
,
3529
3542
.
Jonas
,
T.
,
C.
Marty
, and
J.
Magnusson
,
2009
:
Estimating the snow water equivalent from snow depth measurements in the Swiss Alps.
J. Hydrol.
,
378
,
161
167
.
Kojima
,
K.
,
1966
:
Densification of seasonal snow cover.
Physics of Snow and Ice, Vol. 1, H. Oura, Ed., Institute of Low Temperature Science, 929–951
.
König
,
M.
,
J. G.
Winther
, and
E.
Isaksson
,
2001
:
Measuring snow and glacier ice properties from satellite.
Rev. Geophys.
,
39
,
1
27
.
Koren
,
V.
,
J.
Schaake
,
K.
Mitchell
,
W. Y.
Duan
,
F.
Chen
, and
J. M.
Baker
,
1999
:
A parameterization of snowpack and frozen ground intended for NCEP weather and climate models.
J. Geophys. Res.
,
104
,
19569
19586
.
Lapen
,
D. R.
, and
L. W.
Martz
,
1996
:
An investigation of the spatial association between snow depth and topography in a Prairie agricultural landscape using digital terrain analysis.
J. Hydrol.
,
184
,
277
298
.
Lehning
,
M.
,
I.
Völksch
,
D.
Gustafsson
,
T. A.
Nguyen
,
M.
Stähli
, and
M.
Zappa
,
2006
:
ALPINE3D: A detailed model of mountain surface processes and its application to snow hydrology.
Hydrol. Processes
,
20
,
2111
2128
.
Lemke
,
P.
, and
Coauthors
,
2007
:
Observations: Changes in snow, ice and frozen ground.
Climate Change 2007: The Physical Sciences Basis, S. Solomon et al., Eds., Cambridge University Press, 337–383
.
Liston
,
G. E.
, and
M.
Sturm
,
2002
:
Winter precipitation patterns in arctic Alaska determined from a blowing snow model and snow-depth observations.
J. Hydrometeor.
,
3
,
646
659
.
Liston
,
G. E.
, and
K.
Elder
,
2006
:
A distributed snow-evolution modeling system (SnowModel).
J. Hydrometeor.
,
7
,
1259
1276
.
Liston
,
G. E.
,
R. B.
Haehnel
,
M.
Sturm
,
C. A.
Hiemstra
,
S.
Berezovskaya
, and
R. D.
Tabler
,
2007
:
Simulating complex snow distributions in windy environments using SnowTran-3D.
J. Glaciol.
,
53
,
241
256
.
Lunn
,
D. J.
,
A.
Thomas
,
N. G.
Best
, and
D. J.
Spiegelhalter
,
2000
:
WinBUGS—A Bayesian modeling framework: Concepts, structure, and extensibility.
Stat. Comput.
,
10
,
325
337
.
Martinec
,
I. J.
,
1966
:
Snow cover density changes in an experimental watershed.
Int. Symp. Scientific Aspects of Snow and Ice Avalanches, Davos, Switzerland, International Association of Scientific Hydrology–Commission on Snow and Ice, 29–39
.
McKay
,
G. A.
, and
B. F.
Findlay
,
1971
:
Variation of snow resources with climate and vegetation in Canada.
Proc. 39th Western Snow Conf., Billings, MT, Western Snow Conference, 17–26
.
McKay
,
G. A.
, and
D. M.
Gray
,
1981
:
The distribution of snowcover.
Handbook of Snow, 1st ed. D. M. Gray and E. D. H. Male, Eds., Pergamon Press, 153–190
.
Minoru
,
A.
, and
S.
Hiroshi
,
2006
:
Snow depth derived from airborne LIDAR data and estimation of snow water equivalent volume.
J. Japan Soc. Photogramm. Remote Sens.
,
45
,
24
33
.
Molotch
,
N. P.
,
2009
:
Reconstructing snow water equivalent in the Rio Grande headwaters using remotely sensed snow cover data and a spatially distributed snowmelt model.
Hydrol. Processes
,
23
,
1076
1089
.
doi:10.1002/hyp.7206
.
Molotch
,
N. P.
, and
S. A.
Margulis
,
2008
:
Estimating the distribution of snow water equivalent using remotely sensed snow cover data and a spatially distributed snowmelt model: A multi-resolution, multi-sensor comparison.
Adv. Water Resour.
,
31
,
1503
1514
.
doi:10.1016/j.advwatres.2008.07.017
.
Mote
,
T. L.
,
A. J.
Grundstein
,
D. J.
Leathers
, and
D. A.
Robinson
,
2003
:
A comparison of modeled, remotely sensed and measured snow water equivalent in the northern Great Plains.
Water Resour. Res.
,
39
,
1209
.
doi:10.1029/2002WR001782
.
Meteorological Service of Canada
,
2000
:
Canadian Snow Data.
Climate Research Branch, Environment Canada. CD-ROM
.
Paterson
,
W. S. B.
,
1981
:
The Physics of Glaciers.
2nd ed. Pergamon Press, 380 pp
.
Peck
,
E.
,
T.
Carroll
, and
S.
Vandermark
,
1980
:
Operational aerial snow surveying in the United States.
Hydrol. Sci. Bull.
,
25
,
51
62
.
Peterson
,
N. R.
, and
A. J.
Brown
,
1975
:
Accuracy of snow measurements.
Proc. 43rd Western Snow Conf., Coronado, CA, Western Snow Conference, 1–5
.
Plummer
,
M.
,
N.
Best
,
K.
Cowles
, and
K.
Vines
,
2006
:
CODA: Convergence Diagnosis and Output Analysis for MCMC.
R-News, 6, 7–11
.
Pomeroy
,
J. W.
, and
D. M.
Gray
,
1995
:
Snowcover accumulation, relocation, and management.
NHRI Science Rep. 7, 144 pp
.
R Development Core Team
,
cited
.
2008
:
R: A language and environment for statistical computing.
R Foundation for Statistical Computing, Vienna, Austria. [Available online at http://www.R-project.org]
.
Räisänen
,
J.
,
2008
:
Warmer climate: Less or more snow?
Climate Dyn.
,
30
,
307
319
.
Robinson
,
D. A.
,
1993
:
Hemispheric snow cover from satellites.
Ann. Glaciol.
,
17
,
367
371
.
Rovansek
,
R. J.
,
D. L.
Kane
, and
L. D.
Hinzman
,
1993
:
Improving estimates of snowpack water equivalent using double sampling.
Proc. 50th Eastern Snow Conf., Quebec City, Quebec, Canada, Eastern Snow Conference, 157–163
.
Rutter
,
N.
, and
Coauthors
,
2009
:
An evaluation of forest snow processes models (SnowMIP2).
J. Geophys. Res.
,
114
,
D06111
.
doi:10.1029/2008JD011063
.
Schaffhauser
,
A.
,
M.
Adams
,
R.
Fromm
,
P.
Jorg
,
G.
Luzi
,
L.
Noferini
, and
R.
Sailer
,
2008
:
Remote sensing based retrieval of snow cover properties.
Cold Reg. Sci. Technol.
,
54
,
164
175
.
Sickman
,
J. O.
,
A.
Leydecker
, and
J. M.
Melack
,
2001
:
Nitrogen mass balances and abiotic controls on N retention and yield in high-elevation ecosystems of the Sierra Nevada, California, USA.
Water Resour. Res.
,
37
,
1445
1461
.
Skaugen
,
T.
,
1999
:
Estimating the mean areal snow water equivalent by integration in time and space.
Hydrol. Processes
,
13
, (
12/13
).
2051
2066
.
Spiegelhalter
,
D. J.
,
N. G.
Best
,
B. P.
Carlin
, and
A.
Linde
,
2002
:
Bayesian measures of model complexity and fit.
J. Stat. Soc.
,
64B
,
583
639
.
Staehli
,
M.
, and
D.
Gustafsson
,
2006
:
Long-term investigations of the snow cover in a subalpine semi-forested catchment.
Hydrol. Processes
,
20
,
411
428
.
Steppuhn
,
H.
,
1976
:
Areal water equivalents for prairie snowcovers by centralized sampling.
Proc. 44th Western Snow Conf., Calgary, Alberta, Canada, Western Snow Conference, 63–68
.
Sturm
,
M.
, and
J.
Holmgren
,
1998
:
Differences in compaction behavior of three climate classes of snow.
Ann. Glaciol.
,
26
,
125
130
.
Sturm
,
M.
, and
G.
Liston
,
2003
:
The snow cover on lakes of the Arctic Coastal Plain of Alaska, U.S.A.
J. Glaciol.
,
49
,
370
380
.
Sturm
,
M.
,
J.
Holmgren
, and
G.
Liston
,
1995
:
A seasonal snow cover classification system for local to global applications.
J. Climate
,
8
,
1261
1283
.
Turcan
,
J.
, and
H. S.
Loijens
,
1975
:
Accuracy of snow survey data and errors in snow sampler measurements.
Proc. 32nd Eastern Snow Conf., Manchester, NH, Eastern Snow Conference, 2–11
.
Watson
,
F. G. R.
,
T. N.
Anderson
,
W. B.
Newman
,
S. E.
Alexander
, and
R. A.
Garrott
,
2006
:
Optimal sampling schemes for estimating mean snow water equivalents in stratified heterogeneous landscapes.
J. Hydrol.
,
328
,
432
452
.
Wheeler
,
D. C.
,
D. A.
Hickson
, and
L. A.
Waller
,
2010
:
Assessing local model adequacy in Bayesian hierarchical models using the partitioned deviance information criterion.
Comput. Stat. Data Anal.
,
54
,
1657
1671
.
Wilks
,
D. S.
, and
M.
McKay
,
1996
:
Extreme value statistics for snowpack water equivalent in the northeastern United States using the Cooperative Observer Network.
J. Appl. Meteor.
,
35
,
706
713
.
Work
,
R. A.
,
H. J.
Stockwell
,
T. G.
Freeman
, and
R. T.
Beaumont
,
1965
:
Accuracy of field snow surveys, western United States, including Alaska.
U.S. Army CRREL Tech. Rep. 163, U.S. Army Cold Regions Research and Engineering Laboratory, 43 pp
.

Footnotes

Corresponding author address: Dr. Matthew Sturm, USA-CRREL-Alaska, P.O. Box 35170, Ft. Wainwright, AK 99712. Email: msturm@crrel.usace.army.mil