Observed long-term greening of alpine vegetation—a case study in the French Alps

We combined imagery from multiple sources (MODIS, Landsat-5, 7, 8) with land cover data to test for long-term (1984–2015) greening or browning trends of vegetation in a temperate alpine area, the Ecrins National Park, in the context of recent climate change and domestic grazing practices. We showed that over half (56%) of the Ecrins National Park displayed significant increases in peak normalized difference vegetation index (NDVImax) over the last 16 years (2000–2015). Importantly, the highest proportional increases in NDVImax occurred in rocky habitats at high elevations (> 2500 m a.s.l.). While spatial agreement in the direction of change in NDVImax as detected by MODIS and Landsat was high (76% overlap), correlations between log-response ratio values were of moderate strength (approx. 0.3). In the context of above treeline habitats, we found that proportional increases in NDVImax were higher between 1984 and 2000 than between 2000 and 2015, suggesting a slowing of greening dynamics during the recent decade. The timing of accelerated greening prior to 2000 coincided with a pronounced increase in the amount of snow-free growing degree-days that occurred during the 1980s and 1990s. In the case of grasslands and low-shrub habitats, we did not find evidence for a negative effect of grazing on greening trends, possibly due to the low grazing intensity typically found in the study area. We propose that the emergence of a longer and warmer growing season enabled high-elevation plant communities to produce more biomass, and also allowed for plant colonization of habitats previously characterized by long-lasting snow cover. Increasing plant productivity in an alpine context has potential implications for biodiversity trajectories and for ecosystem services in mountain landscapes. The presented evidence for long-term greening trends in a representative region of the European Alps provides the basis for further research on mechanisms of greening in alpine landscapes.


Introduction
Changes in the productivity of vegetation cover observed by satellite can provide insight into groundlevel shifts in the structure and functioning of terrestrial ecosystems (Kerr andOstrovsky 2003, Pettorelli et al 2005). Potential drivers of spatial and temporal trends in vegetation greening (increased productivity) and browning (decreased productivity), often quantified using the normalized difference vegetation index (NDVI), include climate or land-use change (Jia et al 2003, Metternicht et al 2010, disturbance (Goetz et al 2005) or plant succession dynamics (Walker et al 1995). NDVI provides a remotelysensed proxy of standing biomass and photosynthetic activity (Myneni and Williams 1994), and when integrated over time can generate a number of metrics for the study of vegetation dynamics, including the duration of the growing season, peak greenness, and rates of green-up and senescence (Busetto et al 2010).
In temperate mountain ecosystems such as the European Alps, quantifying potential shifts in NDVI is of interest in order to monitor responses of alpine vegetation to global change and to complement observational field data. Furthermore, shifts in the productivity of mountain vegetation have strong implications for ecosystem services in alpine landscapes, including erosion control through root systems (Meusburger et al 2010) and pasture resource for livestock (Launchbaugh et al 1990).
Greening and browning dynamics in alpine study areas have received substantially less attention than in Arctic and sub-Arctic regions, where multiple field and remote sensing-based studies indicate widespread patterns of recent greening, mostly caused by shrub expansion into tundra (Tape et al 2006, Myers-Smith et al 2011. However, there are a number of reasons to expect temporal trends in NDVI in the context of the European Alps, including recent increases in air temperature above the global average increase (Beniston 2005, Durand et al 2009a, decreases in snow cover duration (Durand et al 2009b, Hantel andHirt-Wielke 2007), glacier retreat and ensuing plant colonization (Cannone et al 2008), tree line rise and forest expansion (Améztegui et al 2010, Carlson et al 2014, observed increases in species richness plant diversity observed on European summits (Pauli et al 2012, Kammer et al 2007 and shifts in alpine land use practices (Gartzia et al 2016).
One possible explanation for the lack of NDVI trend studies in alpine environments is that the 1 km spatial resolution of the widely used Advanced Very High Resolution Radiometer satellite (AVHRR) is poorly suited to quantify change in heterogeneous mountain landscapes (Fontana et al 2008). Currently, there is a need to exploit available time series of moderate to high-resolution satellite imagery in order to improve our understanding of vegetation dynamics in alpine landscapes. The Moderate Resolution Imaging Spectroradiometer on Terra/Acqua platforms (MODIS), which was launched in 1999 and provides 8 day composite NDVI products at 250 m resolution, currently allows for NDVI trend analysis over more than 15 years in temperate mountain study areas. A recent study utilizing MODIS NDVI between 2000 and 2012 demonstrated that inter-annual variability in the productivity of alpine grasslands in the French Alps is largely determined by snow free period, although temporal trends in peak NDVI were not investigated specifically (Choler 2015). Considering that MODIS spectral data have been used for vegetation trend analysis at the global scale (Zhao and Running 2010), and have been validated with respect to ground measurements and other sensors (Fisher and Mustard 2007), trends in NDVI detected with MODIS can be considered to provide reliable measures of change in vegetation cover. In an alpine context, however, the 250 m ground resolution limits the interpretation of NDVI trends with respect to spatially heterogeneous patterns of alpine vegetation (Virtanen and Ek 2014).
With the launch of the Landsat-8 satellite in 2013, the Landsat data program continues to provide the longest available record of earth surface observations at a relatively high ground resolution of 30 m (Roy et al 2014). The higher spatial resolution of Landsat imagery as compared to MODIS allows for the quantification of changes in NDVI closer to the spatial resolution of both (i) heterogeneity of alpine habitats and plant communities (Dedieu et al 2016) and (ii) drivers of vegetation change, for example domestic grazing or other land-use practices. Using Landsat-5 imagery, a study conducted in the Pyrenees between the late-1980s and the late 2000s demonstrated decreased productivity in grasslands characterized by high stocking rates, and increasing productivity in unmanaged alpine grasslands (Gartzia et al 2016). To complement existing research on NDVI dynamics in temperate mountain ecosystems, we believe that it is necessary to consider results from multiple satellite platforms to confirm greening and browning trends in alpine areas, and also to intersect satellite data with independent, ground-level field observations in order to attribute observed changes to specific habitats and vegetation types.
In this observational study, we investigated changes in peak NDVI in a 918 km 2 protected area of the southwestern French Alps. We selected the Ecrins National Park for our study as an emblematic territory representative of vegetation in the European Alps, and also in order to take advantage of spatially continuous field mapping of land cover and domestic grazing activity. To investigate the question of recent greening or browning of vegetation in the Ecrins National Park, we sought to (i) quantify trends in peak NDVI (NDVI max ) between 2000 and 2015; (ii) assess agreement between MODIS and Landsat data in terms of measured NDVI response; (iii) identify which land cover types showed the most change; (iv) quantify the temporal dynamics of shifts in NDVI max in the context of changes in snow cover duration and summer air temperature; and (v) to test for potential influences of domestic grazing pressure on observed NDVI max dynamics. We conclude our study by discussing potential local-scale processes contributing to satellite-detected shifts in alpine plant productivity as well as directions for further research.

Study area and vegetation mapping
The Ecrins National Park is a mountainous area including more than 150 summits over 3000 m a.s.l., with elevations ranging from 700 to 4102 m a.s.l.. Bedrock in the interior of the park consists of granite and gneiss on the highest summits, with zones of sedimentary rock frequently occurring in peripheral areas of the park. Between 1995 and 2000, an exhaustive and fine-scale mapping of habitats within the national park was conducted as part of the DEL-PHINE program at 100 m resolution (Godron and Salomez 1995). The resulting typology of land cover classes utilized in this study is listed in the legend of figure 1(a).  Choler (2015). Level 1T Landsat-5 (Thematic Mapper), Landsat-7 (ETM+, SLC-on) and Landsat-8 (OLI) images were obtained from the Landsat Earth Explorer data portal managed by the United States Geologic Survey (http://earthexplorer.usgs.gov). We downloaded available Landsat images covering the Ecrins National Park for the months of June, July and August for the years 1984, 1985 and 1987 (P1), 2000, 2001 and 2002 (P2) and 2013, 2014 and 2015 (P3; see supplementary figure S1, available at stacks.iop.org/ERL/12/114006/mmedia for the distribution of dates for utilized Landsat-5, 7 and 8 images). We calculated NDVI for each available Landsat date using the formula provided above. Supplementary appendix 1 provides further details on image processing aimed at harmonizing NDVI values between Landsat-5, 7 and 8 satellites.

Image processing and NDVI
To estimate peak photosynthetic activity within the Ecrins National Park, we extracted the highest summer value of NDVI (NDVI max ) for each grid cell and for each study year, as detected by both MODIS and Landsat. We utilized two approaches for estimating temporal change in NDVI max . First, we assessed trends in MODIS NDVI between 2000 and 2015 by applying Theil-Sen single median linear regression (Sen 1968) to each 250 m grid cell using year as the explanatory variable and NDVI max as the response variable. Second, we used the log-response ratio (Hedges et al 1999) to quantify proportional change in NDVI max (LRR NDVI max ) for two different time periods and two spatial resolutions: (i) between P2 (2000-2002) and P3 (2013-2015) at 250 m resolution, to compare change as detected by Landsat and MODIS and (ii) between P1 (1984, 1985, 1987 and P2, and between P2 and P3 at initial Landsat 30 m resolution. NDVI max for a given time period, for example P3, was calculated as follows: where NDVI max2013 is the map of peak NDVI for the year 2013, and so forth. We calculated the log-response ratio of NDVI max for each grid cell and for a given time period, for example between P2 and P3, as follows: where ln is the natural logarithm. To compare LRR NDVI max values between Landsat and MODIS platforms, we (i) mapped LRR NDVI max for both sensors and calculated overlap between positive and negative change in LRR NDVI max and (ii) generated random samples of grid cells stratified by land cover class, and estimated Pearson correlations between MODIS and Landsat LRR values.
We utilized the higher spatial resolution of Landsat to analyze change in NDVI with respect to elevation and land cover. To align layers, the land cover map was resampled to 30 m using the nearest neighbor algorithm, and a digital elevation model at 25 m resolution provided by the French National Institute of Geographic Information was resampled to the same 30 m grid using bilinear interpolation. For further details on NDVI max trend estimation, see supplementary appendix 2.

Climate and domestic grazing data
Annual values of growing degree-days and length of the snow free period between 1959 and 2015 were provided by the SAFRAN-SURFEX/Crocus-MEPRA model chain (S2M) developed by Météo-France for the French Alps (Durand et al 2009a, Durand et al 2009b, Vionnet et al 2012. This non-spatial model provides daily estimates of air temperature and snow cover depth for 300 m elevation bands for massifs within the French Alps. We estimated a combined metric of 'snow-free growing degree-days' for the Oisans massif by summing temperatures above 0 • C in the absence of snow cover. Values presented here represent averages between 1800 and 2500 m a.s.l.. Due to its physical basis and the nonspatial structure of the S2M model chain, we utilized model outputs to provide long-term regional climate trends and not to predict local thermal conditions. Domestic grazing intensity of non-forested habitats (grassland and low-shrub) was quantified based on grazing surveys available for the French Alps (1996-1997-2014SUACI, IRSTEA and CERPAM 2015) that best aligned with P1:P2 and P2:P3 covered by remote sensing data. We estimated summer grazing intensity for 240 pastoral units within the Ecrins National Park (on 15th July) as the coefficient of total livestock units per hectare (LSU ha −1 ). For each temporal period (P1:P2 and P2:P3), we fitted two linear models to test effects of grazing intensity on observed greening dynamics: (i) with initial NDVI max as the explanatory variable and LRR NDVI max as the response variable and (ii) with grazing pressure as the explanatory variable and residuals of the first model as the response variable. Further explanation of grazing data is provided in supplementary appendix 3.

Spatial trends in MODIS and Landsat NDVI max
The distribution of NDVI max linear trend slope values estimated using MODIS NDVI time series had a slight positive skew (mean = 0.0014 yr −1 , max. = 0.018 yr −1 , min. = −0.013 yr −1 ) and values were highly spatially structured ( figure 1(b)). Areas of concentrated negative slope values shown in dark brown in figure 1(b) could generally be attributed to local factors, i.e. (1) a forest fire occurring near the town of Argentière-la-Bessée in 2003; (2a-b) areas of agricultural activity; and (3) a forest that was destroyed by either avalanche or landslide activity in a steep valley, according to aerial photos. Areas of significant positive slope values (shown in green) were more widespread, and the highest slope values occurred in the interior region of the park at mid-to high elevations. South-facing slopes above east-west oriented valleys showed consistent greening trends (4), as well as grassland plateaus (5a-b; figure 1(b)). Overall, 56% of grid cells had a significant positive slope value above 0.001 yr −1 and 16% had a significant negative slope value below −0.001 yr −1 . The remaining 28% of grid cells had non-significant slope values between −0.001 and 0.001 yr −1 . The spatial distributions of log-response ratios of NDVI max (LRR NDVI max ) derived from Landsat and MODIS sensors showed strong agreement (figure 1(c) and (d)), and were spatially coherent with respect to the distribution of slope values shown in figure 1(b

Change in NDVI max in relation to land cover and elevation
Habitats with a mean NDVI max less than 0.75 (nonvegetated cliff and scree, vegetated cliff and scree, grassland and low-shrub) expressed a generally consistent pattern, with the majority of grid cells (60%-65%) exhibiting significant increases in MODIS NDVI max between 2000 and 2015 (figure 2). Above the 0.75 NDVI max threshold (dense scrub, mixed tree/shrub and forest), however, the relationship deteriorated and the number of grid cells showing no significant trend increased sharply (figure 2). It is well known that NDVI tends to saturate when measuring the photosynthetic activity of high leaf area canopies (Gitelson et al 1996) and for this reason, we were unable to quantify change in NDVI max in the case of tree-covered grid cells expressing high initial values of NDVI max . Accordingly, we excluded forest, dense scrub and mixed tree/shrub classes from the remainder of the analysis, and chose to focus on above treeline habitats.
Correlations between Landsat and MODIS LRR NDVI max at 250 m resolution were significant (P < 0.001), with correlation coefficients ranging from 0.18 (non-vegetated cliff and scree) to 0.36 (grasslands; figure 3). While sensors tended to agree on which grid cells showed increases or decreases in NDVI max between 2000 and 2015, the high scatter and fairly weak correlations demonstrated that the magnitude of the measured response varied strongly between MODIS and Landsat (figure 3).
LRR NDVI max for non-vegetated cliff and scree, vegetated cliff and scree, grasslands and low-shrub habitats exhibited a triangular relationship with respect to elevation, with a consistent peak occurring around 2500 m a.s.l. (figure 4). Grasslands showed a second peak in LRR NDVI max around 900 m a.s.l.. Midelevations (between 1000 and 2000 m a.s.l.) were characterized by a consistent majority greening trend across land cover classes. The highest values of LRR NDVI max , occurred in the context of non-vegetated cliff and scree, and were concentrated above 2500 m a.s.l. (figure 4).

Long-term greening dynamics relative to climate and domestic grazing
According to the long-term Landsat record, 67% of above treeline habitats showed consistent increases in LRR NDVI max between 1984 and 2015 (table 1). In comparison, less than two percent of 30 m grid cells in the study area were characterized by long-term browning between 1984 and 2015. Decreases in LRR NDVI max were observed in 23% of habitats from 2000 onward, pointing to the emergence of recent browning trends in certain contexts. While overall only 8% of above treeline habitats showed greening dynamics beginning in 2000, 16% of grid cells in non-vegetated cliff and scree habitats showed increases in LRR NDVI max initiated after 2000 (table 1).
Across land cover classes, greening trends in the Ecrins National Park were stronger between 1984 and 2002 (P1:P2) as compared to between 2000 and 2015 (P2:P3; figure 5). LRR NDVI max mode values for P1:P2 were centered at 0.15, whereas values for P2:P3 were centered at 0.05. Higher values of LRR NDVI max during P1:P2 coincided with a pronounced increase in the quantity of snow-free growing degree-days (GDD) Figure 3. Scatter plots of log-response ratio for NDVI max values calculated using MODIS (x-axis) and Landsat (y-axis), for the 2000 to 2015 period. Panels correspond to land cover classes. For each land cover class, a random of sample of 500 points was generated in order to minimize spatial auto-correlation. 'R' represents the Pearson's product moment correlation coefficient, and ' * * * ' indicates P values < 0.001. above 1800 m a.s.l. during the mid-1980s (figure 6). Figure S2 shows independent time series for the length of the snow free period and for July-August growing degree days, and while we did not disentangle their relative contributions, it appears that change in both snow cover duration and air temperature led to the increase in energy availability shown in figure 6. Greening trends of pastoral units were consistent with those assessed for grasslands and low-shrub habitats overall. Regardless of the period (P1:P2 or P2:P3),   a stronger greening response was observed for the least productive pastoral units (i.e. higher LRR NDVI max for lower initial NDVI max ), suggesting that the most productive pastures were already closer to their potential maximum. This negative linear relationship was significantly more pronounced for P1:P2 than for P2:P3 (P < 0.001, F 3479 = 250.3, R 2 = 0.6102), which is in agreement with the general greening trends observed  for the studied habitats (figure 7(a), table 2). We only found a significant relationship between the residuals of the precedent models (RES LRR) and grazing intensity for P1:P2 (P < 0.1, F 1,90 = 3.297, R 2 = 0.0246), suggesting a slight positive influence of grazing on NDVI max trends. The latter is also supported by the negative skew of RES LRR of non-grazed pastoral units observed also for P1:P2 ( figure 7(b), table 2).

Discussion
Our study combined imagery from multiple sources (MODIS,7,8) with independent field mapping to investigate changes in peak plant productivity in a French high mountain national park over the last 30 years. We found overall concordance between trends detected by Landsat and MODIS platforms (figures 1(c) and (d); figure 3), which supported full exploitation of the Landsat record at 30 m ground resolution. Increases in NDVI max over the last 16 years (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) were concentrated in above-tree line habitats at high elevation (∼2500 m a.s.l.), and the main observed greening dynamic in our study area consisted of increased productivity in grassland, low-shrub and rocky habitats. The long-term  Landsat record indicated ongoing greening trends in over half (67%) of above treeline habitats in the Ecrins National Park, as well as recent greening dynamics initiated since 2000 in previously non-vegetated cliff and scree habitats. We also found that increase in NDVI max was more pronounced between 1984 and 2002 than between 2000 and 2015, suggesting a slowing of greening dynamics during recent years. Finally, greening trends were largely unaffected by summer grazing intensity, although we did observe a slight decrease in greening response in non-grazed areas for the 1984 to 2002 period.

Potential mechanisms and causes of alpine vegetation greening
In the case of alpine grasslands and low-shrub habitats, we hypothesize that increases in NDVI max may have been the result of two possible mechanisms: (i) gradual densification and increases in the height and biomass of herbaceous plant cover (figure S3 (a)−(c)); and (ii) encroachment of low stature shrubs into grassland communities (Dullinger et al 2003). The first mechanism has been visually observed along permanent transects within the Ecrins National Park (Senn 2016), however repeat measures of standing biomass need to be carried out in order to confirm this hypothesis. The second mechanism has also been observed in the study area in the form of local increases in Vaccinium myrtillus (Senn 2011) and Juniperus communis subsp. nana within grassland communities. Expanding shrub cover has also been reported during recent decades in the nearby central Italian Alps at elevations up to 2500 m a.s.l. (Cannone et al 2007). Field observations carried out in the Ecrins National Park suggest that both increased herbaceous biomass and increasing shrub cover likely contributed to the observed NDVI max signal. We consider that decreases in snow cover duration and increases in air temperature since the 1980s (figure 6; Durand et al 2009a, Durand et al 2009b had stronger effects on long-term vegetation dynamics than domestic grazing practices within our study area. Indeed, stronger greening trends observed between 1984 and 2002 than between 2000 and 2015 coincided with a sharp increase in snow-free growing degree-days occurring during the 1980s (figures 5 and 6). Additionally, when grazing intensities remained stable across periods (1984-2002 and 2000-2015) rates of greening were negatively correlated with initial NDVI max , regardless of the magnitude of grazing intensity ( figure 7, table 2).
A lack of influence of grazing intensity on greening trends is contrary to observations in the Spanish Pyrenees by Gartzia et al (2016), where multi-temporal analysis of Landsat-derived vegetation indices showed that grassland biomass and greenness decreased in grazed areas. We explain this discrepancy by the fact that, unlike the Spanish Pyrenees, the southwestern French Alps are subject to relatively low grazing intensities and have undergone less dramatic transformations of land-use practices during recent decades. From a regional-scale standpoint, the pastoral units of the Ecrins National Park are actually among those characterized by the lowest grazing intensity values of all massifs of the French Alps (Dobremez and Borg 2015), despite the existence of local areas of intensive grazing. Pastoral activity has been increasing for several years in the Ecrins massif, however, pointing out a new issue for future vegetation dynamics. While there is no consensus concerning whether the effects of grazing on grassland productivity are positive or negative, we hypothesize that moderate grazing could limit the accumulation of senescent material and improve plant productivity through nutrient inputs, both of which could lead to increased greenness during the growing season. In our study area, a potential positive effect of moderate grazing on greening trends is supported by the observed decrease in greening trend when domestic grazing was totally absent during the 1984-2002 period (figure 7, table 2). Further comparative analysis, however, including a wider range of grazing intensities, is required in order to confirm this hypothesis.
Cliff and scree habitats above 2000 m a.s.l. exhibited strong greening trends during both the 1984-2002 and 2000-2015 time periods (figures 4-5). As illustrated by aerial photographs in figures S3 (d) and (f), steep screes with shrubs and small trees already present appear to have become more densely colonized since the early 1990s, leading to an increase in measured NDVI max signal. It is likely that plant colonization dynamics in steep terrain have contributed to this trend, as tree establishment contributes to slope stability (Meusburger et al 2010) and can ameliorate local environmental conditions (Batllori et al 2009), which in turn further facilitate tree colonization. Field observations suggest that pioneer trees Pinus cembra and Larix decidua, as well as Sorbus spp., Populus tremula, Acer spp. and Rosa spp. on southern slopes, have colonized steep rock faces, screes and grasslands in many areas of the park during recent years and could be contributing to the widespread increase in NDVI max observed in vegetated rocky habitats.
We attribute shifts in NDVI max in the context of low elevation grasslands (<1500 m a.s.l.; figure 4) largely to the effects of local land-use practices, including agriculture, grazing and forestry. Specifically, we consider the increase in NDVI max in montane fields and low elevation grasslands to be linked to forest expansion since 2000 following pasture abandonment. Observational studies elsewhere in the Alps (Carlson et al 2014) and Pyrenees (Améztegui et al 2010) point to the widespread nature of land-use related forest expansion throughout European mountain ranges. In our case, this interpretation is supported by a monitoring program lead by the scientific department of the Ecrins National Park (SOPHIE program;Merhan 2015). Relying on systematic field surveys, similar trends have been described, including a sharp increase in heath and forest cover up to 1998, followed by a slowing expansion of woody vegetation during the last 15 years. Other disturbances affecting the position of the forest-grassland ecotone include fire and logging, however we consider that these drivers are relevant only for small, isolated areas of the Ecrins National Park.
The highest proportional changes in NDVI max between 2000 and 2015 occurred in initially nonvegetated cliff and scree habitats (i.e. mapped in 2000 as scree, talus, ridges and cliffs with less than 10% vegetation cover; figure 4). We propose three external drivers that could be responsible for greening trends in this context. First, glacier retreat during recent decades has been especially pronounced in the southern French Alps (Gardent et al 2014), and as described elsewhere in the Alps (i.e. Raffl et al 2006), has been followed in Ecrins National Park by plant colonization of icefree moraines and exposed bedrock. Observed plant colonization in the study area following glacier retreat consists of forbs or shrubs and small trees, including Salix sp., Acer sp., Alnus alnobetula or Larix decidua. This process is illustrated in figure S3 (g)−(i), where vegetation colonized a recently de-glaciated zone and caused a substantial proportional increase in NDVI max . Second, in addition to vegetation dynamics in periglacial areas, we believe it is possible that the decrease in snow cover duration and increase in growing degree days since the 1980s may have led to the emergence of a brief growing season in high-elevation contexts previously characterized by perennial snowfields. This hypothesis could also help to explain the higher rates of greening observed between 1984 and 2002 as compared to more recent years, given the occurrence of a 'step effect' of sudden increase in air temperature and decrease in snow cover duration occurring in the 1980s followed by more gradual climate change during recent years (figure 6, Durand et al 2009a, Durand et al 2009b. Although we expect that colonization in these contexts likely consists of sparse scree vegetation, with time this process could lead to the emergence of novel snowbed communities. Finally, increased atmospheric deposition of nitrogen has been observed in other temperate alpine study systems (Seastedt et al 2004), and it is possible that enhanced nutrient availability has concurred with observed increases in high-elevation plant productivity. In the context of acidic grasslands, increased nitrogen inputs have been linked to increased biomass and conversion of species-rich communities into tall, species-poor swards (Bobbink et al 1998). Further observational studies are required in an alpine context to improve understanding of interactive effects of nutrient availability and climate parameters on the productivity of alpine plant communities.

Implications of long-term greening dynamics and main conclusions
Increased plant productivity in alpine environments has a number of implications for ecological functioning and land management in protected mountain areas such as the Ecrins National Park. For example, increased biomass in alpine grasslands and in rocky areas could affect the spatial distribution of herbivore habitat, both (i) from the perspective of shepherds concerned by the amount and quality of grazing resource; and (ii) for wild herbivores optimizing foraging strategy based on vegetation quality and availability (e.g. Grignolio et al 2003). However, increased biomass is not synonymous with increased fodder quality, and this difference may have substantial consequences for large herbivores (Corazza et al 2016). A study from the Austrian Alps measured significant herbivory of nival vegetation up to 3000 m a.s.l. (Diemer 1996), which suggests that shifts even in high-elevation plant composition may have implications for other trophic levels, including soil micro-organisms and insects. Also, higher plant cover in alpine and nival habitats may contribute to slope stability and soil water retention capacity, particularly in the context of unstable glacial moraines and steep scree slopes. Decreased erosion would affect local hydrology and runoff (Lehning et al 2006) as well as slope-related risk assessment in alpine watersheds (Descroix and Mathys 2003). Finally, the combination of ongoing glacier retreat and more prevalent plant cover at high elevation has the potential to affect landscape aesthetics and perception, particularly for tourists associating snow and rock with a pristine high-mountain environment and who travel to the mountain areas such as the Ecrins National Park with this landscape image in mind.

Conclusion
This work provides evidence for long-term greening of alpine vegetation in a protected area of the French Alps, especially in above-tree line ecosystems. We hypothesize that increasing air temperatures, decreased snow cover duration and to a lesser extent changes in land-use practices are the main drivers of observed greening dynamics. Although we restricted our study to the Ecrins National Park in order to take advantage of available field data, future studies should expand the spatial scale considered and test for greening and browning trends at the scale of the European Alps. In addition to quantifying extrinsic factors affecting changes in peak productivity, i.e. climate forcing or grazing pressure, further fieldwork is necessary to elucidate vegetation dynamics contributing to the measured greening signal, i.e. shrub encroachment, expansion of herbaceous plant cover, increases in plant height and biomass or colonization by pioneer species. Very promising is the opportunity to combine multi-temporal satellite imagery with repeat field surveys of plant diversity and multi-trophic community composition, in order to explore links between biodiversity and ecosystem functioning in mountain environments.