Changes in the frequency of global high mountain rain-on-snow events due to climate warming

Rain-on-snow (ROS) events can trigger severe floods in mountain regions. There is high uncertainty about how the frequency of ROS events (ROS) and associated floods will change as climate warms. Previous research has found considerable spatial variability in ROS responses to climate change. Detailed global assessments have not been conducted. Here, atmospheric reanalysis data was used to drive a physically based snow hydrology model to simulate the snowpack and the streamflow response to climate warming of a 5.25 km2 virtual basin (VB) applied to different high mountain climates around the world. Results confirm that the sensitivity of ROS to climate warming is highly variable among sites, and also with different elevations, aspects and slopes in each basin. The hydrological model predicts a decrease in the frequency of ROS with warming in 30 out 40 of the VBs analyzed; the rest have increasing ROS. The dominant phase of precipitation, duration of snow cover and average temperature of each basin are the main factors that explain this variation in the sensitivity of ROS to climate warming. Within each basin, the largest decreases in ROS were predicted to be at lower elevations and on slopes with sunward aspects. Although the overall frequency of ROS drops, the hydrological importance of ROS is not expected to decline. Peak streamflows due to ROS are predicted to increase due to more rapid melting from enhanced energy inputs, and warmer snowpacks during future ROS.


Introduction
Rain-on-snow (ROS) is behind of many of the most damaging floods in mountain areas and its estimation requires knowledge not only of the mountain snowpack, but also of rainfall dynamics and alpine micrometeorology (Marks et al 1998, Vionnet et al 2020. The hydrological response of a catchment to ROS involves complex phenomena, as it depends on the specific meteorological and snowpack mass and energy conditions and snow-covered area during the event. This includes enhanced turbulent fluxes driven by temperature, wind and humidity and enhanced longwave irradiance, which are responsible for much of the extra melting energy associated with ROS as shortwave irradiance is reduced by cloudy conditions during ROS (Dadic et al 2013. Such factors also condition the precipitation phase (Jennings et al 2018) and the internal energy, liquid water content and mass of the snowpack at the time of the rainfall (Groisman et al 2006, Würzer andJonas 2018). The energy advected by rainfall itself is usually not the main driver of faster snowmelt during ROS (Marks et al 1998), as rain temperatures are usually not much higher than snowpack temperatures, and rainfall passes relatively quickly through preferential flowpaths in isothermal snowpacks (Leroux et al 2020). The exception is ROS onto cold snowpacks where refreezing and formation of ice layers can add significant energy to the snowpack and slow the release of advected rainfall from the snowpack (Leroux and Pomeroy 2019).
Global warming is changing high mountains and their cryospheric and hydrological components rapidly, and there are concerns about how natural hazards in mountains may evolve as air temperatures increase (Musselman et al 2017, Beniston et al 2018. However, the effect of increasing temperature on the frequency of ROS is difficult to determine. Even without considering the uncertainty of precipitation changes in warming mountain climates, the impacts on mountain hydrology driven only by increasing temperatures have been shown to be highly variable in space, altitude and time (López-Moreno et al 2020). This complexity is explained by the fact that higher temperatures and associated humidities shift the phase of precipitation (Harder and Pomeroy 2014), increasing the rainfall ratio and frequency of winter/spring rainfall, but also lead to a shorter snowcovered season, and thus reducing the period when ROS may occur. The balance of these two factors determines the magnitude and direction of ROS change with climate warming (McCabe et al 2007, Morán-Tejeda et al 2016 and is expected to vary with solar exposure on slopes and elevation and to be impacted by the persistence of alpine snow drifts in summer. Temporal trend analyses of observational data have revealed that different parts of the world, or distinct elevations of the same mountain region, have undergone differing responses in the frequency and intensity of ROS to climate warming. Increases (decreases) in the frequency of ROS have been found in the high (low) elevation zones of the western United States (McCabe et al 2007, Surfleet andTullos 2013) and in the Swiss Alps (Morán-Tejeda et al 2016). Seasonal changes are also reported; Freudiger et al (2014) analyzed ROS in the major basins of Central Europe for the period 1950-2010, and reported that the frequency increased for all elevations in January and February, but decreased in April and May. Large-scale analyses have reported latitudinal patterns in the changing frequency of ROS in the last decade, with more ROS likely at high latitudes, such as the circumpolar regions Fletcher 2007, Ye et al 2008). Sensitivity analyses have therefore reached different conclusions on the changes in mountain ROS frequency with warming, depending on the degree of warming and snowpack properties of the mountain region. López-Moreno et al (2016) reported a continuous increase in ROS in Ny-Ålesund (Svalbard, 79 • N) with a +1 • C to +5 • C change in temperature; while Beniston and Stoffel (2016) showed that floods caused by ROS may increase by 50% in the Swiss Alps with a temperature increase of 2 • C-4 • C, and decrease with temperature increases exceeding 4 • C because of the reduction in snowpack duration.
Despite an increasing interest in the study of ROS, their hydrological implications and their response to climate warming (Musselman et al 2018), there is still a lack of studies comparing the hydrological response of high mountain catchments found in various global mountain climates. This is partially due to the sparse meteorological and hydrological information in many high mountain areas, but also because of the difficulty in comparing processes in catchments that differ in size, hypsometry and land cover (Wayand et al 2015).
This study uses perturbations of a downscaled atmospheric reanalysis dataset used to force a physically based, spatially detailed snow hydrology model of idealized 'virtual basins' to determine the change in ROS and resulting hydrological response in 40 different high mountain areas of the world, representing most of the alpine climatic conditions existing on the planet. Virtual basins (VBs) (Weiler and McDonnell 2004, Armstrong et al 2008, López-Moreno et al 2020 are synthetic drainage basins whose properties reflect the typical spatial organization of basins within the region of interest. Here, a typical alpine headwater basin in a post-glacial mountain landscape was characterized by a VB to represent basins where snow hydrology is important and ROS can occur during seasonal snow cover. The approach permits (a) analysis of the regional sensitivity of ROS to climate warming associated with the counteracting effects of the increasing rainfall ratio and declining snow-covered period, (b) assessing seasonal changes in ROS occurrence as the climate warms, (c) identifying the main predictor variables that explain different regional responses, and (d) quantifying how melt rates during ROS change with warming climate. Exploiting a comparative VB approach, these points are examined under identical conditions of slope/aspect, spatial configuration of topography, land cover, and climate warming, the baseline climate being the only difference among them.

Climatic, snow and hydrological simulations in VBs
A VB comparative methodology was used to ensure that precipitation phase, snowpack dynamics and runoff differences among VBs in different mountain ranges were only due to initial climates. This allowed removal of the impact of other relevant factors for ROS, such as hypsometry, soil characteristics or land cover (Wayand et al 2015). Thus, a 'typical' small high mountain 'alpine' basin of 5.25 km 2 with a 1000 m vertical gradient was chosen to be the VB. In high mountain basins, sparse vegetation, shallow soils, and limited groundwater storage have small influences on hydrology, and wind redistribution by snow, sublimation, solar irradiance and snow-covered area depletion have major effects on hydrology (Fang and Pomeroy 2020). Similar VBs were used previously to analyze the sensitivity of snow accumulation and runoff from melting in 45 high mountain areas of the world (López-Moreno et al 2020). Figure 1 shows a Figure 1. VB study sites in mountain ranges around the world. Colours denote the minimum elevation of each basin. In the bottom-left corner is a representation of the VB and its seven HRUs that were modelled in this study.
representation of the VB, disaggregated into seven hydrological response units (HRUs): (a) a summit area of 0.25 km 2 and 30 • slope angle, (b) a high plateau of 0.5 km 2 and 10 • slope angle, (c) and (d) north and south-facing steep slopes of 0.5 km 2 each, with 25 • slope angles; (e) and (f) north and southfacing moderate slopes of 1.5 km 2 each, with 20 • slope angles; and (g) the mild westerly sloped bottom of the basin at 0.5 km 2 and 10 • slope angle. The soil depth was set to zero at the summit and increased progressively to 50 cm depth at the outlet of the basin. The high plateau and summit were barren, and short (10-15 cm high) meadow grass was the only vegetation included below these HRUs. The VB was 'placed' in different mountain areas of the world under contrasting climatic and snow characteristics (see supplementary ST1 (available online at stacks.iop.org/ERL/ 16/094021/mmedia)). The exact position and elevation of the basin was subjective, with the only requisite being to have a seasonal snowpack that completely melts every year under the current climate. Hence, temperate or low mountain elevations without seasonal snow cover and glaciarized mountain elevations were excluded from this study in order to focus on the impact of ROS on the seasonal snowpack and its hydrology.
Meteorological inputs were obtained from the WFDEI dataset generated in the framework of the WATCH project (www.eu-watch.org) corresponding to a bias-corrected temperature, specific humidity, surface pressure, wind speed, incoming shortwave radiation, and precipitation from ERA-Interim reanalysis for the period 1979-2012 (at 3 h basis), at 0.5 • spatial resolution (Weedon et al 2014). This original resolution is subsequently downscaled to each HRU.
The Cold Regions Hydrological Modelling platform (CRHM) is a flexible, modular, physically based hydrological modelling system that is suitable for snow hydrology (Pomeroy et al 2007, Ellis et al 2010. A flowchart of the different modules used in CRHM for this study is provided in appendix figure S1. Lapse rate gradients for temperature and precipitation (6.5 • C and 50% increase per 1000 m respectively), psychrometric adjustments for atmospheric humidity and precipitation phase and redistribution of wind fields and long-and short-wave radiation according to elevation and topography were used to lapse the input data from the elevation of the WATCH centroid to the elevations of each HRU. This solves the limitation of using an initial 0.5 • resolution of the forcing meteorological data that is insufficient to deal with mountain topographic effects on precipitation, temperature, humidity, radiation and wind fields. This is why the downscaling within CRHM is critical for applying these fields in mountain terrain to the different modules related in the CRHM platform permit to calculate the full range of hydrometeorological processes for each HRU and to aggregate to the basin level hydrological response using the VB. CRHM deployed the psychrometric energy balance method approach (Harder and Pomeroy 2013) to determine precipitation phase; the Prairie Blowing Snow Model (PBSM, Pomeroy and Li 2000) module to calculate blowing snow redistribution and sublimation fluxes; and the Snobal module (Marks et al 1998) to calculate energy balance snowmelt and track the snowpack mass and energy states. Albedo decay is based on the age of snow after the last snowfall, with values ranging between 0.95 and 0.5. CRHM's Evap, Soil and NetRoute modules were used to calculate evapotranspiration, infiltration, soil moisture storage, and subsurface and surface routing (DeBeer and Pomeroy 2017).
The presence of frozen soil and soil depth has a great influence on runoff generation; thus, a typical alpine configuration for our VBs was used. The configuration includes state variables that are responsive to climate warming and ROS events. For instance, the infiltration into frozen soil algorithm takes the soil moisture content from the end of the previous snow-free season to set initial conditions for calculating limited infiltration during the seasonal snowmelt. However, this limited state is adjusted to a restricted state when there is a major mid-winter melt or ROS event (>10 mm), as that can cause a basal ice layer to form at the bottom of the snowpack and restrict infiltration to frozen soils. This restriction of subsequent infiltration is one of the mechanisms by which ROS events can increase runoff dramatically if meltwater is calculated to reach the base of the snowpack during mid-winter. The model can capture this dynamical behaviour and the response of the subsurface hydrology varies with the climate regime and meteorological history of the snow season. The ground surface temperature was estimated using the Radiation-Convection-Conduction approach (Williams et al 2015), and freeze and thaw was estimated using the XG-algorithm dividing the soil into five layers for application of the Stefan Equation (Changwei and Gough 2013).
More details about the configuration of the CRHM model for this study can be found in (López-Moreno et al 2020). The aim of the simulations was not to reproduce the climate, snowpack and runoff for each mountain range exactly, but to ensure that coherent inputs represented the climates of the major snow-dominated mountain headwaters worldwide. This homogenization of the inputs permits a deeper understanding of the influence that climatic characteristics have and will have on ROS since the outputs of the simulations can be directly compared. A CRHM model having an almost identical configuration to the one used in this study was used to satisfactorily reproduce the snowpack and runoff, and to perform a sensitivity analysis to warming, over the same variety of environments considered in this study (López-Moreno et al 2020). CRHM has also been successfully used to analyze the energetic exchanges and flood generation under ROS events  and melt in a wide range of mountain headwater basins-from sub-arctic to cool climate (Rasouli et al 2019).

ROS events identification and sensitivities analysis
Similar to that stated by (Musselman et al 2018), an event was considered ROS when at least 10 mm of daily rainfall fell over a snowpack deeper than 10 cm. The number of ROS days was calculated for each individual HRU, and the average rain and snowmelt for each basin were computed as the areaweighted average of the HRUS. The frequency of ROS at each basin was calculated for the control period  and temperature was progressively increased by +1 • C intervals to +5 • C to each daily value in order to simulate the impact of different magnitudes of climate warming on ROS. The change in the frequency of ROS with warming was calculated for each degree, and the average was considered as the sensitivity of ROS frequency per • C. For the sensitivity analysis, relative humidity was held constant, allowing vapour pressure to rise with T. This influences the longwave radiation, precipitation phase, sublimation, evapotranspiration and snowmelt processes in CRHM.
A linear regression analysis was performed to assess the predictability of the spatial distribution ROS sensitivity. For this, the snow ratio (% percentage of precipitation as snow), snow duration and mean temperature of each VB were used as independent variables for stepwise multiple linear regressions (based on the Akaike Information Criterion, AIC). Other snow and climatic variables were discarded due to lack of explanatory capacity or co-linearity with the aforementioned variables. The adjusted R 2 for each model informed about the quantity of variance in the spatial distribution of ROS sensitivity to climate warming, whereas the beta-coefficients (β) informed about the relative contribution of each variable to the total predictability of the resulting model (Venables and Ripley 2002). R 2 was estimated by comparing our obtained ROS sensitivity at each basin with the predicted values from the linear model after applying a jackknife approach. This resampling technique is especially useful for variance and bias estimation. The estimator is calculated by sequentially deleting a single observation from the sample. The estimator is recomputed until there are n estimates for a sample size of n.
In order to estimate the sensitivity of snowmelt during ROS, the change in snowmelt quantity was calculated between the pairs of ROS that occurred at a given temperature and also under one degree of warming. This procedure was designed to avoid comparing the melt rates for differing numbers of ROS. Finally, the total runoff produced during ROS was compared in order to assess if increasing temperature leads tona rising or declining hydrological response to ROS.

Results
Figure 2(a) shows that the frequency of high mountain ROS under unperturbed conditions (no added warming) varies widely with geographical area. Mountains at mid-latitudes and under oceanic influence generally have more frequent ROS, in contrast with high latitude and continental climate sites. Thus, there are VBs included in this study with more than 10 ROS per year, while in other basins their occurrence is extremely rare during the period 1982-2014. Figure 2(b) shows that there are also noticeable differences in the frequencies of ROS within the basins. The ROS frequency increases with elevation, and there is a lower frequency on south-facing slopes (HRUs 4 and 6), compared to north-facing slopes at the same elevation (HRUs 3 and 5 respectively), likely due to increased snow cover persistence into late spring or summer on north-facing slopes and at higher elevations. Figure 3 shows the response of ROS as the air temperature is increased. Figure 3(a) shows the variability shown by VBs in the percentage of ROS that do not happen due to the disappearance of snow cover due to 1 • C increment, obtained from averaging the observed values from +1 • C to +5 • C. The rest of the events are the increasing ROS due to the change in phase from snowfall to rainfall over snow-covered ground. The balance of decreasing and increasing events delivers the total ROS sensitivity per • C of warming shown in figure 3(b). There are noticeable differences among sites, but in general the number of basins where the number of ROS decreases due to snow disappearance (75% of the total) exceeds that of basins where ROS increases due to a change from solid to liquid precipitation (26%). The mean change for the 40 basins is a decrease of 9% in the number of ROS per • C ( figure 3(b)). In some basins, the decrease in ROS numbers exceeded 20% per • C of temperature increase. One out of four sites considered in this study showed an increase in ROS numbers, an increase that reached 10% per degree in the most extreme cases. Figure 3 also shows that declines in ROS at higher elevations of the basins are moderate compared to the larger decreases at lower elevations (HRUs 5, 6 and 7) where warming causes temperatures to more frequently rise above 0 • C. Aspect also introduces some differences in the sensitivity of HRUs 3-4 compared to HRUs 5-6, with slightly larger declines in ROS for the most irradiated HRUs. In many of the VBs, the ROS frequency increased in the 'high mountain plateau' (HRU2) and the north-facing high-elevation slopes (HRU3) (38% and 41%, respectively). In contrast, the ROS frequency increased in only 21%, 19% and 18% of the basins for the low elevation HRUs 5, 6 and 7, respectively.   : 1979-2012) and for progressive increases in temperature by +1 • C (T1) intervals to +5 • C (T5). Months in VBs located in the Southern hemisphere have been shifted 6 months in order to correspond with cold and warm periods as in the Northern hemisphere. The centre line is the median, the box represents the 25th and 75th percentiles and bars the 10th and 90th percentiles. Figure 4 shows that a warmer climate causes changes not only in the frequency of ROS, but also in their seasonal distribution. In general terms, increasing temperatures lead to a sharp decrease in the importance to total ROS during late spring (mainly May) and a decrease in early winter (December). On the contrary, February, March and April exhibit an increasing frequency of ROS. However, the warmest scenarios (T + 4 and T + 5 • C) show increased variability among basins in March and April, with the frequency for the 25th and 10th percentiles becoming lower with warming. This suggests that the generally observed increase in ROS in late winter and early spring does not occur for most temperature basins where snow cover will have disappeared at this time of year. The median number of ROS remains similar in January for the different warming scenarios, but the 75th and 90th percentiles exhibit a marked increase. Figure 5 shows the relationship between the sensitivity of the ROS frequency to increasing temperature with the snowfall ratio (1-rainfall ratio), average snow cover duration and the mean temperature of each basin during the control period. These variables have been included in a stepwise regression model as statistically significant (p < 0.05) predictors of ROS sensitivity. The most important explanatory variable  is the snowfall ratio with a Beta coefficient, β = 0.59. This suggests that the basins where most precipitation falls as snow are the ones where ROS increases as the climate warms. Conversely, a decrease in ROS is observed in most of the basins where the snowfall ratio is less than 60%, recording the largest decreases for those sites where snowfall ratios are the lowest. Snow cover duration does not exhibit significant correlation with the sensitivity of ROS to climate warming when it is correlated alone ( figure 5(b)); however, it does contribute significantly (p < 0.05), along with the mean temperature of the basin to the explanation of the variance of ROS sensitivity (B = −0.11 and −0.51 respectively) with a sharper decrease in ROS frequency with increasing temperature for shorter duration snow covers and warmer mean temperatures. The three variables (snowfall ratio, snow cover duration and mean temperature) explain 66% of the total variance in the sensitivity of ROS.
The combination of the three predictor variables (snowfall ratio, snow duration and mean temperature during the control period) explains the geographical distribution of different sensitivities shown in figure 6. Thus, the largest decrease of ROS under warmer conditions (<−15% per • C) are found in Mediterranean climate mountains (Central Chile, South Africa, Australia, Morocco and some places around the Mediterranean Sea) where most of the precipitation occurs during the cold season, snowpack duration is rather short, and mean temperatures are among the highest. Changes in ROS for Mediterranean climate mountains are highly influenced by the snowfall ratio during the control period. Large decreases in ROS with warming are also found in northern New Zealand, and the most humid parts of the Himalayas. A more moderate decrease in the number of ROS (−5% to −15% per • C) is predicted in northern Chile and Patagonia, many mountains of southern and central Europe, the Caucasus and Hokkaido (Japan) and some mountains in North America; whilst no significant change or an increase in the number ROS is predicted in continental areas of North America, Anatolia and Himalayas and in the northernmost latitude (Yukon, Northern Quebec, Svalbard and Kamchatka) mountains. This is mainly because most of the precipitation during the snow covered period over these mountains is as snow, and well below the liquid/solid threshold.
Finally, figure 7(a) shows the average change per • C in snowmelt per ROS with increasing temperature (calculated only from pairs of events where ROS occurred when 1 • C increased), and its relationship with the previously quantified sensitivity of ROS frequency to temperature increases. In all basins, snowmelt during ROS increases as the temperature warms, although noticeable differences among basins occur (see boxplot on the right side). Thus, snowmelt increases on average by 16% per • C; basins within the 25th and 75th percentiles range between 12% and 18%, and the most extreme cases to 8% and 27% per • C. Variability in the sensitivity of snowmelt during ROS is positively correlated with the sensitivity of ROS (r 2 = 0.43), suggesting that the largest increases in melt happen in basins with lower snowfall ratio, higher mean temperature and shorter snow cover duration. Figure 7(b) shows the change in total runoff recorded at each basin during ROS per • C of warming, and its relation with sensitivity of ROS. It can be seen how the sensitivity of total runoff generated during ROS is very low and the sign of the change is not consistent with positive and negative values. In most of the cases the magnitude of the change is less than 1% per • C, and it is unrelated to the sensitivity shown by the frequency of ROS.

Discussion
Changes in the frequency of ROS driven solely by temperature increase are spatially complex, with strong differences among different mountain climates of the world. Temperature increases will also cause seasonal shifts in ROS occurrence, with a general tendency to decrease in late spring and to increase in late winter or early spring. In basins with shorter snowpacks, the decrease in ROS may already be evident in March and April. This complexity explains why previous research based on specific geographic areas, mostly from North America and the Alps, has not found consistency in the observed trends and future projections on the occurrence of ROS, as well as on the floods generated from these events (Hock et al 2019). After comparing 40 mountain areas around the world that exhibit a wide variety of mountain climates, the CRHM model predicts that 76% of the basins (30 of them) respond to climate warming with a decreasing frequency of ROS, but there are some sites for which little sensitivity or even an increase of ROS under a warmer climate were predicted. Larger decreases (increases) in ROS frequency occur at sites with lower (higher) snowfall ratio, shorter (longer) snow duration and higher (lower) average temperature. These three factors (snowfall ratio, snow cover duration, and average temperature) explain 66% of the total variance in the sensitivity of ROS events to climate warming. These factors were found to be highly sensitive to climate warming, especially snowfall ratio in those basins with a mean temperature above −8 • C for the period November-June (figure 2 in López-Moreno et al 2020); and they explain the latitudinal patterns found by Ye et al (2008), who reported an increase in ROS in high-latitude mountains in Eurasia, concomitant with a clear warming trend; or Cohen et al (2015), who also detected an increase of ROS at high latitudes in the circumpolar North over the period 1979-2014. In this study, the circumpolar North (>60 • N) of North America, Europe and Asia, as well as continental climate sites in North America, Anatolia-Caucasus and Asia were where ROS will remain constant or increase under warmer conditions. These sites are characterized by a low sensitivity of snowmelt to climate warming ; thus, the decline of snow duration with increasing temperature is attenuated compared to milder and more humid areas. Moreover, in these sites much of the precipitation is snowfall under the current climate, making it possible that many events occurring during the winter may turn from snowfall to rainfall with climate warming. Therefore, the relatively long-lasting snowpack in colder basins promotes increases in rainfall to occur over snow cover, despite climate warming. In contrast, coastal or very humid mountains (New Zealand, eastern and southern Himalayas) and mountains exhibiting Mediterranean climates (Central Chile, Australia, around Mediterranean basin) registered the most intense decreases in the number of ROS with warming. This is a consequence of the milder climates, leading to frequent isothermal snowpacks, enhancing the sensitivity of snow cover duration to climate warming . Also, at these sites, the hydrological response should be faster, since the isothermal snowpack cannot retain and refreeze liquid precipitation falling during the ROS as do snowpacks in colder and dryer conditions (Würzer and Jonas 2018). In addition, the milder the sites, the greater the possibility of rainfall at higher air temperatures, leading to higher melting rates during ROS (Corripio and López-Moreno 2017).
The processes governing the geographical differences in ROS frequency response to warming can explain spatial differences within the same mountain region, and also among the HRUs of the VBs used in this study. Thus, ROS events at high elevation or low insolation aspect HRUs occur onto deeper, colder and longer lasting snowpacks than other HRUs, and so show a small to negligible decrease in ROS frequency with increasing temperature; some basins experience an increase in ROS frequency with warming as spring snowfall shifts to rainfall. The basinscale ROS frequency response to increasing temperature may hide contrasting HRU-scale responses that counterbalance each other. This result is fully coherent with contrasting trends in ROS that have been found at the regional level. This study has shown that the amount of snowmelt during the same ROS but under increased temperature is enhanced in all VBs, especially in those where ROS frequency is expected to decline faster. This is explained by the increased heat content of the snowpack and enhanced sensible and latent heat fluxes, longwave irradiance and advection of energy by the rain under warmer climate . Such an increase in snowmelt rate is expected to often translate into increased peak streamflows during ROS (Musselman et al 2018). For the CRHM-modelled VB, the overall runoff volume produced during ROS is predicted to be unaffected by temperature increase, as a consequence of the general decrease in ROS. Thus, in a warmer climate, it is predicted that in most places there will be fewer ROS floods, but often there may be greater peak streamflows when they happen. This risk is greater for high elevations and cold regions where ROS may still be more frequent with warming, and especially in those areas where climate warming is associated with an increased intensity of precipitation (Prein et al 2016) that has not been considered in this study.
The results presented here are subject to the uncertainty associated with the forcing data (Weedon et al 2014), the way in which the climate is perturbed, the particular threshold value of rainfall used to define a ROS, and the specific configuration (elevation, size, soil characteristics, etc) and the defined land cover of the VBs. For instance, previous studies have already indicated that basin hypsometry and land cover strongly affect the occurrence and hydrological response of ROS (Wayand et al 2015). Downscaling available climate scenarios (e.g. CMIP6 ensembles) to mountain headwater basins introduces substantial uncertainty in temperature and precipitation fields, including orographic and convective effects, which are important for future ROS events. This study uses a delta method approach in order to isolate the role of temperature on the sensitivity of ROS events. However, future research should address changes in compound event characteristics related to snow-drought and heat-drought, but also the opposite: more frequent or more intense rainfall events under sustained warmer conditions or more extreme variable sequencing (AghaKouchak et al 2020). Nonetheless, this is a reliable approach to apply comprehensive snow hydrological models that include all the relevant physical processes involved in alpine snow regimes , and for isolating the role of temperature on changes in ROS for mountains with contrasting climate. These findings are useful for developing a broader picture of the impacts of climate change on hydrological hazards in different high mountain headwaters, and may serve as a baseline for more detailed studies, and with a deeper exploration of the physical processes that drive the ROS sensitivities, in areas where ROS has been identified as particularly sensitive to climate warming and where changes in ROS events may be critical for the management of natural disasters.

Conclusions
This research details and diagnoses the complex hydrological response of ROS events to climate warming. ROS frequency may increase or decrease under a warmer climate depending on the changing balance between more frequent rain events caused by warmer, more humid air, and the shortening snow cover duration with climate warming. A comparison of 40 VBs in mountain areas across the world revealed that 75% (30 of them) are predicted to decrease the frequency of ROS as the climate warms. In addition, the decrease in ROS is higher (lower) within each basin at lower (higher) elevations and for sunny (shadowed) slopes. This is because the snowfall ratio, the length of the snow duration and the average temperature at each site explain 66% of the global variance in the sensitivity of mountain ROS frequency to climate warming. These climatic characteristics also explain if the changes in the frequency of ROS are constant or variable along the +1 • C to +5 • C warming range tested. Generally, the changes will be more constant with temperature in those areas that currently exhibit longer snow cover durations. In all basins, snowmelt volumes during ROS increase with warming, which suggests that peak streamflows generated by ROS will increase as the climate warms, as a consequence of enhanced energy inputs and a decrease in the cold content of the snowpack. However, the general decrease in the frequency of ROS suggests that runoff volumes generated during ROS will not change substantially.

Data availability statement
Reanalysis data was created in the framework of the EU-funded WATCH project (www.eu-watch.org/ data_availability) and can be downloaded from https://rda.ucar.edu/datasets/ds314.2/#!description.
All data that support the findings of this study are included within the article (and any supplementary files).