Change in observed long-term greening across Switzerland – evidence from a three decades NDVI time-series and its relationship with climate and land cover factors

ABSTRACT

Environmental changes are significantly modifying terrestrial vegetation dynamics, with serious consequences on Earth system functioning and provision of ecosystem services. Land conditions are an essential element underpinning global sustainability frameworks, such as the Sustainable Development Goals (SDGs), requiring effective solutions to assess the impacts of changing land conditions induced by various driving forces. At the global scale, long-term increase of vegetation greening has been widely reported notably in seasonally snow-covered ecosystems as a response to warming climate. However, greening trends at the national scale have received less attention, although countries like Switzerland are prone to important changing climate conditions. Hereby, we used a 35-year satellite-derived annual and seasonal time-series of Normalized Difference Vegetation Index (NDVI) to assess vegetation greenness evolution at different spatial and temporal scales across Switzerland and related them to temperature, precipitation, and land cover to investigate possible responses of changing climatic conditions. Results indicate that there is a statistically significant greening trend at the national scale with an NDVI mean increasing slope of 0.6%/year for the 61% significant pixels across Switzerland. In particular, the seasonal mean NDVI shows an important break for winter, autumn and spring seasons starting from 2010, potentially indicating a critical point of changing land conditions. At biogeographical scale, we observed an apparent clustering (Jura-Plateau; Northern-Southern Alps; Eastern-Western Alps) related to landscape characteristics, while forested land cover classes are more responsive to NDVI changes. Finally, the NDVI values are mostly a function of temperature at elevations below the tree line rather than precipitation. The findings suggest that multi-annual and seasonal NDVI can be a valuable indicator to monitor vegetation conditions at different scales and can provide complementary observations for national statistics on the ecological state of vegetation to monitor land affected by changing environmental conditions. This work is aiming at strengthening the insights into the driving factors of vegetation change and supporting monitoring changing land conditions to provide guidance for effective and efficient environmental management and sustainable development policy advice at the national scale.

KEYWORDS:

• Landsat
• climate
• land cover
• NDVI
• anomaly

1. Introduction

Global environmental change and notably climate change are inducing important modifications significantly altering the functioning of the Earth system (Zhu et al., 2016) as well as provision of ecosystem services essential to human well-being (Rumpf et al., 2022). Vegetation is one of the major components of terrestrial ecosystems playing a crucial role in our environment, providing many ecosystem services, such as the habitat for wildlife, regulating the water and carbon cycle, influencing energy balance and absorbing carbon dioxide from the atmosphere (Chu et al., 2019; Yang et al., 2019). To effectively manage and protect natural resources, it is important to have accurate and up-to-date information on vegetation conditions and health.

Climatic factors are predominantly influencing vegetation dynamics and spatial distribution across regions and scales (Piao et al., 2006; Yang et al., 2019). Vegetation dynamics is related to changes in climate conditions due to biophysical responses of plant respiration, photosynthesis, and evapotranspiration (Chu et al., 2019). Vegetation holds about 42% of terrestrial carbon storage and assimilates up to about 20% of annual anthropogenic CO₂ production (Novillo et al., 2019). As greenhouse gas concentrations continue to rise, temperatures are expected to follow a similar trend (Beniston et al., 1994). Therefore, it is fundamental to analyze vegetation dynamics and its response to climate change, and it has been recognized as a major issue for global change in terrestrial ecosystems (Yang et al., 2019).

Over the past two decades, indications of increased vegetation growth and biological activity have been observed in many locations, including Europe, Sahel, parts of Australia, India and China (Brandt et al., 2015; Jiang et al., 2022; Jong et al., 2012; Lu et al., 2021; Lucas et al., 2019; Wu et al., 2019; Zhu et al., 2016). Despite this greening trend, many forested biomes are also showing a decline in vegetation productivity mostly reported in large parts of boreal forests, North America and Southeast Asia (Jong et al., 2012; Murthy & Bagchi, 2018; Nash et al., 2017; Zhu et al., 2016, 2016). This greening trend, mostly driven by temperature, enables plant species growing faster and taller (Rumpf et al., 2022), shows significant spatial and temporal variability (Choler et al., 2021). Therefore, there is a need to better understand the spatio-temporal dynamics of ecosystem’s productivity (Jong et al., 2012).

Switzerland, being covers by 60% of mountains and being in the middle of Europe, is recognized as particularly sensitive to global warming. Indeed, the Swiss territory has experienced a temperature increase twice higher than the Northern Hemisphere average (Gobiet et al., 2014). Changes in precipitation distribution and temperature increase projected by climate scenarios will significantly impact the Swiss vegetation (Vittoz et al., 2013). Mountain ecosystems have undergone rapid changes in response of warming climate, and greening has been detected in the mountains of the Alps (Choler et al., 2021; Thornton et al., 2021).

To monitor and quantify terrestrial ecosystems' responses to climate change at multiple scales, remote-sensing data, have proven to be a useful tool (Cortés et al., 2021; Piao et al., 2015). Using vegetation indices such as the Normalized Difference Vegetation Index (NDVI) as an indicator of chlorophyll abundance correlated to vegetation density, photosynthetic capacities and greenness (Nash et al., 2017; Zhu et al., 2016), satellite Earth Observations (EO) are commonly used for evaluating long-term trends in vegetation cover. Indeed, NDVI being sensitive to phenology (Yang et al., 2019) is considered to have a valuable indicator of aboveground green biomass and photosynthetic activity (Cortés et al., 2021) as well as vegetation growth and coverage change. Spatiotemporal variations of vegetation greenness may be influenced by various factors, such as density and type of vegetation, leaf density, health conditions, climate and land cover change (Carlson et al., 2017), disturbance (Lastovicka et al., 2020) or plant succession dynamics (Knoflach et al., 2021). Integrated over time, this index can provide useful metrics to evaluate vegetation dynamics at various scales (Carlson et al., 2017).

Most of the analysis of trends in vegetation greenness have been performed on large areas (e.g. continents) using coarse spatial resolution sensors, such as the Moderate Resolution Image Spectroradiometer (MODIS) and Advanced Very High-Resolution Radiometer (AVHRR) (Mishra et al., 2015; Wu et al., 2019). However, assessing vegetation greenness evolution is also important at the national scale as vegetation offers a variety of ecosystem services of local importance (Angélica et al., 2022; Chênes et al., 2021). Moreover, in Switzerland, such assessment has never been done; and the use of this coarse spatial resolution sensors is not suitable for a country like Switzerland that is small with a variety of landscapes and mostly covered by mountainous areas that are known as more sensitive to climate change.

Consequently, the aim of this study is to assess the vegetation dynamics, at a higher spatial resolution, in Switzerland from 1984 to 2018 based on NDVI time-series. The objectives are to gain insights on (1) the spatiotemporal patterns of NDVI at national, biogeographical and pixel level; (2) identify which land cover classes are the most responsive to NDVI change; and (3) evaluate the possible relationships between the NDVI time-series and two important climate parameters (temperature and precipitation).

2. Study area

The study area covers all of Switzerland, between the latitudes of 45.82° and 47.81° North and the longitudes of 5.96° and 10.49° East. The total land area is 41,285 km². Switzerland can be divided into the six biogeographical regions (Figure 1), which are based on zones of homogeneity of flora and fauna, as well as hydrologic and topographic features (Gonseth et al., 2001).

Figure 1. Elevation of Switzerland and its six biogeographical regions.

These regions represent a diversity of landscapes and ecological features characterized by different climatic, geological and vegetation conditions. Based on the altitude these regions can be grouped into three major areas: the Alps (North, Central, South) characterized by high elevations and accounting for 60% of the total surface of the country while the Plateau (30%) and the Jura (10%) are with low elevation.

The climate of Switzerland is strongly affected by the Alps and its proximity to the Atlantic Ocean. It is characterized as moderately continental in the plateau, alpine in the mountain regions and more temperate in the Southern Alps. Switzerland has four distinct seasons with varying temperature and precipitation patterns following these different climatic conditions (NCCS, 2018).

3. Methodology and implementation

In this study, to assess vegetation dynamics over time, we used the Normalized Difference Vegetation Index (NDVI) using satellite earth observations data from the Landsat program. NDVI, measuring the state of the vegetation health based on the red (i.e. strongly absorbed) and near-infrared (NIR) (i.e. strongly reflected) bands, is one of the primary measures used for vegetation monitoring. It is calculated as follows (Huang et al., 2020) (1):

(1) $\mathrm{NDVI}=\left(\mathrm{NIR}-\mathrm{RED}\right)/\left(\mathrm{NIR}+\mathrm{RED}\right)$(1)

NDVI values range from − 1 to + 1 depending on the vegetation and land cover type. Negative values correspond to clouds, water, and snow, while positive values indicate roughly vegetation. Values close to zero correspond to bare soils and rocks. When vegetation is affected by dehydration or becomes sick, usually the colors turn from green to brown. This change is particularly significant in NIR light that is increasingly absorbed when vegetation deteriorates. Therefore, NDVI can provide accurate information on the presence of chlorophyll and can be correlated with plant health (Pettorelli et al., 2005).

NDVI production workflow and time-series analysis are summarized in Figure 2. The Swiss Data Cube (Giuliani et al., 2017) has been used to access and process Landsat of 35 years of Analysis Ready Data and generate annual and seasonal time-series, which were further used to detect anomalies and investigate possible correlations with national temperature precipitation and land cover datasets.

Figure 2. Analytical framework of the study.

3.1. Satellite data access and processing

To monitor vegetation greenness in Switzerland, an important challenge relates to the satellite EO data source to be used. Indeed, as mentioned in the introduction, using MODIS or AVHRR is not suitable for such a small country as Switzerland. Therefore, we need higher spatial resolution to measure NDVI to capture representative signal of greenness change. We decided to use Landsat data for two reasons: (1) improved spatial resolution at 30 m (instead of 250 m-500 m-1 km) and (2) the long time-series of observations freely available (e.g. >40 years) with images available every 16 days. Landsat is a series of satellites launched by the United States Geological Survey (USGS) that captures images of the Earth’s surface since 1972 (Wulder et al., 2022). Data collected by Landsat satellites are well suited to track changes in vegetation cover and health over time (Zhou et al., 2023). However, some challenges must be considered to produce reliable and consistent NDVI time-series using Landsat data, such as cloud cover, atmospheric interference, calibration, spectral resolution or data processing.

To tackle the challenges mentioned above, Switzerland has established a satellite EO Analysis Ready Data archive of the entire national territory. The Swiss Data Cube (SDC – http://www.swissdatacube.ch) is a tera-scale analytical cloud-based platform allowing users to access, analyze and visualize more than 40 years of optical (e.g. Sentinel-2; Landsat-5, -7, -8, -9) and radar (e.g. Sentinel-1) EO data (Chatenoux et al., 2021; Honeck et al., 2018; Truckenbrodt et al., 2019). The SDC facilitates national-scale analyses of large volume of spatially aligned and consistently calibrated satellite EO data. It is intended to contribute to the environmental and reporting mandate of the government (Dhu et al., 2019) while at the same time enabling any Swiss institutions to benefit from the wealth of information accessible from satellite EO data. It has already been used to monitor snow cover (Poussin et al., 2023), impacts of droughts on vegetation (Poussin et al., 2021), impacts of climate change on threatened species (Barras et al., 2021), or for SDG monitoring (Giuliani et al., 2020). The archive is updated daily and accounts approximately for 80,000 scenes corresponding to a total volume of 30 TB and more than 3,000 billion observations/pixels. It also recently made available additional official raster datasets from the Swiss government such as Land Use Statistics, Digital Elevation Model, and climate model outputs. The SDC is based on the Open Data Cube (ODC) software stack, an open-source geospatial data management and analysis software designed to handle large volumes of EO data (Dhu et al., 2017; Ferreira et al., 2020; Killough, 2019; Sudmanns et al., 2021, 2022).

The limited temporal resolution of Landsat satellites, combined with the presence of clouds, presents challenges in obtaining cloud-free observations. Therefore, to leverage the high spatial resolution, that is particularly important for a small country like Switzerland with a variety of landscapes on a small surface, we decided to use a multi-date compositing approach. Annual and seasonal NDVI time-series have been computed using Level 2/Collection 1 Landsat-5 and -8 data over the full geographic coverage of Switzerland for the period 1984 to 2018 and made available as a Python script in Jupyter Notebooks (Giuliani et al., 2020) freely and openly available in the SwissEnvEO data repository (Giuliani et al., 2021). Seasonal composites were produced for four periods: (1) December to February (DJF); (2) March to May (MAM); (3) June to August (JJA); and (4) September to November (SON). Winter season has been included in the analysis for two reasons: (1) there can be significant vegetation activity in certain regions of Switzerland, particularly in low altitude areas where temperatures remain relatively mild and snow cover is absent; and (2) in studying long-term trends, considering the winter months can be useful for examining the impacts of climate change on vegetation dynamics. By incorporating these months, we obtain a more comprehensive understanding of the seasonal dynamics of vegetation in Switzerland.

It should be noted that Landsat-7 has not been used because since 2003 it has data gaps induced by a Scan Line Corrector failure (Wulder et al., 2008). Consequently, no data were available for the year 2012. Moreover, to account for cloud perturbations, clouds were filtered out using the cloud cover information included in the Landsat Collection 1 Level-1 Quality Assessment (QA) Band (Dwyer et al., 2018; Ernst et al., 2018). In this analysis, only clear observation pixels labeled as “clear”, “water” or “snow” were used. All the other QA assessment band attributes were considered as “cloud-covered” pixels. The percentage of clear observations per pixel can be determined using the following Equation (2)(2) $\mathrm{Pixel}\mathrm{Quality}\left[\mathrm{\%}\right]=\left(\mathrm{number}\mathrm{of}\mathrm{clear}\mathrm{observations}/\mathrm{total}\mathrm{observations}\right)\times 100$(2) :

(2) $\mathrm{Pixel}\mathrm{Quality}\left[\mathrm{\%}\right]=\left(\mathrm{number}\mathrm{of}\mathrm{clear}\mathrm{observations}/\mathrm{total}\mathrm{observations}\right)\times 100$(2)

Figure 3 shows the percentage of clear observations per pixel across Switzerland for the considered period. The percentage of cloud-free observations per pixel and for the entire period reached up to 35% in some areas in the north of the country and in some valleys (yellow to light blue areas). However, most of all pixels (59%) are only cloud-free in 10% to 20% of all observations. Consequently, only pixels with more than 10% of the clear observations have been considered, and 23% of all pixels (i.e. white areas in Figure 3) had to be excluded for the analysis. Most of these areas are in the Alps that are particularly prone to cloud cover. For the winter season (DJF) around 10% of pixels are excluded, 15% of pixels for the spring season (MAM), 22% of pixels for the summer season (JJA) and 17% of pixels for the autumn season (SON). Likewise, some seasonal periods had to be excluded from the analyses because of cloud cover conditions such as the winter seasons for the year 2003, 2008 and 2013, the 2003 spring season and the 1986 and 2002 autumn seasons.

Figure 3. Percentage of clear observations per pixels across Switzerland for the period 1984–2019 for the annual (top) and seasonal mean (DJF (December/January/February), MAM (March/April/May), JJA (June/July/August), SON(September/October/November)). White zones indicate pixels with less than 10% of clear observations over the study period. The strips are caused by Landsat acquisitions paths that can be overlapped (and may duplicate the amount of data).

In addition of satellite data, we used the land cover dataset obtained from the national Land Use statistics produced by the Federal Office for Statistics (FSO) (https://s.geo.admin.ch/9f4ff25119). The land cover dataset is derived from aerial photography and provides information on the evolution of the Swiss land surface at the hectare scale (i.e. spatial resolution of 100 × 100 m) since 1979 over four survey periods (1979/85, 1992/97, 2004/09 and 2013/18). It will further be used to determine NDVI trends in the different vegetation categories (Table 1) extracted from the 27 land cover classes following the NOLC04 nomenclature for the latest survey campaign 2013/18 (Giuliani et al., 2022; Swiss Federal Statistical Office, 2001). This dataset is static and does not vary annually as explained in: https://www.bfs.admin.ch/bfs/en/home/services/geostat/swiss-federal-statistics-geodata/land-use-cover-suitability/swiss-land-use-statistics.html

Table 1. Vegetation categories following the land cover NOLC04 classification separated in principal domains and basic categories.

Principal domains	Basic categories
10 – Artificial areas	15 – Lawns
	16 – Trees in artificial areas
20 – Grass and herb vegetation	21 – Gras and herb vegetation
30 – Brush vegetation	31 – Shrubs
	32 – Brush meadows
	33 – Short-stem fruit trees
	34 – Vines
	35 – Permanent garden plants and brush crops
40 – Tree vegetation	41 – Closed forest
	42 – Forest edges
	43 – Forest strips
	44 – Open forest
	45 – Brush forest
	46 – Linear woods
	47 – Cluster of trees

3.2. Time-series analysis and statistics

To identify trends, annual and seasonal mean NDVI time-series are assessed through zonal statistics (mean, max, min, standard deviation) across national and biogeographic divisions, as well as relevant land cover classes. These statistics are computed using official national and biogeographic delimitations from the swissBOUNDARIES3D database by swisstopo. Linear regression analysis (de Beurs & Henebry, 2004; Jong et al., 2012) determines NDVI trends over time for different zones (national, biogeographic, land cover). Associated with trend lines, two statistical measures emerge: the correlation coefficient (R²) and p-value. R², ranging from 0 to 1, gauges linear relationship strength, while slope indicates positive/negative correlation. The p-value indicates statistical confidence. A low value (p-value < 0.05) signifies a significant trend, while a higher p-value suggests unrelated predictor and response variable changes.

Additionally, for each NDVI time-series, an autocorrelation (ACF) and partial autocorrelation (PACF) were performed. The autocorrelation of a series indicates that in a time or spatial series, the measurement of an event at time t can be correlated with previous or subsequent measurements. On the other hand, the partial autocorrelation will observe whether a correlation exists between the residuals of effects already explained by the previous lags (Javed et al., 2020, 2021).

To complete the statistical analysis, the Mann–Kendall test will be performed for each annual, regional, and seasonal analysis to determine if there is a greening trend. Unlike simple linear regression, this test does not require that the data be normally distributed or linear. Its null hypothesis is that there is no trend in the data, and the alternative hypothesis is that the data represent a monotonic trend in a time series that may include a seasonal component (Mann, 1945). It is one of the most widely used tests to detect changes in climate time series (Hamed & Ramachandra Rao, 1998) and to determine greening trends (de Beurs & Henebry, 2004; Wessels et al., 2012). Its correlation coefficient (τ) ranges from −1 to 1. The closer it is to 1, the more it indicates a strong positive relationship and inversely. Significance of trends is also tested with a 95% confidence interval level (p-value < 0.05).

3.3. Anomalies

NDVI anomalies were calculated to assess the deviation from the long-term average at a given time. In other words, they are used to identify years that show a particularly high or low greening episodes, based on a specific reference period. These were calculated for each year, season, and decade period based on the standard score (Z-score) using the following formula (3):

(3) $\mathrm{NDVI}{\mathrm{Z}}_{\mathrm{t}}=\left({\mathrm{X}}_{\mathrm{t}}-{\mathrm{X}}_{\mathrm{mean}}\right)/{\mathrm{X}}_{\mathrm{std}}$(3)

X_t represents the pixel value of NDVI for the year, season, or decade of interest. X_mean corresponds to the average NDVI values and X_std is the standard deviation for the reference period 1984–2018.

Anomalies can have positive or negative values depending on whether the data obtained are lower or higher than the long-term average of the entire series. Positive values indicate that greening is higher than the average of the reference period. Conversely, negative values indicate that greening is below the baseline average. These values can be a good indicator of episodes where vegetation greening is very low or very high caused by diverse environmental factors, at the period of year, or during specific climate-related events (e.g. drought).

To facilitate the interpretation, understanding, and visualization of anomalies, the z-score values were reclassified into seven classes ranging from extremely good to extremely bad, following the classification proposed by the World Meteorological Organization (WMO) summarized in Table 2 (Meroni et al., 2019)

Table 2. Classification of NDVI z-score anomaly values based on the methodology of Meroni et al. (2019).

NDVI z-score anomaly
Value	Class
x ≥ 2	Extremely good
1.5 ≤ x < 2	Very good
1 ≤ x < 1.5	Moderately good
−1 < x < 1	Near normal
−1.5 < x ≤ −1	Moderately bad
−2 < x ≤ −1.5	Very bad
x ≤ −2	Extremely bad

3.4. Correlation analysis

Vegetation greening is influenced by different climatic factors, most notably temperature and precipitation (Chu et al., 2019; Gao et al., 2022; Hou et al., 2015; Jong et al., 2012). Therefore, it could be valuable to explore the possible relationship between NDVI and these variables with changing climatic conditions over the study period. Gridded temperature and precipitation datasets have been produced and made available by MeteoSwiss. Temperatures are calculated at 2 m above the surface in degree Celsius (℃) and averaged using daily mean values, while precipitations represent accumulated precipitation, including rainfall and snowfall water equivalent, expressed in millimeters (mm). Both datasets are at 1-km spatial resolution, with seasonal and annual means, covering the entire country over the period 1984–2018 for the temperature and over the period 2005–2018 for the precipitation (MeteoSwiss, 2016a, 2016b).

To ensure that the remote sensing data and the climatic variables have the same spatial resolution, all NDVI data were spatially resampled using the nearest neighbor method at 1 km and realigned to the coordinates of the climate variables. Pixels with less than 10% of clear satellite observations have also been removed from the climatic dataset to have the same number of pixels as the satellite data.

For assessing the linear link between climatic factors and NDVI, Pearson correlation coefficients (r and p-value) gauged annual and seasonal means in Switzerland and its six biogeographic zones over the study period. These coefficients range from −1 to 1, signifying the strength of the linear relationship. A coefficient nearing 1 indicates a robust positive correlation, suggesting when one variable rises, the other does too. Conversely, proximity to −1 suggests a strong negative correlation. Note that this test assumes a normal distribution, confirmed by the Shapiro–Wilk test.

Furthermore, a pixel-by-pixel Spearman correlation analysis between the annual and seasonal NDVI time series was computed to observe the spatial distribution of the relationship over time and between the different climatic variables. This non-parametric test measures the monotonic relationship between two variables.

4. Results

4.1. National and biogeographical NDVI time-series evolution

4.1.1. Annual mean NDVI trend at the national scale

High mean NDVI values for the whole time-period were located in low to medium altitude regions (Jura Mountain, the plain and the Swiss valleys) (Figure 4a). Inversely, smaller NDVI values can be found in higher altitude regions. From the linear regression computed for 77,761,233 pixels across Switzerland over the 34 years, more than 61% of the NDVI pixels have a significant trend at 95% confidence level (Figure 4b). Within these pixels, more than 97% has a positive trend indicating an increase in NDVI over a large part of the Swiss territory for the period. However, significant values in the NDVI trend map are relatively small (i.e. significant trend mean = 0.006) indicating a slightly NDVI change over the period 1984–2018.

Figure 4. (a) overall mean annual NDVI, and (b) annual NDVI significant slope at the 95% significance over the period 1984–2018 across Switzerland.

At the national scale, a positive and significant greening trend is observed for the entire Swiss territory for the period 1984 to 2018 (Figure 5). The correlation coefficient and its p-value (R² = 0.43 and p-value < 0.05) indicate that the strength of the observed relationship between the annual means of NDVI during the observed periods is relatively moderate and significant. The Mann–Kendall test supports the results of the linear regression, i.e. that there is a moderate correlation (tau = 0.49 and p-value < 0.05) between the range of the observed NDVI values and their order in time. Although there is an overall increase in the greening trend over the period 1984 to 2018, the change per year remains small (slope = 1.11 × 10⁻⁵/year). Mean annual NDVI values showed variations ranging from 0.34 in 2002 to 0.57 in 2013 with an overall mean for the whole period is 0.47. Mean seasonal time series (Figure 5 dashed line) described well the seasonal pattern with high NDVI values in summer (mean NDVI of 0.63) and autumn (mean NDVI of 0.46), and low NDVI values in winter (mean NDVI of 0.10) and spring (mean NDVI of 0.30). The ACF and PACF tests reveal that there is no autocorrelation within the NDVI means but pointed out seasonal patterns.

Figure 5. Mean annual NDVI (solid line) and mean seasonal NDVI (dashed line) trend for Switzerland for the period 1984–2018. The regression line is illustrated in red.

To understand the effects of inner-annual variations, the use of the seasonal time-series was considered. The analysis of seasonal NDVI time-series reveals that the greening trend has increased slightly (statistically significant) for each season over the period 1984 to 2018 (Figure 6 and Table 3). However, a stronger trend is detected in summer (R² = 0.31 and p-value < 0.05) as well as in winter (R² = 0.29 and p-value < 0.05) although the correlation coefficients indicate a relatively weak relationship between NDVI, and the period studied (Table 3). The mean NDVI values for summer are generally much higher than for the other seasons and range from 0.44 to 0.71. Logically, the lowest values are found in winter, ranging from 0.00 to 0.29, while in the other seasons, values range from 0.10 to 0.43 for spring and from 0.18 to 0.62 for autumn. The seasonal NDVI means appear to be higher for each season starting in 2010. Indeed, in winter, autumn, and summer, the highest NDVI mean are recorded in 2016 with an NDVI value of 0.29 for the winter season, 0.62 for the autumn season, and 0.71 for the summer season, respectively. In the spring, the highest average was recorded in 2017 (= 0.43).

Figure 6. Mean seasonal NDVI trend for Switzerland (DJF (December/January/February), MAM (March/April/May), JJA (June/July/August), SON(September/October/November)) for the period 1984–2018. Blank values correspond to pixels that do not have enough cloud-free observations.

Table 3. Linear regression coefficient and Mann Kendall coefficient for the mean seasonal NDVI time series for Switzerland for the period 1984–2018 (DJF (December/January/February), MAM (March/April/May), JJA (June/July/August), SON (September/October/November)). The (*) indicate statistically significant correlation at 95% confidence level (p-value < 0.05).

	DJF	MAM	JJA	SON
R^2	0.29	0.23	0.31	0.18
P_value	0.0013*	0.002722*	0.000355*	0. 00965*
Slope	0. 00001136	0.00001091	0.000007851	0.00001110
tau	0.34	0.41	0.35	0.42

4.1.2. Annual mean NDVI trend at biogeographical scale

Regional analyses for the period 1984 to 2018 reveal a positive and significant trend for all biogeographical regions (Figure 7 and Table 4). In 2002, we observe a very low mean NDVI for several regions, especially for the Eastern Central Alps (= 0.03) and the Northern Alps (= 0.29) which reach their lowest NDVI mean values (Table 4). For most of the biogeographic zones, the highest NDVI means are reached from the years after 2010, except for the Northern Alps, which reaches it in 2007 (= 0.58). This indicates that before the 2000s, a generally stable mean NDVI was observed, but it is from the 2000s onwards that a significant increase is observed.

Figure 7. Mean annual NDVI trend for the six biogeographical regions for the period 1984–2018.

Table 4. Linear regression coefficient and Mann Kendall coefficient for the mean annual NDVI time series for the six biogeographical regions. The (*) indicate statistically significant tau correlation at 95% confidence level (p-value < 0.05).

	Jura	Plateau	Northern Alps	Western central Alps	Eastern central Alps	Southern Alps
R^2	0.50	0.40	0.37	0.22	0.19	0.31
P_value	0.00000202*	0.0000334*	0. 000082*	0.00324*	0.0005228*	0.000354*
Slope	0.00001323	0.000009317	0.00001263	0.000009585	0.00001015	0.00001364
tau	0.52	0.49	0.46	0.33	0.37	0.32

In general, the Jura (R² = 0.50, p-value =< 0.05) and the Central Plateau (R² = 0.40, p-value < 0.05) areas show the strongest correlations of the lot with correlation coefficients of medium intensity. This suggests that lower elevation regions seem to be the most responsive to changes in NDVI values during the study period. In these same areas, we observe the highest NDVI means, whether at the minimum, maximum or general level. This is supported by the Mann–Kendall test, which shows relatively high and positive tau coefficient for the Jura (tau = 0.52 and p-value =< 0.05) and for the Central Plateau (tau = 0.49 and p-value < 0.05). Analyses by region revealed three groups of different behaviors by elevation. Indeed, we could notice that the low elevation areas such as the Jura and the Plateau had the same dynamics with higher NDVI values than alpine regions. Also, the Southern and Northern Alps differed from the central Eastern and Western Alps, which had the lowest mean annual NDVI values.

To observe more precisely the variations of the NDVI values, regional analyses for each season over the 35 years were conducted (Figure 8).

Figure 8. Mean seasonal NDVI trends by biogeographical regions for the period 1984–2018 (DJF (December/January/February), MAM (March/April/May), JJA (June/July/August), SON (September/October/November)). Blank values correspond to pixels that do not have sufficient numbers of cloud-free observations.

The Mann-Kendall trend test reveals a positive and significative greening trend in time series except for the eastern central Alps during the autumn season. The seasonal analyses show that the low elevation regions (the Central Plateau and the Jura) have a similar behavior and have for all seasons the highest mean NDVI values. We can also observe two different groups of behavior between the alpine regions. Indeed, the Eastern and the Western Central Alps follow the same dynamics and have the smallest mean NDVI values for all seasons. However, for the winter season, the Western Central Alps region stands out from the rest of the Alps, with the lowest NDVI values. Winter, spring, and autumn seasons present considerable annual and spatial variability, in contrast to the summer season. Since the 2000s, we can observe more events with high NDVI values in all seasons and for all biogeographical regions. Similarly, the majority of the highest mean NDVI values are found after 2010, indicating that the NDVI trend has undergone an increase from the 2000s and especially from 2010. For the winter season, it is interesting to notice that the mean NDVI values are negative between 1984 and 2011 and then become positive. This is in line with previous findings which state that a stronger increase in vegetation productivity is observed from 2010 onwards.

4.1.3. NDVI trend by land cover classes

Results of the linear regressions summarized in Table 5 show that, in general, a positive and significant correlation is observed for most of the land cover classes (p-value < 0.05 for each class). On the other hand, class 32 (Brush meadows) obtains a correlation coefficient too low to indicate that there is a greening trend in this category (R² = 0.16). For the other classes, there is a statistical relationship between the NDVI averages, and the period elapsed within each class indicating a greening trend ranging from low to moderate according to the r-squared intensity. However, classes 41 “Closed forest” (R² = 0.57 and p-value < 0.05), 34 “Vines” (R² = 0.56 and p-value < 0.05) and 44 “Open forest” (R² = 0.4946 and p-value < 0.05) have the highest correlation coefficients of the group, indicating a relatively moderate correlation at a high significance level. In other words, these categories seem to be more responsive to NDVI variations. Moreover, the Mann-Kendall correlation coefficients confirm a strong significant correlation between the rank of the NDVI values and their order between 1984 and 2018 for the three selected classes (see Table 4). Figure 9 shows the evolution of the greening trend for each of the classes and allows us to observe an overall increase of the latter during the period studied.

Figure 9. NDVI annual trend of the land cover class “Open forest” between 1984 to 2018 (slope = 0.00001642; p-value = 0.00000211).

Table 5. Minimum, maximum, and full period (1984–2018) averages by biogeographic zone for each season (summer, fall, winter, spring).

Summer		Mean NDVI
Region	Minimum average	Maximum average	Mean (1984–2018)
Plateau	0.495911 (1984)	0.705233 (2014)	0.641716
Eastern central Alps	0.100609 (2002)	0.70328 (2016)	0.542122
Jura	0.50978 (1984)	0.777003 (2013)	0.698546
Northern Alps	0.459359 (1984)	0.762514 (2016)	0.67113
Southern Alps	0.281476 (2008)	0.755286 (2013)	0.61573
Western central Alps	0.252922 (1984)	0.561585 (2015)	0.44638
Fall		Mean NDVI
Region	Minimum average	Minimum average	Mean (1984–2018)
Plateau	0.227449 (2000)	0.709202 (2013)	0.520063
Eastern central Alps	0.145048 (2000)	0.529749 (2007)	0.339564
Jura	0.021256 (2000)	0.715227 (2016)	0.510779385
Northern Alps	0.199488 (2000)	0.652077 (2016)	0.468775
Southern Alps	0.086188 (2008)	0.641462 (2016)	0.403835
Western central Alps	0.10229 (1998)	0.423621 (2018)	0.260252
Winter		Mean NDVI
Region	Minimum average	Maximum average	Mean (1984–2018)
Plateau	0.024341 (1999)	0.548422 (2014)	0.160238
Eastern central Alps	−0.05841 (1991)	0.124206 (2016)	−0.01764493
Jura	0.001027 (1999)	0.463978 (2014)	0.146457
Northern Alps	0.004408 (2006)	0.299186 (2017)	0.04853
Southern Alps	0.00769 (2009)	0.299677 (2016)	0.054057
Western central Alps	−0.04569 (1987)	0.070245 (2017)	−0.01356073
Spring		Mean NDVI
Region	Minimum average	Maximum average	Mean (1984–2018)
Plateau	0.171140072 (2002)	0.653109456 (2017)	0.48515814
Eastern central Alps	−0.048677218 (1992)	0.147295607 (2017)	0.03633505
Jura	0.11804778 (2005)	0.684873278 (2017)	0.46402739
Northern Alps	0.062479818 (2005)	0.422191889 (2014)	0.23928331
Southern Alps	0.055142098 (1984)	0.418264474 (2013)	0.17921633
Western central Alps	−0.034876464 (1999)	0.151900429 (2018)	0.04370726

4.2. NDVI z-score anomalies

NDVI anomalies were computed at the national scales for annual means for different reference periods to identify greening conditions. Figure 10 shows the NDVI during-score anomalies for the reference period 1984 to 1994. The lowest Z-score is for the year 2002 (i.e. NDVI mean z-score of −1.19 for Switzerland), and the highest for the year 2013 (i.e. NDVI mean z-score of 1.85 for Switzerland). As illustrated in Figure 10, the majority of years with a very positive anomaly compared to the reference period are mostly recorded after 2000 (i.e. dark green and dark blue in Figure 10). Conversely, years 1984, 1992, 2002 and 2005 show unusual low NDVI z-score anomaly values (NDVI mean z-score of −0.58, −0.85, −1.19 and −0.49, respectively) illustrated by red pixels throughout the country classified as extremely bad vegetation conditions. Overall, a significant increase in greening can be observed in Switzerland from the 2000s onwards, with especially high values from 2010 onwards, confirming the previous results indicating that the increase in vegetation greenness has been more consistent over the last decade.

Figure 10. Mean annual NDVI score anomalies for Switzerland for the period 1984–2018. Pixel classes are based on the classification of Meroni et al. (2019). White zones indicate pixels with less than 10% of clear observations over the study period.

NDVI anomalies were calculated at the national and regional scales as well as for relevant land cover classes (Table 6), according to different decades with respect to the reference period from 1984 to 2018 (Table 7). Negative values are generally observed for the periods before 2000 and positive values from 2001 onwards. Significantly higher values are observed for all regions in the year range from 2011 to 2018, again indicating a more consistent increase in the greening trend since 2011 compared to the average of the reference period (1984–2018).

Table 6. Summary table of results from linear regressions and Mann-Kendall tests performed for annual means of NDVI for each land cover class during 1984–2018. Results followed by * indicate a significant trend with 95% confidence interval (p-value < 0.05), ** indicates a highly significant trend (p-value < 0.01) and *** a highly significant trend (p-value < 0.001).

Land Cover classes	slope	p-value	R^2	tau	tau p-value
15 – Lawns	.000007013	.000317 ***	0.3166	0.422		.00046778
16 – Trees in artificial areas	.00000795	.0000329 ***	0.4031	0.462		.00013089
21 – Gras and herb vegetation	.000009111	.000141 ***	0.3488	0.458		.00014758
31 – Shrubs	.00001035	.00238 **	0.2305	0.373		.0020459
32 – Brush meadows	.000008859	.0118 *	0.1567	0.34		.0048528
33 – Short-stem fruit trees	.00001212	.0000636 ***	0.3791	0.497		.00003767
34 – Vines	.00001296	.000000230***	0.5586	0.572		.0000021458
35 – Permanent garden plants and brush crops	.00000931	.0000669***	0.3772	0.483		.000062704
41 – Closed forest	.00001754	.000000132 ***	0.5733	0.551		.0000050068
42 – Forest edges	.00001301	.00000821***	0.4512	0.515		.00001955
43 – Forest strips	.00001231	.0000157***	0.4292	0.497		.00003767
44 – Open forest	.00001642	.00000211 ***	0.4946	0.504		.000029087
45 – Brush forest	.00001251	.001286 **	0.2577	0.383		.0015117
46 – Linear woods	.00001051	.0000275***	0.4096	0.483		.000062704
47 – Cluster of trees	.00001188	.0000323***	0.4038	0.455		.0001663

Table 7. Mean NDVI values and anomalies for Switzerland and the six biogeographical regions by decades. The baseline period is 1984–1994.

	Mean NDVI/NDVI z-score anomalies
Region	1984–1990	1991–2000	2001–2010	2011–2018	1984–1994 (reference period)
Switzerland	0.42/-0.405	0.44/-0.201	0.47/0.034	0.53/0.644	0.46
Jura	0.53/-0.697	0.58/-0.107	0.60/0.099	0.66/0.708	0.59
Plateau	0.53/-0.409	0.56/-0.131	0.56/-0.095	0.63/0.733	0.57
Northern Alps	0.40/-0.480	0.44/-0.168	0.47/0.097	0.52/0.582	0.46
Western central Alps	0.24/0.156	0.24/-0.276	0.28/0.061	0.33/0.505	0.27
Eastern central alps	0.25/-0.196	0.25/-0.336	0.30/0.159	0.34/0.450	0.28
Southern Alps	0.35/-0.279	0.34/-0.424	0.40/0.037	0.50/0.835	0.39

4.3. Correlations between NDVI and climate variables

Annual and seasonal NDVI estimates at the national, biogeographical and pixel scales have been correlated with temperature and precipitation datasets using Pearson’s correlation method for the period 1984 to 2018.

At the national-scale, a moderate positive and significative correlation is observed between annual NDVI means and temperature (r = 0.52 and p-value < 0.05) (Table 8(A)). A positive correlation indicates in this case that the values of NDVI means tend to increase when those of temperatures increase.

Table 8. Correlation results between annual and seasonal NDVI and climate variables ((A) temperature and (B) precipitation) at the national and biogeographical scale for the period 1984–2018. The (*) indicate statistically significant correlation at 95% confidence level (p-value < 0.05).

		Jura	Plateau	Northern Alps	Eastern central Alps	Southern Alps	Western central Alps	Switzerland
(A) Temperature–NDVI
Coefficient of correlation (r)	Annual mean	0.71*	0.57*	0.44*	0.00	0.41*	0.19	0.52*
	Seasonal mean (DJF)	0.47*	0.76*	0.61*	0.16	0.44*	0.19	0.70*
	Seasonal mean (MAM)	0.45*	0.34	0.52*	0.29	0.40*	0.19	0.45*
	Seasonal mean (JJA)	0.49*	−0.41*	0.12	−0.14	−0.08	0.33	0.18
	Seasonal mean (SON)	0.17	0.26	0.48*	0.33	0.42*	0.53*	0.40*
(B) Precipitation–NDVI
Coefficient of correlation (r)	Annual mean	−0.22	−0.20	−0.41	−0.15	−0.41	−0.27	−0.41
	Seasonal mean (DJF)	−0.22	−0.18	0.14	0.43	−0.18	0.10	−0.03
	Seasonal mean (MAM)	−0.62*	−0.63*	−0.37	−0.17	−0.19	−0.31	−0.61*
	Seasonal mean (JJA)	−0.33	−0.46	−0.30	−0.06	−0.40	−0.52	−0.45
	Seasonal mean (SON)	−0.10	−0.18	−0.07	0.09	−0.68*	0.24	−0.29

At a regional scale, the Jura (r = 0.71 and p-value < 0.05) and the Central Plateau (r = 0.57 and p-value < 0.05) regions obtain the strongest NDVI trend indicating strong significant relationship between temperatures and NDVI means for these regions. Conversely, the Western Central Alps (r = 0.19) and Eastern Central Alps (r = 0.00) show a very weak and/or non-significative correlation with temperature (Table 8(A)). In general, it can be observed that the low altitude areas (Jura, Central Plateau) are more strongly correlated with temperatures than the high-altitude areas (Northern Alps, Southern Alps, Western Central Alps, Eastern Central Alps).

The correlations between the seasonal means of NDVI for Switzerland and temperature are positive and significant for most seasons, except for the summer season (r = 0.18 and p-value = 0.327) which shows no significant relationship with the climate variable. Conversely, the winter seasonal NDVI mean for the period of 1984 to 2018 reveals the most significant trend with a very strong and significant correlation (r = 0.70 and p-value < 0.05) between NDVI values and temperatures. In winter, NDVI values in the Central Plateau (r = 0.76 and p-value < 0.05) show a particularly strong correlation with the temperature variable.

Regarding the pixel-by-pixel analyses, for temperature, 45% of all pixels are significant with mainly moderate and positive correlations (99% of the significative pixels) across Switzerland (mean correlation of 0.46*). Illustrated by Figure 10(A), low elevation regions and the Southern Alps are strongly correlated with temperature than the higher elevation regions. As with the previous results, most of the seasons also indicate the presence of a weak and significant correlation, except for the summer season (r = 0.17 and p-value < 0.05).

However, the correlations between annual and seasonal NDVI mean values and precipitation are mostly negative and non-significant at the national and regional scales over the 35 years (Table 8(B)). This tells us that NDVI values tend to decrease when precipitation increases and vice versa. On the other hand, seasonal results of annual NDVI means allow us to observe a negative and significant correlation of moderate intensity (r = −0.61 and p-value < 0.05) in the spring. If we look more precisely at the regional level, we can observe that the low elevation regions (the Central Plateau and the Jura Mountain) are the only regions with a strong negative and significant correlation with precipitation during the spring season (Table 8(B)).

The pixel-wise analysis (Figure 11) highlights a very weak and negative correlation between precipitation and annual mean NDVI values for Switzerland with a mean correlation of −0.14*. Nevertheless, only 6% of the pixel across the country are significative with a few more negative values (59% of the significative pixels). Concerning the seasonal correlations, we observe a moderate and significant negative correlation for the winter season (r = −0.60 and p-value < 0.05) and a weak negative correlation for the spring season (r = −0.30 and p-value < 0.05). In general, Figure 10(B) shows a large spatial variability in the correlation between NDVI means and precipitation indicating that it is difficult to establish a relationship between these two variables even at the regional level.

Figure 11. Pixel-wise correlation between annual mean NDVI, annual mean temperature (a) and annual precipitation (b) over the period 1984–2018. Only significant pixel (correlation at 95% confidence level (p-value < 0.05)) are shown in the figure.

5. Discussion

5.1. NDVI spatiotemporal patterns across Switzerland

Annual, seasonal, and regional findings reveal consistent greening from 1984 to 2018 across Switzerland and its areas, notably in low-altitude zones like Jura and Central Plateau. Contrarily, higher elevations exhibit milder greening. This aligns with the results in Filippa et al. (2019) showing alpine areas above 1800 m have limited vegetation growth due to snow cover. Regionally, three elevation-based behaviour clusters emerge. Low areas like Jura and Plateau exhibit similar dynamics with higher NDVI than alpine regions. Southern/Northern Alps differ from Central Eastern/Western Alps, having lowest NDVI. Nonetheless, the study indicates an earlier snowmelt trend from 2000 to 2018, potentially impacting future vegetation. This greening aligns with climate change-driven long-term greening and rising temperatures in cold/snowy ecosystems (Choler et al., 2021). Additional factors, like extended growing seasons (Vittoz et al., 2013) and higher tree lines (Gehrig-Fasel et al., 2007), contribute to greening in both low- and high-altitude zones.

It can be difficult to observe meaningful trends in vegetation greenness levels across a broad temporal range (Jong et al., 2012), so regional and seasonal analyses allow us to observe variations more accurately in NDVI values and to observe slightly different trends for certain regions during particular seasons. The results of the seasonal means at the national scale show that summer and winter seasons present the stronger greening trend. It is also during these seasons that MeteoSwiss observes the greatest regional and seasonal variations in the warming trend in Switzerland, with temperature increases being greatest in the winter for low-altitude regions and in the summer for alpine regions. In the winter, at the regional scale, it is observed that the Southern Alps and the Eastern Central Alps show the strongest relationship with NDVI during the considered period.

The Southern Alps and Eastern Central Alps display the strongest NDVI relationship in the studied timeframe. The research of Gehrig-Fasel et al. (2007) notes accelerated vegetation growth since 1985 in the Eastern Central Alps, with a doubling forest area on the southern slopes over 120 years (Cioldi et al., 2021). Vittoz et al. (2013) find climate change advances plant phenology and invasive species due to higher temperatures, noting over 20 invasive shrubs/tree species in southern Swiss forests since 1970 due to reduced frost. Decreased snow cover can also spur winter vegetation growth. Winter’s start has shifted by 15 days, with 20–30 days less snow in Swiss Alps. The winter zero temp’s rise from 600 m to 850 m, possibly reaching 1300–1500 m by 2060, reducing snowfall by 30% (NCCS, 2018). This suggests alpine vegetation’s greenness will likely be impacted by diminishing snowpack.

In summer, it is the Jura and western central Alps regions that show a stronger greening trend compared to other biogeographic zones. The Jura is one of the most forested regions in Switzerland (41%) after the Southern Alps (52%), which could explain why the lowland area has one of the highest greening trends. Regarding the Western Central Alps, the phenomenon of earlier snowmelt mentioned above (Filippa et al., 2019) and increasing temperatures in alpine areas in the summer may potentially contribute to the increase in vegetation greenness in these high elevation regions.

All results obtained in this study show that the highest annual and seasonal means of NDVI, both for Switzerland and its biogeographic zones, are predominantly recorded after 2010. This indicates that the trend before the 2000s was stable and that since then there has been a significant increase in vegetation greenness. The produced anomalies maps allow us to visually observe this trend, with higher annual and regional scores especially in the decade from 2011 to 2018. In addition, the maps showing extreme anomalies by season and for all biogeographic areas support the finding of a greater increase in plant greenness from 2010 onwards in every season, and for all regions.

Observations during the exceptional drought years, 2003 and 2018, reveal high NDVI averages, contrasting the study of Jolly et al. (2005) on the 2003 drought’s impact on plant productivity. Jolly’s work suggests lower elevation plant growth decrease due to water stress linked with the 2003 heatwave, while higher elevations might benefit. The 2003 annual average could be high due to no data removal for winter/spring, having lowest NDVI values (averaging 0.10 and 0.30). Another factor is groundwater and evapotranspiration’s influence. Evapotranspiration, along with precipitation, shapes water resource availability for ecosystems. Temperature rise can enhance summer evapotranspiration, potentially leading to drying conditions. Despite Switzerland’s high rainfall (1400 mm/year average), more frequent/intense droughts due to rising temperatures are seen globally and locally. The study of Poussin et al. (2021) on Switzerland’s 35-year drought index detects slight drying trends, mainly at low/medium altitudes and Southern Alps. This trend may decrease vegetation productivity and reverse greening. Remote-sensing offers insights on variables like temperature, albedo, and vegetation cover, enabling evapotranspiration estimation. Swiss data scarcity on this topic is improving but varied collection objectives hinder trend observation (NCCS, 2018).

Climate scenarios predict a greater frequency and intensity of climatic phenomena, such as floods, landslides, or winter storms (Vittoz et al., 2013). These phenomena can potentially impact vegetation dynamics and explain why very low or very high NDVI averages are observed in some years. In the winter of 1999, low-lying areas have the lowest average NDVI recorded during the winter period from 1985 to 2018. This coincides with the Lothar winter storm in December 1999, which devastated forests in the Plateau and Jura regions (Wohlgemuth & Kramer, 2015). Thus, these events could have potentially contributed to the decline in plant greenness seen at this particular time period.

5.2. Response of NDVI to land cover classes

In general, a greening trend is observed throughout Switzerland in terms of land cover classes. The results highlight three vegetation classes with high vegetation greenness during the study period: closed forest, vines, and open forest. Regarding the open and closed forest classes, the increase in forest area by 30% in Switzerland (Gehrig-Fasel et al., 2007) may contribute to the strong greening trend observed in these classes. Indeed, since the last data census in 2009, there has been an increase in forest and wooded areas of more than 2% in Switzerland, and more precisely of 16% and more when one exceeds 2000 m in altitude (Swiss Federal Statistical Office, 2013). The tree line in the mountains is determined by summer temperatures and their increase will therefore imply a raising of this line (Vittoz et al., 2013), which may explain why more vegetation development is observed at higher altitudes.

Pre-2000 winter storms like Viviane and Lothar have led to resilient forest regeneration with diverse species, improving drought and storm resistance compared to monocultures swept away. Elevated temperatures and longer growing seasons further boost plant productivity and greening (Vittoz et al., 2013). Gehrig-Fasel et al. (2007) findings on high elevation forest cover raise two explanations: land abandonment or global warming’s impact. Although initial observations align with land use, climate effects become prominent over the long term. Loran et al.‘s (2017) study shows increased forest area in Switzerland, especially in the Southern Alps due to a shift towards service economy, along with temperature and slope influences. This change leads to land abandonment and less intense management. Yet, recent frequent intense droughts (Piedallu et al., 2019) drive demand for high elevation pasture to escape lower elevation heat. Examining land use evolution amid climate change reveals its impact on greening trends.

Regarding the strong correlation with vines, a decrease in the use of chemical and mechanical weed control between grapevines, and the use of treatment products such as copper and sulfur that control diseases attacking the green leaf surface can potentially explain the increase in the greening trend in this land cover class (Noceto et al., 2020). However, the wine industry is complex and requires the consideration of many external factors such as the various treatments used on this plant class, not to mention the phenological characteristics of the plant itself and would require an analysis of this area by people with expertise in this field.

It is important to emphasize that to characterize NDVI changes and accurate identification resulting from changes in land use and land cover requires the use of dynamic land cover. Nevertheless, Thomas and Giuliani (2023) have reported that Switzerland has experienced subtle, spatially scattered, dynamic, and gradual changes in land cover patterns over the past four decades. They have further highlighted that the capacity to identify driving factors and understand the rate of these changes is constrained by the spatial resolution and temporal update frequency of the ArealStatistik dataset. Consequently, adopting the most recent dataset available, spanning from 2013 to 2018, can be considered as a reasonable choice.

5.3. Response of NDVI to climatic variables

Results of the correlation analyses between NDVI and temperature indicate that there is a relationship between these two variables in Switzerland and for most of its regions during the observed period. This indicates that the increase in NDVI values coincides with the increase in temperature, but there is some variability in terms of regions and seasons. Indeed, the impact of temperature on vegetation will depend on climate types, with generally a positive effect on wet vegetation growth and a negative effect on arid and semi-arid vegetation (Muradyan et al., 2019). For example, the Tibetan Plateau areas and the majority of the Northern Hemisphere show a stronger correlation between NDVI and temperature while plant productivity in the Sahel region and South Africa will tend to be correlated with rainfall (Zhu et al., 2016).

In regional analysis, lower areas like Jura and Plateau exhibit a stronger link between annual NDVI means and temperatures compared to other regions. High-altitude zones like Western and Eastern Central Alps lack temperature correlation due to snow and lower temperatures retarding vegetation growth, unlike low/mid-altitude regions (Filippa et al., 2019). Rumpf et al. (2022) confirms this, noting a 2°C rise affects snow cover less in colder spots than in warmer ones with less snowfall. Nonetheless, their study shows significant vegetation increase in the European Alps, coupled with declining snow cover. This combination could reduce surface albedo, amplifying global warming by increasing solar energy absorption. Thus, even at high altitudes, temperatures are expected to progressively impact alpine vegetation productivity.

Regional results of seasonal averages indicate that a stronger correlation between NDVI and temperature is observed in the winter. The areas of the Central Plateau and the northern slope of the Alps are most affected by this which could be explained by the previously mentioned arguments of an earlier onset of the vegetation growth period due to increasing temperatures for the lower elevation areas (Filippa et al., 2019), and for the alpine areas, a decrease in the snowpack and higher temperatures (Rumpf et al., 2022).

The pixel-by-pixel analysis also indicates a good positive correlation between annual and seasonal means of NDVI for Switzerland and temperatures, except for the summer season which is weaker. Indeed, an increase in temperature can be beneficial for vegetation greening, but up to a certain threshold (Hou et al., 2015). High temperatures coupled with water deficit can reverse the greening trend as mentioned in the study by Piedallu et al. (2019). Corona-Lozada et al. (2019) analyze greening in the European Alps between 2000 and 2020 and find that among the four major heat waves that impacted the Alps during the last 20 years, the only one that did not generate an increase in vegetation productivity but a decrease in vegetation productivity was coupled with a high-water deficit. This could explain why the correlation is weaker in the summer.

In terms of precipitation, regional analysis generally indicates a non-significant negative trend with annual NDVI means across biogeographic zones. Similar conclusions from studies on temperate regions find less correlation between precipitation and NDVI changes compared to temperature (Chu et al., 2019; Hou et al., 2015). Nevertheless, there is a significant negative relationship between spring seasonal NDVI means and precipitation. Lower altitude areas (Jura, Central Plateau) display a stronger negative trend, suggesting increased spring precipitation leads to lower NDVI values. This negative correlation might stem from heightened precipitation and cloud cover reducing solar radiation and temperature, thus hindering photosynthesis (Hou et al., 2015). While lag effects on NDVI-precipitation exist, our study utilizes seasonal and annual NDVI data, with time frames too long for significant impact. The delayed response varies based on factors like plant species, soil moisture, and local climate. Different ecosystems exhibit varied response times to precipitation, from days to weeks or months. In previous work, we tested time-lagged correlation between seasonal NDWI products and climatic variables, observing no significant effect (Poussin et al., 2021).

The pixel-by-pixel analysis again shows that there is a negative and significant correlation between spring seasonal means of NDVI and precipitation. However, unlike the Pearson analysis, the pixel-by-pixel analysis detects a stronger negative relationship between NDVI and precipitation in the winter season. As mentioned earlier, the National Centre for Climate Services (NCCS) observes an increase in heavy precipitation that is more intense and more frequent in Switzerland, with a particular increase in the amount of precipitation in the winter. The naturally lower temperatures during this season, coupled with increased precipitation, could explain the reason why a negative correlation with NDVI is observed during this period.

5.4. Limits and perspectives

Despite a significant Switzerland-wide greening trend, using NDVI and remote-sensing data has limitations. Factors like canopy architecture, leaf water content, soil background, and light can impact index accuracy (Huang et al., 2020, Redowan & Kanan, 2013; Yang et al., 2019). While remotely sensed data offer relevant analysis timescales, sensor errors in geolocation and calibration, atmospheric conditions, climatological data for correction, cloud cover, and regional topography can introduce challenges (Southworth & Muir, 2021). This study excluded multiple years due to unclear observations. Combining earth observation with field data would enhance comprehensive and focused analyses (Honeck et al., 2018).

Furthermore, even if the use of the Landsat covers a long period of time, Sentinel-2, which was launched by the European Space Agency (ESA) in 2015, allows for a better predictive capacity of the models compared to Landsat, due to increased spectral, temporal, and spatial resolutions (Astola et al., 2019).

The NDVI time-series analysis shows greening trends across land cover classes, but this insight is preliminary. Further investigation should examine class changes from 1984 to 2018, assessing their contribution to greening variations. Plants react diversely, identifying vulnerable and positively impacted species are valuable. Therefore, the interpretation of time series results requires caution. Future analysis can incorporate Swiss forest distribution categories (softwood, hardwood, etc.) from FSO (2020) to gauge their response to changing climate and impact on greenness. For instance, Hou et al. (2015) discovers coniferous forests respond more significantly to increasing greening trends.

The results of correlation analyses between NDVI and temperature and precipitation tell us that temperature is more strongly correlated than precipitation. The precipitation data used in this work only allow for an analysis from 2005 to 2018 and this may partially explain the reason why it is more difficult to determine what the relationship is between precipitation and plant greenness. In addition, the pattern of temperature change is clearer than that of precipitation, which is more subject to natural fluctuations, so it is more difficult to see its evolution (NCCS, 2018).

Other variables used by several studies (Novillo et al., 2019; Piedallu et al., 2019; Zhu et al., 2016) such as soil geomorphology, CO₂ concentrations, or quantification of anthropogenic activities could also complement the NDVI trend analysis in Switzerland and potentially explain some of the variations found at the regional level. It is therefore crucial to continue to monitor variations in the rate of vegetation greenness and to develop a precise methodology adapted to the data used and, if possible, to incorporate other explanatory variables to provide a global and accurate picture of the state of vegetation in Switzerland.

Greening trend analysis holds significance, particularly in the context of UN Sustainable Development Goal 15, as it observes crop phenology and land cover shifts (Lucas et al., 2022). Spatial data aids quantifying countries’ global warming efforts, aiding sustainable territorial development. Addressing land degradation is urgent for ecosystem balance, human well-being, and planetary health (Giuliani et al., 2020). Satellite data offers accessible, accurate soil insights, supporting climate change adaptation assessments. It provides clear environmental assessments guiding policy for a sustainable future. WMO states the last 8 years were the world’s warmest, and 2022 is Switzerland’s warmest, sunniest year since 1864. Drought’s impact on beech trees surpasses even Lothar’s damage. Comprehensive analyses across scales are crucial to inform national/regional actions against global warming’s consequences (Brun et al., 2020; Buras et al., 2020).

6. Conclusions

Switzerland with its mid-latitude position, continental climatic conditions and diversity of landscapes is particularly sensitive to climate change. Changing environmental conditions are expected to have important pressures on land ecosystems. Consequently, having means to monitor and assess these changes is crucial for providing ready-to-use information by decision-makers to manage land efficiently and effectively.

We showed that NDVI time-series extracted from 35 years of Landsat data can be used to monitor vegetation greening at different spatial and temporal scales. Even if the national annual mean is useful to obtain an overview of the general condition, disaggregation at lower spatial and temporal scales is needed to gain knowledge on conditions in specific land portions (e.g. biogeographical regions) and to distinguish the contrasting effects at different spatial and temporal scales. Our results suggest that there is a clear and significant greening trend at different scales (e.g. national and landscape) and for all the seasons. We identified 2010 as a possible critical point of changing land conditions showing a clear break in the NDVI trend. Across Switzerland, NDVI appears to be mostly controlled by temperature at elevations below the tree line and forested land cover classes are the most responsive to changing land conditions.

The use of large volume of satellite EO data to generate the necessary measurements to retrieve adequate information on trends is a prerequisite to efficiently monitor environmental variables. To that end, the use of an EO Data Cube (EODC) has been pivotal to effectively manage and process terabytes of ready to analyze data or provide information on land conditions and therefore essential to contribute to measure progresses towards policy outcomes. Ultimately, this work can be seen as an exploratory contribution providing nationwide spatially explicit information on the ecological state of vegetation to monitor land affected by changing environmental conditions using satellite Earth Observations.

Acknowledgements

The views expressed in the paper are those of the authors and do not necessarily reflect the views of the institutions they belong to.

Results of this publication are partly or fully relying on the Swiss Data Cube (http://www.swissdatacube.org), operated and maintained by UN Environment/GRID-Geneva, the University of Geneva, the University of Zurich and the Swiss Federal Institute for Forest, Snow and Landscape Research WSL.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Data availability statement

The NDVI time-series data that support the findings of this study are openly available at (annual mean) https://doi.org/10.26037/yareta:kpmscrogqbdhvjeuev2ydrzk7y and (seasonal mean) https://doi.org/10.26037/yareta:voy277qzczgzbgeiyrcldztuti

Land Use Statistics are available at: https://www.bfs.admin.ch/bfs/de/home/statistiken/raum-umwelt/erhebungen/area.html

National Boundaries and Biogeographical regions can be obtained at: https://www.swisstopo.admin.ch/en/geodata/landscape/boundaries3d.html

Temperature (TabsM and TabsY) and precipitation (RhiresM and RhiresY) dataset from the MeteoSwiss platform – Spatial Climate Analysis data are available for registered users at: https://www.meteoswiss.admin.ch/climate/the-climate-of-switzerland/spatial-climate-analyses.html

Additional information

Notes on contributors

Claire Obuchowicz

Claire Obuchowicz is an MSc student at the University of Geneva/Institute for Environmental Sciences on Biodiversity and Climate Change.

Charlotte Poussin

Charlotte Poussin is a PhD student at the University of Geneva/Institute for Environmental Sciences & GRID-Geneva. Her research topic is on assessing and monitoring the impacts of climate change in Switzerland using time series analyses in the Swiss Data Cube.

Gregory Giuliani

Gregory Giuliani is the Head of the Digital Earth Unit and Swiss Data Cube Project Leader at GRID-Geneva of the United Nations Environment Programme (UNEP) and a Senior Lecturer at the Institute for Environmental Sciences of the University of Geneva. He is a geologist and environmental scientist who specializes in Remote Sensing, Geographical Information Systems (GIS) and Spatial Data Infrastructures (SDI). He also works at GRID-Geneva of the United Nations Environment Programme (UNEP) since 2001, where he was previously the focal point for Spatial Data Infrastructure (SDI) and is currently the Head of the Digital Earth Unit. Dr. Giuliani’s research focuses on Land Change Science and how Earth observations can be used to monitor and assess environmental changes and support sustainable development.