# Intermediate snowpack melt-out dates guarantee the highest seasonal grasslands greening in the Pyrenees

J. Revuelto, D. Gómez, E. Alonso-González, I. Vidaller, F. Rojas-Heredia, C. Deschamps-Berger, J. García-Jiménez, G. Rodríguez-López, J. Sobrino, R. Montorio, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

J. Revuelto, D. Gómez, E. Alonso-González, I. Vidaller, F. Rojas-Heredia, et al., Intermediate snow-pack melt-out dates guarantee the highest seasonal grasslands greening in the Pyrenees. Scientific Reports, 2022, 12 (1), pp.18328. 10.1038/s41598-022-22391-x. hal-04916110

HAL Id: hal-04916110

https://hal.science/hal-04916110v1

Submitted on 28 Jan 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /></figure>

OPEN

# Intermediate snowpack melt-out dates guarantee the highest seasonal grasslands greening in the Pyrenees

J. Revuelto<fnref n="1" />✉, D. Gómez<fnref n="1" />, E. Alonso-González², I. Vidaller<fnref n="1" />, F. Rojas-Heredia<fnref n="1" />, C. Deschamps-Berger<fnref n="1" />, J. García-Jiménez<fnref n="1" />, G. Rodríguez-López³, J. Sobrino⁴, R. Montorio⁵, F. Perez-Cabello⁵ & J. I. López-Moreno<fnref n="1" />

In mountain areas, the phenology and productivity of grassland are closely related to snow dynamics. However, the influence that snow melt timing has on grassland growing still needs further attention for a full understanding, particularly at high spatial resolution. Aiming to reduce this knowledge gap, this work exploits 1 m resolution snow depth and Normalized Difference Vegetation Index observations acquired with an Unmanned Aerial Vehicle at a sub-alpine site in the Pyrenees. During two snow seasons (2019–2020 and 2020–2021), 14 NDVI and 17 snow depth distributions were acquired over 48 ha. Despite the snow dynamics being different in the two seasons, the response of grasslands greening to snow melt-out exhibited a very similar pattern in both. The NDVI temporal evolution in areas with distinct melt-out dates reveals that sectors where the melt-out date occurs in late April or early May (optimum melt-out) reach the maximum vegetation productivity. Zones with an earlier or a later melt-out rarely reach peak NDVI values. The results obtained in this study area, suggest that knowledge about snow depth distribution is not needed to understand NDVI grassland dynamics. The analysis did not reveal a clear link between the spatial variability in snow duration and the diversity and richness of grassland communities within the study area.

Alpine grasslands have exceptional ecological value as they host a very high biodiversity and are the basis of livestock grazing during the summer months<fnref n="1" />⁻³. They occupy large areas in mountain regions⁴ at elevations where a significant snow cover is present during winter and spring⁵. Their growth is controlled by climate variability and snow cover dynamics⁶ as the plants must cope with a growth period limited by low temperatures and snow duration⁷,⁸. The high spatial variability characterizing mountain terrain causes strong variations in snow accumulation patterns over short distances⁹,¹⁰ which translate into snow depth distribution and its persistence¹¹⁻¹⁴. This spatial variability can strongly influence the distribution and phenology of plant communities¹⁵⁻¹⁷, as well as the occurrence of endemic or very rare species with special conservation interest and thus must be included along with other meaningful variables (elevation, soil thickness, slope)⁵,¹⁸ when analyzing communities distribution and its temporal evolution¹⁹,²⁰.

Important shifts in the magnitude and duration of snowpack are occurring due to ongoing climate change and are expected to continue in temperate mountains, with shorter snow seasons and a significant rise in snow line elevation²¹. These changes may have important implications for plant biomass and net above-ground production as a longer growth period may be expected²², but the significance of these changes is still under discussion. At the same time the distribution and abundance of plants may change as temperatures rise²³, and, additionally, precipitation changes are highly uncertain and uneven in mountain areas²⁴⁻²⁶. In recent decades, important changes in plant species distribution and productivity trends have been observed in mountain areas, which have

<footnote n="1">¹Instituto Pirenaico de Ecología, Consejo Superior de Investigaciones Científicas (IPE-CSIC), Zaragoza, Spain. ²Centre d'Etudes Spatiales de la Biosphère, CESBIO, Univ. Toulouse, CNES/CNRS/INRAE/IRD/UPS, Toulouse, France. ³Departamento de Análisis Económico, Universidad de Zaragoza, Zaragoza, Spain. ⁴Facultad de Biología, Universidad Autónoma de Madrid, Madrid, Spain. ⁵Departamento de Geografía y Ordenación del Territorio-Instituto Universitario en Ciencias Ambientales de Aragón (IUCA), Universidad de Zaragoza, Zaragoza, Spain. ☉email: jrevuelto@ipe.csic.es</footnote>
been mainly related to snow cover variations²⁷⁻³⁰, making the comprehension of snow-vegetation interactions highly relevant.

Since the space-borne sensor breakthrough in environmental science, several researchers have paid attention to the relationship between vegetation phenology and snow cover evolution²⁰,³¹⁻³³. These works have shed light on the dependence of alpine vegetation productivity on snow dynamics over large areas and how this vegetation evolves along different temporary windows⁸,¹⁸,³⁴,³⁵. For instance, it is known that variations in the snow-melt date control the greenness onset over extensive mountain regions²². However, a full comprehension of these interactions requires high and very high spatial resolution imagery³¹, for which on-site observations are needed²⁰. In this regard, the combination of on-site observations with space-borne and close-range remote sensing acquisitions has demonstrated good capability for detecting short-term changes with well-established metrics for monitoring plant communities structure evolution and snow dynamics⁴,³³.

In the context of remote sensing imagery, vegetation productivity increases are identified as “greening” and vegetation productivity decreases are referred to as “browning”³⁴,³⁶. In remote sensing studies over extended areas the Normalized Difference Vegetation Index (NDVI)³⁷,³⁸ is considered a good proxy for determining the onset of the growing season and the greenness peak, and many grassland evolution studies are based on the NDVI²⁸,³¹,³⁴,³⁹, obtaining the same conclusions when equivalent indexes are used⁸,⁴⁰. In this study “greening” is used to define the period with a NDVI increase. Over extended areas, the snow melt-out date is usually derived from the Normalized Difference Snow Index (NDSI)⁴¹ from satellite observations, setting a pixel as snow-free/snow-covered depending on a threshold value⁴². Other studies are based on direct snow depth and vegetation height acquisitions, but these typically require on-site observations with automatic weather stations (AWSs) measuring sensor-ground distance⁶. Several studies have investigated the relationship between the snow cover and alpine grassland productivity at high spatial resolution in numerous plant survey plots¹⁵,¹⁶. However, more detailed studies are needed to understand the role of micro-relief characteristics⁴³,⁴⁴. In this sense, Unmanned Aerial Vehicles (UAVs) enable the acquisition of detailed observations of ecological systems not achievable hitherto from satellite acquisitions, allowing the characterization of processes not monitored before⁴⁵.

This study was designed in light of the lack of high spatial resolution studies to understand the influence that snow dynamics exert on alpine grassland growth, while also giving insight into the distribution of plant communities. The main goal is to shed light on the relation between the seasonal evolution of grassland greening and both the snow melt-out date and the snow depth observed the previous winter. The main hypothesis is that the snow disappearance date determines the onset of the grassland greening period, with a primary impact on the maximum NDVI values reached afterwards. As previous studies have noted, the snowpack can exert an influence on plant communities distribution¹⁵,¹⁶. Potential differences in plant species observed over short distances will also be assessed in view of the snow melt timing.

The methodology later detailed relies on UAV observations of the snow depth distribution and the NDVI acquired at Izas Experimental Catchment⁴⁶, in the Central Pyrenees. This site is characterized by a long-lasting snowpack and a heterogeneous topography comprised between 2000 and 2300 m above sea level (masl). Izas Experimental Catchment (Fig. 1) is mostly covered by alpine grassland (about 85%)⁴⁷ dominated by Festuca eskia and Nardus stricta. The high grassland cover that characterizes this site is representative of subalpine areas in the Central Pyrenees⁴⁸⁻⁵⁰. This fact, together with the long persistence of snow above 2000 masl in the Pyrenees⁵¹,⁵², heightens the value of the conclusions obtained in this small study area, which might be extended to other mountain areas in temperate latitudes with similar characteristics.

In the 2019–2020 and 2020–2021 snow seasons, 10 and 7 high spatial resolution (1 m × 1 m grid cell size) snow depth distribution maps covering 48 ha were derived from RGB (red, green, blue) images acquired with a fixed-wing UAV. The spatial distribution of the melt-out date was obtained from these maps. Likewise, over the snow melting period, the NDVI distribution was obtained at the same spatial resolution with a multispectral camera. The NDVI was retrieved for 9 (in 2020) and 5 (in 2021) dates between April and August. With these UAV data (see Table S1 in Supplementary material), it was possible to determine the evolution of grassland greening and relate it to snowpack variables (snow depth accumulated in previous months and melt-out date) and the temporal evolution of the meteorological conditions observed at the AWS installed in the catchment⁴⁶. Figure 1 includes an overview of the spatial data available in the study area.

Moreover, by the end of 2021 melting season (July 2021), seven survey plots distributed across the study area were set up to study vegetation and specifically to analyze its structure (cover and floristic composition) and different diversity indices (Shannon's *H*, Pielou's *J*, Simpson's *D*, alpha, beta, soil coverage), informations which cannot be retrieved by satellite or UAV observations. To the best of our knowledge, this is the first study exploiting very high spatial resolution observations of both snow cover and grassland greening evolution based on UAV acquisitions, and also combining them with plant diversity indexes computed in several vegetation plots.

## Results

**Main meteorological conditions over the study period.** The two snow seasons analyzed had contrasting characteristics, with a short duration and shallow snowpack in 2020–2021 season and an average snow year in terms of duration and accumulation⁴⁶ in 2019–2020 (Fig. S1 in Supplementary material). The first snow season also had milder temperatures than the 2020–2021. Precipitation in spring 2019–2020 was more frequent and intense.

**NDVI response to snow melt-out.** The box plots in Fig. 2 show the temporal evolution of the UAV NDVI acquisitions categorized by snow melt-out date (NDVI observed in the snow-free). These box plots show that grassland greenness onset is closely related to snow disappearance date. In both snow seasons, the NDVI values start to rise in the later UAV snow observations. For the first seasonal NDVI acquisitions, some areas
<figure><img src="image_4.png" /><figcaption>**Figures 1.** Left side panels show the location of Izas Experimental Catchment, an RGB composition obtained for the UAV flight of the 5 August 2020 including the location of plant species survey plots, the validation area of the UAV transported multispectral camera, and the location of the Automatic Weather Station (AWS). In the bottom right side the mask for bare ground and a picture showing an overview of the catchment are included. Right side panels show the distribution of different variables analyzed in this work for 2019–2020 and 2020–2021 snow seasons. The upper right maps show the snow distribution around peak snow accumulation, middle right maps depict the melt-out dates and the bottom right maps show the NDVI distribution by the end of snow melting period. The bare ground mask and an overview of the catchment obtained from the spectrometer validation plots are shown at the bottom left side. All these maps where derived from the UAV acquisitions described in the manuscript. The maps layout, the legends, the scale bars etc., of Figure were generated by J. Revuelto in ArcMAP 10.8 (https://www.esri.com/en-us/arcgis/products/arcgis-desktop/resources). The picture of the Experimental Catchment the 2 was obtained by J. Revuelto.</figcaption></figure>

were snow covered and thus the NDVI was not computed inside. In the following UAV acquisitions the snow of these areas progressively melted-out and thus the NDVI gradually increased in the later observations. Similarly snow free areas in the first UAV acquisitions had overall low NDVI values (average NDVI < 0.3) increasing in the later acquisitions. The variability observed in the initial NDVI values is most likely related to the snow melt-out derived from the recursive snow depth observations. In this way, areas with the same melt-out date depict the snow extent melted between two consecutive UAV acquisitions (from 4 to 20 days among consecutive UAV flights). This hampers determination of the exact onset delay between the snow melt-out date and the greenness onset.

In the 2019–2020 snow season, the higher frequency of UAV acquisitions allows a more precise analysis of NDVI trends. The maximum NDVI values observed for the different melt-out date categories and the greenness onset after the snow disappearance show a rapid NDVI increase for areas with melt-out date around 3 May [day of year (DOY) 124] which reach peak NDVI (~0.7). For areas with later snow disappearance (10–21 June), pasture greening is also quite rapid, but the maximum values reached are smaller (~0.6). The last UAV observation of the NDVI in 2020 showed the early phase of grassland senescent, as NDVI decreased for the entire study area. Despite the lack of more frequent acquisitions and the lack of UAV observations during the 2021 summer, the second snow season shows similar behavior of the NDVI evolution (onset, maximum NDVI and tendency) to that observed in the previous season. In view of the contrasting snow depth evolution and the differences in temperature and precipitation observed in 2019–2020 and 2020–2021 (Fig. S1 in Supplementary material), it is remarkable that, in both seasons, the spatial maximum NDVIs are reached in areas with melt-out date occurring between late April and early May.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Box plots depicting the temporal evolution of the NDVI for the UAV acquisition dates categorized by the melt-out date respectively observed in 2019–2020 (a) and 2020–2021 (b). The boxes show first and third quantiles and the horizontal line inside each box the median. Whiskers of each plot include maximum and minimum values for each category and day.</figcaption></figure>

Overall, areas with the earlier melt-out date do not reach NDVI values as high as those observed in zones where snow disappearance occurs later. This is particularly noticeable for areas where the melt-out date occurs before April (DOY < 105), which have mean NDVI values below 0.65 at peak productivity (observed in July). Similarly, those zones with the latest snow melting do not overtake the maximum mean NDVI values (0.7), reached in areas with melt-out date around DOY 105–124 (15 April–4 May). Interestingly, this NDVI evolution is observed in both snow seasons, illustrating that there is an optimal melt-out date to reach maximum NDVI values. This outcome suggests that the study area can be classified into three main groups in view of the melt-out date: early snow melt-out (DOY < 105), mid snow melt-out (or optimum melt-out, 105 < DOY < 125), and late snow melt-out (DOY > 125). In 2020 the extent covered by these classes is: 30% for early snow melt-out, 22% for mid melt-out and 48% for late melt-out; whereas in 2021 these are respectively 40%, 37% and 23%.

**NDVI response to snow depth.** The box plots in Fig. 3 show the comparison between the maximum NDVI and the maximum snow depth observed in each snow season (seasonal maximum snow depth on each pixel), grouped in 0.25 m snow depth intervals and categorized in the three major melt-out date categories (early, mid and late melt-out DOY). A clear tendency of maximum NDVI values when compared to the maximum snow depths is not observed. The three melt-out date groups exhibit high dispersion for peak NDVI, suggesting that peak snow accumulation (observed between late March and early April) is not related to peak grassland productivity. Only for high snow accumulations (> 2 m) exist a slightly negative tendency between snow depth
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Box plots of maximum snow depth and maximum NDVI values observed in 2020 (a) and 2021 (b). Values have been categorized in the three main melt-out date groups (early melt, mid melt and late melt) and grouped by 0.25 m snow depth intervals. The boxes show first and third quantiles and the horizontal line inside each box the median. Whiskers of each plot include maximum and minimum values for each category.</figcaption></figure>

and NDVI, which is explained by a later melt out date in those locations, tending to prevent high NDVI values. Figure 3 also demonstrates that maximum snow depth cannot be used to derive peak NDVI, as both high and low peak snow depths can later have low NDVI values and thus. Similarly, the correlation between NDVI and snow depth (not shown) was only slightly negative (albeit statistically significant) between observations having an intervening temporal gap (the NDVI and snow depth observations) of two to three weeks. For late NDVI observations (July/August), correlations with previous snow depth observations were close to zero (not shown), showing that the maximum NDVI is not related to the snow depth distribution observed some months before.

**Plant diversity indices; NDVI and topographic characteristics of the plots.** A wide variety of grassland species was observed in the seven survey plots (Supplementary material Table S2); however, in all plots the dominant species is *Festuka eskia*, which shows high and very high plant coverage in most subplots (see plots 4 and 5 in Table S2). Two other species are also present in all plots (Table S2), *Pilosella officinarum* and *Nardus stricta*. The other species observed had smaller plant coverage, with a wide variety of plant composition between the plots. The frequency distribution derived from the plant cover of all species yields a wide variety of indexes describing plant communities. The average species richness (*alpha*) and the plot's species turnover (*beta*) denote that plot 5 had the widest variety of plant species (Fig. 4), while plot 6 had fewest. The other plots had similar
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Box plots showing the dispersion of the topographic variables (a) and grassland communities' indexes (b) observed at the survey plot locations. The boxes show first and third quantiles and the horizontal line inside each box the median). Whiskers include maximum and minimum values for each variable observed in the plots.</figcaption></figure>

species richness. Overall plant cover in the plots was high, and, except for plot 5, which had subplots showing coverage between 50 and 70%, the coverage in all of them ranged between 75 and 100%. *Shannon's H* depicts the marked variety of species observed in plot 5, whereas the species evenness (*Pielou's J*) in subplots 5 and 6 is higher than that observed in the other five plots. Plots 4 and 7 exhibit the highest *Simpson's D* average values, denoting the dominance of a few plant species, whereas in plots 5 and 6 a wider diversity is observed.

Despite the variety in plant species and in plant diversity indices, there was no evident relationship between plots' topographic characteristics and their snow melt-out dates. The topographic characteristics alone (Fig. 4) do not explain the dispersion in plant diversity indices. The slopes of the surveyed plots show that plant coverage is higher for flatter plot locations. Similarly, the higher solar radiation is, the lower *Shannon's H* is. This indicates that, in this study, areas with enhanced solar radiation have fewer species varieties and are dominated to a greater extent by *Festuka eskia*. The melt-out dates, which span about 45 days in 2020 and more than 80 in 2021, do not
show any notable relationship to plant diversity indexes. However, the temporal evolution of the NDVI in both years is closely related to the melt-out date (Fig. S2 in Supplementary Material), in agreement with the NDVI patterns observed for the entire catchment (Fig. 2).

Discussion

This work aimed to shed light on the understanding of snowpack influence on grassland greening in mountain areas over short distances. From previous works in snow-dominated areas<sup>8,22,34,35</sup> this hypothesis was established: the snow disappearance date determines the onset of the grassland greening period, with a primary impact on the maximum NDVI values reached later. The results shown here confirm this hypothesis, demonstrating, with very high spatial resolution observations, the major impact that snow melt-out date has on the temporal evolution of alpine grassland greening. The frequent snow-depth and NDVI observations at Izas Experimental Catchment have also made it possible to confront a question recently raised when investigating shrub growth in snow-dominated areas<sup>44</sup>: "(...) earlier snow melt-out on alpine shrub growth: the sooner the better?". The UAV observations presented here demonstrate that, for alpine grasslands, sooner is not always better<sup>53,54</sup>. Moreover, these acquisitions show that the optimal productivity, in terms of maximum NDVI, is found in areas where melt occurs at intermediate dates (late April–early May, termed mid melt-out or optimum melt-out). Nevertheless, there might be other processes or variables (soil depth and characteristics, terrain curvature, aspect or solar exposure...) contributing to reach peak NDVI values as the dispersion on the NDVI values observed in this study demonstrate. Zones where the melt-out date occurs before (early melt-out) or later than (late melt-out) the optimum melting period show lower maximum NDVI values and thus, it is expected, lesser vegetation productivity. This outcome is significant because many works have already highlighted the major role of snow melt dynamics in grassland greening dynamics<sup>8,22,34,35,55,56</sup> but fewer emphasize the decrease in grassland productivity if the snow disappearance occurs too soon<sup>57-59</sup> or too late. The high spatial resolution of the observations has proven that grassland greening evolution can be highly variable over (very) short distances and it allows a classification in terms of the snow disappearance date as: early snow melt-out, mid snow melt-out, and late snow melt-out.

The discretization of Izas Experimental Catchment on a wide variety of melt-out dates also reveals two greenness trends: one with a moderate NDVI temporal increase which does not reach maximum NDVI values and is observed for early snow-melt areas, and the other with a rapid NDVI increase which can reach maximum NDVI values (for mid melt-out) or not (for late melt-out) but both with a similar temporal tendency. Some works observed that an earlier snow melt implies more productive grassland<sup>55</sup>, but certainly there are other meteorological constraints affecting grassland dynamics. Thus, if snow disappears too soon, the plants are exposed to spring frost which does not enable a longer growing season<sup>57</sup> and thus restricts plant development<sup>60</sup>. Areas with the earliest snow melt match the regions with shallower snowpack<sup>10</sup> at times of peak snow accumulation. These are typically observed in mountain ridges<sup>61,62</sup>, and frequently have less developed and shallow soils<sup>63</sup>. In Izas Experimental Catchment, the shallowest snowpack accumulations at peak accumulation time<sup>14</sup> match areas with very limited soil presence<sup>64</sup> and were masked in the analysis (Fig. 1) because no plants are present here. In these areas, the wind erosion of the snowpack prevents the establishment of a long-lasting snowpack that might protect soil from extreme temperatures. In turn, this affects the soil nutrient cycle, and thus grassland growth and its productivity is also drastically reduced<sup>65</sup>.

The long period covered by the NDVI acquisitions obtained in 2020 allowed determination of the maximum greening, as the last UAV flight in August 2020 already showed a marked NDVI decrease throughout the study area. The impact of snow dynamics in the senescing part is known to be minor as other meteorological variables control this period<sup>8,19,39</sup>, so analogous behavior was expected in 2021 (not observed due to the UAV crashing in July 2021). The grassland greening evolution was observed in two snow seasons with marked snow and climatic differences. Thus, the similar results obtained in both study years allow the conclusions of this work to be extended to years with climatic dissimilarities. Nonetheless, to extract more robust conclusions and catch the high temporal variability that characterize the snowpack accumulation in the Pyrenees it is encouraged to maintain similar observations in this mountain range<sup>66,67</sup>. Similarly, the use of UAV to retrieve high spatial resolution observations has constrains as the limited spatial extent of the study area. In this regard, UAV studies in mountain areas hardly exceed 100 ha<sup>68-70</sup>. This work is based on observations obtained on a small study area where the findings are valid. These might be extended to areas of similar characteristics with equivalent plant coverage and diversity, but would require additional observations to validate them if extended to larger areas.

Another interesting outcome of this research is that the snow depth distribution observed at peak snow accumulation does not explain the peak NDVI distribution. Small correlations between the snow depth and the NDVI were only observed in acquisitions close in time (about 2–3 weeks difference between the snow depth and the later NDVI observation) likely due to the coming snow disappearance in areas with shallow snow depths. This result strongly suggests that NDVI dynamics can be studied in the light of the spatial distribution of the snow melt-out date. The latter variable can nowadays be retrieved for large areas from a wide variety of high-resolution optical satellite sensors (Sentinel, SPOT, MODIS and others). However, a UAV still provides higher spatial and also higher acquisition flexibility than satellite sensors. This is clearly beneficial in areas where clouds limit satellite observations, as a UAV may permit acquisition at the exact times when skies are clear.

Although soil characteristics, humidity content and other variables with notable implications for grassland communities distribution<sup>71,72</sup> have not been investigated here, the major role played by snow melt dynamics in grassland greening evolution has been shown. Other works have found remarkable differences in alpine grassland communities related to snow melt timing<sup>15,73</sup>, but within the seven survey plots established in the Izas Experimental Catchment, no clear relationship between plant communities and snow dynamics or topographic factors was found. The only notable outcome of plant community surveys is that plots with higher solar radiation have less species richness while flatter plots have higher plant cover, which might be connected with the thicker soils
usually found in these areas. This demonstrates that the snowpack is not the unique variable explaining plant communities: although it might have an important role, it must be considered together with other topographic and soil variables⁵,¹⁹. The noteworthy presence in all the plant communities sampled of *Festuca eskia*, *Nardus stricta* and, to a lesser extent, *Festuca nigrescens*, is decisive and relegates the contribution of the rest of the species to the relationship of vegetation greening with the NDVI, but allows a spatial comparison of NDVI evolution. It should be noted that in many plant communities the predominance of biomass is accounted for by a few species (especially of genus *Festuca*, *Carex*, *Trifolium*)⁵, which would allow simplifying the greenish to NDVI ratio estimation by limiting it, in a first approximation, to these species. Moreover, these species dominate plant communities in large alpine and sub-alpine areas of the Pyrenees⁴⁹, what might enable to apply similar analysis to other areas of this mountain range. The short growing season on the alpine floor, is a strong constraint conditioning plants to fast development (growth, flowering and fruiting) after snowmelt in order to complete as soon as possible their life cycle, regardless of their taxonomy, morphology, life form and reproductive biology. This fact suggests a similar temporal behavior towards the greenish of most of the flora.

In view of the expected changes in mountain areas in the coming decades, the snow cover is likely to decrease in both extent and duration²¹. These shifts in spatio-temporal evolution of the snowpack will impact grassland greening dynamics³⁰,³⁶ but will also, if maintained over time, influence plant communities distribution⁷⁴. In some locations, these changes might benefit species migrations from lower to higher elevations, but in other locations earlier snow melt-out in the future could impede the establishment of some species as these might not be protected from atmospheric extremes.

## Conclusions

The results of this work have shown the important role that the snow melt-out date plays in the NDVI evolution in alpine grassland. Peak NDVI values are reached in areas with a melt-out date around late April or early May. Conversely, those zones showing an earlier snow disappearance do not reach NDVI values as high as those obtained in areas of mid melt-out. Similarly, areas with a later snow disappearance do not achieve peak NDVI values. This behavior suggests that there is an optimum melt-out date to reach maximum NDVI values, and that grassland greening is penalized by snow melting that occurs too early or too late. In Izas Experimental Catchment no major differences in grassland communities were related to snow dynamics or NDVI evolution. Nonetheless, there is a predominant species in all plots (*Festuka eskia*) which may have hampered the observation of potential differences in the NDVI evolution of the survey plots. Results also revealed the negligible relationship between maximum NDVI and peak snow depth values, but also the minor contribution made by previous snow depth observations to explaining NDVI. Hence, grassland dynamics studies on mountain areas might only require snow melt-out information, which is easily retrieved from binary snow cover maps from spaceborne, airborne and on-site sensors. Further work for longer periods, over extended mountain areas, and covering a higher altitudinal gradient, is recommended to fully understand grassland dynamics in snow-dominated areas and validate the conclusions obtained here in a small study area.

## Methods

**UAV-based NDVI and snow depth observations.** The high spatial resolution database of both snow depth and NDVI distribution was retrieved in Izas Experimental Catchment by means of a fixed-wing UAV, a *SenseFly eBee-Plus*. Either a RGB camera (S.O.D.A.) or a multispectral camera (*Parrot Sequoia*) was mounted on the UAV and acquired images of the study area during the UAV flights. These images were then processed in *Pix4Dmapper* (version 4.4.12), which is a Structure from Motion (SfM) software that generates orthophotos and 3D point clouds representing the snow/ground surface. The NDVI maps were computed using the red and near-infrared bands⁷⁵,⁷⁶ of the multispectral images. The snow depth was computed by comparing the 3D point clouds with snow presence and a point cloud generated from images retrieved without snow. Following a detailed acquisition protocol and under good illumination conditions, the accuracy of the UAV snow depth observations in this study area is below 0.19 m, as previously demonstrated in Izas Experimental Catchment⁶⁸. By contrast, the UAV-based NDVI observations were not previously evaluated in this study area and, although a thorough methodology was followed, an error assessment for the product generated was needed⁷⁷.

The orthophotos derived from the *Parrot Sequoia* multispectral camera bands and the NDVI for the UAV flight of 5 August 2020 in different sectors of the study area (Fig. 1) were compared to observations retrieved with an Analytic Spectral Devices spectrometer (ASD FieldSpec4 SR), with a bare fore optic (25°). ASD reflectance was calibrated using a white Spectralon panel (Labsphere Inc., North Sutton, NH, USA). Several transects (5 on each validation plot) were established within the validation plots to retrieve the grassland radiometric signal (50 points on each transect). ASD reflectance values were convolved to the spectral response functions of the Sequoia camera for comparison. This comparison showed an RMSE in the NDVI of 0.047 and a linear adjustment between both observation devices (camera Sequoia and ASD spectrometer) with an R² = 0.84 (see Fig. S3 in Supplementary material). These values demonstrate the reliability of the NDVI acquisitions of the UAV exploited here. From the same UAV flight, the bare ground mask was derived based on a threshold value in the blue band of the RGB orthophoto.

The melt-out date computation within the catchment was obtained from the snow depth maps. Those sectors of the study area having snow presence in one UAV acquisition but not in the subsequent flights were set with a melt-out date equal to the last day with snow. This choice was made to ensure an accurate distribution of the last date with snow presence observed within the study area. Thus it would be more accurate to name this variable as "last snow occurrence date", but for simplicity we named it melt-out date.
**Plant communities indexes in survey plots.** The vegetation was sampled in late June 2021. Based on the digital terrain model, seven zones of about 0.1 ha were selected in the study area (Fig. 1), stratified according to their topographic characteristics (elevation, exposure, slope and terrain curvature), which, in turn, determines the snow lifespan and the growing season. Inside each plot, two 10 m perpendicular transects were marked. Along each of these lines ten squares of 0.5 × 0.5 m (subplots), separated by 1 m and placed alternately on both sides of the line, were delimited to record plant cover including that of each plant present within the plot. Inside each subplot all plant species were computed and their coverage determined based on the line-point intercept method<sup>78</sup>. The data obtained from these surveys served to determine species richness and plant coverage, later exploited to determine plant community cover and plant diversity indexes (*Shannon's H*, *Pielou's J*, *Simpson's D*, alpha, beta, soil coverage)<sup>79,80</sup>.

**Data availability**

The database of UAV acquisition obtained in Izas Experimental Catchment will be made freely available in Zenodo repository after publication with this https://doi.org/10.5281/zenodo.7038071.

Received: 6 July 2022; Accepted: 13 October 2022

Published online: 31 October 2022

**References**

1. Battaglini, L., Bovolenta, S., Gusmeroli, F., Salvador, S. & Sturaro, E. Environmental sustainability of alpine livestock farms. *Ital. J. Anim. Sci.* **13**, 3155 (2014).

2. Lavorel, S. *et al*. Historical trajectories in land use pattern and grassland ecosystem services in two European alpine landscapes. *Reg. Environ. Change* **17**, 2251–2264 (2017).

3. Pan, Y., Wu, J. & Xu, Z. Analysis of the tradeoffs between provisioning and regulating services from the perspective of varied share of net primary production in an alpine grassland ecosystem. *Ecol. Complex.* **17**, 79–86 (2014).

4. Rossi, M. *et al*. A comparison of the signal from diverse optical sensors for monitoring alpine grassland dynamics. *Remote Sens.* **11**, 296 (2019).

5. Körner, C. Plant ecology at high elevations. In *Alpine Plant Life: Functional Plant Ecology of High Mountain Ecosystems* (ed. Körner, C.) 1–7 (Springer, 2003). https://doi.org/10.1007/978-3-642-18970-8_1.

6. Jonas, T., Rixen, C., Sturm, M. & Stoeckli, V. How alpine plant growth is linked to snow cover and climate variability. *J. Geophys. Res. Biogeosci.* **113** (2008).

7. Körner, C. Impact of atmospheric changes on high mountain vegetation. In *Mountain Environments in Changing Climates* 155–166 (Routledge, 1994).

8. Choler, P. Growth response of temperate mountain grasslands to inter-annual variations in snow cover duration. *Biogeosciences* **12**, 3885–3897 (2015).

9. Schirmer, M., Wirz, V., Clifton, A. & Lehning, M. Persistence in intra-annual snow depth distribution: 1. Measurements and topographic control. *Water Resour. Res.* **47**, 09516 (2011).

10. Revuelto, J., Jonas, T. & López-Moreno, J.-I. Backward snow depth reconstruction at high spatial resolution based on time-lapse photography. *Hydrol. Process.* **30**, 2976–2990 (2016).

11. López-Moreno, J. I. *et al*. Small scale spatial variability of snow density and depth over complex alpine terrain: Implications for estimating snow water equivalent. *Adv. Water Resour.* **55**, 40–52 (2013).

12. Clark, M. P. *et al*. Representing spatial variability of snow water equivalent in hydrologic and land-surface models: A review. *Water Resour. Res.* **47**, (2011).

13. Wayand, N. E., Hamlet, A. F., Hughes, M., Feld, S. I. & Lundquist, J. D. Intercomparison of meteorological forcing data from empirical and mesoscale model sources in the north fork american river basin in northern sierra Nevada, California. *J. Hydrometeorol.* **14**, 677–699 (2013).

14. Revuelto, J., López-Moreno, J. I., Azorin-Molina, C. & Vicente-Serrano, S. M. Topographic control of snowpack distribution in a small catchment in the central Spanish Pyrenees: Intra- and inter-annual persistence. *Cryosphere* **8**, 1989–2006 (2014).

15. Winkler, D. E., Butz, R. J., Germino, M. J., Reinhardt, K. & Kueppers, L. M. Snowmelt timing regulates community composition, phenology, and physiological performance of alpine plants. *Front. Plant Sci.* (2018).

16. Scherrer, D. & Körner, C. Topographically controlled thermal-habitat differentiation buffers alpine plant diversity against climate warming. *J. Biogeogr.* **38**, 406–416 (2011).

17. Billings, W. D. Arctic and alpine vegetations: Similarities, differences, and susceptibility to disturbance. *Bioscience* **23**, 697–704 (1973).

18. Hua, X., Ohlemüller, R. & Sirguey, P. Differential effects of topography on the timing of the growing season in mountainous grassland ecosystems. *Environ. Adv.* **8**, 100234 (2022).

19. Xie, J. *et al*. Land surface phenology and greenness in Alpine grasslands driven by seasonal snow and meteorological factors. *Sci. Total Environ.* **725**, 138380 (2020).

20. Carlson, B. Z., Choler, P., Renaud, J., Dedieu, J.-P. & Thuiller, W. Modelling snow cover duration improves predictions of functional and taxonomic diversity for alpine plant communities. *Ann. Bot.* **116**, 1023–1034 (2015).

21. Beniston, M. *et al*. The European mountain cryosphere: A review of its current state, trends, and future challenges. *Cryosphere* **12**, 759–794 (2018).

22. Stöckli, R. & Vidale, P. L. European plant phenology and climate as seen in a 20-year AVHRR land-surface parameter dataset. *Int. J. Remote Sens.* **25**, 3303–3330 (2004).

23. Steinbauer, M. J. *et al*. Accelerated increase in plant species richness on mountain summits is linked to warming. *Nature* **556**, 231–234 (2018).

24. Fazeli Farsani, I., Farzaneh, M. R., Besalatpour, A. A., Salehi, M. H. & Faramarzi, M. Assessment of the impact of climate change on spatiotemporal variability of blue and green water resources under CMIP3 and CMIP5 models in a highly mountainous watershed. *Theor. Appl. Climatol.* **136**, 169–184 (2019).

25. Kharin, V. V., Zwiers, F. W., Zhang, X. & Wehner, M. Changes in temperature and precipitation extremes in the CMIP5 ensemble. *Clim. Change* **119**, 345–357 (2013).

26. Engler, R. *et al*. 21st century climate change threatens mountain flora unequally across Europe. *Glob. Change Biol.* **17**, 2330–2341 (2011).

27. Qiao, D. & Wang, N. Relationship between winter snow cover dynamics, climate and spring grassland vegetation phenology in Inner Mongolia, China. *ISPRS Int. J. Geo-Inf.* **8**, 42 (2019).

28. Zong, S. *et al*. Upward range shift of a dominant alpine shrub related to 50 years of snow cover change. *Remote Sens. Environ.* **268**, 112773 (2022).
29. Ernakovich, J. G. et al. Predicted responses of arctic and alpine ecosystems to altered seasonality under climate change. *Glob. Change Biol.* **20**, 3256–3269 (2014).

30. Zheng, J., Gia, J. & Xu, X. Earlier snowmelt predominates advanced spring vegetation greenup in Alaska. *Agric. For. Meteorol.* **315**, 108828 (2022).

31. Dedieu, J.-P. et al. On the importance of high-resolution time series of optical imagery for quantifying the effects of snow cover duration on alpine plant habitat. *Remote Sens.* **8**, 481 (2016).

32. Virtanen, T. & Ek, M. The fragmented nature of tundra landscape. *Int. J. Appl. Earth Obs. Geoinf.* **27**, 4–12 (2014).

33. Fontana, F., Rixen, C., Jonas, T., Aberegg, G. & Wunderle, S. Alpine grassland phenology as seen in AVHRR, VEGETATION, and MODIS NDVI time series—A comparison with in situ measurements. *Sensors* **8**, 2833–2853 (2008).

34. Carlson, B. Z. et al. Observed long-term greening of alpine vegetation—A case study in the French Alps. *Environ. Res. Lett.* **12**, 114006 (2017).

35. Tomaszewska, M. A., Nguyen, L. H. & Henebry, G. M. Land surface phenology in the highland pastures of montane Central Asia: Interactions with snow cover seasonality and terrain characteristics. *Remote Sens. Environ.* **240**, 111675 (2020).

36. Rumpf, S. B. et al. From white to green: Snow cover loss and increased vegetation productivity in the European Alps. *Science* **376**, 1119–1122 (2022).

37. Myneni, R. B. & Williams, D. L. On the relationship between FAPAR and NDVI. *Remote Sens. Environ.* **49**, 200–211 (1994).

38. Pettorelli, N. et al. Using the satellite-derived NDVI to assess ecological responses to environmental change. *Trends Ecol. Evol.* **20**, 503–510 (2005).

39. Asam, S. et al. Relationship between spatiotemporal variations of climate, snow cover and plant phenology over the alps—An earth observation-based analysis. *Remote Sens.* **10**, 1757 (2018).

40. Rossini, M. et al. Remote sensing-based estimation of gross primary production in a subalpine grassland. *Biogeosciences* **9**, 2565–2584 (2012).

41. Dozier, J. Spectral signature of alpine snow cover from the landsat thematic mapper. *Remote Sens. Environ.* **28**, 9–22 (1989).

42. Hall, D. K. & Riggs, G. A. Accuracy assessment of the MODIS snow products. *Hydrol. Process.* **21**, 1534–1547 (2007).

43. Julitta, T. et al. Using digital camera images to analyse snowmelt and phenology of a subalpine grassland. *Agric. For. Meteorol.* **198–199**, 116–125 (2014).

44. Francon, L. et al. Assessing the effects of earlier snow melt-out on alpine shrub growth: The sooner the better?. *Ecol. Ind.* **115**, 106455 (2020).

45. Assmann, J. J., Myers-Smith, I. H., Kerby, J. T., Cunliffe, A. M. & Daskalova, G. N. Drone data reveal heterogeneity in tundra greenness and phenology not captured by satellites. *Environ. Res. Lett.* **15**, 125002 (2020).

46. Revuelto, J. et al. Meteorological and snow distribution data in the Izas Experimental Catchment (Spanish Pyrenees) from 2011 to 2017. *Earth Syst. Sci. Data* **9**, 993–1005 (2017).

47. Nadal Romero, E. et al. Sediment balance in four small catechumen's with different land cover in the Central Pyrenees (Spain). (2009).

48. Gartzia, M., Alados, C. L. & Pérez-Cabello, F. Assessment of the effects of biophysical and anthropogenic factors on woody plant encroachment in dense and sparse mountain grasslands based on remote sensing data. *Progr. Phys. Geogr. Earth Environ.* **38**, 201–217 (2014).

49. Fillat, F., González, R. G., García, D. G., Gómez, D. & Reiné, R. *Pastos del Pirineo*. (Editorial CSIC-CSIC Press, 2008).

50. Gómez-García, D., Ferrández, J. V., Tejero, P. & Font, X. Spatial distribution and environmental analysis of the alpine flora in the Pyrenees. *Pirineos* **172**, e027–e027 (2017).

51. Gascoin, S. et al. A snow cover climatology for the Pyrenees from MODIS snow products. *Hydrol. Earth Syst. Sci.* **19**, 2337–2351 (2015).

52. López-Moreno, J. I. et al. Different sensitivities of snowpacks to warming in Mediterranean climate mountain areas. *Environ. Res. Lett.* **12**, 074006 (2017).

53. Cernusca, A. Standörtliche Variabilität in Mikroklima und Energiehaushalt Alpiner Zwergstrauchbestände. In *Verhandlungen der Gesellschaft für Ökologie Wien 1975: 5. Jahresversammlung vom 22. bis 24. September 1975 in Wien* (ed. Müller, P.) 9–21 (Springer Netherlands, 1976). https://doi.org/10.1007/978-94-015-7168-5-2

54. Cernusca, A. & Seeber, M. C. Canopy structure, microclimate and the energy budget in different alpine plant communities. In *Symposium—British Ecological Society* (1981).

55. Kudo, G., Nordenhäll, U. & Molau, U. Effects of snowmelt timing on leaf traits, leaf production, and shoot growth of alpine plants: Comparisons along a snowmelt gradient in northern Sweden. *Ecoscience* **6**, 439–450 (1999).

56. Baptist, F. & Choler, P. A simulation of the importance of length of growing season and canopy functional properties on the seasonal gross primary production of temperate alpine meadows. *Ann. Bot.* **101**, 549–559 (2008).

57. Baptist, F., Flahaut, C., Streb, P. & Choler, P. No increase in alpine snowbed productivity in response to experimental lengthening of the growing season. *Plant Biol.* **12**, 755–764 (2010).

58. Wipf, S., Rixen, C. & Mulder, C. P. H. Advanced snowmelt causes shift towards positive neighbour interactions in a subarctic tundra community. *Glob. Change Biol.* **12**, 1496–1506 (2006).

59. Sierra-Almeida, A. & Cavieres, L. A. Summer freezing resistance decreased in high-elevation plants exposed to experimental warming in the central Chilean Andes. *Oecologia* **163**, 267–276 (2010).

60. Camarero, J. J., Gutiérrez, E. & Fortin, M.-J. Spatial pattern of subalpine forest-alpine grassland ecotones in the Spanish Central Pyrenees. *For. Ecol. Manag.* **134**, 1–16 (2000).

61. Dadic, R., Mott, R., Lehning, M. & Burlando, P. Parameterization for wind-induced preferential deposition of snow. *Hydrol. Process.* **24**, 1994–2006 (2010).

62. Vionnet, V. et al. Simulation of wind-induced snow transport and sublimation in alpine terrain using a fully coupled snowpack/ atmosphere model. *Cryosphere* **8**, 395–415 (2014).

63. Burns, S. F., Tonkin, P. J. & Thorn, C. E. Soil-geomorphic models and the spatial distribution and development of alpine soils. In *Space and Time in Geomorphology: Binghamton Geomorphology Symposium*, vol. 12 (2020).

64. Lana-Renault, N. et al. Comparative analysis of the response of various land covers to an exceptional rainfall event in the central Spanish Pyrenees, October 2012. *Earth Surf. Proc. Land.* **39**, 581–592 (2014).

65. Freppaz, M., Williams, B. L., Edwards, A. C., Scalenghe, R. & Zanini, E. Simulating soil freeze/thaw cycles typical of winter alpine conditions: Implications for N and P availability. *Appl. Soil. Ecol.* **35**, 247–255 (2007).

66. López-Moreno, J. I. et al. Long-term trends (1958–2017) in snow cover duration and depth in the Pyrenees. *Int. J. Climatol.* **40**, 6122–6136 (2020).

67. López-Moreno, J. I., Vicente-Serrano, S. M. & Lanjeri, S. Mapping snowpack distribution over large areas using GIS and interpolation techniques. *Clim. Res.* **33**, 257–270 (2007).

68. Revuelto, J., López-Moreno, J. I. & Alonso-González, E. Light and shadow in mapping alpine snowpack with unmanned aerial vehicles in the absence of ground control points. *Water Resour. Res.* **57**, e2020WR028980 (2021).

69. Eberhard, L. A. et al. Intercomparison of photogrammetric platforms for spatially continuous snow depth mapping. *Cryosphere* **15**, 69–94 (2021).

70. Harder, P., Schirmer, M., Pomeroy, J. & Helgason, W. Accuracy of snow depth estimation in mountain and prairie environments by an unmanned aerial vehicle. *Cryosphere* **10**, 2559–2571 (2016).
71. Stanton, M. L., Rejmánek, M. & Galen, C. Changes in vegetation and soil fertility along a predictable snowmelt gradient in the mosquito range, Colorado, USA. *Arct. Alp. Res.* **26**, 364–374 (1994).

72. Winkler, D. E., Chapin, K. J. & Kueppers, L. M. Soil moisture mediates alpine life form and community productivity responses to warming. *Ecology* **97**, 1553–1563 (2016).

73. Litaor, M. I., Williams, M. & Seastedt, T. R. Topographic controls on snow distribution, soil moisture, and species diversity of herbaceous alpine vegetation, Niwot Ridge, Colorado. *J. Geophys. Res. Biogeosci.* **113**, (2008).

74. Keller, F., Kienast, F. & Beniston, M. Evidence of response of vegetation to environmental change on high-elevation sites in the Swiss Alps. *Reg. Environ. Change* **1**, 70–77 (2000).

75. Running, S. W. Estimating terrestrial primary productivity by combining remote sensing and ecosystem simulation. In *Remote Sensing of Biosphere Functioning* (eds. Hobbs, R. J. & Mooney, H. A.) 65–86 (Springer, 1990). https://doi.org/10.1007/978-1-4612-3302-2_4.

76. Myneni, R. B., Hall, F. G., Sellers, P. J. & Marshak, A. L. The interpretation of spectral vegetation indexes. *IEEE Trans. Geosci. Remote Sens.* **33**, 481–486 (1995).

77. Huang, S., Tang, L., Hupy, J. P., Wang, Y. & Shao, G. A commentary review on the use of normalized difference vegetation index (NDVI) in the era of popular remote sensing. *J. For. Res.* **32**, 1–6 (2021).

78. Floyd, D. A. & Anderson, J. E. A comparison of three methods for estimating plant cover. *J. Ecol.* **75**, 221–228 (1987).

79. Peet, R. K. The measurement of species diversity. *Annu. Rev. Ecol. Syst.* **5**, 285–307 (1974).

80. Mouillot, D. & Leprêtre, A. A comparison of species diversity estimators. *Res. Popul. Ecol.* **41**, 203–215 (1999).

## Acknowledgements

This work was supported by the research projects of the Spanish Ministry of Economy and Competitiveness projects "El papel de la nieve en la hidrología de la península ibérica y su respuesta a procesos de cambio global-CGL2017-82216-R". J. Revuelto is supported by a "Juan de la Cierva Incorporación" postdoctoral fellow of the Spanish Ministry of Science, Innovation and Universities (Grant IJC2018-036260-I). Authors recognize the important contribution of Antonio Campo in field campaigns for observing plant communities.

## Author contributions

Conceptualization, J.R., J.I.L.M., E.A.G. and D.G.; methodology, J.R., F.R.H.; software, J.R., I.V., F.R.-H., J.G.J., C.D.-B., G.R.L.; validation, J.R., R.M., F.P.C.; formal analysis, J.R., J.S., investigation, J.R., J.S., D.G. J.I.L.M.; resources, J.R. J.I.L.M.; data acquisition, all authors; writing—original draft preparation, J.R.; writing—review and editing, all authors; visualization, J.R.; funding acquisition, J.R., J.I.L.M. All authors have read and agreed to the published version of the manuscript.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary Information** The online version contains supplementary material available at https://doi.org/10.1038/s41598-022-22391-x.

**Correspondence** and requests for materials should be addressed to J.R.

**Reprints and permissions information** is available at www.nature.com/reprints.

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_8.png" /></figure>

**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/.

© The Author(s) 2022