Spring Temperature and Snow Cover Climatology Drive the Advanced Springtime Phenology (1991–2014) in the European Alps

Shifts in phenology are important traces of climate change affecting mountainous ecosystems. We present an analysis of changes in spring phenology using a suite of Earth observation based parameters, that is, start of season (SOS), snow cover extent and meteorological variables from 1991 up to 2012/2014 for the European Alps. Our results show that SOS tends to occur earlier throughout the Alps during this period and spring temperatures have increased in the Eastern Alps. Spring temperatures presented a predominant influence on SOS for both, grasslands and forests across elevations between 500 and 2,200 m asl, while this effect is particularly pronounced in the northeastern Alps. Snow cover duration and snow cover melting days showed secondary impact on SOS. Our research provides a comprehensive observation of spatiotemporal changes in alpine spring vegetation phenology and its driving factors. They improve our understanding of the sensitivity of the European Alps ecosystems to a changing climate.

The average global temperature has been reported rising in recent decades, particularly in mountainous regions . In the European Alps, the average temperature has increased by about 1.1°C over the past 100 years (Böhm et al., 2001). The European Alps were reported to be particularly sensitive to inter-annual variations in climatic drivers such as temperature, precipitation, and solar radiation (Beniston et al., 2003;Gobiet et al., 2014), as well as, depending on elevation, snow cover extent and its depth Xie et al., 2017). Recent climatic warming was reported to have led to considerable variance and shifts in plant phenology across ecosystems (Graven et al., 2013;Parmesan, 2006;Parmesan & Yohe, 2003) and thus may threaten species with synchronized life cycles (Dozier & Painter, 2004;Post et al., 2009). The firm link between phenological changes and seasonal temperature variability suggests that broad-scale and long-term phenology observations could serve as a proxy for detecting and tracing global climate change over space and time (Dannenberg et al., 2018;Myneni et al., 1997;White et al., 1997).
Advances of start of season (SOS) were reported to reflect biological responses to warmer spring temperatures (Delbart et al., 2008;Fu et al., 2014;Linkosalo et al., 2009;Menzel et al., 2006;Parmesan, 2006;Peñuelas & Filella, 2001;Schwartz et al., 2006). More detailed investigations are therefore needed to focus on the detection of winter and spring temperature-induced trends in SOS in mountainous areas (Fontana et al., 2008). In the European Alps, springtime phenology has been well investigated with practice and evidence over the past years (Asam et al., 2018;Cornelius et al., 2013;Fontana et al., 2008;Galvagno et al., 2013;Ide & Oguma, 2013;Oehri et al, 2017Oehri et al, , 2020Pellerin et al., 2012), however, reported findings are mostly focusing on smaller geographical scales or selected aspects in observation. A comprehensive observation of spatiotemporal changes in alpine springtime phenology and its relationship with snow extent and climatologies is still limited in this field. Specifically, the impact and strength of winter and spring temperatures on SOS are not yet well understood across elevations (Xie et al., 2017). Moreover, there are no significant inter-annual trends observed on SOS and its parallel driving seasonal snow and spring meteorological SOS is particularly sensitive to inter-annual variations of climatic factors (Bennie et al., 2018;Fu et al., 2015;Wipf & Rixen, 2010), especially in winter and spring (Euskirchen et al., 2007). Snow cover and snow accumulation are two of the essential environmental factors controlling high-elevation SOS (Badeck et al., 2004;Keller et al., 2005;Wipf & Rixen, 2010). Changes in snow cover have been reported to have reverse implications on SOS (Dorrepaal et al., 2004;Goulden et al., 1998;Jonas et al., 2008;Monson et al., 2006), while snow cover protects vegetation from dry and severe cold conditions and simultaneously supplies moisture (Buckeridge & Grogan, 2008;Desai et al., 2016;Ide & Oguma, 2013;Inouye, 2008;Schimel et al., 2004). Therefore SOS adapted to changes in the duration of snow cover can show distinct responses to snowmelt advancement or delay (Keller & Körner, 2003;Wipf & Rixen, 2010;Wipf et al., 2006). Thus, snow cover seasonality might be the most important control to SOS at high elevations. Changes in seasonal patterns of snow are related to plant photosynthesis and growth in the snow-free season (Galvagno et al., 2013;Rossini et al., 2012), thus influencing ecosystem functioning (Saccone et al., 2013).
Characteristics of snow cover and snow accumulation differ with the variation of geographic factors such as climatic conditions and topography (Dedieu et al., 2014;Gobiet et al., 2014;Xie et al., 2017). Similarly, vegetation distribution is subject to meteorological and geological conditions (Ide & Oguma, 2013). From previous studies, a consistent message has hence emerged: the influence of climatic controls on the SOS is different between climatic regions (Dye & Tucker, 2003;Peng et al., 2010), between vegetation cover types such as grasslands (Cornelius et al., 2013;Zeeman et al., 2017) and forests (Hu et al., 2010;Jönsson et al., 2010), as well as between elevational gradients (Martínez-Berdeja et al., 2019;Trujillo et al., 2012;Walker et al., 2014). Therefore, it is necessary to identify the most sensitive zones to inter-annual variations in climatic controls for the further understanding of management and adaptation strategies responding to climate change. Changes in seasonal snow cover may interact with air temperature and thus affect SOS (Yu et al., 2013). Warming temperatures are likely to cause reduced average snowfall and earlier snow cover melting days in spring (Barnett et al., 2005;Foppa & Seiz, 2012). Together with a reduced snow cover duration, increased air temperature may lead to intensified water stress and ultimately constraining vegetation growth (IPCC, 2014. Meanwhile, these reported factors might lead to an earlier SOS, consequently resulting in advancement and extension of the carbon uptake period (Desai et al., 2016). However, the relative effects of climate and climate-related factors on the variation of SOS at the ecosystem scale across elevational gradients in mountainous regions as the European Alps are still pending to be thoroughly addressed.
Understanding the elevational variation of the seasonal snow and climatic factors influence on SOS and predicting future trajectories of spring phenological shifts are important aspects in mountainous ecological studies. The European Alps were selected as study area because of their strong topographic variability (Auer et al., 2007;Scherrer et al., 2004) with a typical mountainous/alpine climate (Brunetti et al., 2009), as well as their climate warming over the past (Böhm et al., 2001) and expected changes over the next decades (Gobiet et al., 2014). We thus examined the effects of inter-annual variations of the preceding seasonal snow extent and climatic factors on SOS using remote sensing and gridded meteorological data. We focused on four sub-climate regions, grasslands and forested regions, as well as elevation zones. Our aims are to (i) investigate the temporal trends of seasonal snow extent and climatic controls, as well as SOS, over the period 1991-2014, and (ii) test inter-annual association of seasonal snow and climatic controls with SOS between specific subregions, vegetated areas, and elevation.
These subregions are north-west (NW; temperate westerly and Oceanic features), north-east (NE; temperate westerly and continental features), south-east (SE; Mediterranean subtropical and continental features), and south-west (SW; Mediterranean subtropical and Oceanic features). The distribution and variation of natural vegetation types in the European Alps depends on subregion and elevation (Xie et al., 2017). We stratified our analysis using a vegetation cover map at 1-km resolution, generated from imagery of the Advanced Very High Resolution Radiometer (AVHRR) satellites acquired between 1981 and 1994. Satellite data were obtained from the global land cover classification collection of the University of Maryland (Hansen et al., 1998(Hansen et al., , 2000. Spring phenology is different between grasslands and forests (Cornelius et al., 2013;Pellerin et al., 2012;Zeeman et al., 2017); meanwhile, the relationships of spring phenology with snow cover and snow accumulation between high-elevation vegetation types (above 1,800 m asl) are similar (Xie et al., 2017(Xie et al., , 2018. We grouped natural vegetation (NV) into two types, that is, grasslands (represented as GL), including wooded grassland, closed shrubland, open shrubland, and grassland, as well as forest lands (represented as FL), including evergreen needle leaf forest, deciduous broadleaf forest, mixed forest, and woodland (see Text S1 and Figure S1). Grasslands are mainly distributed at elevations above 2,000 m asl, whereas forests are the dominating vegetation type at elevations below 2,000 m asl ( Figures S2 and S3).

Start of Season (SOS)
In this study, we employed SOS metrics from the NASA Making Earth System Data Records for Use in Research Environments (MEaSUREs) Vegetation Index and Phenology (VIP) global data set created from the AVHRR data (1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999) and Moderate Resolution Imaging Spectroradiometer (MODIS) Terra MOD09 (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014). The VIP Vegetation Index (VI) products (i.e., modified Enhanced Vegetation Index [EVI2] and Normalized Difference Vegetation Index [NDVI]) were used to derive SOS (denoted as SOS EVI2 and SOS NDVI respectively) with adjustment of gap-fill missing data and VIP smoothing using time-related moving windows (https://lpdaac.usgs.gov/products/vipphen_ndviv004/, accessed March 2020). SOS EVI2 and SOS NDVI were estimated using a modified Half-Max approach (White et al., 2009), applying a 35% delay in place of the original 50% (White et al., 1997). The 35% delay was found to result in an earlier estimated start of springtime phenology but proved to be more accurate particularly in regions with a slow-developing spring season. Furthermore, we used SOS derived from the Global Inventory Modeling and Mapping Studies (GIMMS 3g ) NDVI (represented as SOS Midpoint ) with adaption of the Midpoint pixel method (Garonna et al., 2014(Garonna et al., , 2016. No smoothing was applied to SOS Midpoint . SOS EVI2 and SOS NDVI are available in gridded format with approximately 5.6-km spatial resolution covering a 24-year period spanning from January 1991 to December 2014. SOS Midpoint is available in a gridded format with approximately 8-km spatial resolution and covering a 22-year period spanning from January 1991 to December 2012. SOS Midpoint data are calculated using the Midpoint pixel approach (Garonna et White et al., 2009), which defines SOS as the day of year at which NDVI first reaches up to half of its inter-annual amplitude. This approach is used to analyze the approximately 8-km/bi-weekly NDVI time series (GIMMS 3g ; Pinzon & Tucker, 2014;Tucker et al., 2005). AVHRR GIMMS NDVI well represents elevation-dependent variation of greenness in mountain areas (Trujillo et al., 2012). The Midpoint pixel approach is well adapted to different vegetation types and various elevation belts in the European Alps (Xie et al., 2017(Xie et al., , 2018. SOS metrics (i.e., SOS EVI2 , SOS NDVI , SOS Midpoint ) are expressed as day of the (calendar) year (DOY).

Snow Cover Seasonality
AVHRR snow data products with approximately 1-km spatial resolution were used to derive snow cover seasonality (SCS) for 24 water years (WY, defined from 01 October to 30 September of the following calendar year) from 1990 to 2014 for the entire European Alps region (Hüsler et al, 2012Metsämäki et al., 2005;Musial et al., 2014). Snow cover melting day (SCMD) was defined as the last date in the WY when a pixel is still snow covered. Snow cover duration (SCD) is the total number of days in the WY when pixels are snow covered (Hüsler et al., 2012. The calibration, preprocessing, and validation of the SCS were described in Hüsler et al., (2014).

Climatologies
Gridded datasets with an approximately 8-km (5-min) spatial resolution of monthly air temperature (Chimani et al., 2013), precipitation and solid precipitation (Chimani et al., 2011) from 1991 to 2014 were collected from the Historical Instrumental Climatological Surface Time Series Of The Greater Alpine Region (HISTALP, http://www.zamg.ac.at/, accessed March 2020). Winter (December, January, and February) and spring (March, April, and May) climatologies were computed by calculating the seasonal average based on the monthly values within a season. They are winter temperature (WT), spring temperature (ST), winter precipitation (WP), and spring precipitation (SP). Snow accumulation (SA) was calculated as the sum of solid precipitation from December to May (i.e., from winter and spring).

Elevation Gradients
Topographic information at a 1 arc-sec scale (∼30 m) obtained from the European Environment Agency (EEA) (http://www.eea.europa.eu/, accessed March 2020) was used to generate a 1-km digital elevation model (DEM) using the Nearest Neighbor approach. The elevational analysis consisted of selecting 200m elevation zones from 100 m up to 3,500 m asl (corresponding to the elevation range of vegetated areas, see Figure S2) within which the mean of each metric was calculated. The zones are iteratively moved by 100-m steps, that is, overlapping the neighboring zones, and thus differentiating well the variations of elevation-dependent distribution of vegetation types (Xie et al., 2017), snow cover seasonality , springtime phenology (Pellerin et al., 2012;Vitasse et al., 2018), and meteorological factors (Xie et al., 2020) across the entire study region.

Statistical Analysis
All data products were resampled to a 1-km resolution using the Nearest Neighbor method to match the grid of the SCMD, SCD, and the vegetation type classification. The annual trend of each metric was calculated by employing linear regression analysis on a pixel-by-pixel basis over time. All fitted slopes were tested for significance using a significance level p < 0.05. Given the effects of other driving variables, partial Spearman's rank correlation (significance defined as p < 0.05) was employed to examine the individual associations between SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) and snow cover seasonality (SCD and SCMD) and climatic metrics (WT, ST, WP, SP, and SA) again based on a pixel-by-pixel basis over time.
The 1-km grids of the statistical results were intersected with the 1-km maps of sub-regions and vegetation types. Elevation gradients with (i) SCD less than 10 days or more than 360 days (i.e., permanent snow areas) (Xie et al., 2017), (ii) less than 1% proportion of each vegetation type, and (iii) a percentage of significant results lower than 5% were filtered out in the regional and elevational analysis. Image data processing was performed with ArcGIS (v10.4, ESRI, USA) and ENVI/IDL (v5.1, EXELIS Inc., McLean, VA, USA), and statistical analysis was performed using R (http://www.r-project.org/, version 3.4.2.).

Regional and Elevational Variation of the Inter-Annual Trends of SOS, SCS, and Climatic Factors Between 1991 and 2014
The average inter-annual trends in SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ), SCS (SCD and SCMD) and climatic factors (WT, ST, WP, SP, and SA) differ within and between subregions and over the investigated period ( Figure 2). All the average inter-annual trend (regression slope) values shown in Figure 2 are significant, and a high number of pixels [% of total] reflects that a significant regression was found on a high proportion of the NV areas.
SOS EVI2 , SOS NDVI , and SOS Midpoint show significant negative trends across the entire Alps. Between 1991-2014, SOS EVI2 shows a significant advance on average by 0.54 days per year (represented as d yr -1 ) over 69.3% of total NV areas (i.e., the corresponding amount of pixels with a significant regression [% of total]) of the Alps, and SOS NDVI presents significant negative trends with an average value of 0.63 d yr −1 over 73.0% of the total NV areas of the Alps. Notably, the negative trends of both SOS EVI2 and SOS NDVI are more pronounced in northern regions than in southern ones. Between 1991 and 2012, SOS Midpoint behaved consistently in XIE ET AL.   Within the measurement period between 1991 and 2014, the inter-annual trends of SCD and SCMD are more pronounced (i.e., with a higher average annual trend and a higher amount of pixels with a significant regression [% of total]) in the southern Alps (SW and SE). SCD presents significant reductions with an average of −1.46 d yr −1 over 10.4% of the total NV areas in the entire Alps but only shows a very few significantly increasing trend above 2,000 m asl in each subregion. SCMD experienced a significant negative trend with an average of 0.41 d yr −1 over 10.8% of the total NV areas across the entire Alps, although with a weak significant trend on only a few NV areas in the northern Alps (NW and NE).
Between 1991 and 2014, ST increased on average by 0.06 degrees Celsius per year (represented as °C yr −1 ) over 73.9% of the total NV areas in NE and, however, by 0.04°C yr −1 only over 8.8% of the total NV areas in SE. The analysis of WT, WP, SP, and SA showed no or only a few areas with a significant temporal trend over 1991-2014, except for SE with an increased trend of WP, SP, and SA. In particular, SA experienced an average increased trend of 3.28 mm per year (represented as mm yr −1 ) over 78.8% of the total NV areas in SE, and of 1.34 mm yr −1 over 18.0% of total NV areas in SW.
Besides, the average values in SOS, SCS, and climatic factors (WT, ST, WP, SP, and SA) are depending on elevation for the entire Alps and their average inter-annual trends are different between the four subregions in both GL and FL (Table S1). Figure 3 shows that the average inter-annual trends of SOS (SOS EVI2 , SOS NDVI , SOS Midpoint ) of both grasslands (GL) and forest lands (FL) are varying with elevation gradient, as well as with elevation-distributed vegetation types ( Figure S2), for the NV areas of the entire Alps. SOS metrics of GL illustrate a slightly higher extent in inter-annual trends and the corresponding amount of NV areas depending on elevation than in the case of FL over the European Alps.
Between 1991 and 2014, both SOS EVI2 and SOS NDVI present an elevation-dependent variation in significant inter-annual advance (more than 0.50 d yr −1 ) with a considerable amount of NV areas (more than 60% of total) and with similar changes in GL and FL at the mid elevations between 1,000 and 2,500 m asl. However, at high elevations above 2,700 m asl, both SOS EVI2 and SOS NDVI show positive trends in GL and FL, with SO-S EVI2 experiencing a pronounced delay in more than 60% of total GL and FL areas. At lower elevations below XIE ET AL.
10.1029/2020JG006150 6 of 15 1,000 m asl, the direction of inter-annual trends of both SOS EVI2 and SOS NDVI is varying in an indecisive way and with a lower amount of NV areas compared with GL and FL areas in the mid elevations.
Meanwhile, between 1991 and 2012, SOS Midpoint advanced with an average of 0.79 d yr −1 over 20.6% of total FL areas, but advanced with an average of only 0.26 d yr −1 over 19.4% of total GL areas. The variation of the average inter-annual trend of SOS is slightly more significant in FL than in GL areas with elevation. The annual trends of SCD, SCMD, and ST also illustrated elevation-dependent variations in both GL and FL areas (Table S1).

Regional and Elevational Variation of the Relationship of SOS with SCS, and Climatic Factors
The average Partial Spearman's correlations of SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) with SCS (SCD and SCMD) and climatic factors (WT, ST, WP, SP, and SA) differ within each subregion and between the four subregions over the investigated period ( Figure 4). All average correlation coefficients r shown in Figure 4 are significant. A higher number of pixels (% of total) suggests that a higher proportion of significant correlations could be found in the NV areas. A high average absolute coefficient r denotes a strong relationship; however, directions and values of the correlation ranges are quite different among metrics. Most important, ST presents the strongest and most widespread correlation with SOS EVI2 , SOS NDVI , and SOS Midpoint for the entire Alps and within each subregion over the study period, followed by SCS (SCD and SCMD) and the other climatic metrics (WT, WP, SP, and SA).
Specifically, the average correlation coefficient r between ST and SOS EVI2 is −0.16 over 18.2% of the total NV areas of the Alps and particularly pronounced in NE (with an average r = −0.23 over 29.1% of the total NV areas). The average r of ST with SOS NDVI is −0.23 over 20.6% of the total NV areas of the Alps and particularly pronounced in NE as well (with an average r = −0.45 over 34.4% of the total NV areas). SOS Midpoint shows a negative correlation with ST (with an average r = −0.57 over 24.5% of the total NV areas) across the Alps, and the correlation is more pronounced in NE (with the average r = −0.60 over 37.3% of the total NV areas), thus being consistent with SOS EVI2 and SOS NDVI .
Compared with ST, SCS (SCD and SCMD) shows the second strongest correlation with SOS EVI2 and SOS NDVI , and these correlations are especially notable in the southern Alps (SW and SE). The correlations of SOS EVI2 and SOS NDVI with other climatic metrics (i.e., WT, WP, SP, and SA) are indifferent and with less proportion of total NV areas (≤10%) within NW, NE and SW, in contrary to SE, where moderate correlations of SOS EVI2 and SOS NDVI with WT, WP, SP, and SA occur on 9.1%-17.5% of the total NV areas.
However, SCS (SCD and SCMD) and the other climatic metrics (i.e., WT, WP, SP, and SA) present less and limited associations with SOS Midpoint because of either weaker correlations (with an average absolute correlation coefficient |r|≤0.45) or fewer number of significant pixels (≤10%) in each subregion, except for NE (with an average r = −0.52 between WT and SOS Midpoint occurring on 12.5% of the total NV areas) and SE (with an average r = −0.51 between WP and SOS Midpoint occurring on 10.7% of the total NV areas).
For both GL and FL, the average correlations of SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) with SCS (SCD and SCMD) and climatic factors (WT, ST, WP, SP, and SA) also differ within each subregion and between the four subregions (Table S2). However, differences in average correlations of SOS metrics with SCS and climatic factors exist between GL and FL areas, but they are not significant.
The Partial Spearman's correlations of SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) with SCS (SCD and SCMD) vary with elevation and vegetation type over the investigated period 1991-2012/2014 ( Figure 5). In this period, the percentage of area with significant correlations of ST with SOS EVI2 and SOS NDVI is higher than for correlations of SCD and SCMD with SOS EVI2 and SOS NDVI (GL between 500 and 2,500 m asl In GL areas at elevations above 2,400 m asl, the positive correlations of SCMD with SOS EVI2 and SOS NDVI and their proportions of the total NV areas with significant correlation are higher than the correlations of ST with SOS EVI2 and SOS NDVI . In FL areas at elevations above 2,800 m asl, SCD and SCMD show more pronounced correlations with SOS EVI2 and SOS NDVI than ST. Above 2,400 m asl in GL and above 2,800 m asl in FL, SCS (SCD and SCMD) shows more pronounced correlations with SOS EVI2 and SOS NDVI than ST.
Notably, the correlations of SCMD with SOS EVI2 and SOS NDVI are positive above 1,100 m asl and negative below 1,100 m asl in both GL and FL areas, while these correlations direction change occurs around 3,000 m asl for FL areas. At some elevations above 3,000 and below 500 m asl, ST also show positive correlations with SOS EVI2 and SOS NDVI in both GL and FL areas. The elevation dependent correlations of SCS (SCD and SCMD) with SOS Midpoint show largely different behavior in comparison to the correlations of SCS with SO-S EVI2 and SOS NDVI .
During the period of 1991-2012, the correlations between ST and SOS Midpoint are stronger than the correlations between SCS (SCD and SCMD) and SOS Midpoint for both GL and FL areas across elevations for the entire study area (Figures 5c and 5f). Negative correlations of ST with SOS Midpoint occur more frequently at elevations above 1,500 m asl for GL and at elevations between 500 and 2,500 m asl for FL.   FL, the elevational variations in average correlations of SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) with SCD, SCMD, and ST also differ between the four subregions ( Figures S4-S7). However, these differences mainly exist between western subregions (NW and SW) and eastern subregions (NE and SE). In NE, compared with SCD and SCMD, ST shows the strongest correlations with SOS EVI2 and SOS NDVI above 1,000 m asl and with SOS Midpoint across all elevations for GL and FL areas.

Regional and Elevational Variation of the Annual Trend of SOS, SCS, and Climatic Factors Between 1991 and 2014
Shifts in the spring timing of plant phenology driven by changing climate have been widely investigated over the past decades (Parmesan & Yohe, 2003). The results from our study point out that start of season (SOS EVI2 , SOS NDVI , and SOS Midpoint ) experienced significant advances in most of the natural vegetation (NV) area, in parallel with increasing spring temperature and significant reduction of snow cover duration (SCD) and advancements of SCMD, in the European Alps over the period 1991-2014.
Specifically, we found that SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) overall advanced significantly during the study period by an average value of 0.54-0.65 d yr −1 over 20.2%-73.0% of the total NV areas of the Alps. The differences in the amount of pixels with a significant trend (% of total) between SOS Midpoint and VIP VI products derived SOS (SOS EVI2 and SOS NDVI ) might be due to the use of time-related moving windows on VIP smoothing and the modified Half-Max method, which strengthens the trend significance of SOS EVI2 and SOS NDVI . This found trend is higher than reported by Defila and Clot (2005) for spring phenological XIE ET AL.
10.1029/2020JG006150 9 of 15  Reduced SCD across the Alps was especially pronounced in the southern Alps (SW and SE) over the study period 1991-2014 ( Figure 2). The decreasing SCD trend observed in our study across the Alps is in agreement with the recent findings that SCD experienced a significant decrease between 2001 and 2014 (Chen et al., 2015) and between 1980 and 2000 (Laternser & Schneebeli, 2003). Meanwhile, we found that SCMD advanced across the entire Alps, particularly in the southern Alps (NW and SW), during the investigated period.
ST showed warming trends by 0.06°C yr −1 over 73.9% of NV areas in NE, however a weaker or non-significant trend was found in the other parts of the Alps over the investigated period. This finding is in accordance with a report stating that the average temperature has increased by about 1.1°C in the European Alps during the 20th century (Böhm et al., 2001), and seasonal warming in spring is expected to vary between +1.2°C until the middle and +2.7°C until the end of the 21st century (Gobiet et al., 2014). We have not found significant WT trends in the Alps over the study period, thus indicating that the European Alps experienced no significant shift in winter temperature over the study period. In contrast to this, winter temperature was reported to have increased by 0.14-0.4°C per decade for the period of 1931-2010 with no significant difference between trends at high and low/medium elevations in the Swiss Alps (Martin Beniston, 2012).
Snow accumulation experienced increasing trends in SE over the study period. However, we have not found any significant trends in the variation of winter precipitation (WP), spring precipitation (SP), and snow accumulation (SA) for the entire European Alps during the same period. These findings might be driven by the limited spatial scale (approximately 8-km resolution) of the selected data that cannot well resolve complex mountain terrain.

Regional and Elevational Variation of the Relationship of SOS with SCS and Climatic Controls
We observed changes of SOS (SOS EVI2 , SOS NDVI , and SOS Midpoint ) across GL and FL areas of the entire Alps between 1991 and 2012/2014 that are significantly correlated with changes in ST, and these effects are more pronounced in the northeastern European Alps. ST is the most significant factor, compared with snow and other climatic factors, to advance spring phenology over the study period 1991-2014. Our results are in agreement with previous reports stating that the spring temperature has a dominating influence on springtime phenology (Chuine et al., 2010;Cornelius et al., 2013;Fu et al., 2014;Martínez-Berdeja et al., 2019;Menzel et al., 2006;Myneni et al., 1997;Peñuelas & Filella, 2001;Piao et al., 2006;Vitasse et al., 2018).
Significantly, our results indicate that the influences of SCD and SCMD on SOS are less important than ST in both GL and FL areas across the European Alps. These findings may result from the fact that SCD and SCMD are themselves influenced by winter and spring temperatures. Chen et al. (2015) reported that dates of green-up were mainly affected by spring air temperatures at locations with less pre-seasonal snowfall and lower elevations, while they were significantly influenced by coupled temperature and precipitation at locations with more pre-seasonal snowfall and higher elevations in the Qinghai-Tibetan Plateau. However, in the European Alps, spring temperature plays a dominant role in the determination of SOS in both GL and FL areas across elevation ranges, rather than pre-seasonal snow cover.
Our study provides comprehensive evidence that the influence of spring temperature on spring phenology is varying with elevation ( Figure 5). The influence of spring temperature on spring phenology in forests (with average r of −0.18 to −0.56 in 17.9%-23.4% NV areas) is similar to that in grasslands (with average r of −0.16 to −0.57 in 18.3%-23.1% NV areas). Our findings are in agreement with previous reports stating that the influence of ST on SOS differs between topographies such as elevation and region (Alikadic et al., 2019;Asam et al., 2018;Pellerin et al., 2012;Vitasse et al., 2018;Xie et al., 2020). The influence of spring temperature on spring phenology is more pronounced in the NE subregion (Figure 4 and Figure S5) and at mid elevations between 1,000 and 2,500 m asl. This finding coincides with reported spring temperature influencing leaf-out dates of trees across elevations in the Swiss Alps (Vitasse et al., 2018).
Compared with ST, our results show that satellite-derived SCD and SCMD are second most correlated with SOS when linking them to spring phenological development across the European Alps over the study period 1991-2014. Nevertheless, our results remain in line with previous findings stating that snow cover seasonality has strong effects on spring phenology in alpine and Arctic regions (Cooper et al., 2011;Dorrepaal et al., 2004;Galvagno et al., 2013;Jonas et al., 2008;Jönsson et al., 2010;Julitta et al., 2014;O Leary et al., 2018;Paudel & Andersen, 2013;Wipf et al., 2009;Wipf & Rixen, 2010;Xie et al., 2017). Besides, we did not find significant correlation between solid precipitation accumulation (i.e., SA) and SOS in the Alps, while snow accumulation can well explain the annual variations of forest greening in the Sierra Nevada mountains over the 1982-2006 period as reported in Trujillo et al. (2012).
In the European Alps, our findings highlight a strong synchronization of the spatiotemporal variations in spring air temperature and the investigated spring land surface phenology. The spring growth of alpine vegetated areas is influenced by many factors; however, our results show that spring phenological events of the Alpine ecosystem are dominated by one factor, namely the spring temperature. Thus, changes in spring temperature maybe an important proxy helping to develop uncertainty measures for the prediction of elevational and temporal variations in the start of vegetation growth.

Limitations and Uncertainties
While the results presented in this study point out some relevant relationships between shifts in climate and vegetation phenology, some methodological limitations have to be taken into account. On the one hand, the results presented here are based on the selection of coarse satellite data and gridded climate metrics. Variations in species composition in each grid cell, mixed signals from multiple canopy types, and the selection of the methodological approach affect the estimation of SOS and snow metrics (Helman, 2018). For instance, the methodological differences in VIP smoothing and the modified Half-Max method applied in the estimation of SOS EVI2 and SOS NDVI , as well as the adapted Midpoint pixel approach in the modeling of SOS Midpoint further confirmed that the annual variations of phenology metrics depend on the selected methods (Garonna et al., 2014). While validation of SOS quality is beyond the scope of this study, VIP products have been reported to exhibit some accuracy discrepancies in SOS estimates between land cover classes and to provide less accurate performance when compared to in situ observations (Berman et al., 2020;Chen & Yang, 2020). Furthermore, the complex conditions of remotely sensed signals face specific challenges in mountainous environments such as cloud coverage (Xie et al., 2017), a strong effect of the AVHRR orbit drifts due to topography, and changing acquisition and illumination geometries partially shifting slopes into the shade in winter and spring. Another strong limitation is the detection of snow in forested areas, especially in spring. It is assumed that snow is much longer present beneath the forest canopy than it can be discovered in optical satellite imagery (Thirel et al., 2012). Snow retrieval accuracy from satellite data and their spatial resolution may also play a role here (Hüsler et al., 2012). Besides, the observed sudden variability of direction and strength in trends and correlations at the highest and the lowest elevations might be caused by the low number of vegetated pixels ( Figures S2 and S3) in the European Alps.
On the other hand, other climatic and environmental factors, such as sunshine duration (Asam et al., 2018;Keller & Körner, 2003), soil temperature (Schimel et al., 2004), and wind exposition still need to be included in the investigation of spring phenology in Alpine regions. For instance, in the Central Alps (2600-3200 m asl) in Austria, Keller and Körner (2003) found that almost half of the high-elevation species are sensitive to photoperiod, which needs to be further addressed in combination with seasonal snow and other climatic factors.
Furthermore, the responses of vegetation communities spring phenology to climate warming maybe complex and nonlinear (Flynn & Wolkovich, 2018). A warming of approximately 0.25°C per decade until 2050, and of roughly 0.36°C per decade until 2100 are expected (Gobiet et al., 2014). Finally, warming directly influences the biodiversity and coexistence patterns of natural systems, which themselves are strongly driven by spring phenology and species interactions (Rudolf, 2019).

Conclusions and Implications
Our study provides two main contributions corresponding with its two main conclusions. First, our results provide empirical evidence that the direction and magnitude of inter-annual trends in springtime land surface phenology differ between climatic subregions and elevations. An overall advance in SOS of 0.54-0.65 d yr −1 was significant for the entire Alps and especially pronounced in the northeastern regions over the years of 1991-2014. Second, our findings reveal a considerable negative relation between spring temperature and springtime land surface phenology on 18.2%-24.5% of the natural vegetated area at both high-elevation grasslands and low-elevation forests across the entire Alps, particularly in northeastern regions. In addition, inter-annual variation of satellite-derived snow cover duration and snow cover melting day play secondary roles in the determination of springtime land surface phenology. Moreover, our findings indicate that other climatic factors have marginal effect on springtime land surface phenology compared to the dominant role of spring temperature, SCD, and SCMD.
Regional and elevational analysis of the direction and magnitude of annual trends of springtime land surface phenology is critical to understand the impacts of climate change on Alpine vegetated ecosystems. While numerous reports on changes in springtime phenology and studies of its response to climatic drivers, seasonal snow, and various environmental factors exist, our understanding of the parallel and combined impacts of these factors on spring phenology of Alpine vegetated ecosystems depending on elevation is still scarce at the landscape scale. Even at a coarse spatial resolution of approximately 8-km scale, the annual trends of spring phenology across mountain landscapes are still strong and significant. Multiple-source data generalization does not weaken this effect, indicating that the European Alps are especially prone to climate changes. Our results suggest that a projected increasing spring temperature results in a continuing advance of springtime phenology in the European Alps, particularly, in the northeastern subregion and at elevations between 1,000 and 2,500 m asl. As alpine vegetated ecosystems may be continuously reshaped through changes in spring phenology, ecosystem studies and corresponding adaptation planning should particularly take into account future elevation-dependent changes in temperature.

Data Availability Statement
The data reported in this study will be available as follows: snow cover seasonality is available through Hüsler et al. (2012Hüsler et al. ( , 2014; spring phenology SOS EVI2 and SOS NDVI are available through MEaSURE (shttps:// lpdaac.usgs.gov/products/vipphen_ndviv004/), and SOS Midpoint is available through Garonna et al. (2016); Climatologies will be available for download from the HISTALP http://www.zamg.ac.at/.