Exploring the Fingerprints of Past Rain-on-Snow Events in a Central Andean Mountain Range Basin Using Satellite Imagery

Abstract

Rain-on-snow (ROS) events can alter nival regimes and increase snowmelt, peak river flow, and reduce water storage. However, detection of ROS events is challenging and only the most intense and obvious cases are identified. Rain is known to reduce snow cover and decrease near-infrared reflectance due to increased grain size. This study explored the fingerprints of ROS events on mountain snowpack with a simple typology that classifies changes in snow reflectance using fifteen years of MODIS imagery, reanalysis, and surface hydrometeorological data. The Maipo River Basin, with strong nival regime and a steep topography, in the western Andean mountain range was selected as a case study. Statistical analysis showed two distinct and opposite responses in the near infrared reflectance distribution of snow-covered pixels after precipitation, consistent with the typology for rain or snow events. For the probable ROS events, the daily maximum and minimum temperature increased in the days preceding the event and subsequently decreased, in some cases followed by a less consistent response in river flow. Although much remains to be studied, this approach can be used to expand historical records and improve modelling and detection schemes.

1. Introduction

1.1. Dynamics of Rain-on-Snow Events

In all mountain basins with a nival regime, temperature decreases with altitude, allowing snow accumulation from precipitation events that occur above the zero-degree isotherm (ZDI). The ZDI can change throughout the year, between storms, and even during individual storms. Changes in storm intensity, wind, or temperature can cause rain to fall in high elevations, modifying existing snow cover [1]. These events are known as rain-on-snow (ROS).

ROS events have complex generation mechanisms, influenced by snow depth, precipitation type (rain or snow), and air temperature [2,3], along with other meteorological variables and the interaction with topography. ROS events are also dependent on the amount and temperature of the raindrops that fall on the snow [4]. Much of the dynamics that precede and trigger a ROS event remain unclear, but their effects have been studied around the world.

Seasonal snow suffers structural changes know as metamorphism, where individual snow grains clump and compact over time, creating larger grains with less air between them. ROS enhances snow metamorphism [5], adds energy to the snowpack [6], increases water content and melting rate [7,8], reduces snow reflectance, and alters the energy balance that maintains low surface temperatures [9]. Even for a deep snowpack, rain forms a more compact and shallow cover, reducing thermal resistivity by up to a factor of eight, thus increasing snowmelt [10]. Short-term hydrological effects include peak-flow surges [11], increased flooding during storms [12,13], and avalanches [14,15].

ROS is common in spring [4,16] where it accelerates melting and sometimes triggers early peak flow and flooding [13]. ROS is less common in winter or fall but can modify snow properties [17] altering local hydrological regimes [13]. The more intense ROS events are widely recognized to dominate much of the flood generation in mountainous and boreal regions [15,18]. Most of the largest floods in British Columbia, Washington, Oregon, and California have been associated with ROS [3,19]. In Germany, ROS has a much larger hydrological influence than rain alone in catchments above 400 meters above mean sea level (m.a.m.s.l.) and are a primary cause of changes in channel morphology due to erosion [20,21,22]. In the Oregon Cascades, 85% of all landslides are associated with snowmelt during rainfall [20]. A combination of rainfall and snow melt has been reported to trigger debris flows in Norway [23]. In the Ésera river in the Pyrenees, Northern Spain, a large flooding was caused by torrential rains from 18 to 19th of June 2013, where snowmelt was estimated to contribute 33% of the total runoff during the main flood and 23% of the peak flow [24]. In the long term, frequent low-intensity ROS can change the distribution of river-flow of a given year, with earlier peak river-flows, less spring-melt, and lower summer flow [7,16,25,26,27]. many studies have been carried out for the more extreme ROS events where large floods occurred [2,16,28]. In the absence of systematic monitoring and records, previous studies have focused on the long-term changes using observation datasets in relation to rising air temperatures as suggestive of occurrence [3,29,30]. Increased occurrence of ROS has been detected in Eurasia [3,31,32] and expected to increase in the Northeast United States [7,16,25,32].

Changes in ROS events from climate change have been explored using GCM (Global Circulation Models) runs and regional climate models for limited areas [33,34,35]. Increasing atmospheric temperatures may decrease the frequency of ROS events in mid-latitude regions [3] by reducing snow cover extent and duration [36]. In mountainous and high-latitude regions like the Andes, ROS events are projected to increase due to increasing rainfall days and occurrence of rain at high-elevations [34]. The mountain range holds a significant area above the ZDI where snow can be accumulated that is vulnerable to ROS events [37]. However, detailed regional scale studies predicting how ROS events will change under future climate are not yet available for most parts of the world that face such risks [36].

For the Western Andean slopes, ROS is known to occur, but no systematic monitoring exists. However, an increase of 127 meters in the ZDI has been observed during the last quarter century [38,39,40,41], likely influencing glacier thinning [42]. By the end of the century, climate change projections anticipate a 10 to 30% decrease in precipitation [41] and a 2 °C increase in annual average temperature [43]. Other changes in precipitation can be expected, such as increased intensity, shorter duration, and the occurrence of early or late precipitation in the season. Prolonged changes in temperature and precipitation can alter the snow regime, with more frequent extreme events like winter ROS and an overall reduction of snow accumulation [13].

Observing and understanding extreme weather events like ROS is fundamental for the future and climate change adaptation. However, the study of ROS events is not that common as it is rarely being detected unless it has significant impact to populated areas. Furthermore, ROS occurs mostly in high-elevation areas where monitoring can be difficult [3,31]. High-elevation weather stations in mountains are scarce, especially for the western Andes. Furthermore, these events can be sudden and sporadic, easily missed by an individual station in the wrong location.

Field measurements of snow properties are sparse or infrequent, often subject to availability of instruments, personnel duty cycles, and constrained by cloud cover or storms. In addition to the difficulties in observation, mountain weather with winter frontal storms is often very complex and dependent on specific interactions between storms, atmosphere, and terrain. Systematic monitoring of ROS from weather stations has been very difficult [44] and in most cases, historical records are mostly anecdotal. Only the most extreme events have been clearly identified by observing the ZDI and the impacts on the resulting river flow [3,17,21,23]. Reconstructing past events from this information can be challenging.

1.2. ROS Monitoring

Weather modelling and forecasting is a fundamental tool for the observation and early warning of extreme events. However, the correct determination of the precipitation phase over changing terrain is critical to identify an ongoing ROS event. Precipitation phase determination schemes rely on weather stations and often assume a constant temperature lapse rate with height and invariant atmospheric conditions [45]. These schemes can also ignore interactions between atmosphere and terrain, warm and cold air masses, or between water droplet and the lower atmosphere as it falls through [46,47,48]. Changes in the ZDI from interactions with terrain and the properties of the melting layer might not be captured in this approach [48]. Other overlooked interactions include: airmass boundaries that introduce a warm air layer over cold air that increases atmospheric melt energy, or earth surface heating or cooling a shallow layer of air allowing rain at cooler temperatures or snow at warmer temperatures [46].

Improving forecasts of rain-on-snow events is hindered by the lack of general knowledge about the character of local historical rain-on-snow events [3]. More complex approaches for weather modelling require a broad understanding of the specific mountain weather, knowledge of the initial hydrometeor phase [26], shape and size [49], and detailed atmospheric information, often not available [26,45]. Considering the importance and need for widespread weather monitoring and knowledge of the complex Andean mountain weather, currently, weather modelling is unlikely to accurately detect or forecast ROS events from strong frontal storms.

1.3. Remote Sensing Applications for ROS

Remote sensing is the most extensive technique to observe snow [44] and can regularly and safely provide maps of snow properties [50]. While optical satellites like the Moderate Resolution Imaging Spectroradiometer (MODIS) and Landsat are frequently used to map snow, to the best of our knowledge they have not been explored to identify ROS events.

The potential use of passive microwave (PMW) sensors to study snow properties has been proven useful and measurements of Temperature Brightness (TB) has been shown to be sensitive to wet snow [51,52,53]. These approaches are based on the detection of change of surface dielectric properties from the presence of water at the air-snow interface, leading to surface scattering in the horizontal polarization and increased emission [54]. Dolant et al. (2016) [55] and Langlois et al. (2017) [56] proposed an empirical detection approach from observed ROS events and ground based PMW measurements of TB to detect a wet snow layer. The research aimed to establish a link between snow properties and signals at 18.7 and 36.5 GHz at both vertical and horizontal polarization [57].

The algorithm was tested for the Advanced Microwave Scanning Radiometer for the Earth observing system (AMSR-E) at 12.5 km pixels. A second testing of this approach was presented for Alaska where ROS events were reconstructed using remote sensing AMSR microwave brightness temperature-based freeze-thaw retrievals (12.5 km resolution) and MODIS snow maps (500m resolution) with a high degree of success [37]. Although very promising, this approach can struggle to differentiate between melt events (wet snow) and ROS events and heavily affected by cloud cover [55]. More robust detection algorithms remain to be explored [51].

Multispectral optical satellites can estimate snow cover extent using the difference between the visible and near-infrared (NIR) using the Normalized Difference Snow Index (NDSI) [58], and also retrieve information of albedo, grain size, liquid water close to the surface, and temperature [50]. While snow reflectance is dependent on grain size, water content, surface roughness, depth, and impurities [59], near-infrared reflectance is mainly affected by grain size [60].

Snow metamorphosis occurs gradually and unidirectionally after a snow fall, causing grain size to increase over time. Increase in grain size can be observed by a near-linear decay in NIR reflectance, from high-reflectance fresh snow to lower-reflectance ice [37,50]. In seasonal snow, this means a potential tracking of new snow from increased NIR reflectance, which can be indicative of ROS. Rain has a similar but more drastic effect on snow, washing away thin layers of snow and forcing wet-snow metamorphosis. This will cause small grains to form large clusters in a matter of hours [10,17]. These clusters behave optically as large grains with low near-infrared reflectance, even after refreezing [50,61]. Observing the reflectance of snow in the near infrared to discern between natural and sudden metamorphism can help to identify the occurrence and extent of ROS events.

The objective of this work is to use near-infrared reflectance from readily available remote sensing data and radiation-reflectance principles [50] and develop a simple typology to discern the likelihood of events from the effects of liquid water on grain size against those of fresh snow. Our central hypothesis is that certain changes in snow extent and status, along with climatic conditions and hydrological responses, consistent with ROS events, can typified to discern the likelihood of ROS.

Although significant efforts are being done to detect events and study its effects and dynamics, we believe a clear dataset and understanding of events are the missing link between detection, modelling, and forecasting. To the best of our knowledge, NIR reflectance has not been explored in the context of ROS events and remote sensing has not been used in the reconstruction of likely past events. This approach has the potential to support and explore the typification of events from historical data, which will contribute to the study and eventual understating of these events for future detection schemes.

2. Materials and Methods

2.1. Study Area

In mountain regions, precipitation can infiltrate the soil, generate outflow, or accumulate as snow until melting, altering temporal water distribution [62]. This accumulation of water and delay in river flow is fundamental for Mediterranean climate regions like the Central Andean mountain range with rainy winters and dry summers when water demands are substantially larger.

The Pacific Westerlies transport most of the water to the Central Andes and regulate precipitation via frontal extratropical systems. Westerlies reach the south of the continent in winter and are diverted north to the central Andes along the mountain range, where precipitation is blocked by the southeast Pacific anticyclone (31°S), creating a Mediterranean climate, and causing aridity further north [63,64]. For the Central Andes, snow is seasonal and peak accumulation occurs in winter, while summer is dry and stable, blocking precipitation and driving snowmelt and river flow [65]. River flow patterns depend mainly on the snowpack characteristics and the energy balance during snowmelt [62,66,67], and the amount and extent of snow helps determine available water supply [68].

The Maipo River basin in central Chile (~33°S, ~70°W), is defined by a strong Mediterranean climate and nival regime, characterized by wet and cold winters and dry hot summers. The Maipo Basin is bounded north by the Aconcagua basin, south by the Cachapoal, east by the Andes and west by the Pacific Ocean. Water flows East to West, beginning in the Andes at 6570 meters above mean sea level (m.a.m.s.l.) and ends in the Pacific Ocean, 250 km away. Like other rivers in the area, the Maipo and its tributaries are short with steep slopes. Maximum precipitation occurs between July and September, with average annual precipitation close to 290 mm [69]. Precipitation amounts are, in general, larger towards south and east. The basin is about 15 thousand square kilometers and is comprised of two main tributaries, the Maipo and the Mapocho Rivers (Figure 1). The Mapocho sub-basin is considerably smaller and with lower attitude (2370 m.a.m.s.l.) than the Maipo sub-basin. The Mapocho sub-basin has less area above the ZDI and therefore less winter snow accumulation. Its fluctuations mostly depend on snowmelt, peak river flow tends to be earlier in spring with a maximum in November, and the extreme low-water flow occurs in the autumn months [63]. This highlights the contrast of nival influence between the sub-basins.

Average annual flow in the Maipo River is 4445 million cubic meters and average annual runoff discharged into the Central Valley, pre-mountain range, amounts to 2620 million cubic meters. River flow in winter is close to 40 cubic meters per second and rises to 170 cubic meters per second in the summer. Average annual flow can vary between 90 and 250 million cubic meters from dry to wet years [63]. Even though there are many factors involved in the precipitation-river flow dynamic, affecting how fast a basin responds to rain/snow events, the Maipo basin is good example of a small basin with a fast and direct response on hydrology. Going from 5500 to 560 m.a.m.s.l. in less than 60 km to the metropolitan region of Santiago and 250 km to the Pacific Ocean, the effects of rain or snow in the upper part of the basin where rain transitions to snow should be easier to discern.

2.2. Datasets

2.2.1. Moderate Resolution Imaging Spectrometer (MODIS) Data

The Moderate Resolution Imaging Spectroradiometer (MODIS, 1999 to present) operates on the National Aeronautics and Space Administration (NASA) Terra and Aqua platforms for a daily repeat cycle with 36 spectral bands ranging from 250 to 1000-meter resolution [70]. Spectral bands between the visible and infrared spectrum are frequently used in snow and sea ice mapping and bands in the short-wave infrared parts of the spectrum (i.e., bands 6 (1.6 mm) and band 7 (2.1 mm)) are adequate for snow/cloud discrimination [71,72,73].

The Terra MOD09A1 Surface Reflectance 8-Day L3 Global 500 m SIN Grid Version 6 product (8-Day composite, L3 processing, Global 500m resolution) provides an estimate of the surface spectral reflectance of Bands 1 through 7 with atmospheric correction to account for aerosols, gasses, and Rayleigh scattering [74]. The 8-day composite period selects the best pixel observed based on cloud and solar zenith criteria [70,74,75]. This product has series of significant improvements that are relevant to the analysis [70]: (1) Improved aerosol retrieval and correction algorithm and (2) refinements to the internal snow, cloud, and cloud shadow detection algorithms.

We selected the MODIS Surface Reflectance 8-Day L3 Global 500m SIN Grid V006 [70] to observe the NIR reflectance data from mostly-snow cover pixels without cloud cover. The main reasons are (i) its temporal frequency and time of acquisition, (ii) NIR band closer to the range of grain size sensitivity (1230–1250 nm), (iii) correction of low-view angle, cloud cover or shadow, aerosol loading and atmospheric gases, (iv) accounts for the possible bidirectional effect of surface roughness in off-nadir angles, specifically relevant for mountain terrain. Even though the 8-day window can potentially hide some events with fresh snow (True Negative), our priority in this analysis is to obtain a clear observation of True Positive events, to set the typology of their characterization. Snow metamorphism would not show any significant effects in the period up to 8 days between MODIS scenes, hence it is unlikely to cause confusion between normal snow metamorphism and the effects of rain, leading to a False Positive. All scenes were downloaded from the United States Geological Survey (USGS) Earth Explorer web portal from 2001 to 2015, providing 678 images.

Snow cover was detected using the NDSI algorithm and processing scheme currently used in NASA’s Earth Observing System (EOS) for the MODIS snow map product [71,73,76]. The NDSI [72,76,77] measures the ratio of reflectance of snow in the visible where is higher, and short-wave infrared (SWIR), where is lower. This algorithm is widely used in the identification of snow- and ice-covered pixels and their separation from cloud cover. Further details can be found in Hall & Riggs (2010) [72].

Pixels with NDSI > 0.4 are classified as snow [72,78,79], and can be assumed to contain snow of depth ≥ 2.5 cm [30,70,80]. NDSI is insensitive to various illumination conditions, atmospheric effects, viewing-geometric distortions, and terrain effects [50,68,71,72]. Making it especially useful for mountainous terrain and shadows.

To focus on seasonal snow and better observe the changes in NIR reflectance, avoiding confusion from patchy snow cover, water, vegetation, dark pixels or cloud cover, a set of pixel discrimination criteria was added to improve the likelihood of observing clear and unobstructed signals from snow:

Elevation and NDVI to exclude vegetation: We used the 30-m ASTER Global Digital Elevation Model (DEM) [81], selecting only pixels higher than 1000 m.a.m.s.l., excluding vegetation and limiting the study area to extents with seasonal snow cover. Normalized Difference Vegetation Index (NDVI) was also used in combination with the Normalized Difference Snow Index (NDSI) to exclude snow covered forest from the observation.

• Pixel fractional snow cover: The NDSI algorithm identifies pixels with over 66% fractional snow cover [50]. A higher threshold was set using the MB3U algorithm (model MB for 3 test sites [68], that uses the NDSI values in a regression to obtain fractional snow cover (MB3U algorithm), excluding pixels with less than 75% snow cover, and avoiding pixels with low snow cover that might reduce reflectance for reasons other than grain size.
• Elimination of dark or cloudy pixels: Pixels with values below 10% in band 4 and 20% in band 6 are identified as dark or cloudy, and thus removed from the analysis [70].
• Decaying seasonal snow: Seasonal snow-cover discrimination excluded data during summer-spring, when snow has the highest state of decay and least amount of snow (leaving 337 usable scenes).

2.2.2. Hydrometeorological Data

Four hydrometeorological parameters were selected to study ROS events along with the MODIS scenes: Minimum and Maximum daily temperature (before and after ROS event), Precipitation (triggering ROS event) and river flow (resulting from ROS event). A high-resolution precipitation and temperature dataset for hydroclimatic research in Chile was selected to use as reference data for interval of 2001–2015 (CR2MET, from the Center for Climate and Resilience Research) [82]. The CR2MET is a spatial dataset that contains precipitation, max and min temperatures for a rectangular grid of 0.05° (~5 km). For the CR2MET dataset, precipitation was based on ERA-Interim atmospheric reanalysis data, statistically downscaled from a ~70 km grid, accounting for local topography, calibrated with local observations of precipitation; the temperature products were built from integration of ground data, reanalysis and MODIS Land Surface Temperature [82].

The spatial precipitation in the original ERA-Interim precipitation product is quite coarse and only a few pixels will cover the snow-covered are of the Maipo basin. Although the final product is statistically downscaled, we cannot confidently assume where a precipitation occurred, only that it occurred in the basin. The same difficulty is found using weather stations, where stations are sparse and rarely found in high elevation areas where winter snow occurs. In both cases, precipitation is not directly observed in the mountain range.

We assume that precipitation occurs in the whole basin, supported with the knowledge of the dynamics of frontal storms in the region that usually cover several kilometers [63]. Values are analyzed visually to verify the occurrence of precipitation at the basin level.

River flow data is collected in situ by river gauges operated by the Chilean Directorate General of Water (DGA) and daily average and monthly average river flow were used. A quality control procedure was applied to remove stations with incomplete records within the period of analysis. All data is accessible from the CR2 Climate Explorer [83] and IDE Chile: Geospatial Data Infrastructure of Chile [84].

2.3. Method

Conceptually, both snow cover and NIR reflectance are expected to increase following snowfall, while an abrupt reduction of both after precipitation could signal the effect of ROS. To be able to detect an event, the response to a precipitation should be opposite to that of fresh snow and significantly different from normal snow metamorphism. The assumptions for this approach are as follows:

• Fresh snow follows a continuous metamorphism in the winter-summer span until the melting process begins.
• Snow metamorphism can be tracked by the observation of NIR reflectance of clear snow-covered pixels.
• The highest levels of metamorphism and lowest NIR reflectance should be observed closest to the snowmelt period, signaling the largest grain size and clumps.
• For the colder months, snow NIR reflectance should remain in relative high levels of reflectance in comparison to those of the snowmelt period due to the frequent addition of fresh snow and the small amount of time for metamorphism.

We compiled daily precipitation record at the basin level and aggregated them for the 8-day interval to match the MODIS scenes. This process was performed with the software Quantum GIS (QGIS) [85] and specific packages from the Geographic Resources Analysis Support System (GRASS) [86]. The MODIS Surface Reflectance products were used to observe the seasonal snow accumulation in the Central Andes for the Maipo River basin, by creating two variables that describe extent and state of snow, respectively:

• Snow-Cover Area (SCA) as the count of pixels with more than 75% snow cover, and
• Near-infrared Reflectance Snow Reflectance (NIR-SR) from Band 5 reflectance pixels (1.0 to 1.3 μm) from the SCA.

A simple typology of likely events was proposed using a 2 × 2 matrix to classify the change in SCA/NIR-SR associated with occurrence of precipitation since the last image (Figure 2).

Snow precipitation (P1S1) is expected to increase the snowpack in extent and state (SCA and NIR-SR, respectively), while rain would reduce both (P1S0). In cases of no precipitation, SCA and NIR-SR can stay unchanged (P0S1) or decrease slightly (P0S0) from the previous image, depending on temperature. P1S0 suggests the occurrence of ROS. To identify and discern the fingerprints of P1S0 from P1S1 events, we study the shapes of the NIR-SR distribution moments (Mean, variance, skewness, and kurtosis), comparing them between the four event types, using non-parametric tests (Mann-Whitney and Kolmogorov-Smirnov [87,88]).

P1S0 events were analyzed to assess changes in daily temperature (max, min), and river flow (immediate effects). Figure 3 describes the implications of each type of event in terms of changes in precipitation, temperature, SCA, NIR-SR, and river-flow.

In Figure 3, the responses on Precipitation, Temperature, SCA, NIR-SR, and River flow are matched for events of snowfall, rain, rain-on-snow, and combinations of events. Plus (+) and minus (−) signs indicate if an increase or decrease is expected from an event (parentheses are used when the response is not certain). Empty cells suggest no response is expected. In snowfall events, precipitation is accompanied by a decrease in temperature, increase in SCA and NIR-SR. Rain events are expected to be matched with a response in river flow. Rain-on-snow events would most likely be related to increases in temperature and reduced SCA and NIR-SR.

There are, however, important sources of error and uncertainty that need to be stated and further explored, including in situ calibration, bidirectional reflectance effects, cloud cover, frequency of observation and other sources of grain size increase. Without in situ grain size measurements and calibrations, we observed the changes on NIR reflectance after precipitation as the most significant response to rain or snow, suggesting the changes in grain size from metamorphism. The influence of surface roughness on the bidirectional reflectance of snow for viewing geometries over rough terrain are mostly negligible for near-nadir but not for off-nadir sensors [89,90,91], Thus a comparison between a series of images needs to account for time of collection. Only the surface snow metamorphism is observable by remote sensing and new snow will reset the observable metamorphism. In the near infrared, the snow is effectively semi-infinite when depth exceeds approximately 3 cm, but different bands can reach different snowpack depth and respond to different grain sizes [50,92], therefore, this approach needs to build upon single band measurements, to eventually explore a multi-band ratio approach.

3. Results

3.1. Likelihood of Rain-on-Snow Events

We used 678 MOD09A1 scenes covering the period 2001 to 2015 for the Maipo basin, to obtain a dataset of SCA, NIR-SR and accumulated precipitation, all in 8-day intervals. Each timestep was classified using the simple typology from Figure 2. For a closer look at the individual analysis of timesteps, Figure 4 shows an example of the arranged datasets to make specific events easy to spot. This figure shows the different SCA/NIR-SR responses to precipitation for the Maipo basin. The Figure is organized as follows:

• X-axis: From left to right, the X-axis is the same for both figures and represents the time steps in an 8-day interval to match the satellite products.
• Bottom figure: the figure uses a heatmap to display NIR-SR in the Y-axis from the lowest reflectance (0) at the bottom and the highest at the top (1). The NIR-SR histogram for each scene is displayed vertically in the heatmap using the color gradient, from white to blue, where the blue means a high frequency of pixels in that reflectance level. The grey pixels are null values where no NIR-SR was found.
• Top figure: two variables are displayed, using two distinct types of graphs: (a) the red line displays SCA as the number of pixels classified as snow, and (b) the blue bars display accumulated precipitation over 2 mm (selecting amounts over 2 mm and above 1000 mamsl to discard low-impact rain).

The patterns of snow accumulation over time can be observed by comparing the shift from one timestep to the next along the X axis. For most cases, precipitation is accompanied by a shift in snow reflectance and cover area (SCA and NIR-SR), consistent with a nival regime. For the year 2012, two events are highlighted where two distinct responses are clearly observed:

Two events are highlighted in Figure 4. The first one, when compared to the previous timestep, shows a normal response for P1S1 with increased high reflectance and snow cover, found in most events of precipitation during winter. However, we found some cases of reduced SCA and NIR-SR after a precipitation event classified as P1S0, highly suggestive of ROS events. The number of P1S0 cases vary between 3 and 10 events a year. However, this is not a clear indication of occurrence since this approach only looks at clear events in an 8-day window. In Figure 5, we take a closer look at the same two scenes from Figure 4 where distinctly different typologies are observed from their resulting snow reflectance distribution (NIR-SR).

The first event, classified as P1S1, is displayed in Figure 5 in a darker color. The second event is displayed in a lighter color and corresponds to the next timestep. The first event displays a normal-looking distribution of NIR-SR ranging between 20% and 40% reflectance, with a symmetric unimodal curve (skew and kurtosis close to zero). The second event is much more skewed, ranging between 10% and 30% reflectance, it has both lower mean and variance. Along with a reduction in SCA, the second event (displayed as P1S0) shows significant changes in the snowpack during a small window of time. These changes are very inconsistent with a snowfall but highly suggestive of ROS.

The mean, variance, skew, and kurtosis for NIR-SR were compared between events, focusing on the differences of snow (P1S1), melting (P0S1 and P0S0), and rain (P1S0). P1S1 events showed and overall higher mean and standard deviation, lower skewness, and kurtosis, the closest to a normal distribution. In P0 events (no precipitation), snow shows either slight changes in moments with no large deviations from fresh conditions (P0S1) or stronger metamorphism (P0S0), often found at the end or beginning of winter. Moments for P1S0 event distribution (likely ROS) show strong changes in NIR-SR, even after precipitation that ought to add fresh snow.

Without assuming a normal distribution for the Mean, Variance, Kurtosis and Skewness of P1S1 and P1S0 events, each distribution moment was compared between types of events using both the Mann-Whitney U and Kolmogorov-Smirnov tests. These are non-parametric tests [49,87] used to determine if two sets of data are significantly different, in this case, P1S1 and P1S0 events (

Table 1. Test of similarity between Snow (P1S1) and rain (P1S0) reflectance distribution moments (W: Mann-Whitney U test; D: Kolmogorov-Smirnov test. * p < 0.0001).

Moment	Statistic	Value
Mean	W D	1647.5 * 0.5695 *
Standard Deviation	W D	2853.5 * 0.365 *
Skew	W D	8647.5 * 0.4682 *
Kurtosis	W D	8336.5 * 0.4297 *

).

High W values for the Mann-Whitney U test and close to 0 values for the D Kolmogorov-Smirnov test indicate the goodness of fit and suggest P1S0 and P1S1 belong to two different types of precipitation. The p-value for all four moments is below 0.0001, rejecting the null hypotheses of P1S1 being equal to P1S0. The results show two significantly different reflectance distributions for precipitation events, one closely related to snow precipitation (P1S1) and the other the likely result of ROS (P1S0). Higher NIR-SR variance and mean, with a lower skew and kurtosis, closely resembling a normal distribution around a reflectance of 30%, are characteristic of P1S1 events. P1S0 events had significantly lower means and variance, with higher skew and kurtosis.

3.2. Characteristics and Patterns

Unlike the previous precipitation data that was aggregated for the 8-day interval, here we used daily precipitation around the 8-day interval to observe changes in minimum and maximum temperature before and after the event. For most P1S0 events in the database, at least one precipitation day inside the 8-day window showed increased temperature leading up to it. We labeled such days as day-0 and tested the changes in daily maximum and minimum temperature on surrounding days for significant differences using t-tests (

Table 2. 5-day interval temperature change test for likely rain-on-snow events (P1S0).

Variable Test	2 Days Prior	1 Day Prior	1 Day Later	2 Days Later
Maximum temperature	−14%	−17%	−16%	−13%
	0.1324	0.003143	0.06251	0.03397
Minimum temperature	−42%	−31%	−26%	−46%
	0.02515	0.04885	0.03224	0.004365

).

Both daily maximum and minimum temperature were observed to have increase briefly, leading up to an event, and subsequently decrease. Minimum temperature shows a higher increase, suggesting a lower temperature amplitude during warmer rain events. A temperature increase, by itself, could directly cause the changes observed in snow for P1S0 events, but this is very unlikely for a snowpack between 4000 and 6000 m.a.m.s.l., at least 2000 meters above the winter ZDI.

River-flow changes were tested using seven gauging stations situated above 1000 m.a.m.s.l. in the Maipo basin, comparing the two major sub-basins (the Mapocho and Maipo). The Mapocho sub-basin is smaller in size and elevation, closer to the ZDI, and further north in the aridity gradient of Central Chile. The effects of the transition to rain between these sub-basins should be different. A four-day interval was used to analyze and test for differences in river-flow between the day prior and the two days following each precipitation event (

Table 3. Test of 4-day interval changes in river-flow for likely rain-on-snow events (P1S0).

Sub-Basin	Station Name	1 Day Prior	1 Day Later	2 Days Later
Mapocho	Rio Mapocho en los Almendros	0.19%	6.9%	3%
		0.9206	0.08302	0.1957
Mapocho	Estero Yerba Loca antes Junta San Francisco	1.26%	1.6%	1.6%
		0.214	0.2015	0.3089
Mapocho	Estero Arrayan en La Montosa	−0.23%	3.27%	0.95%
		0.8631	0.4181	0.7648
Maipo	Rio Colorado antes Junta Rio Maipo	2.6%	3.7%	5.78%
		0.4237	0.2094	0.03362
Maipo	Rio Colorado Antes Junta Rio Olivares	0.67%	0.2%	2.4%
		0.6462	0.9103	0.4772
Maipo	Rio Maipo en El Manzano	−0.06%	3.2%	2.4%
		0.9673	0.1456	0.1569
Maipo	Rio Maipo en Las Hualtatas	0.3%	3.9%	4.4%
		0.7639	0.0327	0.06183
Maipo	Rio Maipo en San Alfonso	−0.6%	2.9%	4.6%
		0.6589	0.05487	0.008764
Maipo	Rio Olivares antes Junta Rio Colorado	−1.5%	1.5%	−0.37%
		0.543	0.562	0.9045
Maipo	Rio Volcan en Queltehues	1.12%	8.17%	8.9%
		0.6968	0.03065	0.06366

).

While no responses were observed for the Mapocho sub-basin, positive responses (albeit non-significant) were found in the Maipo sub-basin. Only three stations (Rio Maipo en las Hualtatas, Rio Maipo en San Alfonso, Rio Volcan en Queltehues), showed significant increases in river-flow after day-0. The results suggest that no immediate general response to river-flow is caused by ROS, but large peak-flows can occur in subsequent days in certain cases. For the Mapocho sub-basin some outliers, whose river flow values were over 300% increase after an event, had to be removed from the analysis because the basin is located at lower elevation and shows pluvial events early and late in the winter. This is not suggestive of ROS all by itself, but when considered in the context of the changes in extent and state of snow presented by the typology in Figure 2, it adds validity the claim that P1S1 and P1S0 events are indeed different from each other.

4. Discussion

Based on the proposed typology [10,50,61,66], between 2001 and 2015 in the Maipo River basins, remote sensing of extent (SCA) and state of snow (NIR-SR) showed two different responses to precipitation. Snowfall was identified by an increase SCA and NIR-SR after a precipitation, with NIR-SR 20 to 40% reflectance, low kurtosis, and skewness. Likely ROS events were identified by a reduction in SCA and NIR-SR after a precipitation, with a NIR-SR below 30%, high skew and kurtosis. This is consistent with the theory of the effect of liquid water on snow [10,17,50,65]. There is a clear difference in the distribution moments between fresh snow and likely ROS, and between likely ROS and normal metamorphism. The distribution moments of likely ROS events are more similar to the snowmelt metamorphism observed at the end of the season.

Moreover, some changes in weather, and hydrology suggested by theory [7,16,27,55] were observed, before and after the occurrence of the selected ROS events. A significant increase in maximum and minimum temperature was observed in the two days prior to most ROS events, followed by a decrease in both. A temperature increase, by itself, could cause the changes observed in snow for P1S0 events, but this is very unlikely for a snowpack between 4000 and 6000 m.a.m.s.l., at least 2000 meters above the winter ZDI. This temperature increase prior to the event suggests the possibility of eventual prediction.

No clear immediate river-flow response was observed. But the more dramatic effects (like river flow surges) commonly associated with ROS, are mostly found in the more extreme and easily observed events. We believe these events more frequent than perceived but less intense, often going unnoticed.

To further discuss the possibility of cumulative effects, a simple comparison was done for the spring-summer river-flow of years with more frequent (likely) ROS events and years with less. This comparison was applied to the sub-basins to observe the difference in response between them. The Mapocho sub-basin is smaller in size and elevation, closer to the ZDI, and further north in the aridity gradient of Central Chile. It is also less dependent on snow accumulation. We examined the effects of high (≥10 P1S0/year) and low (5–6 P1S0/year) occurrence of P1S0 on snowmelt using hydrographs for two stations (Mapocho sub-basin: Los Almendros and Maipo sub-basin: El Manzano; Figure 6).

There are some clear differences in river-flow between years with high and low occurrences of P1S0 events. With a mostly pluvial regime, the Mapocho river has several peak events before summer melting (Figure 6a). The Maipo sub-basin shows a normal-shaped hydrograph with a unimodal peak, more closely related to snowmelt-driven hydrology (Figure 6b). This is likely due to the lower altitude of the Mapocho sub-basin in relation to the ZDI, causing immediate snowmelt with less changes in the later melting period ROS events affect both sub-basins, with higher peak values associated with more ROS events. The Maipo shows a more skewed hydrograph and slightly earlier rise in river-flow with more ROS events. However, the number of detected events is not a clear indication of occurrence since this approach only looks at clear events in an 8-day window and requires a more detailed analysis before exploring cumulative effects of ROS.

There are, however, important sources of error and uncertainty that need to be stated and further explored. The lack of high-elevation weather stations meant that this analysis had to be based on lower altitude data (but over 1000 m.a.m.s.l.). Remote sensing and reanalysis spatial products for precipitation are too coarse and not a good indicator of weather conditions in the Andes. For the analysis, no direct or remote observation of precipitation was done on the mountain range. We build on the assumption that winter frontal storms in central Chile ascend over the Andes, as suggested in the literature [63,93]. We observed snow (SCA) respond to most precipitation records in the low elevation stations, which supports the assumption for frontal storms. This is common problem around the world and any approach that aims at past events needs to overcome the limitation with sturdier analysis. If precipitation did not in fact reach the snowpack in the mountain, P0S1 or P0S0 events can be mistaken for P1S0. However, the abrupt changes observed in those events is significantly different to those of normal winter snow metamorphism. It is unlikely that precipitation detected in a lower elevation station was not involved in the changes observed in the snow further up in the basin. Rain remains the most consistent explanation for the observed rapid metamorphism.

Other sources of uncertainty are the specific sensor limitations like bidirectional reflectance effects, cloud cover, frequency of observation and other sources of grain size increase. The influence of surface roughness on the bidirectional reflectance of snow for viewing geometries over rough terrain are mostly negligible for near-nadir but not for off-nadir sensors [89,90,91]. Thus, a comparison between a series of images needs to account for time of collection. The surface reflectance obtained from the MOD09A1 [70,72,76,77] products can still be affected be affected by atmospheric correction, water, aerosol loading, geometry, vegetation, and ageing of snow [70,74,75]. This study attempted to avoid some sources of error by only observing high elevation without significant vegetation and water bodies, while focusing on areas of seasonal snow with less influence from pollution. The main limitation of the MOD09A1 [70,71,72,77] is the 8-day interval that provides a composite of images that are likely to miss important changes in snow. Within this interval, rain effects can be covered by subsequent snowfall and events can go undetected, leading to False Negatives.

Without in situ grain size measurements and calibrations, we observed the changes on NIR reflectance after precipitation as the most significant response to rain or snow, suggesting the changes in grain size from metamorphism. Only the surface snow metamorphism is observable by remote sensing and new snow will reset the observable metamorphism. Near-infrared bands can reach up to approximately 3 cm, but different bands can reach different snowpack depth and respond to different grain sizes [50,92]. As an initial study it was decided to build upon single-bland measurements, but further studies should explore a sturdier multi-band ratio approach with in-situ hyperspectral calibration.

5. Conclusions

Under normal winter conditions precipitation is expected to contribute fresh snow with a high NIR reflectance that is gradually reduced over time due to snow metamorphism. This gradual reduction in NIR reflectance is reset after each winter snowfall, and only reaches significantly low reflectance levels at the end of winter when snow precipitation stops and temperature increases. However, we detected two different responses to winter precipitation with significantly different reflectance distributions. One of those is inconsistent with what is expected from a winter nival regime, and highly suggestive of the effects of rain on snow [7,16,25,32]. Moreover, most of these events were preceded by an increase in maximum and minimum temperature, followed by a decrease in both. Changes in river flow were less consistent but only expected from intense ROS events where snowmelt is accelerated. River flow surges are commonly associated with ROS in the more extreme and easily observed events. We believe ROS events more frequent than perceived but less intense, often going unnoticed.

The aim, however, is not a detection scheme for all past events but to discern obvious rain from snow, contributing to the historical records of ROS events in poorly monitored areas. This study was focused on True Positives for past events, and a lot of exclusion had to be applied that needs to be addressed to avoid False Positives if this approach can support future detection schemes. Cloud remains a source of confusion if applied for detection [58], since it is causes high reflectance in the NIR, it can mask the viewing of snow and lead to False Negatives. Geometry, cloud cover and aerosol loading might still be a source of confusion that needs to be address in further studies to better define the range of values for NIR-SR [70,74]. Future detection schemes will require in situ validation, daily products and higher resolutions, multi-band indexes and multi-sensor integration, all validated against confirmed ROS events.

Significant efforts are being done to detect events and study its effects and dynamics. However, we believe a clear dataset and understanding of events are the missing link between detection, modelling, and forecasting. To the best of our knowledge, NIR reflectance has not been explored in the context of ROS events and remote sensing has not been used in the reconstruction of likely past events. This research presents a deductive approach to study rain-on-snow by exploring and highlighting likely past rain events on the Andean snowpack. The typology presented can be helpful in discerning the effects of snow from those of rain and natural melt. The responses in snow cover and reflectance to precipitation show that the designated typology and approach are useful to study rain-on-snow events in the absence of observation records, contributing valuable records and datasets to expand understanding. There are clear signs that ROS occur in the study area, even if they do not cause significant impact and can go unnoticed. This approach has the potential to support and explore the typification of events from historical data, which will contribute to the study and eventual understating of these events for future detection schemes and assessing climate change impacts on nival regimes. Clearly, a lot of work remains to be done but a clarity of the fingerprints left by rain on snow open the door for the next steps.