# Evaluating the Performance of Sentinel-1 SAR Derived Snow Depth Retrievals Over the Extratropical Andes Cordillera

N Bulovic, F Johnson, H Lievens, T E Shaw, J Mcphee, S Gascoin, M

Demuzere, N Mcintyre

► To cite this version:

N Bulovic, F Johnson, H Lievens, T E Shaw, J Mcphee, et al.. Evaluating the Performance of Sentinel-1 SAR Derived Snow Depth Retrievals Over the Extratropical Andes Cordillera. Water Resources Research, 2025, 61 (2), pp.e2024WR037766. 10.1029/2024wr037766. hal-05345111

HAL Id: hal-05345111

https://hal.science/hal-05345111v1

Submitted on 3 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Water Resources Research®

## RESEARCH ARTICLE

10.1029/2024WR037766

**Key Points:**

* We derived 4 years of snow depth maps at 1 km² over the Andes mountains with Sentinel-1 SAR and other globally and freely available data

* Snow depth retrievals perform best over sites with negligible evergreen forest cover and in deeper snow

* Snow depth is overestimated for Tundra and Prairie snow class sites, whilst Maritime is unbiased

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

N. Bulovic,
n.bulovic@uq.edu.au

**Citation:**

Bulovic, N., Johnson, F., Lievens, H., Shaw, T. E., McPhee, J., Gascoin, S., et al. (2025). Evaluating the performance of Sentinel-1 SAR derived snow depth retrievals over the extratropical Andes Cordillera. *Water Resources Research*, 61, e2024WR037766. https://doi.org/10.1029/2024WR037766

Received 17 APR 2024
Accepted 18 JAN 2025

**Author Contributions:**

**Conceptualization:** N. Bulovic, F. Johnson, N. McIntyre
**Data curation:** T. E. Shaw, J. McPhee, S. Gascoin
**Formal analysis:** N. Bulovic
**Methodology:** N. Bulovic, H. Lievens
**Software:** H. Lievens, M. Demuzere
**Supervision:** F. Johnson, N. McIntyre
**Visualization:** N. Bulovic
**Writing – original draft:** N. Bulovic
**Writing – review & editing:** F. Johnson, H. Lievens, T. E. Shaw, J. McPhee, S. Gascoin, M. Demuzere, N. McIntyre

## Evaluating the Performance of Sentinel-1 SAR Derived Snow Depth Retrievals Over the Extratropical Andes Cordillera

<figure><img src="image_2.png" /></figure>

N. Bulovic¹ ID, F. Johnson² ID, H. Lievens³,⁴ ID, T. E. Shaw⁵,⁶ ID, J. McPhee⁷,⁸ ID, S. Gascoin⁹ ID, M. Demuzere³,¹⁰ ID, and N. McIntyre¹,¹¹ ID

¹Centre for Water in the Minerals Industry, Sustainable Minerals Institute, The University of Queensland, Brisbane, Qld, Australia, ²Water Research Centre, School of Civil and Environmental Engineering, University of New South Wales, Sydney, NSW, Australia, ³Department of Environment, Ghent University, Ghent, Belgium, ⁴Department of Earth and Environmental Sciences, KU Leuven, Heverlee, Belgium, ⁵Swiss Federal Institute for Forest, Snow and Landscape Research (WSL), Birmensdorf, Switzerland, ⁶Institute of Science and Technology Austria (ISTA), Klosterneuburg, Austria, ⁷Advanced Mining Technology Center, Universidad de Chile, Santiago, Chile, ⁸Department of Civil Engineering, Universidad de Chile, Santiago, Chile, ⁹CESBIO, CNES/CNRS/INRAE/IRD/UPS, Université de Toulouse, Toulouse, France, ¹⁰B-Kode VOF, Ghent, Belgium, ¹¹School of Civil Engineering, The University of Queensland, Brisbane, Qld, Australia

**Abstract** Monitoring and estimating mountain snowpack mass over regional scales is still a challenge because of the inadequacy of observational networks in capturing spatiotemporal variability, and limitations in remotely sensed retrievals. Recent work using C-band synthetic aperture radar (SAR) backscatter data from the Sentinel-1 satellite mission has shown good promise for tracking mountain snow depth over specific northern hemisphere ranges, although the broader potential is still unknown. Here, we extend the new Sentinel-1 based modeling framework beyond the northern hemisphere by only utilizing globally available input data, and evaluate different model parametrization and model performance over the Chilean and Argentine Andes mountains, which contain the largest mountain snowpack in the southern hemisphere. The accuracy of Sentinel-1 snow depth estimates is evaluated against an extensive in situ network available for the region. Satellite-retrieved snow depth is found to have poorer performance across the Andes than observed for northern hemisphere mountain ranges because of greater sensitivity to evergreen forest cover and shallower snowpacks. The algorithm does offer some skill but performance is variable and site-dependent. Algorithm performance is best over regions with limited evergreen forest cover (<15%) and snow depths greater than 0.75 m, although the retrievals over-estimate snow depth across most sites. Systemic errors for specific snow classes and across different snow depths are shown, highlighting specific areas in need of further investigation and development.

## 1. Introduction

Snow is a critical freshwater resource, particularly in mountainous and mid- to high-latitude regions where it is a major component of the water budget (Barnett et al., 2005). Mountain regions commonly behave as natural “water towers” as they accumulate and store snow reserves over winter, which melt and supply freshwater during the warmer seasons (Vivioli et al., 2007). Despite mountain snowpacks accounting for less than a third of the global seasonal snow covered area (Wrzesien et al., 2019), river runoff originating from mountain snowmelt provides water to 1.9 billion people downstream in their respective watersheds (Immerzeel et al., 2020).

The cryosphere of the Chilean and Argentine Andes Cordillera is the largest in the southern hemisphere outside Antarctica (Hammond et al., 2018), and is an important socio-economic and natural resource sustaining eco-systems, agriculture, hydropower, and a globally important mining industry (SERNAGEOMIN, 2021). In the wet southern Patagonian Andes, snowmelt accounts for 26% of annual runoff in Baker River (Krogh et al., 2015); while in the central semi-arid Andes, it contributes between 55% and 95% of annual runoff, depending on annual precipitation, catchment landscape position, and degree of glaciarisation (Ayala et al., 2020; Burger et al., 2019; Masiokas et al., 2006). From 2010 until mid-2023, the central Andes have been affected by a megadrought, characterized by persistent precipitation deficits (Garreaud et al., 2020), depleted snowpacks and substantially reduced river flows (Garreaud et al., 2017). The drought has resulted in reduced water availability and water rights conflicts between users as demand exceeds supply (Rivera et al., 2016)—a scenario that will likely be further exacerbated by climate change (Vicuña et al., 2011).

© 2025. The Author(s).
This is an open access article under the
terms of the Creative Commons
Attribution License, which permits use,
distribution and reproduction in any
medium, provided the original work is
properly cited.
Good medium- to long-term planning and water resource management can alleviate water availability challenges and conflicts, but ideally requires frequent (~7 days), continuous and high resolution (~1 km²) estimates of mountain range snow mass (Bormann et al. 2018), that can adequately capture the high spatial variability in snowpack related to topography. Snow mass is quantified with snow water equivalent (SWE). However, as SWE measurements are more difficult to obtain, snow mass is often approximated by snow depth (SD), assuming that snow density has a smaller influence. One application of high resolution SWE or SD data is through integration into hydrological models to improve the accuracy of seasonal streamflow forecasts (e.g., Li et al., 2019; Shaw, Caro, et al., 2020). Long-term historical data of snow mass can also be used to identify drivers-of (e.g., Cortés & Margulis, 2017) and trends-in (e.g., Pulliainen et al., 2020) snow mass over decadal time frames, which may then be used to assess impacts on future water availability and guide scenario planning. However, estimating snow mass at catchment to regional scales is inherently difficult (Dozier et al., 2016), partially due to mountain SWE heterogeneity caused by variable climatic and topographic drivers (Mendoza et al., 2020). This challenge is compounded as in situ observations are traditionally collected by sparse and unevenly distributed monitoring networks which also struggle to represent local scale heterogeneity of the snowpack. For example, the SNOTEL network is composed of around 900 automated hydroclimatic data collection sites (USDA Natural Resources Conservation Service, 2022), averaging one site per 2,500 km². In comparison, Masiokas et al. (2006) utilized a network averaging one station per 20,000 km² to investigate central Andean snowpack variability.

In an attempt to address this challenge, numerous remote sensing-based approaches have been used to derive spatially continuous estimates of SWE or SD ideally supplemented with snow cover extents. Terrestrial and airborne applications of lidar and photogrammetry have derived SD at high spatial resolutions (~1–50 m) at hillslope (e.g., Mendoza et al., 2020) to catchment scales (e.g., Painter et al., 2016) but acquisitions remain sporadic. In theory, satellite-based platforms such as Pléiades can largely overcome the spatial sampling issues during cloud-free conditions, but the high cost of purchasing data and narrow swath does not allow frequent coverage over large areas (e.g., Deschamps-Berger et al., 2020; Marti et al., 2016). Freely available, widespread and frequent data has predominantly come from passive microwave (PMW) and optical spaceborne missions. PMW data are used to produce quasi-global, daily estimates of SWE (e.g., Pulliainen et al., 2020). However, PMW data applicability over mountainous regions is limited due to the coarse spatial resolution of the sensors (~25 km), and the poor sensor performance in deeper snow (SWE > ~0.25 m.w.e.) commonly associated with mountain environments (Vuyovich et al., 2014). Therefore, mountain SWE is masked in PMW-based products such as GlobSnow (Takala et al., 2011). Optical imagery observes whether snow is present (rather than SD or SWE), and is commonly coupled with in situ data (Yang et al., 2022) or models (Cornwell et al., 2016; Cortés et al., 2016) to indirectly estimate SWE. However, this is not a trivial task for various reasons, including that the relationship between snow covered area and volume is not linear (Shaw, Deschamps-Berger, et al., 2020). Despite the valuable progress and insights gained from these approaches globally and across the Andes (as summarized in Table 1), there still remains a gap in frequent and widespread SD monitoring of Andean mountain snowpacks.

Similar to PMW, radar (active microwave) remote sensing is sensitive to snowpack properties and has all-weather and all-day observing capability, but also shows promise in mapping SD or SWE at spatial resolutions relevant to mountain environments. Radar-based approaches for snow mapping tend to rely either on the “signal phase” or the “microwave scattering” components of the collected data (Tsang et al., 2022). “Signal phase” approaches are related to the travel time of the signal (or path length) which can be correlated to snow properties. Differential SAR interferometry (Tarricone et al., 2023) and co-polar phase difference (Voglimacci-Stephanopoli et al., 2022) are two such approaches used for SWE or SD mapping (see Tsang et al. (2022) for a detailed review). “Microwave scattering” approaches rely on the amplitude of the backscatter response. Most scattering-based snow research has focused on the higher frequency X- to Ku-bands (~8–18 GHz) because of the strong backscatter response to dry snowpack volume (King et al., 2018; Rott et al., 2010; Yueh et al., 2009). Accordingly, various methods to estimate SD or SWE have been developed based on a combination of radiative transfer modeling, empirical and/ or machine learning approaches (e.g., Patil et al., 2020; Zhu et al., 2021). Despite the progress in methods that use X- to Ku-bands, at present the only global and freely available satellite synthetic aperture radar (SAR) observations are in C-band (~5 GHz). After early work found C-band observations to have limited sensitivity to dry snowpack volume (Bernier & Fortin, 1998; Shi & Dozier, 2000), there has been little interest in its application for SD or SWE mapping.

Recently, Lievens et al. (2019) used C-band SAR observations from the European Space Agency and Copernicus Sentinel-1 mission to estimate dry SD over mountainous terrain at a 1 km² resolution on a weekly basis. The
<table><caption>**Table 1**
*Studies Deriving Spatial Estimates of SD or SWE Across the Chilean and Argentine Andes Cordillera*</caption><thead><tr><th>Reference</th><th>Domain</th><th>Extent</th><th>Period</th><th>Approach<sup>a</sup> (M, R, M&R)</th><th>Data</th><th>Res. (m)</th><th>Appl.<sup>b</sup> (N, R)</th></tr></thead><tbody><tr><td>Gascoin et al. (2013)</td><td>29°S, 70°W</td><td>Headwaters</td><td>Winter season 2008</td><td>M—SnowModel</td><td>In situ stations</td><td>90</td><td>—</td></tr><tr><td>Cornwell et al. (2016)</td><td>27°– 38°S</td><td>Region</td><td>2001–2014</td><td>M&R—Couple snow cover data with energy melt model</td><td>MODIS fSCA; in situ stations</td><td>500</td><td>R</td></tr><tr><td>Cortés et al. (2016)</td><td>33°S</td><td>Catchment</td><td>1985–2015</td><td>M&R—Assimilate snow cover data into LSM</td><td>LSM (SSiB3); MERRA; Landsat fSCA</td><td>90</td><td>R</td></tr><tr><td>Cortés and Margulis (2017)</td><td>27°– 37°S</td><td>Region</td><td>1984–2015</td><td>M&R—Assimilate snow cover data into LSM</td><td>LSM (SSiB3); MERRA; Landsat fSCA</td><td>180</td><td>R</td></tr><tr><td>Mernild et al. (2017)</td><td>Full Andes</td><td>Continent</td><td>1979–2014</td><td>M—SnowModel</td><td>MERRA</td><td>4000</td><td>—</td></tr><tr><td>Shaw, Gascoin, et al. (2020)</td><td>34°– 70°S</td><td>Catchment</td><td>September 2017</td><td>R—Photogrammetry</td><td>Pléiades</td><td>4</td><td>N</td></tr><tr><td>Shaw, Deschamps-Berger, et al. (2020)</td><td>34°– 70°S</td><td>Catchment</td><td>September 2017, 2019</td><td>R—Photogrammetry</td><td>Pléiades</td><td>4</td><td>N</td></tr><tr><td>Mendoza et al. (2020)</td><td>31°– 37°S</td><td>Hillslope</td><td>Winter 2018</td><td>R—Lidar</td><td>Ground-based</td><td>1</td><td>—</td></tr></tbody></table>

retrievals were based on the demonstrated response of cross-polarised backscatter to deeper snow, such as found in mountain environments. Comparison against in situ measurements over northern hemisphere mountains showed the algorithm's estimates to be unbiased on average, and to have realistic spatial and temporal patterns, although performance varied between mountain ranges (Lievens et al., 2019). More recently, an updated version of the algorithm developed and optimized for the Alps was assessed over a 2-year period, and found strongest performance at intermediate elevations characterized by lower forest cover and deeper snowpacks (>1 m) (Lievens et al., 2022). Despite the promise of this new approach, the evaluations have been limited to the northern hemisphere, and particularly focused on the Alps and other well-monitored mountain ranges in the USA. There still remain questions about (a) the applicability of the algorithm across other regions with different snow classes, surface attributes and land cover types, which are known to influence radar observations in some circumstances (Feng et al., 2021; King et al., 2018; Yueh et al., 2009), and (b) the drivers of performance more generally.

Accordingly, the aim of this study is to assess the skill of Sentinel-1 snow depth retrievals in the Andes Cordillera over a 4-year period (2017–2020) at a 1 km² resolution. First, temporal correlations between the backscatter data (used in the snow depth retrievals) and an extensive in situ snow depth monitoring network are assessed. The effects of long-term static (location, topography, land cover, climate, and snow class) and dynamic (snow persistence, snow depth range, and retrieval frequency) attributes on temporal correlations are considered to identify factors influencing the algorithm performance. Importantly, the in situ monitoring network was specifically compiled to maximize observational coverage as opposed to focusing only on the highest quality stations or stations with data for the full 4 year evaluation period. Thereafter, some of the retrieval input data are updated so that all ancillary data sets are globally available, rather than being northern hemisphere (or regionally) specific. Lastly, we evaluate the performance of the original algorithm (Lievens et al., 2019), and discuss the drivers of performance and the implications for broader scale applications.

## 2. Study Area

The SD algorithm is evaluated over the central-to-southern Andes of Chile and Argentina, running in a north-south direction (Figure 1). The tropical Andes are out of scope of this study because in that region snow is mostly limited to the highest mountain peaks and glaciers where there are minimal snowpack observations.

<footnote><sup>a</sup>Approach abbreviations: M—snow model; R—remote sensing; M&R—remote sensing integrated into snow model. <sup>b</sup>Application abbreviations: N—near-real time; R—retrospective.</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Study area incorporating the Chilean and Argentine Andes: (a) Elevation (m a.s.l.) from SRTM DEM; (b) Seasonal snow classes (Sturm & Liston, 2021): ice (I), Tundra (T), Boreal forest (BF), Montane forest (MF), Prairie (P), and Maritime (M); (c) Percentage days of snow cover from 2017 to 2020 calculated from VNP10A1F satellite snow cover product; (d) % Forest cover (Tuanmu & Jetz, 2014). Location of in situ monitoring stations shown by white dots, inset locations are shown in subplot (a).</figcaption></figure>

The study area encompasses multiple climatic zones and a large range of elevations, topographic relief, vegetation types and cover. Snow across the region is classified into one of six snow classes (related to average snow depth, stratigraphy, temperature and snow-grain characteristics) (Sturm & Liston, 2021), or by the snow cover persistence (Figures 1b and 1c). In the desert Andes (<23°S), snow cover is intermittent and limited to the highest elevations and peaks. The largest proportion of seasonal and permanent snow occurs between 28° and 37°S (Cornwell et al., 2016; Cortés & Margulis, 2017) due to the winter-dominant precipitation, high elevations (Figure 1a) and low temperatures. Around two-thirds of the snow mass occurs on the western Andean slopes (Cortés & Margulis, 2017), and the region is dominated by two snow classes—Tundra (higher elevation, colder) and Prairie (intermediate elevation, warmer)—both of which typically have shallow, wind-affected snowpacks. South of 40°S, precipitation occurs throughout the year, mountain peaks drop in elevation and so does the snow line. The Maritime snow class, associated with deep, wet snow, predominantly occurs here, although detailed mapping of snow cover is limited due to the very high cloud cover (Saavedra et al., 2017). Vegetation patterns
broadly reflect the climatic variability (Figure 1d). Forest cover is limited in the arid regions and at high elevations; in contrast it is widespread across the lower elevation, wet, southern Andes (Huerta et al., 2019).

This study focuses on a 4 year period from 2017 to 2020. The extratropical Andes were under unprecedented drought conditions from 2010 to 2023 (Garreaud et al., 2017), with 2019 being the second driest year on record (Shaw, Deschamps-Berger, et al., 2020). Regionally, reduced snow cover and higher snowline elevations have been observed (Saavedra et al., 2018), whilst 2019 snow depths and volumes were substantially diminished in comparison to 2017 over a high, headwater catchment (Shaw, Deschamps-Berger, et al., 2020).

# 3. Data and Methods

## 3.1. Snow Depth Modeling

Synthetic aperture radar (SAR) is an active (rather than passive) form of remote sensing, where a sensor emits microwave energy of a specific frequency and polarization toward a surface and measures the return signal, commonly referred to as backscatter. Backscatter from a dry snow covered surface is predominantly the result of a combination of (a) scattering from the ground surface attenuated by the snowpack, and (b) volume scattering from the snowpack (Rott et al., 2010). For wet snow, radiative transfer is dominated by strong absorption of the signal (because of the high dielectric losses within water) (Mätzler, 1987) and therefore backscatter is substantially reduced.

The relative importance of each scattering mechanism for dry snow depends on SAR sensor frequency and polarization. For Sentinel-1 (C-band frequency of 5.4 GHz), dual polarimetric backscatter data is acquired in vertical-vertical (VV) and vertical-horizontal (VH) polarisation, also known as co- and cross-polarised backscatter, respectively. The empirical framework of Lievens et al. (2019) for snow depth retrieval using Sentinel-1 backscatter observations is based on a number of assumptions about the scattering mechanisms of C-band SAR over dry snow covered terrain. Co-polarised backscatter ($\sigma_{vv}^0$) shows little sensitivity to dry snow and the signal is predominantly due to ground surface scattering (depending on ground surface roughness, temperature and moisture content). In comparison, cross-polarised backscatter ($\sigma_{vh}^0$) shows increased sensitivity to snow volume scattering for two primary reasons (Lievens et al., 2022). First, as the radar signal travels through a dry snowpack it depolarizes and thus increases scattering in cross-polarisation. The depolarization occurs for multiple reasons, including scattering on anisotropic snow crystals (and clusters of crystals) within a dense snowpack and multiple scattering on snow layer and snow-ground interfaces. Consequently, an increase in snow depth is expected to result in an increase in cross-polarised backscatter due to enhanced opportunity for depolarization. Second, the strength of ground surface scattering in cross-polarisation is less than that in co-polarisation. Changes in ground surface properties over time are expected to have a similar relative influence in both co- and cross-polarised backscatter, and so the scattering contribution of snow attributes can be better extracted by combining observations in both co- and cross-polarisation ($\sigma_{vv}^0$ and $\sigma_{vh}^0$). Here, only certain snow attributes (particularly snow volume) are assumed to be important in influencing the backscatter signal, while others—such as snow crystal shape, snow anisotropy and snow metamorphism—are not explicitly considered. This is a practical but necessary limitation for large-scale regional studies, and the implications are discussed in Section 5.

Figure 2 provides a summarized depiction of the snow depth retrieval algorithm (Lievens et al., 2019), and the necessary input data and pre-processing requirements. The algorithm relates snow depth variability over time to changes in successive backscatter observations, utilizing three main information sources: (a) Pre-processed and normalized Sentinel-1 data, used to produce a merged single time series of Sentinel-1 co- and cross-polarised backscatter data ($\sigma^0$) at each pixel (Figures 2a and 2b). Mean orbital bias correction was used to normalize Sentinel-1 $\sigma^0$ values across different orbits to account for viewing geometry/incidence angle differences (as detailed in Text S1.1 of Supporting Information S1). This $\sigma^0$ normalization removes biases in the long-term mean backscatter between different orbits, but does not resolve differences that occur due to variations in footprint size (Small, 2011) or different scattering mechanisms at different angles. The implications on results are further discussed in Section 5; (b) Snow cover data, derived from the optically based VNP10A1F data set, used to define snowy periods (Figure 2c); and (c) Land-cover data for correcting or masking snow depth estimates (Figure 2d). For brevity, only the snow depth retrievals are detailed below whilst data background and pre-processing information is provided in Texts S1.1 to S1.3 of Supporting Information S1. The algorithm was applied to derive snow depths for four austral winters (2017–2020).
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Flowchart of data processing steps for the snow depth retrievals. (a) Sentinel-1 data acquisition and standard processing; (b) Backscatter normalization for incidence angle; (c) Snow cover gap filling and standardization; (d) Ancillary data acquisition and standardization. Processing steps for parts (a) and (b) are detailed in Text S1.1 in Supporting Information S1, part (c) is detailed in Text S1.2 in Supporting Information S1, and part (d) in Text S1.3 in Supporting Information S1.</figcaption></figure>

Using the daily normalized co- and cross-polarised backscatter ($\sigma_{vv}^0$ and $\sigma_{vh}^0$), the first step was to quantify the cross-polarisation ratio ($\sigma_{cr}^0$) defined as a weighted difference in dB scale as follows:

$$ \sigma_{cr}^0 = A \cdot \sigma_{vv}^0 - \sigma_{vv}^0 \quad (1) $$

where $\sigma_{cr}^0$, $\sigma_{vv}^0$ and $\sigma_{vh}^0$ are in dB, and A is a dimensionless constant used to weight $\sigma_{vv}^0$ values. The $\sigma_{cr}^0$ data were temporally smoothed to reduce the impact of outliers from propagating in the retrievals and to minimize the impact of observational noise intrinsic to backscatter data. The temporally smoothed $\sigma_{cr}^0$ was equal to the mean of $\sigma_{cr}^0$ at the current, prior and posterior timesteps $t$, $t_{pri}$, and $t_{post}$, respectively, where $\sigma_{cr}^0(t)$ was obtained directly from Equation 1. The posterior estimate $\sigma_{cr}^0(t_{post})$ was the weighted average of the cross-polarisation ratios over the subsequent 12 days (timesteps $t + 1$ to $t + 12$):

$$ \sigma_{cr}^0(t_{post}) = \frac{1}{\sum_{j=t+1}^{t+12} w_j} \cdot \sum_{j=t+1}^{t+12} w_j \cdot \sigma_{cr}^0(j) \quad (2) $$

where the weights $w_j$ are the inverse distance in time from $t$ such that $w_j = \frac{1}{j-t}$. This normally involved two to three observations of the cross-polarisation ratio because most of the region had an overpass frequency of 4–5 days. No prior estimate was available for $t = 1$, but thereafter $\sigma_{cr}^0(t_{pri})$ would be equal to the final, smoothed $\sigma_{cr}^0$ derived for the previous timestep where $t$ was equal to $t_{pri}$. Next, the change in the cross-polarisation ratio ($\Delta\sigma_{cr}^0$) at timestep $t$ and point $i$ was calculated between successive timesteps $t$ and $t_{pri}$:
$$ \Delta\sigma_{cr}^{0}(i,t) = \sigma_{cr}^{0}(i,t) - \sigma_{cr}^{0}(i,t_{pri}) \quad (3) $$

From here the snow index (SI) was estimated at timestep t and point i:

$$ SI(i,t) = SI(i,t_{pri}) + \Delta\sigma_{cr}^{0}(i,t) \quad (4) $$

SI was set to 0 in instances where there was no snow cover (SC(i,t) = 0) (as depicted in Figure 2c) or if SI became negative. Finally, SI was linearly rescaled and converted into snow depth (SD) in meters, although SD values were marked NA (rather than 0) for time steps with negative SI to limit inconsistencies with the snow cover data:

$$ SD(i,t) = \frac{C}{1 - B.FC(i)} \cdot SI(i,t), \quad (5) $$

where FC is evergreen forest cover fraction (as shown in Figure 2d and detailed in Text S1.3 in Supporting Information S1) with values ranging between 0 and 1. FC was used to rescale the SD retrievals to correct for signal attenuation due to forest cover. B (dimensionless) and C (m/dB) are constant scaling parameters. The same parameter values derived by Lievens et al. (2019) were used here, namely: A = 1, B = 0.6, and C = 1.1 m/dB. These parameters were derived with data from numerous northern hemisphere mountain ranges, and thus should be broadly representative of a quasi-global snowpack, unlike parameter sets specifically derived for an individual mountain range (Lievens et al., 2022) or region (Hoppinen et al., 2024). The potential implications of using the Lievens et al. (2019) parameter set, as opposed to deriving one for the Andes, is further discussed in Section 5.3.

Note that snow depth retrievals during wet snow conditions were masked with a wet snow detection algorithm. This is because (a) backscatter from wet snow is strongly decreased as the radar signal is dominated by strong absorption (Marin et al., 2020; Nagler et al., 2016), and (b) the snow depth retrieval algorithm was specifically developed and optimized for dry snow conditions (Lievens et al., 2019). The first step in the wet snow detection algorithm was to calculate the change in mean $\sigma_{v,h}^{0}$ at point i for 12 days prior to and 12 days after (and including) timestep t:

$$ \Delta\sigma_{v,h}^{0}(i,t) = \frac{1}{n_j} \sum_{j=t-1}^{t-12} \sigma_{v,h}^{0}(i,j) - \frac{1}{n_k} \sum_{k=t}^{t+11} \sigma_{v,h}^{0}(i,k) \quad (6) $$

where $n_j$ and $n_k$ are the number of observations in the 12 days prior to and post timestep t. Timestep t was flagged as wet if snow cover was present (SC(i,t) = 1) and $\Delta\sigma_{v,h}^{0}$ was greater than a threshold value set to 1 dB. The applied threshold was lower than in other studies (e.g., 2 dB; Nagler et al., 2016) because wet snow was determined here using a moving-window as opposed to differences in $\sigma_{v,h}^{0}$ relative to a reference snow-free period. If point i was classified as wet, the wet snow classification would persist until snow was no longer present at point i. To prevent points from being classified as wet snow early in the season when snow cover can be highly dynamic, the wet snow detection algorithm was only enabled from September onwards to capture the dominant spring snowmelt period.

## 3.2. In Situ Snow Measurements

### 3.2.1. Data and Data Pre-Processing

The Andes mountains have historically suffered and continue to suffer from a lack of continuous snow data (Mendoza et al., 2020), with monitoring points sparsely and unevenly distributed across the region (Masiokas et al., 2020). Hence, a concerted effort was made to collate all potential in situ data sets during the analysis period to maximize the amount of validation data across both the Chilean and Argentine Andes. Continuous in situ point measurements of snow depth (SD) and (or) snow water equivalent (SWE) were obtained from a variety of sources, including the Chilean Dirección General de Aguas (DGA), Argentine government agency data stored on the Sistema Nacional de Información Hídrica (SNIH, https://snih.hidricosargentina.gob.ar/), Universidad de
Chile, and Chilean-owned mining company Codelco. Similar data supplied by mining companies have been used previously for remote sensing product evaluations in the Andes (e.g., Bulovic et al., 2020).

The initial step in the quality control process was to identify stations with high quality data, and remove those stations or individual snow years that were of inadequate quality from the record. First, time series data from each station were plotted and the quality was visually inspected for missing data and obvious issues with sensor performance. These issues included unrealistic snowpack evolution over time, a non-static reference measurement point elevation, and large inconsistencies between (a) SD and SWE measurements, and (b) automatic and manual measurements (where both were available). Depending on the extent of the data quality issues, either individual snow years or the complete station were removed. Next the accuracy of station coordinates was manually assessed by cross-checking the timing of snow cover with PlanetScope and Sentinel-2 imagery, using the Planet online image archives (Planet Team, 2017). From an initial 42 stations, a total of 12 were completely removed (11 had poor quality data, one had questionable coordinates), while five stations lost at least 1 year of data due to quality issues. A total of 30 in situ stations (equivalent to 77 water years) were determined to be of suitable quality: 10 from Argentina and 20 from Chile (14 from DGA, three each from Universidad de Chile and Codelco). In situ monitoring station locations are shown in Figure 1, while their site, data, and snow class attributes are detailed in Table 2. The stations are distributed between 29 and 44°S, the mean elevation is 3,026 m a. s.l. (ranging between 1,165. and 4,459 m a.s.l.), and they cover the seasonal snow classes Tundra, Prairie, Maritime, and Montane Forest (Sturm & Liston, 2021).

Next, all data from the 30 in situ stations were quality-controlled and the timestep was standardized to produce a final data set of daily snow depth. As the data were collected from multiple sources, the quality control steps had to be customized by source to account for the differences in data attributes (e.g., temporal resolution, variable measured, etc.) and the level of post-processing applied by the data owner. Here, some of the main steps, attributes and differences by data source are discussed, while the exact steps applied by data source are detailed in Supporting Information S1 (Table S4). SNIH provided daily post-processed SWE data for the snow season (1 May–1 December), and as such the main steps were to first extend the time series to the entire year (by adding zeros), and then convert SWE measurements to SD. On the other hand, all Chilean SD and SWE data were sub-daily and provided in a raw format that typically required outlier and reference measurement point corrections, definition of snow/no-snow periods, and aggregation to a daily timestep.

SWE measurements were converted to SD based on four different approaches of increasing complexity. The simplest approach (model 1) used a constant snow density of 370 kg/m³—the median value measured across numerous Andean snow pits (Cornwell et al., 2016), some of which were co-located with the in situ monitoring stations (see Table 2). The other three approaches (models 2–4) utilized the empirical models of Hill et al. (2019), which were developed and fitted against an extensive network of North American stations covering a wide range of snow types. Each of these models is a more complex version of the prior—namely: the first models the relationship between SWE and SD with a power law (model 2); the second incorporates separate snow accumulation/ablation phases (model 3); whilst the third also incorporates climatic and day of year attributes (model 4). For more details on all models refer to Text S2.1 in Supporting Information S1. Model performance were evaluated on the DGA in situ stations with coincident SWE and SD values (Table 2). Based on model performance, a subset of models were selected, SD time series were derived and the impact of model choice on various results evaluated. Furthermore, using the same DGA in situ stations, additional analyses of some results (given in Section 4.2) were performed to test sensitivities to data type (i.e., SWE or SD).

### 3.2.2. Site, Snowpack, and Sentinel-1 Acquisition Attributes

For the purpose of understanding physical controls on algorithm performance, data on long-term and dynamic physical attributes were prepared. Twelve long-term physical attributes were extracted at each in situ site (Table 3). The attributes are grouped by location, topography, climatology, and land cover. Location and topographic attributes were selected as they are known to have an impact on SD distribution (Cornwell et al., 2016; Cortés & Margulis, 2017; Shaw, Gascoin, et al., 2020). Location information was directly obtained from the station metadata. Topographic variables were derived from a 3-arc second Shuttle Radar Topography Mission (SRTM) Digital Elevation Model (DEM) that was resampled to 1 km to match the snow depth retrievals. Climatology attributes were included because of their influence on snow quantity and atmospheric conditions, which govern snowpack characteristics (e.g., stratigraphy, microstructure, texture, and density; Sturm & Liston, 2021) potentially relevant
<table><caption>**Table 2**

In Situ Monitoring Station Attributes Organized by Country and Data Source</caption><thead><tr><th rowspan="2">ID</th><th rowspan="2">Station</th><th rowspan="2">Lat.</th><th rowspan="2">Long.</th><th rowspan="2">Elev. (m a.s.l.)</th><th rowspan="2">Yrs. of Data</th><th rowspan="2">Data type<sup>a</sup></th><th colspan="3">Mean snow characteristics</th></tr><tr><th>Snow class</th><th>Peak depth (m)</th><th>Days of snow cover</th></tr></thead><tbody><tr><td colspan="10"><strong>Argentina (Source: SNIH)</strong></td></tr><tr><td>A01</td><td>Cerro Negro</td><td>-29.89</td><td>-69.56</td><td>4172</td><td>4</td><td>SWE</td><td>Tundra</td><td>0.29</td><td>109</td></tr><tr><td>A02</td><td>Paso de Agua Negra</td><td>-30.19</td><td>-69.80</td><td>4459</td><td>1</td><td>SWE</td><td>Tundra</td><td>0.55</td><td>136</td></tr><tr><td>A03</td><td>Toscas</td><td>-33.16</td><td>-69.89</td><td>3000</td><td>4</td><td>SWE</td><td>Prairie</td><td>0.30</td><td>115</td></tr><tr><td>A04</td><td>Laguna Diamante</td><td>-34.20</td><td>-69.70</td><td>3300</td><td>4</td><td>SWE</td><td>Prairie</td><td>0.40</td><td>147</td></tr><tr><td>A05</td><td>Laguna Atuel</td><td>-34.51</td><td>-70.05</td><td>3420</td><td>4</td><td>SWE</td><td>Tundra</td><td>1.76</td><td>210</td></tr><tr><td>A06</td><td>Las Leñas</td><td>-35.14</td><td>-70.09</td><td>2545</td><td>3</td><td>SWE</td><td>Prairie</td><td>0.38</td><td>112</td></tr><tr><td>A07</td><td>Valle Hermoso</td><td>-35.14</td><td>-70.20</td><td>2250</td><td>4</td><td>SWE</td><td>Prairie</td><td>1.51</td><td>154</td></tr><tr><td>A08</td><td>Paso Pehuenche</td><td>-35.98</td><td>-70.39</td><td>2545</td><td>4</td><td>SWE</td><td>Prairie</td><td>1.46</td><td>177</td></tr><tr><td>A09</td><td>La Hoya</td><td>-42.83</td><td>-71.25</td><td>1530</td><td>1</td><td>SWE</td><td>Prairie</td><td>0.64</td><td>150</td></tr><tr><td>A10</td><td>Lago Gaucho</td><td>-43.83</td><td>-71.41</td><td>1165</td><td>4</td><td>SWE</td><td>Montane Forest</td><td>0.38</td><td>124</td></tr><tr><td colspan="10"><strong>Chile (Source: DGA)</strong></td></tr><tr><td>C01</td><td>Quebrada Larga</td><td>-30.73</td><td>-70.29</td><td>3500</td><td>2</td><td>Both (C)</td><td>Tundra</td><td>0.89</td><td>139</td></tr><tr><td>C02</td><td>Cerro Vega Negra*</td><td>-30.90</td><td>-70.52</td><td>3600</td><td>3</td><td>Both (I)</td><td>Tundra</td><td>1.01</td><td>146</td></tr><tr><td>C03</td><td>El Soldado</td><td>-32.01</td><td>-70.32</td><td>3290</td><td>1</td><td>SWE</td><td>Prairie</td><td>1.29</td><td>105</td></tr><tr><td>C04</td><td>Nacimiento Del Sobrante</td><td>-32.19</td><td>-70.48</td><td>3380</td><td>2</td><td>SD</td><td>Prairie</td><td>0.59</td><td>113</td></tr><tr><td>C05</td><td>Portillo*</td><td>-32.84</td><td>-70.11</td><td>3000</td><td>4</td><td>Both (C)</td><td>Prairie</td><td>0.96</td><td>131</td></tr><tr><td>C06</td><td>Glaciar Olivares Gamma</td><td>-33.16</td><td>-70.15</td><td>3628</td><td>1</td><td>SD</td><td>Tundra</td><td>0.76</td><td>108</td></tr><tr><td>C07</td><td>Valle Olivares</td><td>-33.19</td><td>-70.11</td><td>2786</td><td>1</td><td>SD</td><td>Prairie</td><td>0.52</td><td>115</td></tr><tr><td>C08</td><td>Laguna Negra*</td><td>-33.67</td><td>-70.11</td><td>2780</td><td>4</td><td>Both (C)</td><td>Prairie</td><td>0.94</td><td>130</td></tr><tr><td>C09</td><td>Las Melosas</td><td>-33.86</td><td>-70.27</td><td>3320</td><td>2</td><td>Both (I)</td><td>Tundra</td><td>1.18</td><td>143</td></tr><tr><td>C10</td><td>Termas Del Flaco</td><td>-34.89</td><td>-70.33</td><td>2650</td><td>4</td><td>Both (C)</td><td>Maritime</td><td>1.45</td><td>157</td></tr><tr><td>C11</td><td>Lo Aguirre</td><td>-35.97</td><td>-70.57</td><td>2000</td><td>4</td><td>Both (C)</td><td>Maritime</td><td>1.52</td><td>95</td></tr><tr><td>C12</td><td>Los Condores</td><td>-36.01</td><td>-72.55</td><td>2441</td><td>1</td><td>Both (C)</td><td>Maritime</td><td>3.69</td><td>206</td></tr><tr><td>C13</td><td>Volcan Chillan</td><td>-36.89</td><td>-71.41</td><td>2060</td><td>3</td><td>Both (C)</td><td>Maritime</td><td>1.08</td><td>159</td></tr><tr><td>C14</td><td>Alto Mallines</td><td>-37.16</td><td>-71.24</td><td>1770</td><td>4</td><td>Both (C)</td><td>Maritime</td><td>2.12</td><td>161</td></tr><tr><td colspan="10"><strong>Chile (Source: Codelco)</strong></td></tr><tr><td>C15</td><td>Laguna Angela</td><td>-33.08</td><td>-70.27</td><td>3575</td><td>1</td><td>SD</td><td>Prairie</td><td>1.56</td><td>206</td></tr><tr><td>C16</td><td>Laguna Barrosa</td><td>-33.11</td><td>-70.23</td><td>3765</td><td>1</td><td>SD</td><td>Tundra</td><td>0.70</td><td>177</td></tr><tr><td>C17</td><td>Machu Pichu</td><td>-33.17</td><td>-70.26</td><td>4250</td><td>1</td><td>SD</td><td>Tundra</td><td>0.95</td><td>208</td></tr><tr><td colspan="10"><strong>Chile (Source: Universidad de Chile)</strong></td></tr><tr><td>C18</td><td>PQN10</td><td>-33.31</td><td>-70.25</td><td>3570</td><td>2</td><td>SD</td><td>Tundra</td><td>0.57</td><td>150</td></tr><tr><td>C19</td><td>PQN07</td><td>-33.31</td><td>-70.26</td><td>3483</td><td>2</td><td>SD</td><td>Tundra</td><td>0.57</td><td>138</td></tr><tr><td>C20</td><td>PQN06</td><td>-33.32</td><td>-70.26</td><td>3556</td><td>2</td><td>SD</td><td>Tundra</td><td>0.93</td><td>167</td></tr></tbody></table>

Note. Snow class was extracted from Sturm and Liston (2021) while other characteristics were derived from monitoring data. * marks stations with snow pit density measurements (Cornwell et al., 2016). <sup>a</sup>Abbreviations for “Data Type” attribute “Both”: C—Coincident SWE & SD data; I—Independent SWE & SD data.

to radar observations. Land cover is known to influence radar observations (e.g., Meroni et al., 2021; Veloso et al., 2017) and the Sentinel-1 snow depth retrievals (Lievens et al., 2022). Three land cover categories were considered to try and distinguish the influence of evergreen forests, general vegetation and barren ground.

A total of five dynamic snowpack attributes were derived at each in situ site (Table 3) for each year of monitoring (i.e., annually), and also a mean of the individual years was calculated as a representative value for the full
<table><caption>**Table 3**
List of Physical Site Characteristics and Dynamic Snowpack Attributes Extracted or Derived for Each In Situ Monitoring Station</caption><thead><tr><th>Attribute</th><th>Definition</th><th>Unit</th></tr></thead><tbody><tr><td><strong>Location</strong></td><td></td><td></td></tr><tr><td>Latitude</td><td>Latitude of station (WGS84)</td><td>°S</td></tr><tr><td>Longitude</td><td>Longitude of station (WGS84)</td><td>°W</td></tr><tr><td><strong>Topography</strong></td><td></td><td></td></tr><tr><td>Elevation</td><td>Mean elevation of 1 × 1 km grid derived using SRTM DEM</td><td>m a.s.l.</td></tr><tr><td>Slope</td><td>Hillslope angle of 1 × 1 km grid derived using SRTM DEM</td><td>°</td></tr><tr><td>Northness</td><td>Cosine of the aspect of 1 × 1 km grid derived using SRTM DEM</td><td>—</td></tr><tr><td>Relief</td><td>Topographic relief within 1 × 1 km grid measured as the difference in elevation range derived using SRTM DEM</td><td>m</td></tr><tr><td><strong>Climate and snow classifications</strong></td><td></td><td></td></tr><tr><td>Climate type</td><td>Regional climate types based on 30 Köppen-Geiger climate classes (Peel et al., 2007)</td><td>—</td></tr><tr><td>Snow class</td><td>Global seasonal snow classification into one of seven categories (Sturm &amp; Liston, 2021)</td><td>—</td></tr><tr><td>Wind orientation</td><td>Leeward (east)/windward (west) relative to Andes Cordillera</td><td>—</td></tr><tr><td><strong>Land cover</strong></td><td></td><td></td></tr><tr><td>Evergreen forest cover fraction</td><td>Fraction of 1 × 1 km grid covered in evergreen forest (Tuanmu &amp; Jetz, 2014)</td><td>%</td></tr><tr><td>Vegetative cover fraction</td><td>Fraction of 1 × 1 km grid covered in vegetation (Tuanmu &amp; Jetz, 2014)</td><td>%</td></tr><tr><td>Barren land fraction</td><td>Fraction of 1 × 1 km grid with barren land (Tuanmu &amp; Jetz, 2014)</td><td>%</td></tr><tr><td><strong>Snowpack attributes (Mean and annual)</strong></td><td></td><td></td></tr><tr><td>Site representativeness</td><td>Classification of how similar (or representative) are the in situ measurements to the surrounding 1 × 1 km retrieval grid based on PlanetScope imagery (Planet Team, 2017). Sites manually assigned into one of three classes: low, medium or high</td><td>—</td></tr><tr><td>Number of Sentinel-1 observations</td><td>Number of Sentinel-1 observations during snow season per station</td><td>—</td></tr><tr><td>Maximum SD</td><td>Maximum observed (in situ) snow depth per station</td><td>m</td></tr><tr><td>Mean SD</td><td>Average observed (in situ) snow depth per station</td><td>m</td></tr><tr><td>Number snow days (in situ)</td><td>Number of days with snow cover based on in situ observations</td><td>Days</td></tr><tr><td>Number snow days (remote sensing)</td><td>Number of days that 1 × 1 km grid is classified to have snow cover based on VNP10A1F remotely sensed data set (details in Text S1.2 of Supporting Information S1)</td><td>Days</td></tr></tbody></table>

Note. Physical site attributes are grouped by location, topography, climate and snow classification, and land cover and are assumed to be static over the monitoring period.

monitoring period. In situ monitoring data were used to calculate the annual maximum snow depth, annual mean snow depth and the annual number of days with snow cover based on the in situ data and the VNP10A1F remotely sensed snow cover product. Sentinel-1 sampling (overpass) frequency during periods of snow cover was used to calculate the number of available data points for each snowy season to explore whether it has an influence on the retrieval algorithm, particularly as the retrievals are based on changes between successive timesteps. A “site representativeness” category was assigned to each site to capture how representative (or similar) the point measurements were to the 1 × 1 km retrieval grid. Site representativeness was visually assessed with PlanetScope (Planet Team, 2017) and manually assigned as “low,” “medium” or “high” based on (a) the consistency in snow cover timing between the in situ data and satellite imagery, and (b) the heterogeneity of snow cover over the 1 × 1 km grid. For instance, a site was deemed to have “low” representativeness if the snow cover was highly dynamic in the surrounding grid, and thus the measurements were capturing highly localized snow dynamics. In contrast, sites were categorized to have “high” representativeness if the snow cover was relatively homogeneous across most of the grid and the snow cover timing was well matched between the measurements and PlanetScope imagery.
### 3.3. Evaluation

The evaluation of the results is divided into two main components. The first part describes how the relationships (and influencing factors) between in situ snow depth and four backscatter coefficients (raw coefficients—$\sigma_{vv}^0$, $\sigma_{vh}^0$, $\sigma_{cr}^0$; and temporally smoothed $\sigma_{cr}^0$ as defined in Section 3.1) were assessed. This provides insights into what input signals and site attributes may or may not be affecting variability in performance of the algorithm. The second part details the approach used to evaluate the Sentinel-1 snow depth retrievals. Due to the strong influence of wet snow on backscatter coefficients, all assessments focus on dry snow periods only whilst Sentinel-1 observations flagged as “wet snow” by the wet snow mask (detailed in Section 3.1) are excluded.

#### 3.3.1. Backscatter Coefficient Correlations

For each monitoring site, temporal correlations ($R_t$) were calculated between time series of in situ snow depth measurements and each of the backscatter coefficients in co-polarisation ($\sigma_{vv}^0$), cross-polarisation ($\sigma_{vh}^0$) and the raw and temporally smoothed cross-polarisation ratio ($\sigma_{cr}^0$). The Pearson correlation coefficient was used to calculate $R_t$ using the full monitoring period ($R_{t,full}$) and also annually per site ($R_{t,ann}$) as follows:

$$R_{t, tp, pol} = \frac{\sum_{i=1}^{n} (\overline{SD}_{obs,i} - \overline{SD}_{obs}) \cdot (\sigma_{pol,i}^0 - \overline{\sigma}_{pol}^0)}{\sqrt{\sum_{i=1}^{n} (\overline{SD}_{obs,i} - \overline{SD}_{obs})^2 \cdot \sum_{i=1}^{n} (\sigma_{pol,i}^0 - \overline{\sigma}_{pol}^0)^2}} \quad (7)$$

where tp refers to the time period (full or annual), n is the number of co-incident Sentinel-1 and snow depth observations (when in situ snow depth > 0 m and the snow is dry), $SD_{obs,i}$ and $\sigma_{pol,i}^0$ are the in situ station snow depths and backscatter coefficients (for polarization pol), respectively, at timestep i, and $\overline{SD}_{obs}$ and $\overline{\sigma}_{pol}^0$ are the mean snow depths and backscatter coefficients. pol takes one of four values (vv, vh, raw cr, temporally smoothed cr). It is worth noting that $R_{t,full}$ values are also somewhat influenced by inter-annual variability in backscatter values prior to the start of the snow season, due to differences in soil properties and vegetation conditions.

Regional relationships between $R_t$ and monitoring site attributes (Table 3) were explored to determine whether any (and which) factors may influence backscatter coefficient response in the different polarizations across the study area. Again, the Pearson correlation coefficient was used to calculate the regional correlations ($R_r$) as follows:

$$R_{t, pol} = \frac{\sum_{a=1}^{s} (A_a - \bar{A}) \cdot (R_{t, pol, a} - \overline{R_{t, pol}})}{\sqrt{\sum_{a=1}^{s} (A_a - \bar{A})^2 \cdot \sum_{a=1}^{s} (R_{t, pol, a} - \overline{R_{t, pol}})^2}} \quad (8)$$

where s is the number of monitoring sites (or the total number of snow years for annually varying attributes), $A_a$ and $R_{t, pol, a}$ are the monitoring site attributes and temporal correlations (for polarization pol), respectively, for site a, and $\overline{A}$ and $\overline{R_{t, pol}}$ are the mean site attributes and temporal correlations across all sites (or snow years). Note that temporal correlations for the full monitoring period ($R_{t, full}$) were used when assessing relationships with long-term static (or mean) site attributes (including location, topography, snow classification and land cover), while annual correlations ($R_{t, ann}$) were used with annually varying snowpack attributes (such as mean annual snow depth or number of days with snow cover).

#### 3.3.2. Snow Depth Retrievals

The snow depth estimates ($SD_{mod}$) were evaluated against in situ data ($SD_{obs}$). Algorithm performance was evaluated on a per-site basis using five performance statistics: Pearson correlation coefficient ($R$), mean absolute error (MAE), normalized mean absolute error (nMAE), bias and normalized bias (nBias), defined as follows:

$$R = \frac{\sum_{i=1}^{n} (SD_{mod,i} - \overline{SD_{mod}}) \cdot (SD_{obs,i} - \overline{SD_{obs}})}{\sqrt{\sum_{i=1}^{n} (SD_{mod,i} - \overline{SD_{mod}})^2 \cdot \sum_{i=1}^{n} (SD_{obs,i} - \overline{SD_{obs}})^2}} \quad (9)$$
$$MAE = \frac{1}{n} \sum_{i=1}^{n} |SD_{\text{mod},i} - SD_{\text{obs},i}| \qquad (10)$$

$$nMAE = \frac{\frac{1}{n} \sum_{i=1}^{n} |SD_{\text{mod},i} - SD_{\text{obs},i}|}{SD_{\text{obs}}} \qquad (11)$$

$$Bias = \frac{1}{n} \sum_{i=1}^{n} (SD_{\text{mod},i} - SD_{\text{obs},i}) \qquad (12)$$

$$nBias = \frac{\frac{1}{n} \sum_{i=1}^{n} (SD_{\text{mod},i} - SD_{\text{obs},i})}{SD_{\text{obs}}} \qquad (13)$$

Here, n is the number of observations, SD<sub>mod,i</sub> and SD<sub>obs,i</sub> are the Sentinel-1 retrieval and in situ station snow depths, respectively, at timestep i, and SD<sub>mod</sub> and SD<sub>obs</sub> are the mean snow depths. All observations marked as wet snow were masked and excluded from the analysis.

Statistics were quantified for each site using the full time series and then also per water year. Statistics based on the full time series were used to evaluate algorithm performance. Relationships between the performance statistics and long-term physical and dynamic snowpack attributes were also assessed to identify potential drivers of performance.

## 4. Results

The results are presented below in three main components. In situ snowpack and Sentinel-1 observational data are described and explored in Section 4.1 to understand behavior and potential collinearities in site attributes (as listed in Table 3). In Section 4.2, temporal correlations between the in situ snow depth data and the Sentinel-1 backscatter coefficients (including raw σ<sup>0</sup><sub>vv</sub>, σ<sup>0</sup><sub>vh</sub>, σ<sup>0</sup><sub>cr</sub>, and temporally smoothed σ<sup>0</sup><sub>cr</sub>) are presented and possible factors influencing relationships are explored. Finally, Section 4.3 focuses on the actual snow depth retrievals and their evaluation.

It should be noted that all the results presented below are based on the most complex SWE to SD conversion approach (i.e., model 4) for the monitoring stations only collecting SWE data (see Section 3.2.1). Since models 1 and 4 were found to have the best and similar performance—as detailed in Text S2.2 in Supporting Information S1—relevant model 1 results are also reported below as necessary and included in Supporting Information S1.

### 4.1. In Situ Snowpack and Remote Sensing Snow Data Attributes

Figure 3 presents annually derived in situ snow attributes and Sentinel-1 overpass frequency during periods of snow cover at the monitoring sites. Mean snowpack depth across the monitoring sites was greatest in 2020, followed by 2017, 2018 and lastly 2019 (median depths of 0.87, 0.48, 0.38, and 0.25 m, respectively); while the number of Sentinel-1 observations increased on a year-to-year basis nearly doubling from a median of 21 in 2017 to a peak of 37 in 2020 (Figure 3a). Surprisingly, snow cover did not follow the same annual trend as mean snow depth. The site-averaged number of snow cover days was highest in 2017, lowest in 2019, while 2018 and 2020 had the same number of snow days (median number of snow days of 157, 117, 148, and 148, respectively).

Trends in annual days of snow cover were the same when derived with in situ observations or remote sensing (VNP10A1F), although VNP10A1F had a larger range in values. The reason for the differences between the snow depth and snow cover day trends is that the 2017 snow season started 20 days earlier, on average, than the 2020 season. Figure 3d shows that mean snow depth and number of snow days are positively correlated but with substantial variability particularly for more persistent snowpacks. For instance, across sites with 200 days of snow cover, mean snow depth can vary from 0.5 m to greater than 2 m. Assessment of year-to-year variation in snowpack attributes at individual monitoring sites suggests that site snow depth is consistently more variable than number of days of snow cover (Table S5 in Supporting Information S1). For example, at the majority of sites the annual number of days with snow cover seldom varies by more than 25% of the mean 2017–2020 value, while annual variability in mean snow depth is commonly greater than 60%. Similar observations have been made in an
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Annually derived attributes at the monitoring sites (Table 3) including mean snow depth, number of snow days based on in situ observations and VNP10A1F, and number of Sentinel-1 observations during periods of snow cover by (a) year, (b) snow class, (c) site representativeness of the sub-grid topographical variability of the SD retrieval grid, and (d) scatter plots colored by snow class.</figcaption></figure>

Andean headwater catchment where relatively small changes in snow cover equated to substantially larger differences in snow volume (Shaw, Deschamps-Berger, et al., 2020).

Annually derived snowpack attributes by seasonal snow class and site representativeness are shown in Figures 3b and 3c. Overall, the Maritime snow class has the highest mean snowpack depth (median SD value ~0.69 m across all sites and individual years) with high variation (ranges from 0.25 to 2.26 m), and also a consistently persistent snowpack (median number of in situ snow days ~159 days, range of 116–206 days). In comparison to Maritime, Tundra and Prairie have shallow mean snowpack depths (medians of 0.32 and 0.40 m), a similar number of days of snow cover but greater variability (medians of 150 and 139 days, ranges of 90–224 and 88–206 days, respectively). For example, the longest observed snow cover of 224 days was in 2017 for station Laguna Atuel (A05, Tundra snow class, Table 2) with a relatively shallow mean snow depth of 0.75 m. Montane Forest snow class only has one site and is therefore excluded from the analysis. Monitoring sites that do not represent well the SD retrieval grids (i.e., low representativeness) are characterized by shallow snowpacks (mean snow depth <
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Time series data for (left) Paso Pehuenche (ID: A08, Argentina) and (right) Laguna Negra (ID: C08, Chile) stations (Table 2) of (a) in situ snow depth (m) measurements, backscatter coefficients in (b) VV polarisation ($\sigma_{vv}^0$), (c) VH polarisation ($\sigma_{vh}^0$), (d) cross-polarisation ratio ($\sigma_{cr}^0$), where black dots show the raw backscatter values (as defined in Equation 1 in Section 3.1) whilst the red dots show the temporally smoothed data. In all plots, dots show actual values while general trends are shown by smoothed lines. Periods of snow cover based on VNP10A1F are shaded in gray, and the y-axis extents differ between the stations.</figcaption></figure>

0.5 m) and the shortest number of days with snow cover, and the opposite is observed for sites with high representativeness. The most obvious sampling bias is undersampling of Maritime compared to other snow classes.

## 4.2. Assessment of Sentinel-1 Backscatter Coefficients

### 4.2.1. Backscatter Coefficient Time Series at Sample Sites

Figure 4 shows 4 years of time series data of in situ snow depth measurements and backscatter coefficients at monitoring stations Paso Pehuenche in Argentina and Laguna Negra in Chile (IDs A08 and C08, respectively, Table 2). The two monitoring stations were selected as examples of the full set of results because of the general
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Cumulative distribution functions of temporal correlation ($R_t$) between in situ measurements and raw $\sigma^0$ (dB) observations for individual stations in VV and VH polarizations, and the raw and temporally smoothed cross-polarisation ratio. Range in $R_t$ values are on the x-axis, whilst y-axis shows the respective probability that a value (X) is less than or equal to $R_t$. Correlations are calculated for the full monitoring period when in situ snow depth > 0 m and are also masked for wet snow. See Figure S4 in Supporting Information S1 for the same SWE to SD model 1 results.</figcaption></figure>

similarity in their site attributes, including that both sites have a temperate, dry, and warm summer climate (Köppen-Geiger class Csb) (Peel et al., 2007), a Prairie snow class (Sturm & Liston, 2021), and are rather bare with limited vegetation (~60% and 50% barren land) (Tuanmu & Jetz, 2014) and nearly no tree cover (Table 2). However, they have notably different backscatter responses (Figures 4b–4e), which highlights some of the challenges that affect the algorithm performance generally. A detailed analysis of the difference in responses at these sites (see detail in Text S3 in Supporting Information S1) led to the following observations.

The $\sigma_{vh}^0$ response to snow cover is quite variable between Paso Pehuenche and Laguna Negra (Figure 4c). Based on visual inspection, $\sigma_{vh}^0$ observations at Paso Pehuenche are very sensitive to the presence of snow, with values consistently increasing between 2 and 3 dB from the start of the snow season to the pre-melt peak depth, irrespective of the interannual variability in peak snow depth (0.6–2 m). Observations coinciding with wet snow are clearly visible and characterized by substantial reductions in backscatter, commonly dropping below –25 dB. On the other hand, at Laguna Negra $\sigma_{vh}^0$ is substantially more variable from year-to-year and less sensitive to the presence of snow. Yet, the magnitude of change in $\sigma_{vh}^0$ from the start of the season to pre-melt peak depth appears to be linked to snowpack depth as the smallest (largest) changes in $\sigma_{vh}^0$ were observed in the driest (wettest) years (2019 and 2020, respectively). Although, the magnitude of change in $\sigma_{vh}^0$ does not appear to be linearly related to snow depth. One reason for the differences in $\sigma_{vh}^0$ response between the two stations may be that the snow-free ground surface at Paso Pehuenche has very little backscatter naturally (–22 dB as opposed to—17 dB at Laguna Negra), and so it ends up being highly sensitive to change in snow depth. Although, some of the variation in signal may also be due to

processing artifacts, including remaining impacts of the viewing geometry that could not be accounted for. At Laguna Negra, the decrease in $\sigma_{vv}^0$ with snow onset indicates a few possibilities, including: early wet snowfall, an initial freezing of soil or differences in dominant scattering mechanisms, with attenuation of ground scatter being larger than the snow volume scatter contribution at shallower depths. This could be influenced by Laguna Negra's shallower snow (1.49 vs. 0.94 m) and shorter snow cover period (177 vs. 130 days annually).

During periods of snow cover, observations of $\sigma_{vv}^0$ at Paso Pehuenche have far more variability between years than $\sigma_{vh}^0$ and no consistent relation to snow depth (Figure 4b). In comparison to basal properties, $\sigma_{vv}^0$ generally has low sensitivity to snow (Bernier & Fortin, 1998; Shi & Dozier, 2000; Strozzi & Matzler, 1998) and so the findings at Paso Pehuenche are not unexpected. The potential influence of soil moisture and temperature could not be evaluated due to the low reliability of 1 km² Noah-MP land surface model (Niu et al., 2011) simulations. Overall, the $\sigma_{vv}^0$ response was very similar to $\sigma_{vh}^0$ at Laguna Negra.

Cross-polarisation ratio backscatter—the primary input to the snow depth retrievals—is notably different between the two sites (Figure 4d). In the instance of Laguna Negra, the noise in $\sigma_{vv}^0$ and $\sigma_{vh}^0$ observations is reduced by merging both observations into a cross-polarisation ratio, and the variability in $\sigma_{cr}^0$ is further reduced through temporal smoothing (as shown by the red dots), elucidating the relationship between $\sigma_{cr}^0$ with snow depth, particularly in 2017 and 2020. In contrast, there is notably more variation in raw $\sigma_{cr}^0$ for Paso Pehuenche, where a response to snow depth is hard to distinguish, although the variability in $\sigma_{cr}^0$ is noticeably reduced by temporal smoothing. Orbital differences in incidence angles were smaller at Paso Pehuenche than Laguna Negra, and so are not a likely reason for the greater noise in cross-polarisation ratio backscatter observed for the former.

### 4.2.2. Variability in Correlations Between Backscatter Coefficients and In Situ Snow Depth

The distribution of temporal correlations $R_t$ between in situ measurements and $\sigma^0$ observations across the individual sites is shown in Figure 5. Overall, $\sigma_{vv}^0$ tends to have the lowest correlations with in situ measurements,
followed by $\sigma_{vh}^0$, $\sigma_{cr}^0$, and finally temporally smoothed $\sigma_{cr}^0$ which has highest correlations. The median values ($P = 0.5$) of $R_{t,vv}$, $R_{t,vh}$, $R_{t,cr}$, and smoothed $R_{t,cr}$ are -0.06, 0.17, 0.32, 0.38, respectively. This general pattern holds across most exceedance probabilities and individual sites, however there are a few exceptions. $R_{t,vv}$ was on occasion greater than $R_{t,vh}$, but only where both correlations with snow depth were negative. The few sites where $R_{t,vh}$ was greater than $R_{t,cr}$ coincided with sites with the greatest $R_{t,vv}$ values. The impacts of smoothing $R_{t,cr}$ differed depending on whether raw $R_{t,cr}$ was greater than or less than zero, where smoothing generally increased $R_{t,cr}$ for the former and decreased it for the latter. Furthermore, the gains to $R_{t,cr}$ by smoothing tended to be greater across sites with top 50% of $R_{t,cr}$ values (~0.3).

The fact that $R_{t,cr}$ values (raw and smoothed) are generally higher than $R_{t,vv}$ and $R_{t,vh}$ justifies the use of the cross-polarisation ratio for snow depth retrievals (as opposed to directly using the $\sigma_{vh}^0$ time series data), and indicates that snowy surfaces over the study region have a similar backscatter response to those observed in the northern hemisphere. The general enhancements in $R_{t,cr}$ due to temporal smoothing highlight the value of further reducing differences in viewing geometries and minimizing outlier impacts. It should be noted that these results are the same even when SWE was converted to SD assuming a constant density (Figure S4 in Supporting Information S1). However, based on an analysis of all stations with coincident SD and SWE observations (Figures S5 and S6 in Supporting Information S1), it is likely that $R_t$ values (especially $R_{t,cr}$) are slightly underestimated for some sites (Table 2 - data type “SWE”), and thus the reported ranges in $R_t$ are slightly conservative.

The spatial distribution of temporal correlations $R_{t,vv}$, $R_{t,vh}$, and temporally smoothed $R_{t,cr}$ are shown in Figure 6. All three temporal correlations have similar spatial patterns, where positive correlations are predominantly clustered in the central region of the study area (between latitudes of 32° and 36°S), while sites with negative values tends to occur in the northern or mid-southern extents of the study area. Monitoring stations Paso de Agua Negra (A02, in the north), and Volcan Chillan and Alto Mallines (located around 37°S; C13 and C14, Table 2) consistently have the most negative temporal correlations.

### 4.2.3. Site and Snowpack Drivers of Temporal Correlation Variability

Possible influence of physical site characteristics and dynamic snowpack attributes on the four temporal correlations ($R_{t,vv}$, $R_{t,vh}$, $R_{t,cr}$, and temporally smoothed $R_{t,cr}$) were assessed by deriving regional correlations ($R_{r,vv}$, $R_{r,vh}$, $R_{r,cr}$, and smoothed $R_{r,cr}$) separately between each of the site attributes and the four $R_t$ values (Table 4). Boxplots of $R_t$ were used to analyse the influence of categorical variables (Figure S7 in Supporting Information S1).

The relationships between $R_t$ and site attributes were generally low, and consistent irrespective of the SWE to SD conversion approach (see Table S6 in Supporting Information S1). The relationships were particularly weak when $\sigma_{vv}^0$ was used to calculate $R_t$; whilst across the other three source backscatter, the same site attributes generally had the highest $R_r$ values. The results are predominantly focused on $R_{r,cr}$ from hereon in, as smoothed $R_{r,cr}$ was typically within 0.03 points of the raw value. The strongest regional correlations were with land cover variables (particularly evergreen forest fraction $R_{r,cr} \sim -0.47$ and barren fraction $R_{r,cr} \sim 0.38$) suggesting that land cover is the strongest indicator of performance. Scatter plots show that $R_{t,cr}$ is inversely related to evergreen forest fraction (Figure 7a), indicating that the temporal correlations between $\sigma_{cr}^0$ and snow depth is greater when the land is less forested (or more bare) likely due to less attenuation and interference of the snowpack backscatter signal by vegetation, or it is less dominated by depolarization from vegetation compared to depolarization from snow. Interaction between vegetation and backscatter is well-known and has been observed in other regions (Lievens et al., 2019, 2022); although degradation in $R_{t,cr}$ across the study sites occurs at notably low levels of forest cover (~15%). This may be partly attributable to seasonal snow class where $R_{t,cr}$ are lower for sites classified as Maritime than Prairie despite having a deeper mean snowpack and higher site representativeness. Similarly, the large spread in $R_{t,cr}$ for sites with no forest cover may be associated with the Tundra snow class. Longitude and northness were the other two non-snowpack attributes with $R_{r,cr}$ greater than 0.2 (Table 4). Further examination of scatter plots (not shown here), highlights that regional correlations with northness and longitude are actually influenced by collinearities with evergreen forest fraction. For instance, longitude $R_{r,cr}$ is reduced to -0.02 once sites with an evergreen forest fraction greater than 15% are excluded.
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Maps of (a) the entire study area and (b) inset (included for improved visibility of high-density region) of temporal correlations *R*<sub>t,vv</sub>, *R*<sub>t,vh</sub>, and temporally smoothed *R*<sub>t,cr</sub> at monitoring stations. Raw *R*<sub>t,cr</sub> is not shown because of similarity to smoothed values. Correlations were calculated for the full monitoring period when in situ snow depth > 0 m and are masked for wet snow.</figcaption></figure>

All snowpack attributes had regional correlations (*R*<sub>r,cr</sub>) greater than 0.2 for at least one of the evaluation time periods, that is, derived across the full monitoring period or annually (Table 4). Overall, highest regional correlations were for the frequency of Sentinel-1 observations, followed by the snow depth attributes, and lastly number of days of snow cover attributes. Scatter plots of annually derived temporal correlation *R*<sub>t,cr</sub> against these attributes are used to further examine the relations (Figure 7b). *R*<sub>t,cr</sub> and mean snow depth are positively correlated (*R*<sub>t,cr</sub> = 0.24, Table 4). For shallow mean snow depths (less than 0.25 m), there is a large, relatively uniform spread between both negative and positive *R*<sub>t,cr</sub> values (*R*<sub>t,cr</sub> ~ -0.4–0.6). *R*<sub>t,cr</sub> tends to become more positive than
<caption>**Table 4**</caption>

Regional Correlations *R*<sub>r</sub> Derived Between Individual Temporal Correlations (*R*<sub>t,yy</sub>, *R*<sub>t,yh</sub>, Raw and Smoothed *R*<sub>t,cr</sub>) and Physical Site Characteristics or Dynamic Snowpack Attributes

<table><thead><tr><th rowspan="2">Name</th><th rowspan="2">Values</th><th colspan="4">Regional correlations (R<sub>r</sub>)</th></tr><tr><th>VV</th><th>VH</th><th>CR</th><th>CR (smooth)</th></tr></thead><tbody><tr><td colspan="6"><strong>Location</strong></td></tr><tr><td>Latitude</td><td>-33.3 (-43.8--29.9)</td><td>0.02</td><td>0.11</td><td>0.15</td><td>0.13</td></tr><tr><td>Longitude</td><td>-70.3 (-71.4--69.6)</td><td>0.18</td><td>0.25</td><td>0.32</td><td>0.30</td></tr><tr><td colspan="6"><strong>Topography</strong></td></tr><tr><td>Elevation</td><td>3288 (1105-4607)</td><td>0.14</td><td>0.15</td><td>0.14</td><td>0.15</td></tr><tr><td>Slope</td><td>6 (1-19)</td><td>0.26</td><td>0.08</td><td>-0.09</td><td>0.03</td></tr><tr><td>Northness</td><td>-0.17 (-1-0.99)</td><td>0.26</td><td>0.36</td><td>0.31</td><td>0.30</td></tr><tr><td>Relief</td><td>230 (12-687)</td><td>-0.06</td><td>-0.03</td><td>-0.01</td><td>0.04</td></tr><tr><td colspan="6"><strong>Land cover</strong></td></tr><tr><td>Evergreen forest cover fraction</td><td>4 (0-30)</td><td>-0.13</td><td>-0.37</td><td><strong>-0.47</strong></td><td><strong>-0.43</strong></td></tr><tr><td>Vegetative cover fraction</td><td>29 (0-95)</td><td>-0.17</td><td>-0.38</td><td>-0.36</td><td>-0.34</td></tr><tr><td>Barren land fraction</td><td>71 (5-100)</td><td>0.20</td><td><strong>0.40</strong></td><td>0.38</td><td>0.36</td></tr><tr><td colspan="6"><strong>Snowpack attributes</strong></td></tr><tr><td>Number of Sentinel-1 observations</td><td>A: 27 (10-60)<br>F: 27 (11-47)</td><td>-0.03<br>(0.17)</td><td>0.19<br>(0.25)</td><td>0.34<br>(0.3)</td><td>0.37<br>(0.34)</td></tr><tr><td>Maximum SD</td><td>A: 0.99 (0.11-3.69)<br>F: 0.94 (0.18-3.69)</td><td>-0.04<br>(-0.09)</td><td>0.14<br>(0.19)</td><td>0.23<br>(0.25)</td><td>0.16<br>(0.18)</td></tr><tr><td>Mean SD</td><td>A: 0.41 (0.05-2.26)<br>F: 0.45 (0.08-2.26)</td><td>0.04<br>(-0.03)</td><td>0.22<br>(0.26)</td><td>0.24<br>(0.29)</td><td>0.17<br>(0.23)</td></tr><tr><td>Number snow days (in situ)</td><td>A: 144 (88-224)<br>F: 144 (105-210)</td><td>0.01<br>(0.16)</td><td>0.10<br>(0.24)</td><td>0.15<br>(0.23)</td><td>0.18<br>(0.24)</td></tr><tr><td>Number snow days (remote sensing)</td><td>A: 151 (46-227)<br>F: 149 (63-216)</td><td>0.04<br>(0.18)</td><td>0.11<br>(0.25)</td><td>0.14<br>(0.22)</td><td>0.17<br>(0.23)</td></tr></tbody></table>

Note. For the snowpack attributes, two *R*<sub>r</sub> values were calculated using: annually derived values, as well as site-based means representing the full monitoring period (denoted as “A” and “F,” respectively). For each attribute, median values (and ranges) are shown. Attributes with *R*<sub>r</sub> ≥ |± 0.2| are italicized, while *R*<sub>r</sub> ≥ |± 0.4| are in bold. See Table S6 in Supporting Information S1 for the SWE to SD model 1 results.

negative for mean snow depths up to around 0.75 m; after which variability increases and *R*<sub>t,cr</sub> correlations range approximately between -0.6 and 0.8. In this instance, most sites with negative *R*<sub>t,cr</sub> have greater evergreen forest cover. Exclusion of sites with an evergreen forest fraction greater than 15% increases *R*<sub>t,cr</sub> for annually derived maximum and mean SD from 0.23 to 0.44 and 0.24 to 0.41, respectively. Distributions of snow classes and individual stations were checked but neither had an obvious impact on spread in *R*<sub>t,cr</sub> and mean snow depths. Despite deeper mean snow depths generally being associated with longer periods of snow cover across the monitoring stations (Figure 3), the relationship between *R*<sub>t,cr</sub> and both annually derived number of snow days is weak (in situ observations and VNP10A1F *R*<sub>t,cr</sub> = 0.15 and 0.14, respectively, Table 4). However, *R*<sub>t,cr</sub> values based on the full monitoring period are notably larger than the annually derived counterparts (and even more so when limited to sites with evergreen forest fraction below 15%), suggesting that climatologic data on snow cover persistence may be somewhat indicative of performance. Finally, *R*<sub>t,cr</sub> and Sentinel-1 observation frequency have the highest correlation of all the dynamic attributes (*R*<sub>t,cr</sub> =0.34, Table 4). This may be indicative of the potential importance of recurring Sentinel-1 observations in capturing temporal variability in snow depth, or that larger numbers of observations may facilitate better orbit bias correction and backscatter normalization.
<figure><img src="image_9.png" /><figcaption>**Figure 7.** (a) Scatter plots of *R*<sub>t,cr</sub> against evergreen forest fraction with points colored by (left) seasonal snow class, (center) site representativeness, and (right) mean snow depth in meters for the full monitoring period. Seasonal snow class abbreviations are the same as in Figure 1; (b) Scatter plots of annually derived temporal correlations *R*<sub>t,cr</sub> against mean snow depth, number of snow days based on in situ observations and VNP10A1F, and number of Sentinel-1 observations during periods of snow cover with points colored by evergreen forest fraction.</figcaption></figure>

In comparison to raw *R*<sub>t,cr</sub>, temporally smoothed *R*<sub>t,cr</sub> resulted in slightly greater regional correlations for most snowpack attributes (Table 4, last column). Maximum and mean SD were the two exceptions, where temporal smoothing resulted in reductions of 0.07 compared to raw *R*<sub>t,cr</sub> values. Considering that smoothed *R*<sub>t,cr</sub> values at sites were generally greater than the raw values (Section 4.2.1), this indicates that temporal smoothing reduced sources of noise—such as differences in orbital viewing geometries—so that site snow depth became less of a determining factor than prior to smoothing.

## 4.3. Evaluation of Sentinel-1 Snow Depth Estimates

### 4.3.1. Modeled Snow Depth

Maps and time series of the Sentinel-1 snow depth retrievals are shown in Figure 8. Snow depth estimates for the month of July are specifically shown in the maps, because snowpacks are widely accumulating mass at this time of year across the study region, but widespread snow melt has not commenced yet. The snow depth distribution across the region for July 2020 (Figure 8a) looks reasonable, with deeper snow observed on hilltops, at higher elevations, and in the latitude band of 34–36°S, broadly aligning with the findings of Cortés and Margulis (2017). Across most of the region, these estimates are based on eight or 10 acquisition dates (not shown in the figure), but do drop to as little as two dates in the south. Here, the snow depth retrievals do exhibit some artifacts due to the infrequent satellite acquisitions, as is highlighted by the “snow-free” diamond in the south.
<figure><img src="image_10.png" /><figcaption>**Figure 8.** (a) Map of Sentinel-1 snow depth across the study area for July 2020. Map insets and monitoring locations detailed in (b, c) and (d, e) are marked in red, respectively. Maps of mean July snow depth for two snow years (2019, 2020) are shown for (b) inset 1, and (c) inset 2. Glacier extents are masked in yellow. Snow depth time series of in situ measurements and Sentinel-1 retrievals (dry snow only) are shown for (d) Laguna Negra, and (e) Paso Pehuenche. Days with snow cover based on VNP10A1F are marked in gray.</figcaption></figure>

Differences in snow depth between July 2019 and July 2020 are evident for both insets and algorithm versions (Figures 8b and 8c). The 2020 snow season was far wetter than 2019, resulting in deeper snowpacks (see Section 4.1) and this is also reflected in the retrievals (Figure 8: B1 vs. B2, C1 vs. C2). The time series data show that modeled snow depth tends to be overestimated but this does vary between sites.

### 4.3.2. Sentinel-1 Snow Depth Performance

Figure 9 shows the cumulative distribution functions of correlation R, mean absolute error (MAE), normalized mean absolute error (nMAE), bias and normalized bias (nBias) derived across the in situ monitoring stations. The performance metrics were calculated for both the full time series and periods with snow cover only (based on the remotely sensed snow cover data set VNP10A1F detailed in Text S1.2 in Supporting Information S1). Mean correlation is 0.39 for periods of snow cover only (Figure 9a) but there is a large spread in correlations between sites with values varying uniformly from -0.17-0.94. Compared to the smoothed temporal correlations $R_{t,cr}$ in Section 4.2.2, the snow depth algorithm correlations are broadly similar but not identical for two main reasons: (a) Periods of snow cover are no longer based on in situ observations but rather remote sensing, meaning there may be shifts in snow cover duration and/or timing; (b) Some dates with Sentinel-1 observations may not result in a valid snow depth estimate because of a negative snow index (Section 3.1). The main differences in correlation occur for the poorest performing sites, where R tends to be greater than $R_{t,cr}$. For the full time series, mean correlations increased to 0.73 (and rarely dropped below 0.5 for any site), highlighting the substantial influence snow free
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Cumulative distribution functions of five performance metrics (a) temporal Pearson correlation (R), (b) mean absolute error (MAE), (c) normalized mean absolute error (nMAE) (d) bias, and (e) normalized bias. Values were derived for periods with full time series ("ALL": black line) and periods with snow cover only ("SC": red line). And masked for wet snow. The mean and median (in brackets) are in the bottom right corner for each metric. See Figure S8 in Supporting Information S1 for the SWE to SD model 1 results.</figcaption></figure>

periods have in inflating performance. These findings were consistent when model 1 was used to convert SWE to SD (Figure S8 in Supporting Information S1).

For periods with snow cover, modeled snow depth bias ranges from -0.66 to 1.95 m across monitoring stations (Figure 9d), with a median bias of 0.1 m meaning that modeled SD is more likely to be overestimated than underestimated. Normalized bias is relatively high (Figure 9e) - only half the sites have nBias within ±0.5 (50%), whilst nBias is greater than 1 (100%) across eight sites. Similarly, nMAE is relatively high across all sites (Figure 9c) where the minimum and median values are 0.25 and 0.7, respectively, but nMAE is especially large across two sites exceeding 6.6. The normalized data shows that differences between both metrics derived across
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Scatter plots of (a) correlation, (b) bias and (c) normalized bias for the modeled snow depth against evergreen forest fraction with points colored by seasonal snow class. Metrics are derived per station using the full monitoring period (2017–2020) and periods with snow cover only.</figcaption></figure>

the two time periods (“SC” vs. “ALL”) are actually negligible. The findings were consistent when model 1 was used to convert SWE to SD (Figure S8 in Supporting Information S1). However, the model errors were slightly larger—especially for normalized metrics—likely because snow depth derived with model 1 tends to be smaller than that derived using model 4 (as shown in Figure 9).

### 4.3.3. Key Drivers of Sentinel-1 Snow Depth Performance

The influence of different long-term (static) and annually derived factors on Sentinel-1 snow depth algorithm performance are shown in Figures 10 and 11, respectively. Only results for the correlation ($R$), bias and normalized bias (nBias) metrics are presented, while MAE (nMAE) is strongly correlated to bias (nBias) and was accordingly excluded. Rather than assessing the potential influence of all site attributes (Table 3) again, the factors shown in Figures 10 and 11 were selected based on the findings in the previous section (Section 4.2).

Over the full monitoring period, correlations between observed and modeled snow depth generally got poorer as evergreen forest fraction increased at the monitoring stations (Figure 10a), consistent with Section 4.2.3 findings. No stations with evergreen forest cover greater than 15% had correlations exceeding 0.5. The negative relationship between correlation and evergreen forest cover was not observed across all sites with very low evergreen forest fraction (less than ~5%) because of the large range in values (demarcated as “low EF” in Figure 10a). This variability may be related to snow conditions associated with the Prairie and Tundra snow classes—both of which are shallow and occur at high elevations. The snowpack at monitoring stations with very low evergreen forest cover may have more diverse snow stratigraphy and metamorphism, due to a lack of sheltering by trees and higher exposure to winds, potentially resulting in wind slabs which are common for both snow classes (Sturm & Liston, 2021).

On an annual basis, Figure 11 (top) highlights that correlations ($R$) between observed and modeled snow depth are influenced by site mean snow depth, whereby $R$ is generally greater for deeper snowpacks. Similarly, annual $R$ values are also positively correlated with number of snow days and the number of Sentinel-1 data points. Here, annual number of snow days appears to act as a surrogate variable for snow depth, which by extension also typically reflects snow persistence (e.g., intermittent vs. seasonal snow). This is sensible, considering that seasonal snowpacks will have a longer period of continuous snow cover than intermittent snowpacks, resulting in a continuous period of snowpack accumulation (and thus deeper snow), as opposed to intermittent snowpacks that have a mix of accumulation and ablation. These findings are broadly consistent with Section 4.2.3, other than that evergreen forest fraction has a lesser influence on modeled snow depth, whilst snowpack attributes become more important.

Evergreen forest fraction has no obvious impact on modeled snow depth bias once normalized (Figures 10b and 10c). However, nBias appears to be dependent on snow class, where Maritime sites are generally unbiased, Prairie sites have a positive bias, whilst largest outliers occur for Tundra sites. These findings suggest that the northern
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Scatter plots of annually derived (top) correlation, (middle) bias and (bottom) normalized bias against mean snow depth, number of snow days based on in situ observations and VNP10A1F, and number of Sentinel-1 observations. Correlation points are colored by evergreen forest fraction, whilst bias is colored by snow class. Metrics are derived only using periods with snow cover.</figcaption></figure>

hemisphere sites used in optimizing the snow depth algorithm may be more similar to the Maritime snow class. Figure 11 (center, bottom) suggests that mean annual snow depths at a monitoring station have a major influence on year-to-year model bias and nBias. Snow depth is generally overestimated over shallow snowpacks ($<$~0.5 m) but is underestimated for deeper snowpacks, although nBias indicates that the underestimation is proportional to the snow depth. This suggests that the backscatter response is not linearly related to snow depth as assumed in the algorithm and may be influenced by other snow and/or ground properties, such as anisotropy, freeze-thaw cycles, etc (Yueh et al., 2009). Variability in bias may be compounded by point-area mismatches between the monitoring stations and Sentinel-1 retrieval grid.

# 5. Discussion

## 5.1. Applicability of Sentinel-1 Snow Depth Retrievals Over the Andes

At a given location, Sentinel-1 snow depth algorithm utility in the Andes was principally related to land cover, followed by snow class, snowpack attributes and Sentinel-1 observation frequency (Table 4). Consistent with PMW snowpack retrievals (e.g., Smith & Bookhagen, 2016; Vuyovich et al., 2014), best performance was observed over barren land or areas with no more than 15% evergreen forest cover, suggesting that the Andes have greater sensitivity to evergreen forests than observed in the Alps (Lievens et al., 2022). Importantly, this
relationship may not just be indicative of the influence of evergreen forest on retrievals, but also capture the
differences in snowpack attributes commonly associated with forest-cover due to wind-interception, partitioning
of longwave and shortwave radiation, etc. (refer to Tennant et al. (2017) for vegetation impacts on snow).

The influence of station snow depth and snow cover on Sentinel-1 algorithm performance became evident when excluding sites with higher evergreen forest fraction (Section 4.2.3) and considering individual snow seasons, because of the substantial inter-annual variability of snow cover properties (Section 4.1). Generally, correlation between observed SD and modeled SD (or R<sub>t,cr</sub>) was greater for deeper and more persistent snowpacks. The positive relationship with annual number of snow days is useful as it may be a valuable surrogate for snow depth and determining where the Sentinel-1 retrievals may be more accurate. However, there is still a lot of variability in correlations across individual snow depth ranges and sites, indicating there are unaccounted influencing factors yet to be incorporated into the Sentinel-1 snow depth model. This potentially includes (a) different snow attributes including liquid water in the snowpack (Painter et al., 2016), snow microstructure and stratigraphy, freeze-thaw cycles (Yueh et al., 2009); (b) basal surface properties such as ground surface cover (Feng et al., 2021), freezing/thawing of soil, and (c) influences of Sentinel-1 data preprocessing and viewing geometry. For instance, recent findings from the Alps show snow volume scattering has stronger contributions to backscatter at higher—as opposed to lower—local incidence angles, due to the longer travel path through the snow (Jans et al., 2025). Additionally, some of the σ⁰ temporal variability would be as a result of observation footprint size differences between orbits, although the temporal smoothing did reduce some of the orbital variability and improve σ⁰<sub>cr</sub> correlations to snow depth (see Section 4.2.2).

Sentinel-1 snow depth algorithm evaluation showed systematic errors in bias related to snow class (Figure 10c) and mean annual snow depth (Figure 11—bottom). Normalized SD tends to be overestimated for shallower depths (<0.5 m) and Tundra and Prairie snow classes, whilst deeper snowpacks (>1 m) are slightly underestimated. The non-linear response in backscatter to snow depth across the Andes may be due to a number of factors, including signal saturation at shallower depths, snowpack moistening toward peak depth in the snow season, or algorithm sub-optimality.

This study indicates that there are some substantial differences in Sentinel-1 snow depth retrieval accuracy be-
tween the Andes and northern hemisphere mountain ranges, with notably poorer algorithm performance across
the former. Mean correlation, MAE, and bias for Northern Hemisphere mountains was reported to be 0.65,
0.31 m, and -0.02 m, respectively, in comparison to 0.41, 0.51, and 0.21 m for the Andes. The weaker per-
formance may be due to a number of factors. Lievens et al. (2022) found that snow depth-backscatter correlations
were poorest for shallow snow while algorithm performance was best for depths exceeding 1 m. Comparatively,
snowpacks were relatively shallow across the Andes, which was likely conflated by this analysis being performed
during an ongoing drought, with 2018 and 2019 being extremely dry years (Shaw, Deschamps-Berger,
et al., 2020). Another likely compounding factor is the higher Sentinel-1 observation frequency over the primary
northern hemisphere monitoring sites (European Space Agency, 2023), as this had an important influence on
temporal correlations in the Andes. The Andean monitoring stations used in this study are likely far higher than
the dominant cluster of Alpine (mean elevation of 3,026 and 1,395 m), Scandinavian and North American sta-
tions, and therefore experience different weather systems with different snow classes and snow characteristics.
Snow is likely to be drier at higher elevations, which impacts radar volume scattering behavior.

An important implication of this study is that it evaluated the suitability of a new remote sensing technique for regional snow depth retrievals across the Andes, extending knowledge on remote sensing applications beyond the northern hemisphere. To the best knowledge of the authors, it is one of the only continuous and free approaches for monitoring mountain snowpack depths at quasi-global scales. Unlike the prior implementations of the algorithm (Hoppinen et al., 2024; Lievens et al., 2019, 2022), all data sets used in this study are globally available, ensuring consistency in the retrievals and enabling easy implementation of the workflow to other regions. Finally, the outputs and findings of this work can (and should) be used to inform locations of intensive data collection and modeling campaigns for improved understanding of Sentinel-1 backscatter response and the advancement of Sentinel-1 snow depth retrievals across the Andes.
## 5.2. Limitations of Sentinel-1 Approach

Poorest performance in the Sentinel-1 snow depth algorithm were observed in regions of shallow snow and evergreen forest cover. Therefore the approach is likely not suited to different parts of the Andes including (a) the northern Chilean Andes characterized by shallow and intermittent snow, and (b) regions south of 42°S where evergreen forests are abundant and the “wet” Maritime snow class dominates (Figure 1). However, very little data was available to evaluate algorithm performance across the southern Andes and therefore point (b) could not be substantiated.

For algorithm evaluation, this study utilized an opportunistic snow depth data set, in the sense that existing data across the region was compiled as opposed to being strategically collected. However, opportunistic data collection did result in (a) inconsistencies in the measured variable between stations (SD or SWE), (b) a mismatch in model (grid) and evaluation (point) scale, and (c) the use of a relatively small data set by global standards that was not able to definitively cover all snow and land characteristics. Even so, the impact of the limitations were minimized as practically possible. For instance for point (a) multiple SWE to SD conversion approaches were tested; however the findings were consistent irrespective of whether the simplest or most complex conversion approach was applied (i.e., assuming a constant snow density based on regional snow pit measurements (Cornwell et al., 2016) or accounting for spatiotemporal variability in snow densities (Hill et al., 2019)). While in the instance of point (b), differences in point-grid scale were somewhat alleviated by the analysis predominantly focusing on temporal correlation, which are less susceptible to the monitoring location. Due to the relatively small evaluation data set (point c), a simple approach was used to assess the influence of the likely main drivers on algorithm performance. However, it is acknowledged that assessments could benefit from more advanced multivariate approaches and to explore the influence of additional variables on performance, particularly as longer evaluation data sets become available.

Finally, snow cover data is an important source of uncertainty in the Sentinel-1 snow depth retrieval algorithm. The algorithm uses a separate binary data set to determine whether snow is, or is not, present on the ground for each pixel; and then either derives snow depth from the Sentinel-1 data or sets it to zero. However, it is not uncommon for a pixel to have both a mix of snow and bare land. However, as this is not currently accounted for in the algorithm, the “mixed” pixels need to be reclassified into snow or bare land. A snow cover data set had to be derived for this study (Section 3.1 and Text S1.2 in Supporting Information S1) because there are currently no binary snow cover products that also extend into the southern hemisphere. As a part of this process, assumptions had to be made on what constituted a “snowy pixel,” in terms of NDSI and pixel count applied when averaging. Multiple thresholds were tested to select the best fit based on monitoring station data.

## 5.3. Recommendations and Future Potential

This study is one of the first that applies SAR data to snow monitoring in the Andes, and is also one of the few that is retrieving regional snow depths from Sentinel-1 observations. Despite the “Northern Hemisphere” algorithm performing well across some sites, it is likely improved performance can be achieved through algorithm enhancements specifically targeted at the Andes. For instance, it may be advantageous to shift away from the current practice of using a “global” or “regional” parameter set to exploring options for targeted parameter optimization based on important site attributes, such as focusing on snow class due to its influence on model bias in the Andes (Figure 1c). This should not only reduce mean errors across all sites, but reduce the large spread in biases observed between sites. High-resolution, spatial snow depth data sets of Andean catchments and/or hillslopes should also be used in conjunction with ancillary data (e.g., basal conditions, local incidence angles) to explore spatial patterns in backscatter change (to pre-snow conditions) and elucidate important influencing factors. Existing snow depth data from Pleiades (e.g., Shaw, Deschamps-Berger, et al., 2020; Shaw, Gascoin, et al., 2020) or terrestrial Lidar (e.g., Mendoza et al., 2020) could provide a good starting point. Furthermore, detailed site-specific models capturing soil hydrology (e.g., Niu et al., 2011) and snowpack dynamics (e.g., Vionnet et al., 2012) should be utilized to understand the difference in response between seemingly similar sites (e.g., Paso Pehuenche as opposed to Laguna Atuel) or even to unravel annual differences in response at certain stations. This should be complemented with detailed assessment of data processing impacts on backscatter time series trends, including the approach implemented in Lievens et al. (2022). Lessons from these analyses would ideally be used to improve algorithm implementations and performance across the Andes.
A number of limitations were identified in the previous section that would be worth addressing in future work. Ongoing measurements across the monitoring stations will be helpful in reducing uncertainties, especially by collecting additional data representative of above-average rainfall conditions, to compensate for the below average snowpack depths experienced during the study monitoring period due to ongoing drought. Limitations in uneven sampling across different snow classes and snow depths by the monitoring stations used in this study could be addressed by using a spatially distributed evaluation data set, such as that of Cortés et al. (2016); Cortés and Margulis (2017). Unfortunately, this data does not currently extend into the Sentinel-1 period and would be influenced by other uncertainties. However, it is the best available spatially resolved data set for the Andes at present.

Remote-sensing has vast potential in helping understand the global water cycle and to monitor resources regionally. Thus, it is worth considering (a) what future potential do the Sentinel-1 snow depth retrievals have for global snowpack monitoring, and (b) how does this integrate with other current/future satellite missions. This study showed that two important influencing factors on snow depth-backscatter correlations were the number of Sentinel-1 observations and mean annual snow depth (Section 4.2). The influence of number of Sentinel-1 observations is relevant because Sentinel-1 data observation scenarios vary by region (European Space Agency, 2023). For instance, observation frequency over Australian snow fields is some of the lowest globally (European Space Agency, 2023). Furthermore, global Sentinel-1 data acquisition rates have been negatively impacted since December 2021, when satellite Sentinel-1B went out of operation. This should be rectified with the launch of two new satellites (Sentinel-1C and Sentinel-1D), extending the Sentinel-1 satellite mission beyond 2030. Another point of consideration into the future is the impact of climate change on snowpacks. In regions such as central Chile, snowpack depths are likely to become shallower as temperatures are projected to increase while precipitation decreases (Salazar et al., 2024); whilst elsewhere climate change may result in increased precipitation or a change in precipitation phase shifting from snow to rain at lower elevations. Hence, Sentinel-1 retrieval suitability may be impacted in some regions because of the strong dependence on snow depth and wet snow.

In addition to C-band, future satellite missions should also incorporate X- and Ku-band SAR because of the better sensitivity to shallower snowpacks (such as the intermittent snowpacks currently found in the northern study area extents). In the meantime, a potential pathway for satellite snow depth retrievals could involve leveraging data from multiple sources to (a) further optimize the Sentinel-1 retrievals with site-specific parameters for instance, and (b) develop a merged snow depth product. Possible sources include altimetry-based snow depths (e.g., ICESat-2; Deschamps-Berger et al., 2023) and snow depths derived using the signal interferometry (e.g., Marshall et al., 2021) from future L-band InSAR satellite missions with freely available data, such as ROSE-L and NISAR.

## 6. Conclusions

In this study, C-band synthetic aperture radar data from the Sentinel-1 satellite mission were used to estimate snow depth over the Chilean and Argentine Andes mountains. A recent approach applied to the northern hemisphere (Lievens et al., 2019) was generalized and extended to the Andes by only using globally available remote sensing data instead of regional data sources of the original approach. The utility of the original northern hemisphere algorithm version (Lievens et al., 2019) was evaluated across 30 Andean snow depth monitoring stations over four winters (2017–2020). First temporal correlations between snow depth observations and raw backscatter coefficients in different polarizations were assessed to determine important factors affecting backscatter temporal variability and by extension performance of the snow depth model. Next, correlation, bias and MAE of snow depth estimates were evaluated. In summary, the main findings of this study are:

1. Overall, Sentinel-1 retrievals have some potential over the Andes despite performance being poorer than in the northern hemisphere. Model skill is site-specific and highly variable across the study region.

2. Evergreen forest cover is the principal factor affecting temporal correlations between observed snow depth and Sentinel-1 backscatter data, although the importance of mean snow depth and snow cover persistence increases across sites with lower evergreen forest cover.

3. Correlations between observed and estimated snow depth are better over regions with limited evergreen forest cover (<15%) and for snow depths greater than 0.75 m.
4. Systemic overestimation of modeled snow depth are observed for shallow snowpacks and Tundra and Prairie snow classes. Unbiased snow depth estimates are associated with the Maritime snow class, characterized by deeper, wetter snow.

To enhance the utility of Sentinel-1 SAR snow depth retrievals over the Andes, future in-depth studies looking at spatial variability in backscatter are needed to investigate algorithm limitations, data processing and local incidence angle implications, and elucidate basal-impacts on backscatter. This work should also be complemented by time-series modeling of snow and soil properties at two similar sites with differing backscatter responses. Finally, other more detailed site-sensitive methods for model parameterization should be explored to improve snow depth retrievals over individual locations. This is foreseen as further research.

## Acronyms

<table><tr><td>DEM</td><td>Digital elevation model</td></tr><tr><td>DGA</td><td>Dirección General de Aguas, Chile</td></tr><tr><td>M</td><td>Maritime</td></tr><tr><td>MAE</td><td>Mean absolute error</td></tr><tr><td>MF</td><td>Montane forest</td></tr><tr><td>P</td><td>Prairie</td></tr><tr><td>PMW</td><td>Passive microwave</td></tr><tr><td>SAR</td><td>Synthetic aperture radar</td></tr><tr><td>SD</td><td>Snow depth</td></tr><tr><td>SNIH</td><td>Sistema Nacional de Información Hídrica, Argentina</td></tr><tr><td>SRTM</td><td>Shuttle Radar Topography Mission</td></tr><tr><td>SWE</td><td>Snow water equivalent</td></tr><tr><td>T</td><td>Tundra</td></tr><tr><td>VH</td><td>Vertical-horizontal</td></tr><tr><td>VV</td><td>Vertical-vertical</td></tr></table>

## Data Availability Statement

The in situ station data is bound by data sharing agreements from the respective owners, but can be obtained directly from the source as follows: Argentine government agency data can be freely accessed on the Sistema Nacional de Información Hídrica at https://snih.hidricosargentina.gob.ar/; Chilean government data needs to be requested through the "Chilean Transparency Law" from the Chilean Dirección General de Aguas at https://www.portaltransparencia.cl/PortalPdT/ley-de-transparencia, whilst the Andina mine data can be requested directly from CODELCO (https://www.codelco.com/). Access to the Planet image collections was granted via the PlanetLabs research and education initiative and an account can be requested at https://www.planet.com/. Sentinel-1 data are available from Google Earth Engine (image collection name: COPERNICUS/S1_GRD_FLOAT). VNP10A1F data are available from the National Snow and Ice Data Center at https://nsidc.org/data/vnp10a1f/versions/1. All other ancillary data used here is referenced, wherein data access details are provided.

## References

Ayala, A., Farías-Barahona, D., Huss, M., Pellicciotti, F., McPhee, J., & Farinotti, D. (2020). Glacier runoff variations since 1955 in the Maipo River basin, in the semiarid Andes of central Chile. *The Cryosphere*, 14(6), 2005–2027. https://doi.org/10.5194/TC-14-2005-2020

Barnett, T. P., Adam, J. C., & Lettenmaier, D. P. (2005). Potential impacts of a warming climate on water availability in snow-dominated regions. *Nature*, 438(7066), 303–309. https://doi.org/10.1038/nature04141

## Acknowledgments

This research was supported by the University of Queensland's PhD scholarship program, the Australian Research Council under the Future Fellowship program (Project ID: FT140100977), and the Sustainable Minerals Institute International Centre of Excellence (Chile). Fiona Johnson is supported by a UNSW Scientia Funding and ARC Training Centre in Data Analytics for Resources and Environments (Grant IC190100031). The authors also thank Liliana Pagliero, Maxi Viale and Rodrigo Correa for their support with obtaining the DGA, SNIH, and Codelco data sets, and the PlanetLabs research and education initiative for free imagery. Open access publishing facilitated by The University of Queensland, as part of the Wiley - The University of Queensland agreement via the Council of Australian University Librarians.
Bernier, M., & Fortin, J. P. (1998). The potential of times series of C-band SAR data to monitor dry and shallow snow cover. *IEEE Transactions on Geoscience and Remote Sensing, 36*(1), 226-243. https://doi.org/10.1109/36.655332

Bormann, K. J., Brown, R. D., Derksen, C., & Painter, T. H. (2018). Estimating snow-cover trends from space. *Nature Climate Change, 8*(11), 924-928. https://doi.org/10.1038/s41558-018-0318-3

Bulovic, N., McIntyre, N., & Johnson, F. (2020). Evaluation of IMERG V05B 30-min rainfall estimates over the high-elevation tropical Andes mountains. *Journal of Hydrometeorology, 21*(12), 2875-2892. https://doi.org/10.1175/JHM-D-20-0114.1

Burger, F., Ayala, A., Farias, D., Shaw, T. E., MacDonell, S., Brock, B., et al. (2019). Interannual variability in glacier contribution to runoff from a high-elevation Andean catchment: Understanding the role of debris cover in glacier hydrology. *Hydrological Processes, 33*(2), 214-229. https://doi.org/10.1002/HYP.13354

Cornwell, E., Molotch, N. P., & McPhee, J. (2016). Spatio-temporal variability of snow water equivalent in the extra-tropical Andes Cordillera from distributed energy balance modeling and remotely sensed snow cover. *Hydrology and Earth System Sciences, 20*(1), 411-430. https://doi.org/10.5194/hess-20-411-2016

Cortés, G., Girotto, M., & Margulis, S. (2016). Snow process estimation over the extratropical Andes using a data assimilation framework integrating MERRA data and Landsat imagery. *Water Resources Research, 52*(4), 2582-2600. https://doi.org/10.1002/2015WR018376

Cortés, G., & Margulis, S. (2017). Impacts of El Niño and La Niña on interannual snow accumulation in the Andes: Results from a high-resolution 31 year reanalysis. *Geophysical Research Letters, 44*(13), 6859-6867. https://doi.org/10.1002/2017GL073826

Deschamps-Berger, C., Gascoin, S., Berthier, E., Deems, J., Gutmann, E., Dehecq, A., et al. (2020). Snow depth mapping from stereo satellite imagery in mountainous terrain: Evaluation using airborne laser-scanning data. *The Cryosphere, 14*(9), 2925-2940. https://doi.org/10.5194/tc-14-2925-2020

Deschamps-Berger, C., Gascoin, S., Shean, D., Besso, H., Guiot, A., & López-Moreno, J. I. (2023). Evaluation of snow depth retrievals from ICESat-2 using airborne laser-scanning data. *The Cryosphere, 17*(7), 2779-2792. https://doi.org/10.5194/TC-17-2779-2023

Dozier, J., Bair, E. H., & Davis, R. E. (2016). Estimating the spatial distribution of snow water equivalent in the world's mountains. *WIREs Water, 3*(3), 461-474. https://doi.org/10.1002/wat2.1140

European Space Agency. (2023). Sentinel-1 observation scenario. Retrieved from https://sentinels.copernicus.eu/web/sentinel/missions/sentinel-1/observation-scenario

Feng, T., Hao, X., Wang, J., Li, H., & Zhang, J. (2021). Quantitative evaluation of the soil signal effect on the correlation between Sentinel-1 cross ratio and snow depth. *Remote Sensing, 13*(22), 4691. https://doi.org/10.3390/RS13224691

Garreaud, R. D., Alvarez-Garreton, C., Barichivich, J., Boisier, J. P., Christie, D., Galleguillos, M., et al. (2017). The 2010-2015 megadrought in central Chile: Impacts on regional hydroclimate and vegetation. *Hydrology and Earth System Sciences, 21*(12), 6307-6327. https://doi.org/10.5194/hess-21-6307-2017

Garreaud, R. D., Boisier, J. P., Rondanelli, R., Montecinos, A., Sepúlveda, H. H., & Veloso-Aguila, D. (2020). The Central Chile Mega Drought (2010-2018): A climate dynamics perspective. *International Journal of Climatology, 40*(1), 421-439. https://doi.org/10.1002/JOC.6219

Gascoin, S., Lhermitte, S., Kinnard, C., Bortels, K., & Liston, G. E. (2013). Wind effects on snow cover in Pascua-Lama, Dry Andes of Chile. *Advances in Water Resources, 55*, 25-39. https://doi.org/10.1016/j.ADVWATRES.2012.11.013

Hammond, J. C., Saavedra, F. A., & Kampf, S. K. (2018). Global snow zone maps and trends in snow persistence 2001-2016. *International Journal of Climatology, 38*(12), 4369-4383. https://doi.org/10.1002/JOC.5674

Hill, D. F., Burakowski, E. A., Crumley, R. L., Keon, J., Michelle Hu, J., Arendt, A. A., et al. (2019). Converting snow depth to snow water equivalent using climatological variables. *The Cryosphere, 13*(7), 1767-1784. https://doi.org/10.5194/TC-13-1767-2019

Hoppinen, Z., Palomaki, R. T., Brencher, G., Dunmire, D., Gagliano, E., Marziliano, A., et al. (2024). Evaluating snow depth retrievals from Sentinel-1 volume scattering over NASA snowwex sites. *EGUSphere, 2024*, 1-35. https://doi.org/10.5194/tc-18-5407-2024

Huerta, M. L., Molotch, N. P., & McPhee, J. (2019). Snowfall interception in a deciduous Nothofagus forest and implications for spatial snowpack distribution. *Hydrological Processes, 33*(13), 1818-1834. https://doi.org/10.1002/hyp.13439

Immerzeel, W. W., Lutz, A. F., Andrade, M., Bahl, A., Biemans, H., Bolch, T., et al. (2020). Importance and vulnerability of the world's water towers. *Nature, 577*(7790), 364-369. https://doi.org/10.1038/s41586-019-1822-y

Jans, J.-F., Beernaert, E., De Breuck, M., Brangers, I., Dunmire, D., De Lannoy, G., & Lievens, H. (2025). Sensitivity of Sentinel-1 C-band SAR backscatter, polarimetry and interferometry to snow accumulation in the Alps. *Remote Sensing of Environment, 316*, 114477. https://doi.org/10.1016/j.rse.2024.114477

King, J., Derksen, C., Toose, P., Langlois, A., Larsen, C., Lemmetyinen, J., et al. (2018). The influence of snow microstructure on dual-frequency radar measurements in a tundra environment. *Remote Sensing of Environment, 215*, 242-254. https://doi.org/10.1016/j.rse.2018.05.028

Krogh, S. A., Pomeroy, J. W., & McPhee, J. (2015). Physically based mountain hydrological modeling using reanalysis data in Patagonia. *Journal of Hydrometeorology, 16*(1), 172-193. https://doi.org/10.1175/JHM-D-13-0178.1

Li, D., Lettenmaier, D. P., Margulis, S. A., & Andreadis, K. (2019). The value of accurate high-resolution and spatially continuous snow information to streamflow forecasts. *Journal of Hydrometeorology, 20*(4), 731-749. https://doi.org/10.1175/JHM-D-18-0210.1

Lievens, H., Brangers, I., Marshall, H.-P., Jonas, T., Olefs, M., & De Lannoy, G. (2022). Sentinel-1 snow depth retrieval at sub-kilometer resolution over the European Alps. *The Cryosphere, 16*(1), 159-177. https://doi.org/10.5194/TC-16-159-2022

Lievens, H., Demuzere, M., Marshall, H.-P., Reichle, R. H., Brucker, L., Brangers, I., et al. (2019). Snow depth variability in the Northern Hemisphere mountains observed from space. *Nature Communications, 10*(1), 4629. https://doi.org/10.1038/s41467-019-12566-y

Marin, C., Bertoldi, G., Premier, V., Callegari, M., Brida, C., Hürkamp, K., et al. (2020). Use of Sentinel-1 radar observations to evaluate snowmelt dynamics in alpine regions. *The Cryosphere, 14*(3), 935-956. https://doi.org/10.5194/tc-14-935-2020

Marshall, H. P., Deeb, E., Forster, R., Vuyovich, C., Elder, K., Hiemstra, C., & Lund, J. (2021). L-band InSAR depth retrieval during the NASA SnowEx 2020 campaign: Grand Mesa, Colorado. In *International Geoscience and Remote Sensing Symposium (IGARSS), 2021-July* (pp. 625-627). https://doi.org/10.1109/IGARSS47720.2021.9553852

Marti, R., Gascoin, S., Berthier, E., De Pinel, M., Houet, T., & Laffly, D. (2016). Mapping snow depth in open alpine terrain from stereo satellite imagery. *The Cryosphere, 10*(4), 1361-1380. https://doi.org/10.5194/TC-10-1361-2016

Masiokas, M. H., Rabatel, A., Rivera, A., Ruiz, L., Pitte, P., Ceballos, J. L., et al. (2020). A review of the current state and recent changes of the Andean cryosphere. *Frontiers in Earth Science, 8*, 99. https://doi.org/10.3389/FEART.2020.00099

Masiokas, M. H., Villalba, R., Luckman, B. H., Le Quesne, C., & Aravena, J. C. (2006). Snowpack variations in the central Andes of Argentina and Chile, 1951-2005: Large-scale atmospheric influences and implications for water resources in the region. *Journal of Climate, 19*(24), 6334-6352. https://doi.org/10.1175/JCLI3969.1

Mätzler, C. (1987). Applications of the interaction of microwaves with the natural snow cover. *Remote Sensing Reviews, 2*(2), 259-387. https://doi.org/10.1080/02757258709532086
Mendoza, P. A., Shaw, T. E., McPhee, J., Musselman, K. N., Revuelto, J., & McDonell, S. (2020). Spatial distribution and scaling properties of lidar-derived snow depth in the extratropical Andes. *Water Resources Research*, 56(12), e2020WR028480. https://doi.org/10.1029/2020WR028480

Mernild, S., Liston, G., Hiemstra, C., Malmros, J., Yde, J. C., & McPhee, J. (2017). The Andes Cordillera. Part I: Snow distribution, properties, and trends (1979–2014). *International Journal of Climatology*, 37(4), 1680–1698. https://doi.org/10.1002/joc.4804

Meroni, M., D'Andrimont, R., Vrieling, A., Fasbender, D., Lemoine, G., Rembold, F., et al. (2021). Comparing land surface phenology of major European crops as derived from SAR and multispectral data of Sentinel-1 and 2. *Remote Sensing of Environment*, 253, 112232. https://doi.org/10.1016/J.RSE.2020.112232

Nagler, T., Rott, H., Ripper, E., Bippus, G., & Hetzenecker, M. (2016). Advancements for snowmelt monitoring by means of Sentinel-1 SAR. *Remote Sensing*, 8(4), 348. https://doi.org/10.3390/RS8040348

Niu, G.-Y., Yang, Z.-L., Mitchell, K. E., Chen, F., Ek, M. B., Barlage, M., et al. (2011). The community Noah land surface model with multi-parameterization options (Noah-MP): 1. Model description and evaluation with local-scale measurements. *Journal of Geophysical Research*, 116(D12), 12109. https://doi.org/10.1029/2010JD015139

Painter, T. H., Berisford, D. F., Boardman, J. W., Bormann, K. J., Deems, J. S., Gehrke, F., et al. (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. https://doi.org/10.1016/j.rse.2016.06.018

Patil, A., Singh, G., Rudiger, C., Mohanty, S., & Kumar, S., & Snehmani (2020). A novel approach for the snow water equivalent retrieval using X-band polarimetric synthetic aperture radar data. *IEEE Transactions on Geoscience and Remote Sensing*, 59(5), 3753–3763. https://doi.org/10.1109/tgrs.2020.3016527

Peel, M. C., Finlayson, B. L., & McMahon, T. A. (2007). Updated world map of the Köppen-Geiger climate classification. *Hydrology and Earth System Sciences*, 11(5), 1633–1644. https://doi.org/10.5194/hess-11-1633-2007

Planet Team. (2017). Planet application program interface: In space for life on Earth. Retrieved from https://api.planet.com/

Pulliainen, J., Luojus, K., Derksen, C., Mudryk, L., Lemmetyinen, J., Salminen, M., et al. (2020). Patterns and trends of Northern Hemisphere snow mass from 1980 to 2018. *Nature*, 581(7808), 294–298. https://doi.org/10.1038/s41586-020-2258-0

Rivera, D., Godoy-Faúndez, A., Lillo, M., Alvez, A., Delgado, V., Gonzalo-Martín, C., et al. (2016). Legal disputes as a proxy for regional conflicts over water rights in Chile. *Journal of Hydrology*, 535, 36–45. https://doi.org/10.1016/J.HYDROL.2016.01.057

Rott, H., Yueh, S. H., Cline, D. W., Duguay, C., Essery, R., Haas, C., et al. (2010). Cold regions hydrology high-resolution observatory for snow and cold land processes. *Proceedings of the IEEE*, 98(5), 752–765. https://doi.org/10.1109/JPROC.2009.2038947

Saavedra, F. A., Kampf, S. K., Fassnacht, S. R., & Sibold, J. S. (2017). A snow climatology of the Andes Mountains from MODIS snow cover data. *International Journal of Climatology*, 37(3), 1526–1539. https://doi.org/10.1002/joc.4795

Saavedra, F. A., Kampf, S. K., Fassnacht, S. R., & Sibold, J. S. (2018). Changes in Andes snow cover from MODIS data, 2000–2016. *The Cryosphere*, 12(3), 1027–1046. https://doi.org/10.5194/tc-12-1027-2018

Salazar, Á., Thatcher, M., Goubanova, K., Bernal, P., Gutiérrez, J., & Squeo, F. (2024). CMIP6 precipitation and temperature projections for Chile. *Climate Dynamics*, 62(3), 2475–2498. https://doi.org/10.1007/s00382-023-07034-9

Servicio Nacional de Geología y Minería (SERNAGEOMIN). (2021). *Anuario de la Minería de Chile* (Tech. Rep.). Servicio Nacional de Geología y Minería.

Shaw, T. E., Caro, A., Mendoza, P., Ayala, Á., Pellicciotti, F., Gascoin, S., & McPhee, J. (2020). The utility of optical satellite winter snow depths for initializing a glacio-hydrological model of a high elevation, Andean catchment. *Water Resources Research*, 56(8), e2020WR027188. https://doi.org/10.1029/2020WR027188

Shaw, T. E., Deschamps-Berger, C., Gascoin, S., & McPhee, J. (2020). Monitoring spatial and temporal differences in Andean snow depth derived from satellite tri-stereo photogrammetry. *Frontiers in Earth Science*, 8, 579142. https://doi.org/10.3389/feart.2020.579142

Shaw, T. E., Gascoin, S., Mendoza, P. A., Pellicciotti, F., & McPhee, J. (2020). Snow depth patterns in a high mountain Andean catchment from satellite optical tri-stereoscopic remote sensing. *Water Resources Research*, 56(2), e2019WR024880. https://doi.org/10.1029/2019WR024880

Shi, J., & Dozier, J. (2000). Estimation of snow water equivalence using SIR-C/X-SAR, Part II: Inferring snow depth and particle size. *IEEE Transactions on Geoscience and Remote Sensing*, 38(6), 2475–2488. https://doi.org/10.1109/36.885196

Small, D. (2011). Flattening gamma: Radiometric terrain correction for SAR imagery. *IEEE Transactions on Geoscience and Remote Sensing*, 49(8), 3081–3093. https://doi.org/10.1109/TGRS.2011.2120616

Smith, T., & Bookhagen, B. (2016). Assessing uncertainty and sensor biases in passive microwave data across High Mountain Asia. *Remote Sensing of Environment*, 181, 174–185. https://doi.org/10.1016/j.rse.2016.03.037

Strozzi, T., & Matzler, C. (1998). Backscattering measurements of Alpine snowcovers at 5.3 and 35 GHz. *IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing*, 36(3), 838–848. https://doi.org/10.1109/36.673677

Sturm, M., & Liston, G. E. (2021). Revisiting the global seasonal snow classification: An updated dataset for earth system applications. *Journal of Hydrometeorology*, 22(11), 2917–2938. https://doi.org/10.1175/JHM-D-21-0070.1

Takala, M., Luojus, K., Pulliainen, J., Derksen, C., Lemmetyinen, J., Kärnä, J. P., et al. (2011). Estimating northern hemisphere snow water equivalent for climate research through assimilation of space-borne radiometer data and ground-based measurements. *Remote Sensing of Environment*, 115(12), 3517–3529. https://doi.org/10.1016/J.RSE.2011.08.014

Tarricone, J., Webb, R. W., Marshall, H.-P., Nolin, A. W., & Meyer, F. J. (2023). Estimating snow accumulation and ablation with L-band interferometric synthetic aperture radar (InSAR). *The Cryosphere*, 17(5), 1997–2019. https://doi.org/10.5194/tc-17-1997-2023

Tennant, C. J., Harpold, A. A., Lohse, K. A., Godsey, S. E., Crosby, B. T., Larsen, L. G., et al. (2017). Regional sensitivities of seasonal snowpack to elevation, aspect, and vegetation cover in western North America. *Water Resources Research*, 53(8), 6908–6926. https://doi.org/10.1002/2016WR019374

Tsang, L., Durand, M., Derksen, C., Barros, A. P., Kang, D. H., Lievens, H., et al. (2022). Review article: Global monitoring of snow water equivalent using high-frequency radar remote sensing. *The Cryosphere*, 16(9), 3531–3573. https://doi.org/10.5194/TC-16-3531-2022

Tuanmu, M. N., & Jetz, W. (2014). A global 1-km consensus land-cover product for biodiversity and ecosystem modelling. *Global Ecology and Biogeography*, 23(9), 1031–1045. https://doi.org/10.1111/geb.12182

USDA Natural Resources Conservation Service. (2022). SNOpack TELEmetry Network (SNOTEL). Retrieved from https://data.nal.usda.gov/dataset/snowpack-telemetry-network-snotel

Veloso, A., Mermoz, S., Bouvet, A., Le Toan, T., Planells, M., Dejoux, J. F., & Ceschia, E. (2017). Understanding the temporal behavior of crops using Sentinel-1 and Sentinel-2-like data for agricultural applications. *Remote Sensing of Environment*, 199, 415–426. https://doi.org/10.1016/J.RSE.2017.07.015

Vicuña, S., Garreaud, R. D., & McPhee, J. (2011). Climate change impacts on the hydrology of a snowmelt driven basin in semiarid Chile. *Climatic Change*, 105(3–4), 469–488. https://doi.org/10.1007/s10584-010-9888-4
Vionnet, V., Brun, E., Morin, S., Boone, A., Faroux, S., Le Moigne, P., et al. (2012). The detailed snowpack scheme Crocus and its implementation in SURFEX v7.2. *Geoscientific Model Development*, 5(3), 773–791. https://doi.org/10.5194/GMD-5-773-2012

Viviroli, D., Dürr, H. H., Messerli, B., Meybeck, M., & Weingartner, R. (2007). Mountains of the world, water towers for humanity: Typology, mapping, and global significance. *Water Resources Research*, 43(7), W07447. https://doi.org/10.1029/2006WR005653

Voglimacci-Stephanopoli, J., Wendleder, A., Lantuit, H., Langlois, A., Stettner, S., Schmitt, A., et al. (2022). Potential of X-band polarimetric synthetic aperture radar co-polar phase difference for arctic snow depth estimation. *The Cryosphere*, 16(6), 2163–2181. https://doi.org/10.5194/tc-16-2163-2022

Vuyovich, C. M., Jacobs, J. M., & Daly, S. F. (2014). Comparison of passive microwave and modeled estimates of total watershed SWE in the continental United States. *Water Resources Research*, 50(11), 9088–9102. https://doi.org/10.1002/2013WR014734

Wrzesien, M. L., Pavelsky, T. M., Durand, M. T., Dozier, J., & Lundquist, J. D. (2019). Characterizing biases in mountain snow accumulation from global datasets. *Water Resources Research*, 55(11), 9873–9891. https://doi.org/10.1029/2019WR025350

Yang, K., Musselman, K. N., Rittger, K., Margulis, S. A., Painter, T. H., & Molotch, N. P. (2022). Combining ground-based and remotely sensed snow data in a linear regression model for real-time estimation of snow water equivalent. *Advances in Water Resources*, 160, 104075. https://doi.org/10.1016/J.ADVWATRES.2021.104075

Yueh, S. H., Dinardo, S. J., Akgiray, A., West, R., Cline, D. W., & Elder, K. (2009). Airborne Ku-band polarimetric radar remote sensing of terrestrial snow cover. *IEEE Transactions on Geoscience and Remote Sensing*, 47(10), 3347–3364. https://doi.org/10.1109/TGRS.2009.2022945

Zhu, J., Tan, S., Tsang, L., Kang, D., & Kim, E. (2021). Snow water equivalent retrieval using active and passive microwave observations. *Water Resources Research*, 57(7), e2020WR027563. https://doi.org/10.1029/2020WR027563

## References From the Supporting Information

Friedl, M., & Sulla-Menashe, D. (2019). MCD12Q1 MODIS/Terra+Aqua Land Cover Type Yearly L3 Global 500m SIN Grid V006 Land Cover Type 1. NASA EOSDIS Land Processes DAAC. https://doi.org/10.5067/MODIS/MCD12Q1.006

Gascoin, S., Hagolle, O., Huc, M., Jarlan, L., Dejoux, J. F., Szczypta, C., et al. (2015). A snow cover climatology for the Pyrenees from MODIS snow products. *Hydrology and Earth System Sciences*, 19(5), 2337–2351. https://doi.org/10.5194/hess-19-2337-2015

Hall, D. K., Riggs, G. A., DiGirolamo, N. E., & Román, M. O. (2019). Evaluation of MODIS and VIIRS cloud-gap-filled snow-cover products for production of an Earth science data record. *Hydrology and Earth System Sciences*, 23(12), 5227–5241. https://doi.org/10.5194/hess-23-5227-2019

Härer, S., Bernhardt, M., Siebers, M., & Schulz, K. (2018). On the need for a time- and location-dependent estimation of the NDSI threshold value for reducing existing uncertainties in snow cover maps at different scales. *The Cryosphere*, 12(5), 1629–1642. https://doi.org/10.5194/tc-12-1629-2018

RGI Consortium. (2017). *Randolph Glacier Inventory—A dataset of global glacier outlines, Version 6.0*. NSIDC: National Snow and Ice Data Center. https://doi.org/10.7265/4m1f-gd79

Riggs, G. A., Hall, D. K., & Román, M. O. (2019). VIIRS/NPP CGF Snow Cover Daily L3 Global 375m SIN Grid, Version 1. NASA National Snow and Ice Data Center Distributed Active Archive Center. https://doi.org/10.5067/VIIRS/VNP10A1F.001