Abstract

The annual consistency of spatial patterns of snow accumulation and melt suggests that the evolution of these patterns, known as depletion curves, is useful for estimating basin water content and runoff prediction. Theoretical snow cover depletion curves are used in models to parameterize fractional snow-covered area (fSCA) based on modeled estimates of snow accumulation and snowmelt. Directly measuring the spatio-temporal snow distribution, characterization of depletion curves, and understanding how they vary across mountainous landscapes was not possible until the recent U.S. National Aeronautics and Space Administration (NASA) Airborne Snow Observatory (ASO). Herein, for the first time, high-resolution spatio-temporal snow depth information from the ASO is used to derive observation-based snow cover depletion curves across physiographic gradients by estimating the slope of the fSCA–snow depth relationship (i.e. depletion slope). The depletion slope reveals important insights into snow processes as it is strongly related to snow depth variability (r2 = 0.58). Regression tree analysis between observed depletion slopes and physiography, particularly vegetation height and terrain roughness, displays clear nonlinear dynamics and explains 31% of the variance in depletion slope. This unique observation-based analysis of snow cover depletion curves has implications for energy and water flux calculations across many earth system models.

HIGHLIGHTS

  • Relationships between snow depth and fSCA (i.e. depletion slope) were robust over the 4 years of study.

  • Significant spatial variability in depletion slope is well correlated with snow depth variability.

  • Increased vegetation height and decreased terrain roughness were associated with more homogeneous snowpacks and lower depletion slopes.

INTRODUCTION

The spatial distribution of snow water equivalent (SWE) is affected by numerous processes at multiple scales including, but not limited to, orography, wind and avalanche redistribution, and ablation dynamics driven by the land-surface energy balance. Inherently, these factors are strongly influenced by topography and vegetation and therefore result in annually repeating spatial patterns (Kirnbauer & Bloschl 1994; König & Sturm 1998; Deems et al. 2008; Sturm & Wagner 2010). Our ability to characterize these persistent patterns of snow cover and their relationships with snow accumulation have direct relevance for forecasting water supplies (Potts 1937; Hannaford et al. 1979; Kolberg & Gottschalk 2010) and estimating land–atmosphere energy exchange and atmospheric processes (Flanner & Zender 2006). In addition, patterns of snow distribution have relevance for understanding the distribution of soil moisture (Harpold & Molotch 2015), ecosystem response to water availability (Trujillo et al. 2012), and wintertime animal movement and habitat (Parker et al. 1984; Bhattacharyya et al. 2014). Depletion curves that link the snow-covered area (SCA) with SWE are a frequently used method for accurately estimating the contributing area to snowmelt runoff (Anderson 1973; Luce et al. 1999) and for characterizing snow-albedo feedbacks within land-surface and general circulation models (Niu & Yang 2007; Swenson & Lawrence 2012).

The existing literature regarding depletion curves typically relates SWE and SCA using empirical relationships or probability distribution functions for modeling applications. These cover a variety of spatial scales from point-based derivations over hillslopes (Anderson 1973; Luce et al. 1999; Luce & Tarboton 2004), physically based gridded models applied to the watershed scale (Clark et al. 2011), and regional to global-scale land-surface models (Niu & Yang 2007; Swenson & Lawrence 2012; Driscoll et al. 2017). However, the scale at which snow heterogeneity is consistent is ambiguous, prompting statistical parameterizations of depletion curves based on the coefficient of variation (CV) of SWE (Liston 2004). Estimates of the CV of SWE and snow depth for different environments have been provided in the literature based on field studies (Liston 2004; Clark et al. 2011); however, these are broad classifications that may not hold true at horizontal scales larger than 100 m, in mixed landscapes, or for elevation differences larger than 200 m (Essery & Pomeroy 2004; Clark et al. 2011). Importantly, satellite-based methods of SWE estimation have limited utility resolving the CV of SWE and snow depth, given that the accuracy of SWE estimates is highly sensitive to SWE spatial variability (Vander Jagt et al. 2013).

Multiple studies have investigated sub-grid depletion curves of SWE based on terrain characteristics (Donald et al. 1995; Niu & Yang 2007; Swenson & Lawrence 2012; Helbig et al. 2015; Fassnacht et al. 2016; Driscoll et al. 2017). In particular, Niu & Yang (2007) and Swenson & Lawrence (2012) observed that snow cover depletion patterns across the USA varied substantially between cells with differing sub-grid variances in elevation. Accordingly, they developed scale-variant SWE depletion curve parameterizations based on sub-grid elevation variability. Topographic controls on SWE–SCA relationships exhibit significant interannual consistency (Fassnacht et al. 2016; Driscoll et al. 2017), suggesting that the shape of depletion curves is largely controlled by the terrain that affects both snow accumulation and ablation.

Few studies have explicitly investigated the relationship between snow depth and SCA, perhaps because snow depth has been perceived as less relevant than SWE for estimating water supply. Helbig et al. (2015) parameterized a snow depth depletion curve with the sub-grid standard deviation of snow depth at peak accumulation. Similarly, Egli et al. (2012) indicated that estimates of maximum snow depth distribution were critical to model the evolution of basin mean snow depth and SCA throughout the melt season. López-Moreno et al. (2014) observed a consistent spatial pattern of the CV of snow depth across multiple lidar acquisitions in time, and Egli & Jonas (2009) found that variation in year-to-year snow depth depletion was largely controlled by the standard deviation of snow depth at peak snow depth. Other works have also related the standard deviation of snow depth to the variability of the underlying terrain (Marchand & Killingtveit 2005; López-Moreno et al. 2014; Helbig et al. 2015). These studies suggest that topographic variability that affects the mean and standard deviation of snow depth will influence the relationship between snow depth and SCA.

These previous studies were largely based on theory, modeling results, or relatively sparse measurements of the snowpack in time or space because regular, spatially extensive measurements such as those of the U.S. National Aeronautics and Space Administration (NASA) Airborne Snow Observatory (ASO) to assess depletion curve characteristics have not been available until the last 7 years. Hence, we still lack a process-level understanding of the spatial variability in depletion curve characteristics and how individual terrain elements at scales less than 100 m affect this behavior.

The NASA ASO dataset offers a new opportunity to observe and investigate spatial differences in depletion curves that result from topographic influences on snow accumulation and ablation processes (Painter et al. 2016). The ASO dataset provides high spatial resolution time series of lidar-derived snow depth and hyperspectral measurements of snow properties over multiple years from which relationships between SCA and snow depth can be determined. ASO also produces a SWE dataset that combines the lidar-derived snow depth with modeled density, but we focus on snow depth depletion curves in this study to keep the analysis anchored by direct measurements. However, improvements in characterization of snow depth depletion curve variability are relevant to the estimation of basin SWE since snow depth varies an order of magnitude more than snow density (Mizukami & Perica 2008).

Our objective is to gain insight into the processes by which snow and terrain interact to produce differences in SCA. We do not attempt to parameterize a universal depletion curve but rather to assess the influence of terrain characteristics on depletion curve shape. Specifically, we ask: what are the primary physiographic controls on the characteristics of snow depth depletion curves?

METHODS

Site description

This study was conducted on data collected by the ASO in the Tuolumne River basin in the Sierra Nevada Mountains in CA, USA (Figure 1). The basin area is 1,175 km2 with land cover that is 48% vegetation, 50% rock/talus, 2% open water, and << 1% permanent snow/ice. The elevation range of 1,127–3,997 m encompasses four distinct ecological zones: the lower montane forest; upper montane forest; a sub-alpine mixed-meadow/coniferous forest zone; and an alpine zone with large granitic features, talus slopes, and limited herbaceous vegetation (NPS 2016).

Figure 1

Landcover at 500 m resolution with 250 m elevation contours and location map of the Tuolumne basin, CA, USA. Landcover was determined based on a spectral reflectance classification by ASO. Elevation contours were extracted from an ASO lidar-derived DEM. Additionally, water and ice were delineated from the National Hydrography Dataset (U.S. Geological Survey 2016).

Figure 1

Landcover at 500 m resolution with 250 m elevation contours and location map of the Tuolumne basin, CA, USA. Landcover was determined based on a spectral reflectance classification by ASO. Elevation contours were extracted from an ASO lidar-derived DEM. Additionally, water and ice were delineated from the National Hydrography Dataset (U.S. Geological Survey 2016).

The climate is characterized as Mediterranean with relatively mild winter temperatures and dry, warm summers. The snow season typically lasts from November to March with more than 70% of the yearly precipitation falling as snow. A long-term record of SWE since 1967 (snow depth is only available since 2004) from 17 snow pillows within 20 km of the Tuolumne Basin recorded a mean peak SWE of 0.8 m with a range of 0.3–1.5 m. The four study years were characterized by below average snowpack with 2014 and 2015 characterized as a severe dry snow drought (Harpold et al. 2017), experiencing only 36 and 35% of the climatological mean peak SWE, respectively. The years 2013 and 2016 also experienced below average snowpack conditions, but less severely, with the data from the snow pillows reporting 62 and 71% of climatological peak SWE (data available from http://cdec.water.ca.gov/). Runoff from the watershed flows into the Hetch Hetchy Reservoir, the water supply for the City of San Francisco and other Bay Area municipalities.

Data sources

We use the ASO dataset which is derived from a paired scanning lidar and imaging spectrometer (Painter et al. 2016). We use the ASO 3 m resolution snow-free digital elevation model (DEM), 3 m vegetation height, and 3 m resolution snow depth datasets. The ASO flew approximately weekly from near peak SWE through the snowmelt season, with six flights in 2013, ten in 2014, eight in 2015, and eight in 2016. Painter et al. (2016) report a mean absolute vertical error of 8 cm and a bias of <1 cm at 3 m spatial resolution when compared with manually measured snow depths. Further details about ASO and data processing can be found in Painter et al. (2016).

The 3 m resolution gridded snow depth data are aggregated to two new 495 m (nominally 500 m) resolution datasets, containing the mean and standard deviation of snow depth for each 500 m grid cell, respectively. Similarly, fractional SCA (fSCA) maps at 500-m resolution are derived from the 3 m gridded snow depth data based on a binary snow-free/snow-on schema. Pixels at the edge of the domain were not included if the 500 m grid cell extended past the boundary of the basin and pixels containing water bodies were also removed. We chose 500 m resolution for the analyses herein because it is the scale at which daily satellite observations of fSCA exist from the MODerate-resolution Imaging Spectroradiometer (MODIS); hence, our results are theoretically transferable to applications using MODIS data. We consider the ASO snow-free classification to be robust, because the original 3 m gridded snow depth data were masked for snow-free areas based on spectral data from the spectrometer (Painter et al. 2016). ASO SCA has also compared well against SCA derived from 0.3 m resolution World View imagery, with Recall, Precision, Accuracy, and F statistics all >0.99 (Bair et al. 2016).

Deriving depletion curves

We derived depletion curves for each 500 m grid cell through time using the 33 flights available from the ASO. For each grid cell, we analyzed 33 snow depth data points versus 33 fSCA data points and characterized the depletion curve with a bilinear regression between snow depth (dependent variable) and fSCA (independent variable). Theoretical considerations would motivate a curvilinear fit because fSCA can, by definition, only increase to a value of one, whereas maximum snow depth is unbounded (Liston 2004). However, we observed distinctly linear relationships between snow depth and fSCA for low fSCA and snow depth values in each grid cell, which is consistent with previous studies (Swenson & Lawrence 2012). As stated above, the objective of this work is not to parameterize a universal depletion curve but rather to evaluate the relationships between topographic characteristics and snow depletion dynamics. Hence, we quantify the depletion curve for each 500 m pixel based on a segmented bilinear regression, resulting in two slopes and a breakpoint value between the two linear segments.

We tested the statistical relationship for a significant breakpoint (p < 0.1) in the bilinear regression based on a Davies' Test as implemented in the R package ‘segmented’ (Muggeo 2003, 2020). The coordinate of the breakpoint represents the fSCA and snow depth at which the fitted slope of the depletion curve changes (Figure 2(a)–2(c)). Grid cells without a significant breakpoint were fit with an ordinary least squares linear regression and the slope tested for significance. Grid cells with no breakpoint but with significant slopes less than 1.6 were included in our analyses under the assumption that these slopes reflect a similar process as the first slope of the bilinear regression but did not accumulate snow beyond the theoretical breakpoint (Figure 2(d)). Visual inspection of the data indicated that significant linear slopes greater than 1.6 typically represented grid cells that did not fully deplete within the ASO observation period and the x-intercept for these linear models were often much greater than zero (Supplementary Figure S1).

Figure 2

(a–d) Example depletion curves from four of the 2,596 grid cells from the Tuolumne. These were chosen because they represent a range of depletion slopes and nicely illustrate the observed depletion dynamics. All breakpoints and slopes are statistically significant (p < 0.1). Differences from 0 of the y-intercept are not statistically significant (p < 0.05). (e) The distribution of variance explained by the models in each of the 2,596 grid cells.

Figure 2

(a–d) Example depletion curves from four of the 2,596 grid cells from the Tuolumne. These were chosen because they represent a range of depletion slopes and nicely illustrate the observed depletion dynamics. All breakpoints and slopes are statistically significant (p < 0.1). Differences from 0 of the y-intercept are not statistically significant (p < 0.05). (e) The distribution of variance explained by the models in each of the 2,596 grid cells.

Consequently, we examined the topographic controls on the first slope of the bilinear regression and slope values of the ordinary linear regression below values of 1.6, and we denote these as the ‘depletion slope’ – directly defined by the breakpoint coordinates and a y-intercept of zero. Significant y-intercept differences from zero (p < 0.05) are rare, and we consider this to be a function of vertical uncertainty in the measurements. The upper slope of the bilinear regression was found to be uncorrelated with topography. This is consistent with the concept that the upper portions of the curves describe dynamics of snow depths greater than the sub-grid topographic variability. Thus, they are not explicitly addressed herein.

Relating depletion slope to topography

We analyze 2,596 grid cells (500 m) for relationships between topographic variables and the depletion slope. Focusing on the depletion slope is particularly revealing with regard to the compounded result from various snow distribution processes. Fundamentally, steeper depletion slopes describe a relatively heterogeneous snow depth distribution where snowpacks covering rough terrain result in areas of relatively deep snow accumulations. As such, for an equivalent decrease in snow depth, the decline in fSCA for grid cells with these deep pockets of snow is relatively small compared with grid cells with shallower, more uniform snow distributions. Intuitively, these different pockets of snow are influenced by the underlying and surrounding topography.

We use regression tree analysis (Breiman 1984) to identify potential physiographic controls (independent variables) on the depletion slope (dependent variable). We consider a range of physiographic variables that have demonstrated physical or statistical relationships with snow depth previously under the assumption that the snow depth–fSCA relationship is determined both by the magnitude and variability of snow depth. We investigate three different scales of control including grid cell scale (500 m) physiography, inter-grid cell effects, and sub-grid variability within each grid cell (Table 1).

First, the pixel-scale physiographic variables represent potential physical process controls on the magnitude of snow depth: elevation (elev), slope (zness), aspect (northness; eastness), and vegetation height (vegheight). Second, the inter-pixel topographic variables, topographic position index (tpi) and relief, quantify the topographic variability within a 3 × 3 pixel window around each 500 m grid cell and also represent potential controls on the magnitude of snow depth. Third, variables representing sub-grid cell scale variability describe terrain roughness based on elevation, slope, and curvature within each 500 m grid cell and represent potential controls on the variability of snow depth. Further information on the pixel-scale and inter-pixel-scale variables are provided in Table 1 and the references therein. Further detail on the sub-grid cell variables is included below as these variables are less prevalent in the literature.

Table 1

Variables used in a regression tree to explain depletion curve slopes

Variable nameShort nameResolutionDefinition/UnitsSource
Grid cell scale physiography 
Elevation elev 495 m Elevation from the summer flight; meters Elder et al. (1991, 1998
Slope zness 495 m Sin(slope); ranges 0–1; dimensionless Balk & Elder (2000); Erxleben et al. (2002); Fassnacht et al. (2003)  
Northness northness 495 m Cos(aspect); ranges 0–1; dimensionless Balk & Elder (2000); Erxleben et al. (2002); Fassnacht et al. (2003)  
Eastness eastness 495 m Sin(aspect); ranges 0–1; dimensionless Balk & Elder (2000); Erxleben et al. (2002); Fassnacht et al. (2003)  
Vegetation height vegheight 495 m Mean vegetation height measured by ASO from the summer flight; meters Deems et al. (2006); Trujillo et al. (2007); Painter et al. (2016)  
Inter-grid cell variability 
Topographic position index tpi 495 m from a 3 × 3 pixel window around each pixel Elevation difference of a pixel from the mean of the surrounding pixels; meters Revuelto et al. (2014); GDAL (2015)  
Relief relief 495 m from a 3 × 3 pixel window around each pixel Maximum elevation difference between any two of the surrounding pixels; meters Bookhagen & Burbank (2010); GDAL (2015)  
Sub-grid variability 
Standard deviation of elevation stdelev 495 m from 45 m Standard deviation of elevation from 45 m DEM; meters Marchand & Killingveit (2005)  
Standard deviation of slope stdslope 495 m from 45 m Standard deviation of slope from 45 m DEM; shown to detect changes in slope at multiple scales; radians Grohmann et al. (2011); Fassnacht et al. (2016)  
Standard deviation of maximum curvature stdmaxcurv 495 m from 45 m Standard deviation of maximum curvature from 45 m DEM; maximum curvature from SAGA GIS; dimensionless Marchand & Killingveit (2005); López-Moreno et al. (2014); Conrad et al. (2015)  
Terrain variability Δhw 495 m from 3 m Terrain variability metric defined in Equation (1) in Helbig et al. (2015); computed from the 3 m DEM; dimensionless Helbig et al. (2015)  
Variable nameShort nameResolutionDefinition/UnitsSource
Grid cell scale physiography 
Elevation elev 495 m Elevation from the summer flight; meters Elder et al. (1991, 1998
Slope zness 495 m Sin(slope); ranges 0–1; dimensionless Balk & Elder (2000); Erxleben et al. (2002); Fassnacht et al. (2003)  
Northness northness 495 m Cos(aspect); ranges 0–1; dimensionless Balk & Elder (2000); Erxleben et al. (2002); Fassnacht et al. (2003)  
Eastness eastness 495 m Sin(aspect); ranges 0–1; dimensionless Balk & Elder (2000); Erxleben et al. (2002); Fassnacht et al. (2003)  
Vegetation height vegheight 495 m Mean vegetation height measured by ASO from the summer flight; meters Deems et al. (2006); Trujillo et al. (2007); Painter et al. (2016)  
Inter-grid cell variability 
Topographic position index tpi 495 m from a 3 × 3 pixel window around each pixel Elevation difference of a pixel from the mean of the surrounding pixels; meters Revuelto et al. (2014); GDAL (2015)  
Relief relief 495 m from a 3 × 3 pixel window around each pixel Maximum elevation difference between any two of the surrounding pixels; meters Bookhagen & Burbank (2010); GDAL (2015)  
Sub-grid variability 
Standard deviation of elevation stdelev 495 m from 45 m Standard deviation of elevation from 45 m DEM; meters Marchand & Killingveit (2005)  
Standard deviation of slope stdslope 495 m from 45 m Standard deviation of slope from 45 m DEM; shown to detect changes in slope at multiple scales; radians Grohmann et al. (2011); Fassnacht et al. (2016)  
Standard deviation of maximum curvature stdmaxcurv 495 m from 45 m Standard deviation of maximum curvature from 45 m DEM; maximum curvature from SAGA GIS; dimensionless Marchand & Killingveit (2005); López-Moreno et al. (2014); Conrad et al. (2015)  
Terrain variability Δhw 495 m from 3 m Terrain variability metric defined in Equation (1) in Helbig et al. (2015); computed from the 3 m DEM; dimensionless Helbig et al. (2015)  

Variables were chosen based on their existing use in the literature; sources listed are by no means exhaustive.

We evaluated two scales of sub-grid cell variability at 3 and 45 m resolutions. At the finer, 3 m scale, we derive the ‘Δhw’ parameter, which represents the ratio of the typical height of a terrain feature to the typical width of a terrain feature. Larger Δhw values indicate a greater relative protrusion of the terrain. We calculated Δhw for each 500 m grid cell from the 3 m ASO DEM after the following equation (Helbig et al. 2015):
formula
(1)
where ∂xz and ∂yz are the first partial derivatives of elevation, i.e. the orthogonal slope components (Neteler et al. 2012).

At the coarser scale, we quantify the standard deviation of elevation (stdelev), the standard deviation of slope (stdslope), and the standard deviation of maximum curvature (stdmaxcurv). Previous studies have shown the length scale at which snow distribution and depletion processes change to be between 6 and 40 m (Deems et al. 2006; Trujillo et al. 2007; Schirmer et al. 2011), so we use a 45 m resolution DEM in order to explicitly treat variation native to processes operative at, or greater than, this scale and the 3 m DEM to capture the processes at finer scales.

Regression trees have commonly been used to estimate the spatial distribution of snow (e.g. Elder et al. 1998; Balk & Elder 2000; Molotch et al. 2005). A regression tree can model nonlinear interactions between the topography and depletion slope. We grew regression trees until additional splits no longer improved the cross-validated r2 by 0.001. We subsequently pruned the regression trees to an optimal complexity by selecting the number of splits corresponding to within 1 standard error of the minimum mean squared error (Breiman 1984; R Core Team 2015; Therneau et al. 2015).

RESULTS AND DISCUSSION

Depletion curves

Relationships between snow depth and fSCA are remarkably consistent in all 4 years despite significant interannual variability in weather patterns, i.e. points from different years plot along the same depletion curve in Figure 2. On average, the model for each grid cell explained 81% of the variance in snow depth (Figure 2(e)). Moreover, 84% of depletion slopes fall within the 95% confidence interval of the mean cross-validated depletion slope (where four depletion slopes were calculated with the data for each year iteratively removed). This suggests that the snow depth–fSCA relationships presented here are robust even in a year with an otherwise anomalous snowpack.

The range of the 2,596 depletion slopes analyzed with respect to topography is 0.20–3.06 with a mean of 0.73 (Figure 3(a)). The spatial pattern of grid cells that had statistically significant depletion slopes corresponds to higher elevation grid cells (mean 2,822 m with a range of 1,276–3,607 m) of the basin (Figure 3(a) and 3(b); see Figure 1 for elevation contours). At lower elevations, shallow snow depths and a lack of clear depletion dynamics – likely driven by the drought conditions during the study period – rendered the regressions statistically insignificant (white, Figure 3(a)). The highest depletion slopes (darkest colors, Figure 3(a)) roughly correspond with deeper snow and occur along the southern and northern basin divides. The lowest values (lightest colors, Figure 3(a)) are clustered in the southeast.

Figure 3

The significant depletion slope values (a), the fraction of grid cells with significant depletion slopes by elevation (b), and temporally averaged standard deviations of snow depth (c). The value of each 500 m grid cell was computed based on the standard deviations of the 3 m snow depth observations within each 500 m grid cell on each date. The map values are divided into 20% quantiles, so that each color represents the same number of observations. 250 m elevation contours are shown for context. (d) Depletion slopes versus average standard deviations of snow depth with the best-fit linear relationship. Please refer to the online version of this paper to see this figure in color: http://dx.doi.10.2166/nh.2020.267.

Figure 3

The significant depletion slope values (a), the fraction of grid cells with significant depletion slopes by elevation (b), and temporally averaged standard deviations of snow depth (c). The value of each 500 m grid cell was computed based on the standard deviations of the 3 m snow depth observations within each 500 m grid cell on each date. The map values are divided into 20% quantiles, so that each color represents the same number of observations. 250 m elevation contours are shown for context. (d) Depletion slopes versus average standard deviations of snow depth with the best-fit linear relationship. Please refer to the online version of this paper to see this figure in color: http://dx.doi.10.2166/nh.2020.267.

The temporally averaged standard deviations of snow depth (Figure 3(c)) exhibit a similar spatial pattern as the depletion slope (Figure 3(a)). A positive relationship between depletion slope and snow depth standard deviation exists across the basin with r2 = 0.58; p < 0.01 (Figure 3(d)). This supports our interpretation of the depletion slope that grid cells with steeper depletion slopes have more heterogeneous snowpacks.

Topographic controls on depletion curves

None of the topographic variables has a strong individual relationship with the depletion slope (Table 2).

Table 2

The single variable Kendall rank correlation (two-sided) and associated p-value between each of the physiographic predictor variables and the depletion slope

Topographic variableKendall rank correlationp-value
stdmaxcurv 0.27 < 0.01 
Δhw 0.26 < 0.01 
stdslope 0.25 < 0.01 
vegheight − 0.23 < 0.01 
stdelev 0.21 < 0.01 
elev 0.20 < 0.01 
northness 0.19 < 0.01 
relief 0.13 < 0.01 
tpi 0.13 < 0.01 
zness 0.09 < 0.01 
eastness 0.01 0.41 
Topographic variableKendall rank correlationp-value
stdmaxcurv 0.27 < 0.01 
Δhw 0.26 < 0.01 
stdslope 0.25 < 0.01 
vegheight − 0.23 < 0.01 
stdelev 0.21 < 0.01 
elev 0.20 < 0.01 
northness 0.19 < 0.01 
relief 0.13 < 0.01 
tpi 0.13 < 0.01 
zness 0.09 < 0.01 
eastness 0.01 0.41 

The highest magnitude correlation with the depletion slope is displayed by the standard deviation of maximum terrain curvature (r = 0.27) followed by the Δh/Δw parameter (r = 0.26) and the standard deviation of terrain slope (r = 0.25). The prevalence of three different sub-pixel topographic variables indicates that sub-pixel terrain variability exerts substantial influence on the depletion slope. This finding agrees well with previous work, showing snow depth variability is well correlated with topographic variability, in part due to wind redistribution (Winstral et al. 2002; Erickson et al. 2005). All correlations are statistically significant except for eastness, possibly due to gradients in precipitation shadowing from north to south in the Tuolumne Basin.

The regression tree analysis explains 31% (95% CI: ±0.2%) of the variance in the depletion slope (Figure 4). Vegetation height explains the most variance with the depletion slope, with taller vegetation corresponding with smaller depletion slopes and shorter to no vegetation corresponding with larger depletion slopes. This implies that taller vegetation leads to more homogeneous snowpacks and shorter vegetation leads to more heterogeneous snowpacks, which agree with the previous work showing that areas with forests have more homogenous snow distribution due to reduced wind redistribution (Hiemstra et al. 2002; López-Moreno & Latron 2008). We interpret this split in the regression tree as indicative of differences in accumulation and depletion dynamics above and below the timberline. We note that northness, elev, and multiple sub-grid cell variables are used for subsequent splits both above and below the tree line. A key difference between the regression tree results and the single variable correlation presented in Table 2 is the nonlinearity inherent in the regression tree results. This is apparent with the appearance of elev, stdmaxcurv, and northness multiple times at different levels of the tree. They explain a higher variance than other variables for individual subsets of the data.

Figure 4

The resultant regression tree for predicting depletion slope using the physiographic variables. The boxes are shaded light to dark indicating low to high depletion slope, respectively. The top number in each box is the depletion slope estimated by the regression tree. The percent in each box is the percent of data in each bin. The value under the box is the value of the independent variable at which the model split the data.

Figure 4

The resultant regression tree for predicting depletion slope using the physiographic variables. The boxes are shaded light to dark indicating low to high depletion slope, respectively. The top number in each box is the depletion slope estimated by the regression tree. The percent in each box is the percent of data in each bin. The value under the box is the value of the independent variable at which the model split the data.

We observe that increased terrain roughness at sub-grid scales – represented by the stdmaxcurv, stdelev, and Δh/Δw variables – correlates with larger depletion slopes at multiple nodes in the regression tree. Lower sub-grid terrain variability is associated with lower depletion slope values, consistent with expectations of relatively homogenous snow distribution. This finding is intuitive, given that these areas have shallower terrain depressions in which to trap snow. In the latter stages of depletion, decreases in snow depth in these areas will yield relatively large decreases in SCA because large portions of bare ground will be exposed nearly simultaneously across the grid cell. We also tested sub-grid standard deviation of vegetation height as another measure of sub-grid variability that could potentially be important near the tree line, but found mean vegetation height to explain more variance.

Northness occurs twice in the regression tree, and in each case, larger northness values are associated with increased depletion slopes, i.e. north-facing grid cells exhibit more heterogeneous snow distribution which is congruent with Lehning et al. (2011). This may be a result of preferential snow accumulation on northeast-facing slopes due to preferential deposition of precipitation and snow drift (Lehning et al. 2008; Dadic et al. 2010). Moreover, flow separation across ridge-tops can produce cornices which can exhibit particularly heterogeneous snow depths (Anderton et al. 2004). We also considered increased small-scale terrain variability, e.g. boulder fields, as a confounding factor but found no correlation between Δh/Δw and northness.

Higher elevations are associated with steeper depletion slopes and therefore a more heterogeneous snowpack. Further work is warranted to explore the physical processes linking the depletion slope and elevation, but several possibilities exist. Increased snow heterogeneity at higher elevations can be produced by more prevalent wind redistribution effects, i.e. at higher elevations, increased snowfall and wind speeds and lower-density snow result in greater wind redistribution relative to lower elevations. In addition, ridgelines that exhibit cornices are typically at the higher elevations. Conversely, shallow depletion slopes would result from a less-distinct accumulation season, e.g. at the lowest elevations, where snowfall results in broad, shallow coverage but quickly melts away. In addition, we might also expect grid cells with high elevation and high terrain curvature to have steep depletion slopes as a consequence of avalanche redistribution which leaves narrow areas of high snow depth.

Implications

The relationship between snow depth and fSCA appears to be consistent across 4 years during the depletion period, suggesting that remotely sensed fSCA contains a large amount of information regarding snow distribution. Moreover, we show, for the first time, that repeat lidar scans of snow depth distribution are a novel data source for developing fSCA–snow depth relationships over large areas (i.e. >1,000 km2). The establishment of these relationships from explicit observations provides a framework for improving depletion curve parameterizations in hydrologic, land-surface, and climate models with the expansion of the ASO program to other mountain ranges and other efforts to explicitly map snow depth and fSCA (e.g. Nolan et al. 2015). Incorporating physically representative empirical factors that capture accumulation and ablation dynamics over multiple scales could result in improved spatial representation of snow-free and SCAs in models thus improving energy balance calculations.

The findings presented in this paper are conceptually consistent with a global classification system of SWE variability used for parameterizing probability distribution depletion curves (Sturm et al. 1995; Liston 2004). It is encouraging that depletion curves developed from high-resolution spatio-temporal data similarly describe more heterogeneous snowpacks in areas with greater topographic variability and greater wind speeds.

We see a significantly larger range in CV of snow depth than previously reported with the mean observed sub-grid cell CV (3 m spatial resolution snow depth within each 500 m grid cell) ranging from 2.4 to 34.0 across all dates. This is compared with field values of CV of snow depth for mountainous regions of 0.1–3 (Clark et al. 2011). We consider this inconsistency to be a function of scale. The 3 m spatial resolution snow depth dataset from ASO samples the complete range of terrain and snow depths compared with the studies reviewed in Clark et al. (2011), and thus captures snow depth extremes.

The results of our regression tree analyses are specific to the Tuolumne basin for two very dry years and two dry to average years. Future analyses of repeat lidar measurements across time and space are warranted to determine the spatio-temporal transferability of our results. In this regard, the regression tree analysis misses the extremes of the range of depletion slope values, which is not surprising given an r2 value of 0.31. There is considerable variability in the observed depletion slopes (Figure 3(a)) that may have resulted from local processes that our predictor variables do not capture. This suggests that there were either too few grid cells with very high depletion slopes for a meaningful statistical relationship at the extremes of the observed parameters, or that a key variable is missing from our analysis. Furthermore, grid cells without a breakpoint but with linear slopes greater than 1.6 were excluded from the topographic analysis (see Methods and Supplementary Information for details); ASO flights did not extend late enough into the year to observe the lowest values of fSCA and depth for these grid cells.

Herein, we show that snow depletion curve shape is inherently linked to snow depth distribution. Hence, further analysis of remotely sensed fSCA depletion may improve understanding of the processes controlling snow distribution. In particular, the fSCA and snow depth components of the depletion slope may be analogous to the scale break or fractal dimension reported in previous studies that investigated fractal scaling in snow depth distributions (Deems et al. 2006, 2008; Trujillo et al. 2009; Schirmer & Lehning 2011). Further work aimed at identifying relationships between topography and the snow depth value corresponding to the depletion curve breakpoint could provide a means to parameterize depletion curve shape within hydrologic models. This could provide useful information not only about the depletion slope during ablation, but also about the capacity of the terrain to trap snow in concavities.

We also performed the same analysis with SWE from ASO to link to the previous literature but did not find a consistent bilinear relationship. We might expect differences in the shapes of depletion curves derived from SWE rather than snow depth because density strongly increases early in the ablation season, thus offsetting decreases in snow depth. Consistent density during melt as fSCA is decreasing may then contribute to a more dynamic relationship between SWE and fSCA. ASO SWE estimates use field-measurement constrained, model-based snow density information with density uncertainties of 12–30 kg m−3 (Painter et al. 2016). Future work evaluating SWE depletion with respect to topography is warranted. Given the extensive use of SWE depletion curves in hydrologic modeling, improved linkage between snow depth depletion curves and SWE depletion curves would also prove valuable.

CONCLUSIONS

We derived individual depletion curves for 500 m grid cells using a lidar-derived, high-resolution, spatio-temporal dataset of observed snow depth and analyzed the effect of physiography on these depletion curves. We found that relationships between snow depth and fSCA (i.e. depletion slope) were robust over the 4 years of study with a mean r2 value of 0.81. Additionally, 84% of presented depletion slopes fall within the 95% confidence interval of the mean cross-validated depletion slope. These depletion slopes exhibited significant variability across the watershed. In this context, we show that a positive relationship with r2 = 0.58; p < 0.01 between depletion slopes and snow depth variability exists across the basin. We also show that sub-pixel and pixel-scale terrain variables explain 31% of the spatial variability in the depletion slope. In particular, increased vegetation height and decreased sub-pixel terrain variability were associated with more homogeneous snowpacks and lower depletion slopes. These results illustrate that repeat, distributed snow depth measurements such as those from the NASA ASO can provide insights to the influence of topography on the evolution of snow distribution. Such understanding has important implications for developing parameterizations of snow cover depletion curves across physiographic gradients. Given that parameterizations of snow cover depletion underpin sub-grid representation of energy and water fluxes across a range of earth system models, the results and approach presented herein has potentially broad applicability.

ACKNOWLEDGEMENTS

Data are available from https://doi.org/10.17605/OSF.IO/A96BP. D.S. was supported by a NASA Earth and Space Science Fellowship (NNX14AL27H). N.P.M. was supported by NASA Grants NNX17AF50G and 80NSSC17K0071, and USDA grant 2012-67003-19802. T.H.P. acknowledges the NASA Terrestrial Hydrology Program, NASA Applied Sciences, and the California Department of Water Resources. Part of this work was performed at the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. Special thanks to contributors to R, the tidyverse and raster packages, and GDAL (Bivand et al. 2015; GDAL/OGR contributors 2015; Hijmans 2015; R Core Team 2015; Wickham et al. 2019).

DATA AVAILABILITY STATEMENT

All relevant data are available from an online repository or repositories (https://doi.org/10.17605/OSF.IO/A96BP).

REFERENCES

Anderson
E.
1973
Snow Accumulation and Ablation Model
.
U.S. Department of Commerce, Silver Spring
,
MD
.
Anderton
S. P.
White
S. M.
Alvera
B.
2004
Evaluation of spatial variability in snow water equivalent for a high mountain catchment
.
Hydrological Processes
18
(
3
),
435
453
.
doi:10.1002/hyp.1319
.
Bair
E. H.
Rittger
K.
Davis
R. E.
Painter
T. H.
Dozier
J.
2016
Validating reconstruction of snow water equivalent in California's Sierra Nevada using measurements from the NASA Airborne Snow Observatory
.
Water Resources Research
52
(
11
),
8437
8460
.
doi:10.1002/2016WR018704
.
Balk
B.
Elder
E.
2000
Combining binary decision tree and geostatistical methods to estimate snow distribution in a mountain watershed
.
Water Resources Research
36
(
1
),
13
26
.
doi:10.1029/1999WR900251
.
Bhattacharyya
S.
Adhikari
B. S.
Rawat
G. S.
2014
Influence of snow, food, and rock cover on Royle's Pika abundance in Western Himalaya
.
Arctic, Antarctic, and Alpine Research
46
(
3
),
558
567
.
doi:10.1657/1938-4246-46.3.558
.
Bivand
R.
Keitt
T.
Rowlingson
B.
2015
rgdal: Bindings for the Geospatial Data Abstraction Library.
.
Bookhagen
B.
Burbank
D. W.
2010
Toward a complete Himalayan hydrological budget: spatiotemporal distribution of snowmelt and rainfall and their impact on river discharge
.
Journal of Geophysical Research
115
(
F3
),
F03019
.
doi:10.1029/2009JF001426
.
Breiman
L.
1984
Classification and Regression Trees
, 1st edn.
Routledge
,
New York
.
Clark
M. P.
Hendrikx
J.
Slater
A. G.
Kavetski
D.
Anderson
B.
Cullen
N. J.
Kerr
T.
Örn Hreinsson
E.
Woods
R. A.
2011
Representing spatial variability of snow water equivalent in hydrologic and land-surface models: a review
.
Water Resources Research
47
(
7
).
doi:10.1029/2011WR010745
.
Conrad
O.
Bechtel
B.
Bock
M.
Dietrich
H.
Fischer
E.
Gerlitz
L.
Wehberg
J.
Wichmann
V.
Boehner
J.
2015
System for Automated Geoscientific Analyses (SAGA) v. 2.1.4
.
Geoscientific Model Development
8
(
7
),
1991
2007
.
doi:10.5194/gmd-8-1991-2015
.
Dadic
R.
Mott
R.
Lehning
M.
Burlando
P.
2010
Wind influence on snow depth distribution and accumulation over glaciers
.
Journal of Geophysical Research
115
(
F1
).
doi:10.1029/2009JF001261
.
Deems
J. S.
Fassnacht
S. R.
Elder
K. J.
2006
Fractal distribution of snow depth from lidar data
.
Journal of Hydrometeorology
7
(
2
),
285
297
.
doi:10.1175/JHM487.1
.
Deems
J. S.
Fassnacht
S. R.
Elder
K. J.
2008
Interannual consistency in fractal snow depth patterns at two Colorado mountain sites
.
Journal of Hydrometeorology
9
(
5
),
977
988
.
doi:10.1175/2008JHM901.1
.
Donald
J. R.
Soulis
E. D.
Kouwen
N.
Pietroniro
A.
1995
A land cover-based snow cover representation for distributed hydrologic-models
.
Water Resources Research
31
(
4
),
995
1009
.
doi:10.1029/94WR02973
.
Driscoll
J. M.
Hay
L. E.
Bock
A. R.
2017
Spatiotemporal variability of snow depletion curves derived from SNODAS for the conterminous United States, 2004–2013
.
JAWRA Journal of the American Water Resources Association
53
(
3
),
655
666
.
doi:10.1111/1752-1688.12520
.
Egli
L.
Jonas
T.
2009
Hysteretic dynamics of seasonal snow depth distribution in the Swiss Alps
.
Geophysical Research Letters
36
(
2
).
doi:10.1029/2008GL035545
.
Egli
L.
Jonas
T.
Grünewald
T.
Schirmer
M.
Burlando
P.
2012
Dynamics of snow ablation in a small Alpine catchment observed by repeated terrestrial laser scans
.
Hydrological Processes
26
(
10
),
1574
1585
.
doi:10.1002/hyp.8244
.
Elder
K.
Dozier
J.
Michaelsen
J.
1991
Snow accumulation and distribution in an Alpine Watershed
.
Water Resources Research
27
(
7
),
1541
1552
.
doi:10.1029/91WR00506
.
Elder
K.
Rosenthal
W.
Davis
R. E.
1998
Estimating the spatial distribution of snow water equivalence in a Montane watershed
.
Hydrological Processes
12
(
10–11
),
1793
1808
.
doi:10.1002/(SICI)1099-1085(199808/09)12:10/11 < 1793::AID-HYP695 > 3.3.CO;2-B
.
Erickson
T. A.
Williams
M. W.
Winstral
A.
2005
Persistence of topographic controls on the spatial distribution of snow in rugged mountain terrain, Colorado, United States
.
Water Resources Research
41
(
4
),
W04014
.
doi:10.1029/2003WR002973
.
Erxleben
J.
Elder
E.
Davis
R.
2002
Comparison of spatial interpolation methods for estimating snow distribution in the Colorado Rocky Mountains
.
Hydrological Processes
16
(
18
),
3627
3649
.
doi:10.1002/hyp.1239
.
Essery
R.
Pomeroy
J.
2004
Implications of spatial distributions of snow mass and melt rate for snow-cover depletion: theoretical considerations
.
Annals of Glaciology
38
(
1
),
261
265
.
doi:10.3189/172756404781815275
.
Fassnacht
S. R.
Dressler
K. A.
Bales
R. C.
2003
Snow water equivalent interpolation for the Colorado River Basin from snow telemetry (SNOTEL) data
.
Water Resources Research
39
(
8
),
1208
.
doi:10.1029/2002WR001512
.
Fassnacht
S. R.
Sexstone
G. A.
Kashipazha
A. H.
Jasinski
M. F.
Kampf
S. K.
Von Thaden
B. C.
López Moreno
J. I.
2016
Deriving snow-cover depletion curves for different spatial scales from remote sensing and snow telemetry data
.
Hydrological Processes
30
(
11
),
1708
1717
.
doi:10.1002/hyp.10730
.
Flanner
M. G.
Zender
C. S.
2006
Linking snowpack microphysics and albedo evolution
.
Journal of Geophysical Research: Atmospheres
111
(
D12
),
D12208
.
doi:10.1029/2005JD006834
.
GDAL/OGR contributors
2015
GDAL/OGR Geospatial Data Abstraction Software Library. Open Source Geospatial Foundation. http://gdal.org
.
Grohmann
C. H.
Smith
M. J.
Riccomini
C.
2011
Multiscale analysis of topographic surface roughness in the Midland Valley, Scotland
. In
IEEE Transactions on Geoscience and Remote Sensing
. Vol.
49
(
1
), pp.
1200
1213
.
doi:10.1109/TGRS.2010.2053546
.
Hannaford
J. F.
Hall
R. L.
Brown
A. J.
1979
Application of snowcovered area to runoff forecasting in the southern Sierra Nevada
. In
Western Snow Conference
,
Sparks, Nevada
.
Harpold
A. A.
Molotch
N. P.
2015
Sensitivity of soil water availability to changing snowmelt timing in the western U.S
.
Geophysical Research Letters
42
,
8011
8020
.
doi:10.1002/2015GL065855
.
Harpold
A.
Dettinger
M.
Rajagopal
S.
2017
Defining snow drought and why it matters
.
EOS – Earth & Space Science News
98
.
doi:10.1029/2017EO068775
.
Helbig
N.
van Herwijnen
A.
Magnusson
J.
Jonas
T.
2015
Fractional snow-covered area parameterization over complex topography
.
Hydrology and Earth System Sciences
19
(
3
),
1339
1351
.
doi:10.5194/hess-19-1339-2015
.
Hiemstra
C. A.
Liston
G. E.
Reiners
W. A.
2002
Snow redistribution by wind and interactions with vegetation at Upper Treeline in the Medicine Bow Mountains, Wyoming, U.S.A
.
Arctic, Antarctic, and Alpine Research
34
(
3
),
262
273
.
doi:10.2307/1552483?ref = no-x-route:13870eb3f86249a29540a7419b9312fd
.
Hijmans
R. J.
2015
raster: Geographic Data Analysis and Modeling
.
Kirnbauer
R.
Bloschl
G.
1994
Wie ähnlich sind Ausaperungsmuster von Jahr zu Jahr (How similar are snow cover patterns from year to year)?
Deutsche Wässerkundliche Mitteilungen
37
(
5/6
),
113
121
.
Kolberg
S.
Gottschalk
L.
2010
Interannual stability of grid cell snow depletion curves as estimated from MODIS images
.
Water Resources Research
46
(
11
),
279
.
doi:10.1029/2008WR007617
.
König
M.
Sturm
M.
1998
Mapping snow distribution in the Alaskan Arctic using aerial photography and topographic relationships
.
Water Resources Research
34
(
12
),
3471
3483
.
doi:10.1029/98WR02514
.
Lehning
M.
Löwe
H.
Ryser
M.
Raderschall
N.
2008
Inhomogeneous precipitation distribution and snow transport in steep terrain
.
Water Resources Research
44
(
7
),
123
.
doi:10.1029/2007WR006545
.
Lehning
M.
Grünewald
T.
Schirmer
M.
2011
Mountain snow distribution governed by an altitudinal gradient and terrain roughness
.
Geophysical Research Letters
38
(
19
),
1
5
.
doi:10.1029/2011gl048927
.
Liston
G. E.
2004
Representing subgrid snow cover heterogeneities in regional and global models
.
Journal of Climate
17
(
6
),
1381
1397
.
doi:10.1175/1520-0442(2004)017 < 1381:RSSCHI > 2.0.CO;2
.
López-Moreno
J. I.
Latron
J. B. P.
2008
Influence of canopy density on snow distribution in a temperate mountain range
.
Hydrological Processes
22
(
1
),
117
126
.
doi:10.1002/hyp.6572
.
López-Moreno
J. I.
Revuelto
J.
Fassnacht
S. R.
Azorín-Molina
C.
Vicente-Serrano
S. M.
Morán-Tejeda
E.
Sexstone
G. A.
2014
Snowpack variability across various spatio-temporal resolutions
.
Hydrological Processes
29
(
6
),
1213
1224
.
doi:10.1002/hyp.10245
.
Luce
C. H.
Tarboton
D. G.
2004
The application of depletion curves for parameterization of subgrid variability of snow
.
Hydrological Processes
18
(
8
),
1409
1422
.
doi:10.1002/hyp.1420
.
Luce
C. H.
Tarboton
D. G.
Cooley
K. R.
1999
Sub-grid parameterization of snow distribution for an energy and mass balance snow cover model
.
Hydrological Processes
13
(
12–13
),
1921
1933
.
doi:10.1002/(SICI)1099-1085(199909)13:12/13 < 1921::AID-HYP867 > 3.0.CO;2-S
.
Marchand
W. D.
Killingtveit
A.
2005
Statistical probability distribution of snow depth at the model sub-grid cell spatial scale
.
Hydrological Processes
19
(
2
),
355
369
.
doi:10.1002/hyp.5543
.
Mizukami
N.
Perica
S.
2008
Spatiotemporal characteristics of snowpack density in the mountainous regions of the Western United States
.
Journal of Hydrometeorology
9
(
6
),
1416
1426
.
doi:10.1175/2008JHM981.1
.
Muggeo
V. M. R.
2003
Estimating regression models with unknown break-points
.
Statistics in Medicine
22
(
19
),
3055
3071
.
doi:10.1002/sim.1545
.
Muggeo
V. M. R.
2020
segmented: Regression Models with Break-Points/Change-Points Estimation. http://cran.r-project.org/package=segmented
.
Neteler
M.
Bowman
M. H.
Landa
M.
Metz
M.
2012
GRASS GIS: a multi-purpose Open Source GIS
.
Environmental Modelling & Software
31
,
124
130
.
doi:10.1016/j.envsoft.2011.11.014
.
Niu
G.-Y.
Yang
Z.-L.
2007
An observation-based formulation of snow cover fraction and its evaluation over large North American river basins
.
Journal of Geophysical Research
112
(
D21
).
doi:10.1029/2007JD008674
.
Nolan
M.
Larsen
C. F.
Sturm
M.
2015
Mapping snow-depth from manned-aircraft on landscape scales at centimeter resolution using structure-from-motion photogrammetry
.
The Cryosphere
9
(
1
),
1445
1463
.
doi:10.5194/tc-9-1445-2015
.
NPS
2016
Plants – Yosemite National Park
.
Available from: https://www.nps.gov/yose/learn/nature/plants.htm (accessed 5 October 2016).
Painter
T. H.
Berisford
D. F.
Boardman
J. W.
Bormann
K. J.
Deems
J. S.
Gehrke
F.
Hedrick
A.
Joyce
M.
Laidlaw
R.
Marks
D.
Mattmann
C.
McGurk
B.
Ramirez
P.
Richardson
M.
McKenzie Skiles
S.
Seidel
F. C.
Winstralf
A.
2016
The Airborne Snow Observatory: fusion of scanning lidar, imaging spectrometer, and physically-based modeling for mapping snow water equivalent and snow albedo
.
Remote Sensing of Environment
184
,
139
152
.
doi:10.1016/j.rse.2016.06.018
.
Parker
K. L.
Robbins
C. T.
Hanley
T. A.
1984
Energy expenditures for locomotion by Mule Deer and Elk
.
The Journal of Wildlife Management
48
(
2
),
474
416
.
doi:10.2307/3801180
.
Potts
H. L.
1937
Snow surveys and runoff forecasting from photographs
. In
Proceedings of the 5th Annual Western Interstate Snow Survey Conference
.
R Core Team
2015
R: A Language and Environment for Statistical Computing
.
Vienna
,
Austria
.
Revuelto
J.
López-Moreno
J. I.
Azorín-Molina
C.
Vicente-Serrano
S. M.
2014
Topographic control of snowpack distribution in a small catchment in the central Spanish Pyrenees: intra- and inter-annual persistence
.
The Cryosphere
8
(
5
),
1989
2006
.
doi:10.5194/tc-8-1989-2014
.
Schirmer
M.
Lehning
M.
2011
Persistence in intra-annual snow depth distribution: 2. Fractal analysis of snow depth development
.
Water Resources Research
47
(
9
),
W09517
.
doi:10.1029/2010WR009429
.
Schirmer
M.
Wirz
V.
Clifton
A.
Lehning
M.
2011
Persistence in intra-annual snow depth distribution: 1. Measurements and topographic control
.
Water Resources Research
47
(
9
),
W09516
.
doi:10.1029/2010WR009426
.
Sturm
M.
Wagner
A. M.
2010
Using repeated patterns in snow distribution modeling: an Arctic example
.
Water Resources Research
46
(
12
).
doi:10.1029/2010WR009434
.
Sturm
M.
Holmgren
J.
Liston
G. E.
1995
A seasonal snow cover classification system for local to global applications
.
Journal of Climate
8
,
1261
1283
.
Swenson
S. C.
Lawrence
D. M.
2012
A new fractional snow-covered area parameterization for the community land model and its effect on the surface energy balance
.
Journal of Geophysical Research
117
(
D21
).
doi:10.1029/2012JD018178
.
Therneau
T.
Atkinson
B.
Ripley
B.
2015
rpart: Recursive Partitioning and Regression Trees
. .
Trujillo
E.
Ramirez
J. A.
Elder
K. J.
2007
Topographic, meteorologic, and canopy controls on the scaling characteristics of the spatial distribution of snow depth fields
.
Water Resources Research
43
(
7
).
doi:10.1029/2006WR005317
.
Trujillo
E.
Ramirez
J. A.
Elder
K. J.
2009
Scaling properties and spatial organization of snow depth fields in sub-alpine forest and alpine tundra
.
Hydrological Processes
23
(
11
),
1575
1590
.
doi:10.1002/hyp.7270
.
Trujillo
E.
Molotch
N. P.
Goulden
M. L.
Kelly
A. E.
Bales
R. C.
2012
Elevation-dependent influence of snow accumulation on forest greening
.
Nature Geoscience
5
(
10
),
705
709
.
doi:10.1038/ngeo1571
.
U.S. Geological Survey
2016
National Hydrography Dataset
.
NHDPlus High Resolution
.
Vander Jagt
B.
Durand
M. T.
Margulis
S. A.
Molotch
N. P.
Kim
E. J.
2013
The effect of spatial variability on the sensitivity of passive microwave measurements to snow water equivalent
.
Remote Sensing of Environment
136
,
163
179
.
doi: 10.1016/j.rse.2013.05.002
.
Wickham
H.
Averick
M.
Bryan
J.
Chang
W.
D'Agostino McGowan
L.
François
R.
Grolemund
G.
Hayes
A.
Henry
L.
Hester
J.
Kuhn
M.
Pedersen
T. L.
Miller
E.
Milton Bache
S.
Müller
K.
Ooms
J.
Robinson
D.
Paige Seidel
D.
Spinu
V.
Takahashi
K.
Vaughan
D.
Wilke
C.
Woo
K.
Yutani
H.
2019
Welcome to the Tidyverse
.
JOSS
4
(
43
),
1686
1686
.
doi:10.21105/joss.01686.
Winstral
A.
Elder
K.
Davis
R. E.
2002
Spatial snow modeling of wind-redistributed snow using terrain-based parameters
.
Journal of Hydrometeorology
3
(
5
),
524
538
.
doi:10.1175/1525-7541(2002)003 < 0524:SSMOWR > 2.0.CO;2
.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY 4.0), which permits copying, adaptation and redistribution, provided the original work is properly cited (http://creativecommons.org/licenses/by/4.0/).

Supplementary data