Twentieth century temperature and snow cover changes in the French Alps

Changes in snow cover associated with the warming of the French Alps greatly influence social-ecological systems through their impact on water resources, mountain ecosystems, economic activities, and glacier mass balance. In this study, we investigated trends in snow cover and temperature over the twentieth century using climate model and reanalysis data. The evolution of temperature, precipitation and snow cover in the European Alps has been simulated with the Modèle Atmospherique Régional (MAR) applied with a 7-km horizontal resolution and driven by ERA-20C (1902-2010) and ERA5 (1981–2018) reanalyses data. Snow cover duration and snow water equivalent (SWE) simulated with MAR are compared to the SAFRAN - SURFEX-ISBA-Crocus - MEPRA meteorological and snow cover reanalysis (S2M) data across the French Alps (1958–2018) and in situ glacier mass balance measurements. MAR outputs provide a realistic distribution of SWE and snow cover duration as a function of elevation in the French Alps. Large disagreements are found between the datasets in terms of absolute warming trends over the second part of the twentieth century. MAR and S2M trends are in relatively good agreement for the decrease in snow cover duration, with higher decreases at low elevation (∼\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\sim $\end{document} 5–10%/decade). Consistent with other studies, the highest warming rates in MAR occur at low elevations (< 1000 m a.s.l) in winter, whereas they are found at high elevations (> 2000 m a.s.l) in summer. In spring, warming trends show a maximum at intermediate elevations (1500 to 1800 m). Our results suggest that higher warming at these elevations is mostly linked to the snow-albedo feedback in spring and summer caused by the disappearance of snow cover at higher elevation during these seasons. This work has evidenced that depending on the season and the period considered, enhanced warming at higher elevations may or may not be found. Additional analysis in a physically comprehensive way and more high-quality dataset, especially at high elevations, are still required to better constrain and quantify climate change impacts in the Alps and its relation to elevation.


Introduction
Since the late ninteenth century, the Greater Alpine Region (GAR,(4)(5)(6)(7)(8)(9)(10)(11)(12)(13)(14)(15)(16)(17)(18)(19) • N) as a whole has seen a temperature increase of 1.2 • C, which is twice the Northern Hemispheric average (Auer et al. 2007;Brunetti et al. 2009). It is often claimed that the rate of warming is amplified with elevation (elevation-dependent warming) and that most high-mountain regions of the world are experiencing faster changes than surrounding lowlands (Pepin et al. 2015). Pepin et al. (2015) identified physical mechanisms contributing to modified warming rates with elevation, and stated that all of them should lead to enhanced warming with height or at a critical elevation. The most commonly cited mechanism is the snow-albedo feedback (e.g., Pepin and Lundquist (2008) and Scherrer et al. (2012)). However, these results have been contradicted by several other studies (Auer et al. 2007;Isotta et al. 2019;Hock et al. 2019) indicating that enhanced warming may or may not be found at higher elevations depending on the period, the season or the region considered. In the Swiss Alps, winter (Begert and Frei 2018;Scherrer 2020) and yearly (Philipona 2013) temperature have been increasing more rapidly at lower elevations (< 1000 m a.s.l). In the Mont Blanc massif, reconstructed temperatures over the last century from Dôme du Goûter (4300 m a.s.l, Gilbert et al. (2012)) were found to increase at similar rates than nearby low elevation areas (Vincent et al. 2020) .
Over the last few decades, one of the most iconic and noticeable consequence of anthropogenic climate change at high elevations has been the decrease of snow cover and the mass loss of glaciers. Meanwhile, snow water equivalent (Marty et al. 2017) and snow cover duration were consistently found to decrease in each part of the Alps Valt and Cianfarra 2010), with sharper decreases at low elevations. This decrease is caused more by an earlier melt in spring than by a later autumn onset (Klein et al. 2016;Matiu et al. 2021). The retreat of glaciers has been reported in numerous studies (e.g., Beniston et al. (2018) and Hock et al. (2019)). During the twentieth century, the loss of glaciers ice volume was estimated to reach about 50%, with accelerated mass loss since the late 1980s (e.g., Huss (2012) and Bolibar et al. (2020)). The ongoing shrinking of glaciers is expected to continue in the coming decades with the disappearance of about 90% of alpine glaciers by 2100 (Huss 2012;Zekollari et al. 2019;Vincent et al. 2019). These decreases, associated with changes in temperature and precipitation, have been and will increasingly impact many components of the alpine social-ecological system: -Agriculture (Fuhrer et al. 2014), ecology and local water availability , plant phenology, and wildlife (Mignatti et al. 2012). -Water availability downstream and hydropower production (Hanzer et al. 2018). -Winter tourism and the ski industry (Spandre et al. 2019;Berard-Chenu et al. 2020). -Natural hazards such as rock falls and debris flows (Ravanel et al. 2013;Phillips et al. 2017), avalanches (Castebrunet et al. 2014;Eckert et al. 2013), and extreme precipitation and rain-on-snow events (Beniston and Stoffel 2016; Scherrer et al. 2016).
In this context, it is crucial to identify and quantify temperature and snow cover trends as a function of elevation and to understand associated feedbacks. More than elsewhere, assessing climatic trends using in situ measurements in high-mountain regions is hindered by the lack of high-quality long-term records at high elevations.
Moreover, many stations at intermediate to high elevations are located on the slopes of deep valleys where temperature trends might be largely influenced by local topographic effects such as temperature inversions (Beniston et al. 1997), while stations on isolated peaks that are more representative of free atmospheric conditions are generally underrepresented. Similar issues arise for reanalysis data: trends are strongly influenced by the quality, number and distribution of observational data assimilated. As an example, Scherrer (2020) showed that winter-time warming trends at high elevation in the Swiss Alps are regularly overestimated in commonly used global or regional reanalyses data. The resolution of global reanalysis are also too coarse to assess climate trends in high-mountain regions.
With their finer resolution, regional climate models (RCMs) allow a better representation of the interactions between the atmosphere and the surface, particularly in moutainous areas where the complex topography and the snow cover strongly affect local-scale climate. RCMs are therefore valuable tools to investigate climatic trends in mountainous areas even though simulated trends are strongly influenced by the large scale model used as boundary condition.
In this study, we take advantage of a relatively highresolution (7 km) experiment of the RCM "Modèle atmosphérique régional (MAR)" driven by different climate reanalysis datasets to investigate changes in air temperature and snow cover in the French Alps. The model set-up and the other datasets used in this study are described in "Data and method". In "Evaluation and comparison", we evaluate and compare MAR simulations, high-resolution snow and climate reanalyses data, in situ stake measurements for late-spring snow water equivalent and homogenized weather station data. Trends over the twentieth century, with emphasis on the 1959-2010 period common to most datasets are presented in "Climatic trends". The physical processes that likely affect the relation of the trends with elevation are discussed in "Discussion" along with a comparison with previous studies.

The Modèle Atmosphérique Régional (MAR) model and set-up
The Modèle Atmosphérique Régional (MAR, Gallée and Schayes (1994), Gallée (1995), and Gallée et al. (1996)) RCM is a hydrostatic, primitive equation, limited-area model with constant sigma coordinates on the vertical axis. MAR has been developped for polar regions (e.g., Gallée et al. (2013) and Fettweis et al. (2017)) and has also been used over mountainous areas: e.g., the Himalayas  and Patagonia (Collao 2018). One asset of using MAR in high mountain regions is its surface scheme SISVAT (Ridder and Schayes 1997;Gallée et al. 2001) which includes a multi-layer snow cover model (Brun et al. 1992;Gallée and Duynkerke 1997) with prognostic equations for snow density, temperature, water content, and albedo.
In this study, MAR has been laterally driven at the border of the integration domain with the ERA-20C reanalysis (Poli et al. 2016), this simulation is referred to as MAR-ERA-20C hereafter. ERA-20C is one of the first reanalysis to cover the whole twentieth century , and has been used to drive MAR over western Europe (Wyard et al. 2018). Although this reanalysis is less accurate than more recent ones as it only assimilates surface pressure and surface marine winds observations, this dataset is expected to be consistent trough time and the risk to introduce spurious trends is limited (Fettweis et al. 2017). Nevertheless, in order to asses the dependency of the trends simulated with MAR to the atmospheric reanalysis used as boundary conditions, we drove MAR with the more recent 35-km horizontal resolution ERA5 reanalysis (Hersbach et al. 2020) to produce a second experiment over the 1981-2018 period (MAR-ERA5 hereafter). A first MAR-ERA-20C simulation has been recently used to investigate precipitation trends over the Alps (Ménégoz et al. 2020). Here, we performed this experiment with identical horizontal (7 km) and vertical (24 levels) resolution using a more recent version of the model (MARv3.9.10, http://mar.cnrs.fr/index.php) while saving additional variables needed to study the snow cover such as snow water equivalent (SWE) and surface albedo. Moreover, the integration domain (1.5-18.5 • E to 41.5-49.5 • N) has been extended to the South in order to better simulate low-pressure systems arriving from the Mediterranean Sea, as well as to the East, to include the easternmost parts of the Alps (e.g., Julian Alps). The domain considered in this study is shown in Fig. 1.
In both model and reanalysis data, the snow cover duration is computed using a 20 kg m −2 (or mm w.e.) snow water equivalent minimum threshold for snow covered days for hydrologic years starting on the 1st August and ending on the 31st July. In the MAR experiments, the snow cover ranges from values close to zero in the plains to almost 365 days per year at highest elevation areas (Fig. 1).

The S2M snow and climate reanalysis
The S2M ) snow and climate reanalysis (SAFRAN -SURFEX/ISBA-Crocus -MEPRA) dataset is available on the period 1958-2018. The SAFRAN reanalysis ) assimilates in situ meteorological observations to downscale large-scale reanalysis or numerical weather prediction models first-guess. Meteorological variables (wind, temperature, precipitation, radiation, etc.) are then aggregated by 300-m elevation bands at the scale of mountain massifs and used to drive the land surface scheme SURFEX/ISBA-Crocus. This scheme simulates snow cover with a multi-layer snow cover model based on prognostic equations for the evolution of physical snow properties (Brun et al. 1992;Vionnet et al. 2012). Each sub-massif, which has an area of about 500-1000 km 2 , is considered to have horizontally homogeneous meteorological conditions as the emphasis is put on the relation of meteorological conditions and snow cover properties with elevation. Here, we used S2M data representative of horizontal surfaces.

The SPAtialisation en Zone de Montagne (SPAZM) climate reanalysis
The SPAZM (SPAtialisation en Zone de Montagne, Gottardi (2009) and Gottardi et al. (2012)) dataset is a 1-km resolution gridded product for daily minimum and maximum temperature and precipitation. It covers the French Alps over the period 1948-2012. SPAZM is built on the interpolation of temperature and rain gauge measurements from a dense network and combines observations from Météo-France and the French electricity company (EDF). The mean temperature from SPAZM is estimated as the mean between minimum (TN) and maximum (TX) temperature (TN+TX/2).

Glacier mass balance measurements
The French National Observation Service GLACIOCLIM (https://glacioclim.osug.fr/) has monitored the evolution of surface mass balance (SMB) of six glaciers located in the French Alps since 1993 (e.g., Six and Vincent (2014)). The observational dataset includes winter and summer mass balances. The winter-time snow accumulation is measured in spring (late April/early May) by drilling snow cores and measuring snow density along the glaciers, while summer mass balance is estimated in late summer (early September) by measuring the emergence of stakes set in the ice.
In this study, we compare the observed winter mass balance (in meter of water equivalent, m w.e.) with the accumulation of snow simulated with MAR experiment and S2M reanalysis data for three glaciers : Mer de Glace and Argentière located in the Mont Blanc massif, and Gebroulaz located in the Vanoise massif. For each glacier and each year, we select the MAR grid points and the S2M values from the corresponding Mont Blanc and Vanoise massifs for the days corresponding to the dates of the field campaign on each glacier (weighted mean as a function of the number of observation per day). In terms of elevation resolution, there are, as an example, 9 corresponding MAR grid points covering the Mont Blanc massif with surface elevation ranging from 1800 to 3000 m a.s.l.

Homogenized in situ temperature measurements
Minimum (TN) and maximum (TX) homogenized temperature series from Météo-France data cover the 1959-2012 period and are available at a monthly time step. The data was homogenized using the HOMER method (Mestre et al. 2013). In this study, we use a set of 19 stations located in the French Alps or within the nearby vicinity, with only one station above 1500 m a.s.l, for which both TN and TX are available. For these data, mean temperature is also approximated as the average between TN and TX. More information about these data and the complete list of stations used with their location and elevation are given in the supplementary materials.

Evaluation and comparison
The MAR-ERA-20C and MAR-ERA5 simulations have been evaluated and compared with S2M and SPAZM for winter and summer mean temperature and with S2M for snow cover duration over 1981-2010. Results and figures are presented in details in the supplementary material (Section B1).
Overall, a larger altitudinal temperature gradient is found in MAR experiments (6 to 7 • C/km) with respect to S2M and SPAZM (4 to 6 • C/km). In winter, MAR-ERA-20C is warmer (∼2.5 • C) at lowest elevations and colder at highest elevations (around 3000 m a.s.l.), with fairly good agreement at intermediate elevations. In summer, MAR-ERA-20C agrees better with the other datasets, especially S2M, with no major biases at low elevations. The MAR-ERA5 simulation is generally 0.5 to 0.8 • C warmer than MAR-ERA-20C. Overestimation of temperature and precipitation altitudinal gradient is a typical RCM bias (Kotlarski et al. 2015).
As a result of warmer winter temperature and a higher rain to total precipitation ratio (even at elevations without significant temperature bias), MAR simulations systematically underestimate (∼ 10 to 60%) snow cover duration (compared to S2M) at low elevation (< 1500 m a.s.l). MAR-ERA-20C and S2M are in fairly good agreement for snow cover duration at intermediate elevations (1800 to 2700 m a.s.l) while MAR-ERA-20C shows longer snow cover duration at its highest elevations. MAR-ERA-20C and S2M agree on total precipitation rates at low elevation, while the annual total precipitation increases with elevation is twice as fast for MAR-ERA-20C compared to S2M (not shown). Due to larger warm biases at low elevations, MAR-ERA5 has a reduced agreement with S2M snow cover duration.
Comparisons with late-spring mass balance (in m w.e.) from glaciers in situ measurements (GLACIOCLIM) are shown in Fig. 2. The method used to generate this comparison is described in "Glacier mass balance measurements". Stake estimation has measurement uncertainty of ∼0.2 m w.e. and show large spatial variability (typically ∼30%) due to local processes such as blowing snow and avalanches (Dadic et al. 2010). As a result, snow accumulation is generally higher on glacier's flat surfaces than on nearby slopes. Therefore, these comparisons should be considered with caution, since these processes are not explicitly accounted for in MAR and S2M.
In the Mont Blanc massif (Mer de Glace and Argentière glaciers), MAR-ERA-20C and S2M show a fairly good agreement with the observations in the ablation area (< 2700 m). Higher up, S2M seems to underestimate the mass balance, while MAR-ERA-20C slightly overestimates it, although the lack of MAR grid points around 3000 m a.s.l. and beyond hampers a proper evaluation. This pattern (underestimation in S2M and overestimation in MAR-ERA-20C) is more visible in the comparison at Gebroulaz glacier (Vanoise massif). Perpetual snow (365 days of snow cover) is almost reached at the highest MAR grid point at 3000 m a.s.l (see supplementary material), a result in accordance with the actual equilibrium line observed around 2900 m a.s.l for glaciers of the Mont Blanc massif (Fig. 2).

Climatic trends
Annual mean temperature trends over periods ranging from 20 years to the longest time series available are shown in Fig. 3a for all the French Alps. Trends are computed using Theil-Sen slopes and significance is assessed using Mann-Kendall test (95% confidence interval).

Warming trends
The warming trend found in MAR-ERA-20C simulation over the whole twentieth century (0.8 • C/century) is somewhat lower than for previous studies (e.g., Auer et al. (2007) and Vincent et al. (2020)). For more recent period , large discrepancies in terms of mean temperature trends are found between the datasets. Nevertheless, MAR-ERA-20C, SPAZM and S2M agree on accelerating warming trends during the second half of the twentieth century. MAR-ERA-20C trends are consistent with those from the SPAZM reanalysis, while MAR-ERA5 is in better agreement with S2M, except for low elevation areas (below 1000 m). S2M trends below 1000 m are spurious : many negative and non significant trends are  found for annual or winter (see Fig. 4a) mean temperature, even over long periods .
Over 1981-2010, trends simulated with MAR-ERA-20C and MAR-ERA5 are very different, suggesting that MAR trends are largely driven by the large-scale model used as boundary conditions. For shorter periods (20 years), climatic trends are affected by the large inter-annual and decadal variability of the Alpine climate (and European climate in general), which has experienced a warming hiatus during the early 2000s (Meehl et al. 2014). However, many of the trends found for the time series extending up to 2018 (MAR-ERA5 and S2M) are significant as a result of the return of fast warming rates occurring over the last 8 to 10 years. All datasets agree on faster summer warming compared to winter warming.

Relationship with elevation
Focusing on a potential dependency between warming and the elevation, MAR-ERA-20C shows larger warming at low elevations over the twentieth century , especially in winter (see supplementary material, Table S1). Considering the 1959-2010 period, the relationship between warming trends and the elevation as well as its seasonality can be seen in Fig. 4a and b. Larger winter warming is systematically found at low elevations (below 1000 m a.s.l) in MAR-ERA-20C while enhanced warming at high elevations is found in summer. S2M shows generally enhanced warming at intermediate elevations (1200 to 1800 m a.s.l.) in winter, a result consistent with . Winter temperature trends also show exceptionally large spatial variability in S2M. In summer, the elevation of the maximum warming is shifted upwards and is more pronounced.
SPAZM shows a similar and constant slight increase of warming with elevation independently from the period or season considered, which somehow questions the reliability of this product to assess the dependency between warming rates and the surface elevation.  Fig. 3b and c. Regarding snow cover duration trends, there is surprisingly less disagreement between the data sets. When long periods are considered, all datasets agree on stronger negative trends at lower elevations. Negative trends also intensify over the second half of the twentieth century. Regarding shorter and more recent periods, trends at intermediate elevations are affected by the large interannual and decadal variability of the snow cover in the European Alps, and accordingly, positive but not significant trends may occur. At high elevation, trends are slightly

Seasonal trends and relationship to elevation
In this section, we investigate more in detail warming trends relationship with the elevation for the four seasons. We focus on the results with MAR-ERA-20C as this simulation allows to assess the interactions between the atmosphere and the snow-pack over a long period  for which local meteorological observations are available.
Warming trends in the French Alps over 1959-2010 simulated with MAR-ERA-20C and observed in Météo-France station are shown in Fig. 5 for the four seasons. trends as function of the elevation for model grid points located in the Northern and Southern French Alps are shown in Fig. 6a. In winter (DJF), the warming trend is maximal below 1000 m a.s.l and decreases with the elevation higher up (large negative correlation). The agreement between MAR-ERA-20C and station data is rather poor in winter. However, even at similar elevation station data trends shows high spatial heterogeneity (Fig. 6a). In spring (MAM), warming is maximum at elevation close to 1500 m in the Northern and around 1800 m in the Southern Alps. In summer (JJA), warming trends show few variations for elevations up to 2000 m a.s.l, and then increase sharply with the elevation higher up (large positive correlation). In autumn, temperature trends show few variations with the elevation (slight decrease with height, no correlation found). There a is good agreement between MAR-ERA-20C and homogenized station trends for spring trough autumn. Unfortunately, the lack of meteorological stations above 1500 m a.s.l hampers from discussing and confirming the relationship between warming trend and the elevation found in MAR-ERA-20C.
The seasonal trends for net absorbed shortwave radiation (SWnet) as a function of the elevation are shown in Fig. 6b. SWnet changes as a function of the elevation show high similarity to those found for temperature in spring and summer, with a positive correlation (r ∼ 0.73 to 0.92). SWnet patterns during these two seasons are largely driven by decreases in surface albedo (not shown), themselves caused by an earlier onset of snow cover decrease at intermediate elevation in spring and high elevation in summer. Elevations below 2000 m a.s.l were already snowfree in summer at the beginning of the period. In winter, a maximum increase in SWnet occurs at intermediate elevation (1300 to 1500 m a.s.l) that does not exactly match with the maximum warming. No significant correlation is found between mean temperature and SWnet trends during this season (r from 0.24 to 0.54), as well as in autumn (r ∼ 0).

Comparison with previous studies
Warming rates simulated with MAR-ERA-20C over the twentieth century (0.6 to 0.8 • C/century) are lower than the warming observed at nearby Swiss stations, (around 1.0 to 1.3 • C/century, Begert and Frei (2018) and Isotta et al. (2019)), or estimated by Vincent et al. (2020). ERA-20C reanalysis is based on the assimilation of surface pressure and maritime surface winds. The number of assimilated observations increased significantly between 1900 and 1950, which impacts the reliability of the atmospheric trends in the North-Atlantic region during the first half of the twentieth century (Bloomfield et al. 2018;Wohland et al. 2019).
Conversely, simulated trends (∼ 0.35 to 0.45 • C/decade) in MAR-ERA-20C over the 1959-2010 period are in better agreement with studies carried out over similar periods (e.g., Ceppi et al. (2012) and Lejeune et al. (2019)). Trends found in MAR-ERA5 experiment over 1981-2010 are lower than the one simulated with MAR-ERA-20C. However, simulated winter trends are consistent the observations in the Swiss Alps described in Scherrer (2020) with positive (∼ 0.2 • C/decade) yet insignificant trends found over 1981-2017. Nevertheless, comparisons of trends simulated with MAR and trends found in other studies is complex due to the differences in periods or areas studied.
The large difference in temperature trends over the common period 1981-2010 found with MAR-ERA5 (0.14 • C/decade) and MAR-ERA-20C (0.41 • C / decade) is surprising. We investigated trends directly in the largescale forcing fields of ERA5 and ERA-20C over their common period (see supplementary materials). We found smaller warming trends near the surface in spring and autumn in ERA5 compared to ERA-20C, and much smaller trends in lower and mid-troposphere in ERA5 in all seasons but winter. Near-surfaces trends in MAR-ERA5 are substantially different than those in its largescale driving fields ERA5, reflecting more the reduced warming trends in the upper levels in ERA5, and are possibly influenced by the different resolution of the surface processes. Such differences in near-surface warming trends between regional simulations and large-scale driving fields are not found for MAR-ERA-20C. Finally, the common period  used to perform this comparison is short for computing climatic trends and this result should be considered with caution.
Large disagreements between the reanalysis datasets S2M and SPAZM are also worth questioning. Negative trends in winter (−0.40 • C/decade) and annual (−0.10 • C/decade) temperature found in S2M at low elevations (< 1000 m a.s.l), even when considering a long period  are spurious and unreported in other studies. This is likely due to the low number of observations at low elevation assimilated in S2M. Therefore, S2M trends at low elevations are not further discussed in the remainder of the paper and should be considered with the greatest caution in future studies. Furthermore, the unreliable trends below 1000 m a.s.l. in S2M hinders us from properly discussing the evolution of warming with elevation using this dataset. Higher up, the large spatial variability of S2M trends at similar elevations (see Fig. 4a and b) suggests that trends in this product are very sensitive to the number and quality of observations assimilated. Temporal trends in S2M can suffer from the fact that the observation network changed over time (addition or discontinuation of observation stations). The same issues could be present in the trends derived from SPAZM or ERA5 reanalysis as the assimilation of temperature observations for the reconstruction of long-term trends is complex, especially Blue triangle (orange diamonds) represent Northern (Southern) Alps grid points. SWnet is estimated using shortwave downward radiation (SWD) and surface albedo (SWnet = SWD*(1-albedo)). Significant trends in MAR-ERA-20C (Météo-France station) are shown with a black dot (black countour). Spearman correlation coefficient between trends and elevation in MAR-ERA-20C simulation are shown in the legend in mountainous terrains. The quantification of climatic trends relation to elevation is best performed using spatiotemporally consistent (gridded) products that assimilate homogenized series of temperature observations such as done in Isotta et al. (2019). In this study, comparison with trends in MAR-ERA-20C simulation with those from homogenized station data from Météo-France (Figs. 6a and 5) suggest that trends in MAR-ERA-20C are quite realistic (except for a reduced agreement in winter). The facts that trends in MAR-ERA-20C are also in fairly good agreement with trends over the last 50 years found in the nearby Swiss Alps (Isotta et al. 2019) gives confidence in the trends simulated for this period.
Trends in the French Alps over the twentieth century found in this study with MAR-ERA-20C do not suggest a warming rate twice as fast as the Northern Hemisphere average. Auer et al. (2007) insisted on the fact that it is the Greater Alpine Region as whole that has been warming twice as fast as the hemispheric average over the twentieth century, while many studies have suggested that within the GAR itself temperatures are not rising faster at high elevations compared to nearby lowlands (Brunetti et al. 2009;Auer et al. 2007;Vincent et al. 2020). Moreover, the GAR average warming rates is influenced by the warming of peripheral regions, such as the Mediterranean regions (e.g., the Po Plain), which are experiencing much faster warming and drying, especially in the summer as a result of changing lapse rates (Kröner et al. 2017).
The small differences in warming rates between low and high elevation is explained by contrasting seasonal patterns. In this regard, trends simulated with MAR over the second part of the twentieth century agree with other studies that also found higher warming at low elevations (< 1000 m a.s.l) in winter (Begert and Frei 2018;Scherrer 2020), while Begert and Frei (2018) also found larger warming at high elevation in summer at Swiss stations. The fact that such seasonal patterns are found for both MAR-ERA5 and MAR-ERA-20C (Fig. 6) for the second part of the century gives confidence in this model result. Unfortunately, the lack of homogenized station data above 1500 m a.s.l in the French Alps hampers from evaluating the dependency between warming trends and elevation found in MAR simulations. Moreover, the lack of MAR grid points above 3000 m a.s.l at the resolution used prevent us from drawing conclusions on warming trends at very high elevations.
Trends in snow cover duration estimated from MAR-ERA-20C or S2M (Fig. 3) agree with previous studies reporting higher decreases at low and intermediate elevations Valt and Cianfarra 2010). The fact that we find more negative trends for shorter periods centered around the mid-1980s, with slight (non significant) increases for periods afterwards is consistent with a steplike decrease in the late 1980s followed by less significant trends afterwards also reported in these studies.

Snow-albedo feedback and other processes
The analysis of MAR-ERA-20C outputs for temperature and net shortwave radiation trends and their relation to elevation suggest a large effect of the snow-albedo feedback in spring and summer, which is not found in winter and autumn. These results are consistent with Kotlarski et al. (2015) that found a strong impact of the snow-albedo feedback in spring and summer at similar elevations to those found in this study. However, a proper quantification of the snow-albedo feedback as well as other feedbacks (e.g., lapse-rate feedback, surface temperature feedback) requires more complex or run-time diagnostics (e.g., radiative kernels, Soden et al. (2008)) that are beyond the scope of this study. The capability of climate models to represent realistic snow-albedo feedback has been linked to their capacity to represent the albedo of snow and vegetation (Thackeray et al. 2018). In our MAR RCM simulations, there is likely room for improvement in the latter as no prognostic or seasonal evolution of the vegetation cover is accounted for.
The fact that MAR simulation matches observed higher warming at low elevations in winter is interesting. Philipona (2013) attributed this higher warming at low elevation in Switzerland to higher increase in incoming shortwave radiation as a result of the reduction of aerosol particle emissions in Europe since the 1980s (also known as solar brightening). MAR RCM does not have interactive and explicit representation of the evolution of the aerosol concentration (Wyard et al. 2018) and is therefore not able to represent changes in concentration in alpine valleys in winter time. This suggests that other physical processes are also possibly determinant in this low elevation winter warming and correctly reproduced by the MAR RCM.

Application and perspectives
In this section, we briefly discuss the relevance and associated perspectives of MAR RCM simulations (and S2M reanalysis) within the framework of interdisciplinary initiatives. Many impact assessment studies rely on high quality and/or high resolution climate data. Prior to the launch of MAR simulations, interviews have been carried out with researchers from largely different background in order to identify essentials climate variables that should be saved and climate indices that should be calculated. Thanks to the relatively high-spatial resolution of MAR simulations and to the emphasis on the relation of meteorological conditions with altitude, MAR and S2M are complementary tools for the study of climate-related processes occurring at much higher spatial resolution that the whole French Alps scale used in this study. An example of snow cover duration for the Mont Blanc, Haute-Maurienne and Oisans massifs is shown in the supplementary materials. These three mountain ranges respectively surround the Arve, Maurienne and Romanche valleys.
Furthermore, MAR data are also available over the whole European Alps, and could be used for studies at this scale once they are evaluated outside of the French Alps. Because it reproduces physically consistent interactions between the atmosphere and the snow-pack, MAR RCM is an adequate tool for downscaling climate projections produced with global climate models such as those participating to the ongoing CMIP6 exercise. Similarly, recent 20CRv3 reanalysis (1836( , Slivinski et al. (2019) could be used to drive MAR RCM and extend the reconstruction of the Alpine climate back to the 1850s.

Summary and conclusion
In this study, we evaluated the MAR RCM driven by two different sets of climate reanalysis data for the representation of temperature and snow cover in the French Alps with focus on the representation of the elevation gradient. MAR's slight overestimation of the temperature decreases with the elevation is typical of what is found for many RCMs in high-mountain regions. The model is nevertheless able to reasonably represent the increase in snow cover duration and snow water equivalent with altitude, even though an underestimation of snow cover duration at low elevation (< 1500 m a.s.l) is found.
The warming trends found over the French Alps during the second part of the century show large disagreement between the reanalyses data and model simulations analyzed in this study. This highlights (i) the large remaining uncertainties on climatic trends in high-mountain regions, even in a region relatively well covered by observations such as the European Alps, and (ii) the need for better integration and coordination between modelling, reanalyses and observational studies.
Concerning the dependency of warming trends to elevation, MAR-ERA-20C and MAR-ERA5 agree with previous findings in the European Alps. When annual mean temperature is considered, little variation of warming rates is generally found. This finding hides large differences in seasonal patterns, with larger warming at low elevations (< 1000 m a.s.l) in winter and increasingly larger warming rates at high elevations (> 2000 m a.s.l) in summer. In spring, a maximum in warming is found at intermediate elevations (1500-1800 m a.s.l), where no relation with altitude is found for autumn warming. An investigation of MAR model outputs suggests that the snow-albedo feedback is likely to explain the relation between warming trends and altitude in spring and summer, while this is not the case in autumn and winter. However, more investigations are needed to better quantify the snow-albedo feedback in MAR model as well as other feedback processes. In the light of these results, we propose that elevation-dependent warming should not be considered as increasing warming rates with elevation but rather as differences in warming rates depending on the season, the elevation and/or the period considered.
Because the MAR RCM was shown to reproduce a realistic evolution of the snow cover with elevation and provides physically consistent atmosphere and snow cover interactions, the model, at the resolution used in this study, is a relevant tool for downscaling climate projections over the European Alps. Assessing the evolution of the relation between warming trends and elevation in the context of decaying snow cover in a warmer climate will be crucial for impact assessment studies as well as adaptation and mitigation strategies.

Acknowledgements
The authors thanks the providers of the observational datasets : Frédéric Gottardi and the French Electric Compagny (EDF) for providing SPAZM reanalyses data, Matthieu Vernay (CNRM, CEN) and Béatrice Dubuisson (Météo-France, DCSC) for Météo-France homogenized station data. We also thanks Matthieu Lafaysse and Marie Dumont for helpful discussions. The figures have been produced with the python package basemap (https://matplotlib. org/basemap/, last access: 20 December 2020). The authors sincerely thanks Jai Chowdhry Beeman and Nathan Maier for editing the English of this manuscript. Finally, we thank the two referees for their useful and constructive comments and the editor for the efficient management of the review process.
Author contribution All authors contributed to the design to the study. XF provided ERA5 and ERA-20C files to drive MAR RCM. DS and CV provided the data and helpful recommendations for comparisons with GLACIOCLIM in situ data. JB ran the MAR experiments, produced the figures and wrote the article and other authors contributed with suggested changes and comment.
Funding As part of the project CDP TRAJECTORIES, this work is funded by the French National Research Agency in the framework of the "Investissements d'avenir" programme (ANR-15-IDEX-02). We sincerely thank CDP TRAJECTORIES for its support. The authors also thank the GRICAD project (https://gricad.univ-grenoble-alpes.fr/, last access: 15 December 2020) for providing computer time for the simulations presented in this paper. The authors also thank the "Institut du Développement et des ressources en Informatique scientifique" (IDRIS, CNRS, project no. A0080101523).

Conflict of interest The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons. org/licenses/by/4.0/.