Exploring the Use of Sentinel-2 Data to Monitor Heterogeneous E ﬀ ects of Contextual Drought and Heatwaves on Mediterranean Forests

: The use of satellite data to detect forest areas impacted by extreme events, such as droughts, heatwaves, or ﬁres is largely documented, however, the use of these data to identify the heterogeneity of the forests’ response to determine ﬁne scale spatially irregular damage is less explored. This paper evaluates the health status of forests in southern Italy a ﬀ ected by adverse climate conditions during the hot and dry summer of 2017, using Sentinel-2 images (10m) and in situ data. Our analysis shows that the post-event—NDVI (Normalized Di ﬀ erence Vegetation Index) decrease, observed in ﬁve experimental sites, well accounts for the heterogeneity of the local response to the climate event evaluated in situ through the Mannerucci and the Raunkiaer methods. As a result, Sentinel-2 data can be e ﬀ ectively integrated with biological information from ﬁeld surveys to introduce continuity in the estimation of climate change impacts even in very heterogeneous areas whose details could not be captured by lower resolution observations. This integration appears to be a successful strategy in the study of the relationships between the climate and forests from a dynamical perspective.


Introduction
The health of forests has been under increased scrutiny over the last few decades and climate change is a major contributor [1,2]. In particular, the stochasticity, frequency and magnitude of extreme weather events outside the range of normal values or rare within a statistical range at a particular place and time of year [3] have increased on a global scale [4][5][6], and their damaging impacts are expected to increase in the future [7]. These negative trends are mainly increasing in Europe [8], especially in the Mediterranean regions, where the climate has been warmer and drier [9,10] since the end of the last century. Mediterranean areas, commonly considered as rather resistant to drought event [11], have been recently characterized by more frequent and severe drought episodes which have already significantly impacted forests [12,13], triggering forest decay and mortality episodes [14][15][16][17]. For this reason, the beginning of the new millennium has been identified as a general period of Mediterranean forest decline related to drought [18]. Especially when drought is contextual to a prolonged period of abnormally hot weather (heatwave), the effects on forests can be destructive.
Under limited availability of soil water, the first plant response is to close stomata, which decreases the CO 2 inflow into the leaves and determines the consequent reduction of structural carbohydrates [19].
Thus, forest ecosystem responses to drought and their adaptation depend on the drought intensity and duration [20,21]. Moreover, some tree species are more sensitive to drought-triggered growth decline and damage [21][22][23][24] and sometimes the response of forests to drought can be exacerbated by local environmental and site conditions [24][25][26][27][28]. Although, several studies have highlighted a diverse vulnerability of vegetation to drought damage, a full knowledge of the response of Mediterranean forest ecosystems, characterized by environments with high heterogeneity, is still lacking. In short, there are many reasons why the local response of vegetation can be highly heterogeneous, and the assessment of such heterogeneity is a challenging task. In particular, local vegetation cover can be composed of different plant types (trees, shrubs, herbaceous species) which add further complexity to the total forest response. The most accurate evaluations are commonly carried out through data collected in situ, however, the high costs of field surveys do not allow the continuous collection of data in situ, especially for large areas.
Nowadays, remote sensing technologies are an essential tool to study all types of anthropic (agricultural lands, e.g., see [29,30]; seminatural (grassland, e.g., see [31])), and natural covers (shrublands, forests, prairies, e.g., [32][33][34]). Satellite remote sensing signals retrieve indirect indicators of forest health which are useful in the detection or monitoring of the stress phenomena due, for example, to drought, parasite attacks, soil pollution [35][36][37][38] however, they provide no direct measures of forest state variables such as height or diameter of trees.
The use of satellite-derived NDVI (Normalized Difference Vegetation Index) to identify areas impacted by climate events has been widely discussed in literature (e.g., [39]). A decrease in the NDVI is expected in areas where external stress causes a decline in photosynthetic activity and damages occur, e.g., after severe fires, the NDVI time series offers helps to identify impacted zones [40][41][42][43][44].
In southern Japan, for example, the NDVI evaluation before and after a typhoon and its negative correlation with electric conductivity indicated that the reduction of NDVI values is attributable to increased salinity from flooding, proving the NDVI as an indicator of salinity damage to vegetation [45]. Moreover, a variation in NDVI was found to be the most effective in predicting total floc soil phosphorus values, due to its capacity to capture small variations in chlorophyll associated with soil total phosphorus levels in periphyton, especially in aquatic/non-impacted areas [46].
The integration of in situ and remote sensing data can improve the identification of vegetation processes (i.e., senescence, photosynthesis, evapotranspiration). Temporal, thematic, and spatial aspects for operational implementation of satellite remote-sensing data in forestry are described well in Pause et al. [47].
The main drawback of these technologies is that remote indicators combine information from spatial units (pixels), which generally cover several square meters. Low resolution removes local details and the average information obtained from these observations could not reflect differences in the response of different areas. The response of forests to adverse meteorological events is generally irregular and the damage can be non-uniform even in small areas, particularly due to the heterogeneity of the vegetation cover. Although the link between the NDVI and the amount of photosynthetic active biomass is well established and many studies testify as to its benefits in the detection of areas damaged by external stresses [39], the ability of satellite observations to provide actual quantitative estimates of the damage at the local scale cannot be assumed directly as satellite sensors retrieve spatial averaged responses.
Very high spatial resolution (a few meters), images such as WorldView, RapidEye, etc., which have been widely used for forest research [48,49] are often very expensive and unsuitable for continuous monitoring.
The launch of the new Sentinel-2 satellite, one of the core missions of the Copernicus Earth Observation Programme of the European Union, provided opportunities for new research activities into the assessment and monitoring of forests. Sentinel 2, freely available, is characterized by a wide swath (290 km), a short revisit time (5 days with 2 satellites under cloud-free conditions which results in 2-3 days at mid-latitudes), 13 spectral bands (4 bands at 10 m of spatial resolution; 6 bands at 20 m of spatial resolution, and 3 bands at 60 m of spatial resolution) at a wider scale [50] that guarantee consistent time series minimizing any artifacts introduced by atmospheric variability and allowing for the systematic monitoring of the land surface [20]. The resolution (10 m) of this sensor could be sufficient to retain most of the details of the local response of forests to external stresses accounting both for damaged crowns and understory vegetation.
This explorative work analyzes the response of forest areas to negative and prolonged meteorological conditions, in five case studies investigating the relationships between the decrease in Sentinel-2 NDVI observed in the post-event period and in situ observations concerning both the vitality of single trees and the structural layers of the vegetation species which are present in the forest communities.
Our analysis aims at evaluating the consistency between the drop in NDVI and the average vitality of the plants in different vegetation layers (trees, shrubs, herbaceous plants). Plants are firstly classified according to the Raunkieur [51] biological spectrum (percent of species belonging to each life-form in flora) and their vitality is then classified according to the Mannerucci method [52].
The final goal is to evaluate the use of Sentinel 2 data to follow the evolution over time of Mediterranean forests under the impact of climate change by retaining information on the local heterogeneity.

Study Area
The study areas are located in Basilicata, southern Italy, 46.8% of its territory is mountainous, 45.2% is hilly and only 8% is flat ( Figure 1). In general, natural forested areas are found at high elevation, whereas man managed areas are mostly found in low-level areas. The five sites fall within areas negatively impacted by 2017 summer heatwave and represent the major forest vegetation types covering the inner part of Basilicata (mixed broadleaves stands with dominant Quercus and Carpinus genus species) within the sub-montane belt.
Gole del Marmo Platano and Costa dell'Angelo are located in the Marmo Platano-Melandro basin (municipalities of Vietri di Potenza and Savoia di Lucania) an enchanting landscape of forests and mountains of about 87,700 hectares located in the north-western part of the region.
Beech (Fagus sylvatica L.) inhabits at high altitudes, whilst oak species (Quercus cerris L. and Quercus pubescens Willd.), Carpinus orientalis Mill., and Fraxinus ornus L. characterize the sub-montane belt. Abies alba Mill. and other minor conifer species are also present in Costa dell'Angelo due to reforestation activities in the 1970s and 1980s.
The remaining sites are located in the Gallipoli Cognato-Piccole Dolomiti Lucane Park which covers over 27,000 hectares in the heart of Basilicata and includes the municipalities of Accettura, Calciano, and Oliveto Lucano in the province of Matera, and Castelmezzano and Pietrapertosa in the province of Potenza, all areas with a typically Mediterranean climate [53]. The park includes several tree species such as downy oak (Q. pubescens, Quercus ilex L.), oriental hornbeam (C. orientalis, black hornbeam (Ostrya carpinifolia Scop.), manna ash (F. ornus, maples (Acer campestre L., Acer opalus Mill., Acer monspessulanum L. ), hollies (Ilex aquifolium L.), and Turkey oak (Q. scerris ) which are the most widespread species.
According to the bioclimatic classification [54,55], the study areas belong to the Mediterranean pluviseasonal oceanic-weak semicontinental macro-bioclimate, lower subhumid ombrotype and upper meso-mediterranean termotype belt. However, local factors such as altitude and aspect determine the transitions towards the drier and warmer types.
The ecoregional classification includes the study areas into the southern Tyrrhenian section, Lucanian subsection [56]. Additional information on the field survey areas is reported in Table 1.

Materials and Methods
The approach we followed to evaluate the coherence between in-situ and Sentinel-2 satellite data is based on the use of NDVI data as proxy for the average photosynthetic activity at the pixel level (10 m for Sentinel-2 data). Three datasets (meteorological, remote, and in situ) were analyzed.

Meteorological Data
Meteorological data (total monthly precipitation; maximum, minimum, and average monthly temperatures) for 2016-2018 from 4 stations were provided by the Basilicata Italian Department for Civil Protection and were selected on the basis of their proximity to the study areas: 1.
Balvano station is closest to the Gole del Marmo Platano site; 2.
Sant'Angelo Le Fratte -Fiume Melandro station is closest to the Costa dell'Angelo site; 3.
Albano di Lucania station is closest to the Castelmezzano and Pietrapertosa sites; 4.
San Mauro Forte station is closest to the Palazzo site.

Satellite Data
The satellite data selected were based on the phenological phases of the forest stands analyzed: (i) early stage (April/May), (ii) middle stage (July), and (iii) late stage (September/October) of vegetation activity. We processed images acquired in 2017 and before (2016) and after (2018) the extreme climatic event.
Nine Sentinel-2 images with tile 33TWE were retrieved from the European Space Agency's Scientific Data Hub (https://scihub.copernicus.eu/) as standard Level-1C (Top of Atmosphere reflectance) ( Table 2). Data were selected on the basis of the low percentage of cloud cover in the investigated areas (completely free of clouds in the test-sites). Diffuse cloudiness in the study areas limited the availability of images in the corresponding dates for all the three years. Table 2. List of the processed Sentinel-2 images and corresponding cloud cover percentage. Satellite image processing was fundamentally based on the analysis of the NDVI, which is a proxy for plant photosynthetic activity [57,58] and a valuable tool at different spatial scales for the Mediterranean ecosystems investigated [53,[59][60][61]. Several studies have established the efficiency of the NDVI as an indicator for anomalies associated with climatic changes [39,62,63].
The responsiveness of the NDVI to precipitation changes is well established at various spatial and temporal scales [27].
Drought condition effects are described by a deviation of the NDVI from their typical values, usually estimated over extended periods. Based on the concept that healthy vegetation reflects more near-infrared (NIR) light and absorbs more RED light compared to other wavelengths, the NDVI is calculated using the following formula: where ρ NIR and ρ RED correspond to the reflectances at the NIR and RED wavelengths, respectively [57]. For each acquired image, we calculated the NDVI using bands 8 (near-infrared channel) and 4 (red channel) at 10 m spatial resolution. Successively, we calculated the average NDVI for each test-area by considering a buffer of 2 × 2 pixels centered on the GPS coordinates (see Table 1).
Finally, in order to capture potential differences in the phenological behavior of the vegetation we compared 2017 NDVI values with corresponding 2016 (pre-event) and 2018 (post-event) values.

In Situ Measurements
In situ measurements were collected between September and November 2017 (post-event period) in the study areas described in Section 2 in order to evaluate the impacts on forests caused by heat and drought in the summer of 2017.
The sampling sites were randomly chosen as a result of on-site observations in different forest typologies in the Lucanian Apennine [65] which showed heat/drought stress.
The major problems were found in deciduous forests located in an altitude range of 500-800 m a.s.l. They were identified by the decline symptoms of the tree and shrub layers highlighted mainly by shoot dieback, withering and early leaf loss compared to the normal seasonal pattern. In the Central Basilicata territory, the dormancy of deciduous forests occurs from November with a drop in temperature.
In September 2017 the deciduous forests already fully manifested this condition, determined by the simultaneous occurrence of prolonged drought and particularly high temperatures during the spring-summer period.

The Mannerucci Method
For each site, a stand area of about 200 m 2 was surveyed, identified by two geographic coordinates (center of the area). For each site, a field-form was completed with the description of: (a) number of trees; (b) tree species; (c) diameter of trees measured at 1.30 m height using a tree calliper; (d) height of the trees measured using a vertex ipsometer; (e) dendrotype to identify the origin of the trees (agamic or gamic); (f) tree vitality.
These measurements are also useful to perform the calculation of the basal area and volume, in order to evaluate the biomass present and to monitor, with subsequent measurements, variations in growth and therefore wood biomass.
The evaluation of the vitality of the vegetation was performed using the Mannerucci method [52] based on visual and qualitative analysis, assigning a vitality class, i.e., a numeric value in the range 1 (healthy plant) to 6 (dead plants), to each observed plant (Figure 2).

Floristic-Vegetational Sampling and Plant Functional Type Classification
The analysis of the biological spectra of the floristic-vegetational surveys and of the vulnerability classes of the species present provides a summary of the consequences of the climatic event on the structure and functional recovery capacity of the forest communities examined.
Available data for vegetational analysis (see Supplementary Materials, Table S1) were collected according to the Braun-Blanquet method [66]. Frequency and vulnerability classes were detected for each species according to the Mannerucci method (see Section 3.3.1). The average of the vulnerability class of the species observed on the site was reported to quantify the structural and functional damage to the forest community. In order to highlight and evaluate the synecological damage and vulnerability, the structural composition of the plants communities is indicated and the functional type based on the life forms frequency is taken into consideration for each species according to the Raunkiaer method [51]. The plant functional types are considered a possible framework for predicting ecosystem response to human-induced changes such as global warming [67].
The life forms method was based on the observation of the different protection of the vegetative buds as a strategy of adaptation to environmental stress (drought or cold) and allows for the evaluation of each plant and in the plant community of the potential resistance to periods of adverse climatic conditions [68,69].
Plant life-form usually displays an obvious relationship between the plants and the environmental factors [70][71][72].
In Raunkiaer's system [51], the more pronounced the unfavorable season the more protected the renewing buds. The effectiveness of the protection depends on the position of the perennial buds with respect to the surface of the soil during the adverse seasons (e.g., perules in the aerial buds of woody plants such as fanerophytes and chamaephyta; underground stems such as rhizomes, tubers and bulbs in geophytes; seed in therophytes).
For each site, the analysis of the functional types was performed by calculating the biological spectrum both on the percentage of the life forms and on their average frequency.
Life form for each species is derived from Pignatti [73]. Vascular plants nomenclature followed Bartolucci et al. [74].  The similarity between thermal curves measured in different study sites reflects the non local character of this large scale climatic event that can be classified as a heatwave. This warm anomaly reached its peak in August 2017, with monthly maximum values of about 40 • C in all sites examined, more than 5 • C higher than the values measured in other years. If we consider that the elevation of the sites is higher than 500 m, this temperature appears to be really high [60]. The positive anomaly occurs between early spring and early autumn which, on the contrary, show the lowest values observed in these periods.

Analysis of Meteorological Data
As for rainfalls, 2017 was drier than 2016 and 2018 for most of the year and in particular during the heatwave (Figure 4). Severe drops in precipitation from 2016 to 2017 were observed (90% in the San Mauro Forte meteorological station). Altogether, the contextual occurrence of temperature and rainfall anomalies in summer, which is already the warmest and driest period of the year, defines the conditions observed between June and August 2017 as a climate event potentially harmful to forest health.

Analysis of Sentinel-2 NDVI
The sequence of the NDVI values estimated for the surveyed areas ( Figure 5) summarizes the phenological period April to October observed by Sentinel-2 during 2016, 2017, and 2018. This sequence provides a proxy for the photosynthetic activity before, during, and after the climate event. The most striking feature which emerges from the comparison of the plots in Figure 5 is the low level of photosynthetic activity in October 2017, after the event, which is the lowest estimated in this month over the three years. According to the NDVI values, at the end of the meteorological critical period a decrease in photosynthetic activity is detected in all sites.
The 2016 spring values are particularly high, as largely testified in literature (e.g., [39]) while the values in August are similar for all sites, except for Castelmezzano where the 2017 spring values were already low and the decrease in vegetation activity in October leads to an NDVI value that is typical of less densely vegetated areas.
Nevertheless, the magnitude of the decrease is site-specific, in other words each site experienced different negative impacts, according to satellite data.
Finally, we can hypothesize vegetation recovery as Figure 5 shows similar NDVI trends for the 2016 and 2018 summer-autumn periods for all sites and the difference between NDVI values for the 2018 autumn/summer period is around zero (Table 3).

Analysis of In Situ Data
The analysis of in situ data aims to provide a quantitative indicator of the damage to vegetation produced by drought and a heatwave for comparison with the decrease of NDVI in October.
Trees form the structure of the forest with an evergreen canopy tree layer (height about 15 m; coverage 70%) with lower deciduous broadleaf trees and shrub layers, which showed the greatest damage (Figure 7). Among these plants, the species in the most critical conditions were C. orientalis (94% of this species belonging to class 4 of the Mannerucci method) and especially A. alba, whose young trees and seedlings were dead. The remaining plants were class 1 (40%), 2 (3%), 3 (5%), and 6 (9%) (Figure 6).
No damage was observed on the woody plants in the undergrowth (Pn, Plian) or on the herbaceous plants (Geophytes and Emicriptophytes; see Supplementary Materials, Table S1-Relevè number 2).

Palazzo Site
This site is currently used as pasture and the soil showedmarked erosion. The vegetation cover was discontinuous: the few trees (about 50% cover) and the scattered bushes of Q. pubescens, Spartium junceum L. and C. monogyna suffered significant damaged from the climatic event (Vulnerability class > 3) (Figures 6 and 7). Two different tree species are present in this area: 27 are Q. pubescens specimen, only 1 F. ornus specimen. These species are classified as 3 (50%) and 4 (50%) according to the Mannerucci method ( Figure 6).
Stress levels were intense between July and August 2017 (Figures 3 and 4) and the shrubs and herbaceous species with epigean vegetative buds (Hemicriptophytes) were mainly affected (Figure 7 and Supplementary Materials, Table S1-Relevè number 3).

Comparative Analysis
The comparative analysis between satellite and in situ data was carried out considering the percentage of damaged plants (from class 2/3 to 6 evaluated for each vegetation layer). Figure 8 shows the distribution by class of the tree layers independently of the species. Costa dell'Angelo is the site with less average damage, whereas the Palazzo site was the most vulnerable site with all specimens belonging to class 3 and 4.
The percentage of damaged trees (classes 2-6) appears to be proportional to the decrease in the NDVI between July and October 2017 (R 2 = 0.98; Figure 9). As such, the decrease in NDVI accounts for the average impact of the climate event on the forest, Costa dell'Angelo is paradigmatic. Despite the presence of many damaged plants and some dead trees (class 6), the high percentage of unaffected trees (40%) limits the decrease of the NDVI in autumn 2017. Similarly, in the Castelmezzano site, the presence of about 20% of unaffected plants produces a decrease in the NDVI which is about twice that observed in Costa dell'Angelo. Finally, for the other three sites, where vulnerability levels are mainly distributed between classes 2, 3, and 4, the decrease of the NDVI is doubled.   Table 3).
The analysis of the remaining vegetation layers is consistent with Figure 9 and indeed completes the explanation of the results. In fact, in Costa dell'Angelo, the healthy condition of the dominant evergreen trees and of the undergrowth and herbaceous plants explains the modest decrease of the NDVI in autumn 2017.
The NDVI recorded very low values in Castelmezzano ( Figure 5) as the soil was almost devoid of vegetation due to the widespread death of the leaf systems due to the severe and prolonged drought during September and October 2017. Moreover, in 2018, the spring rainfalls and the short period of summer drought favored the vegetative growth of less damaged plants such as Olea europea L., P. terebinthus, and F. Ornus, as highlighted by the NDVI values thereby favoring the full 2018 recovery. The decrease in the NDVI in Gole del Marmo-Platano can be mainly ascribed to tree and shrub layers, whereas the floristic richness of herbaceous layers has allowed the community to renew itself in a short time and to return in 2018 to slightly higher NDVI values than those in 2016 ( Figure 5; Table 3).
The NDVI in Palazzo showed a remarkable recovery in 2018 due to the vegetative growth of the little damaged woody species (e.g., F.ornus,) and to the heliophilous herbaceous plants with hypogean vegetative buds (Geophytes; e.g., Asphodelus microcarpus Salzm. & Viv. ) not eaten by cattle ( Figure 5; Table 3).
Finally, the climatic trend over the years in Pietrapertosa is similar to Castelmezzano (Figures 3 and 4), however, the microclimate conditions and the remarkable floristic richness of the undergrowth facilitated efficient and effective recovery of the forest vegetation. This was confirmed by the NDVI values which had already returned to the level close to that of 2016 by autumn 2018, before the climatic stress ( Figure 5).

Conclusions
This study investigated the consistency between in-situ data and Sentinel-2 imagery to assess the health status of forests affected by adverse climate conditions (rainfall deficit and high temperatures) in 2017, classified as one of the driest and warmest years since 1800 (https://www.ncdc.noaa.gov/sotc/ national/201713;http://www.isac.cnr.it/climstor/climate_news.html).
The use of Sentinel-2 for a continuous monitoring which accounts for local damage heterogeneity was evaluated, beyond the mere identification of affected areas which is usually investigated in the current literature. Such a continuous monitoring is fundamental to follow the health status of Mediterranean forests, which are among the most vulnerable to climate risk, not only because of the expected increase of extreme events but also for possible changes in the rainfall regimes [75] that can alter the ecosystem phenological cycles. Our analyses focused on in situ and remote data separately to make independent assessments of the forest status. Overall, field surveys showed different local scale response of forests in the different sites studied. The worst case was Gole del Marmo Platano, where a high percentage of trees and shrubs was damaged. On the other hand, Costa dell'Angelo was characterized by good post-event conditions of the dominant evergreen trees and of the undergrowth and herbaceous plants. This diversity is ascribable both to different vegetation and to different magnitude of the extreme meteorological event. The 2017 spring-summer period in Gole del Marmo Platano was the driest in the area, whereas the warmest conditions were observed in the Palazzo site where, in fact, all the tree crowns were damaged.
Although drought and heatwaves are large-scale events, local differences due, e.g., to elevation and aspect can mitigate or exacerbate their impact on vegetation. In spite of these differences, the overall response of the Mediterranean forests analyzed in this work seems to be rather efficient as only a low percentage of trees, only in one site, died and the recovery was rather fast in all the area investigated.
Our results showed a good consistency between the different assessments including a separate evaluation of the health status of each vegetation layer (trees, shrubs, and herbaceous plants). The role of tree cover in determining the drop in the NDVI was striking as it is proportional to the number of damaged trees estimated using the Mannerucci method. The consistency of the contribution of the other vegetation layers was also impressive and helps to interpret both the response of the forest ecosystems and their recovery favored by the understory vegetation.
These findings suggest that the Sentinel-2 NDVI is a reliable index for the assessment of damage to vegetation at the pixel level (10 m) and shows its suitability for reporting on different levels of damage; this is an encouraging result in view of the climatic changes witnessed in recent years.
Moreover, the availability of a free and robust Sentinel-2 database and the frequent revisit time of its constellation makes these data especially suited to the continuous monitoring of forests after extreme events as they are able to provide information on large forest areas while preserving local heterogeneity. This information is useful to study the complex forest climate dynamics and to support decision makers in the management of Mediterranean forests which are particularly vulnerable to climate change.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-445X/9/9/325/s1, Figure S1: Monthly rainfalls values estimated in the ground meteorological stations from January 2016 to October 2018, Figure S2: Monthly mean temperature values estimated in the ground meteorological stations from January 2016 to October 2018, Table S1: Floristic-vegetational surveys with vulnerability class for each species.