Analysis of topographic controls on depletion curves derived from airborne lidar snow depth data

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 ﬁ rst 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 ( r 2 ¼ 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 ﬂ ux calculations across many earth system models.


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 influ- 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. ).
Multiple studies have investigated sub-grid depletion curves of SWE based on terrain characteristics (Donald et Helbig et al. ). 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. ).
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 ).
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?

Site description
This study was conducted on data collected by the ASO in the Tuolumne River basin in the Sierra Nevada Mountains

Data sources
We use the ASO dataset which is derived from a paired scanning lidar and imaging spectrometer (Painter et al.

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 ). 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 ). 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 , ). 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).
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 ) 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  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. (1) where ∂ x z and ∂ y z are the first partial derivatives of elevation, i.e. the orthogonal slope components (Neteler et al. ).
At the coarser scale, we quantify the standard deviation of elevation (stdelev), the standard deviation of slope 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 have a strong individual relationship with the depletion slope ( Table 2).
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   Uncorrected Proof accumulation and depletion dynamics above and below 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.
We observe that increased terrain roughness at sub-grid scalesrepresented by the stdmaxcurv, stdelev, and Δh/Δw variablescorrelates 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. ). 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 redistri-

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 r 2 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 r 2 ¼ 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.