On the Impacts of Historical and Future Climate Changes to the Sustainability of the Main Sardinian Forests

: The Mediterranean Basin is affected by climate changes that may have negative effects on forests. This study aimed to evaluate the ability of 17 forests located in the Island of Sardinia to resist or adapt to the past and future climate. Sardinia is experiencing a decreasing anthropic pressure on forests, but drought-triggered dieback in trees was recently observed and conﬁrmed by the analysis of 20 years of satellite tree-cover data (MOD44B). Signiﬁcant negative trends in yearly tree cover have affected the broad-leaved vegetation, while signiﬁcative positive trends were found in the bushy sclerophyllous vegetation. Vegetation behavior resulted in being related to the mean annual precipitation (MAP); for MAP < 700 mm, we found a decline in the tall broad-leaved stands and an increase in the short ones, and the opposite was found for bushy sclerophyllous vegetations. In forests with MAP > 700 mm, both stands are stable, regardless of the growing trends in the vapor-pressure deﬁcit (VPD) and temperature. No signiﬁcative correlation between bushy sclerophyllous tree cover and the climate drivers was found, while broad-leaved tree cover is positively related to MAP1990–2019 and negatively related to the growing annual VPD. We modeled those relationships, and then we used them to coarsely predict the effects of twelve future scenarios (derived from HADGEM2-AO (CMIP5) and HadGEM3-GC31-LL (CMIP6) models) on forest tree covers. All scenarios show an annual VPD increase, and the higher its increase, the higher the trees-cover loss. The future changes in precipitation were contrasting. SC6, in line with past precipitation trends, predicts a further drop in the mean annual precipitation ( − 7.6%), which would correspond to an average 2.1-times-greater reduction in the tree cover ( − 16.09%). The future changes in precipitation for CMIP6 scenarios agree on a precipitation reduction in the range of − 3.4% (SC7) to − 14.29% (S12). However, although the reduction in precipitation predicted in SC12 is almost double that predicted in SC6, the consequent average reduction in TC is comparable and stands at − 16%. On the contrary, SC2 predicts a turnaround with an abrupt increase of precipitation (+21.5%) in the upcoming years, with a reduction in the number of forests in water-limited areas and an increase in the percentage of tree cover in almost all forests.


Introduction
In recent years, a growing awareness of climate changes and related consequences has developed, and many countries are promoting mitigation and adaptation actions. However, despite this, climate changes are still taking place at an increasing rate across the globe [1,2] In the last few decades, Mediterranean countries, for example, have experienced increasing drought frequency and severity [3][4][5][6][7][8]. These findings are in line with the comparative analysis by Giorgi (2006) [9] that places the Mediterranean Basin among the most responsive regions to global climate change [10]. Future climate predictions are not encouraging, since they foresee a further air-temperature increase, particularly in the summer season, and a precipitation decrease, especially in winter [11,12].
Differences in temperature or precipitation patterns have produced a direct impact on the distribution of the terrestrial ecoregions and biomes [13]. The Mediterranean topic are focused on Africa, America, and Australia [38,39,[41][42][43], while analogous research in Mediterranean regions, with a similar precipitation gradient, has not yet been seen.
Trees do not all respond in the same way to changing climates and environments, and some species are naturally resistant and resilient to climate-influenced disturbances, whereas others are sensitive and vulnerable [44][45][46]. The importance of past climate-change conditions in forest ecosystems is still under investigation [28,44,45,47], and accumulating evidence indicates that warmer temperatures and the increasing vapor-pressure deficit (VPD) are raising the mortality risks for forests during drought [48] and are negatively impacting vegetation growth [49]. Recently, Gazol and Camarero (2022) analyzed the compound effect of high summer VPD and low soil moisture on defoliation and tree mortality across European forests, highlighting that the steady increase in VPD is raising the risk of increased drought-triggered mortality rates of forests in Southern and Eastern Europe particularly.
Nowadays, long time series of optical remote-sensor observations that allow the monitoring of tree-cover dynamics are starting to be available [50,51], although not as extended as the hydrological database, and can be an appropriate support for climate-change studies. Although the observations of the Advanced Very High Resolution Radiometer (AVHRR) optical sensor started early (from 1984), they are at a coarse spatial resolution (1.1 km) that is inappropriate for heterogeneous Mediterranean ecosystems. Instead, the observations of the MODIS sensor on the Aqua and Terra satellite platforms are robust and still extended enough (from 2000), and with a more appropriate spatial resolution (250 m) [9,52,53]. In particular, a MODIS product, the MOD44B Version 6 Vegetation Continuous Fields-Collection 6 provides robust estimates of yearly tree-cover dynamics [43,[54][55][56][57][58][59] In the last decades, Sardinia has experienced long droughts that impacted the forest cover and tree density in several forests (see Supplementary Figures S2 and S3), showing alarming vulnerable conditions and raising the question of the actual and future hydrologic sustainability of the Sardinian forest ecosystems. Given this, and considering that it is well documented that MAP has a downward trend [60,61], while air temperature has a positive trend [27], in the Mediterranean Basin, we investigate the effects of climate change on the hydrological sustainability of the Sardinian forests. The specific objectives are to (1) quantify past trends in the amount, type, and distribution of tree cover in the main Sardinian forests; (2) evaluate the historical trends of the precipitation, air temperature, and vapor-pressure deficit in the Sardinian forests; and (3) evaluate the impacts of past and also future changing climate conditions on the Sardinian forest tree cover. The outcome of this study is expected to identify the key climate factors controlling tree growth and mortality in the main Sardinian forests, as well as the responses of different tree species, and therefore provide suggestions for forest management practice.

Study Area
The island of Sardinia is a region of Italy located between 38 • 510 N and 41 • 150 N latitude and 8 • 80 E and 9 • 500 E longitude, with a surface of about 24,100 km 2 , and is second in size only to Sicily among the islands of the Mediterranean Basin ( Figure 1). Sardinia has an ancient geoformation and is not earthquake prone. Its granitic basement rocks date, in fact, from the Palaeozoic and the pre-Palaeozoic Era, and they are overlain by discontinuous sedimentary and volcanic units deposited in marine or continental environments between the Permo-Carboniferous and the Quaternary [23,62]. Due to long erosion processes, the island's highlands, formed of granite, schist, trachyte, basalt, sandstone, and dolomite limestone, range between 300 m and 1000 m. The highest peak is Punta La Marmora (1834 m), part of the Gennargentu mountain range in the center of the island. Although the island's ranges and plateau are separated by wide alluvial valleys and flatlands (e.g., Campidano, Nurra, etc.), and although the island has an average elevation of 334 m, the slopes are so steep that most of the island is considered mountainous [62]. The climate is maritime Mediterranean, with mild, rainy winters contrasting with hot, dry summers. Two macrobioclimates (Mediterranean pluviseasonal oceanic and Temperate oceanic), one macrobioclimatic variant (Sub-Mediterranean), and four classes of continentality (from weak semihyperoceanic to weak semicontinental) can be found in Sardinia [63].
Although the island's ranges and plateau are separated by wide alluvial valleys and flatlands (e.g., Campidano, Nurra, etc.), and although the island has an average elevation of 334 m, the slopes are so steep that most of the island is considered mountainous [62]. The climate is maritime Mediterranean, with mild, rainy winters contrasting with hot, dry summers. Two macrobioclimates (Mediterranean pluviseasonal oceanic and Temperate oceanic), one macrobioclimatic variant (Sub-Mediterranean), and four classes of continentality (from weak semihyperoceanic to weak semicontinental) can be found in Sardinia [63].  Broad-leaved forest (BLF) Vegetation formation composed principally of trees, including shrub and bush understory, where broad-leaved species predominate. This class is applicable for mature forests of natural or anthropogenic origin. In Sardinia, it mainly includes evergreen broad-leaved woodlands composed of sclerophyllous trees (mainly Quercus Ilex, Quercus suber, Quercus Rotondifolia). Broad-leaved area must cover at least 75% of the forest tree component.

3.1.2.
Coniferous forest (CF) Woodlands consisting mainly of trees, but also of bushes and shrubs, in which coniferous forest species dominate. This class is applicable for mature coniferous (needle-leaved) forests of natural or anthropogenic origin pure or mixed stands of fir (Abies), pine (Pinus), spruce (Picea), cedar (Cedrus), cypress (Cupressus), juniper (Juniperus), yew (Taxus), and species, among others. The coniferous surface must constitute at least 75% of the forest tree component.

Tree Cover of the Main Sardinian Forests
Sardinia's wildlands were selected based on the Corine Land Cover (CLC) classification produced within the framework of the Copernicus Land Monitoring Service, referring to the land cover/land use status of years 2018 and 2000 (the source dataset is described and accessible from https://land.copernicus.eu/pan-european/corine-land-cover (accessed on 1 September 2022)). This database relies on standard methodology and nomenclature with the following base parameters: 44 classes in the hierarchical 3-level CLC nomenclature, a minimum mapping unit for status layers of 25  land-cover types that have not been significantly modified by human activity, and they exclude urbanized and agriculture lands. Table 1 describes the main properties of the three  land-cover classes. Tree cover was estimated from satellite information, using a product of the widely used optical MODIS sensor, the MOD44B Version 6 Vegetation Continuous Fields-Collection 6 [42,[53][54][55][56][57][58]. MOD44B provides a yearly global representation of the percentage of three ground-cover components at a 250 m spatial resolution: tree-canopy coverage (TC), nontree vegetation (NTV), and non-vegetated (bare) coverage (BS). Moreover, trees are defined as woody vegetation higher than 5 m. Data for the 2000-2019 period were used. CLC information for the years 2000 and 2018 was resampled at the 250 m spatial resolution of the MODIS product. In this study, we did not consider any pixels falling into one of the following categories: (i) nonhomogeneous land cover, (ii) land-cover change between CLC 2000 and CLC 2018, and (iii) affected by forest fires between 2000 and 2020. We selected 17 forests (Figure 1) that have been stable since 1932 [64] and included mainly tall trees with a sparse understory layer of Q. ilex subsp. Ilex, dense forests with a garrigue or maquis formation, and open woodlands of Q. ilex subsp. ballota (frequently described as savannah-like ecosystems, with low canopy cover, accompanied by grass or shrub understory). Table 2 summarizes the main properties of the forests under investigation. The forest areas range from 11.03 km 2 of Su Sartu (F16) to the 447.5 km 2 of the Monti del Gennargentu (F3), which is also the forest at the highest average altitude (1147 m asl).

Meteorologic Data
The rain-gauge network over Sardinia consists of 414 stations [23], and data are available for almost one century (1922-2019). We further selected 179 rain gages with more than 80 complete years of data. In more detail, 51.4% of them have ≥80 complete years of data, 86% have ≥70 complete years of data, 91.62% have ≥60 complete years of data, and 100% have more than 59 years of data.
Daily maximum and minimum air temperatures at 112 stations are available for the period 1981-2010 (39 years), but only a few of them have remained available in the last ten years, for which tree-cover data are instead available. To overcome the reduced spatial availability and density of the observed temperature data, and to extend the data set to 2019, the ERA5 monthly reanalysis dataset for the 1980-2019 period for both air temperature and dew point was used. We successfully tested the accuracy of the reanalysis dataset (correlation coefficient > 0.62, RMSE < 1.14 • C) in 6 of the observing Sardinian stations, with almost-complete datasets (Carbonia C.ra Flumentepido, Tertenia, Senorbi, Macomer, Villa Verde, Alghero, Italy). VPD was calculated by using the ERA5 monthly near-surface air temperature (T) and dew point (Td), using the following equation [65]: where c 1 = 0.611 kPa, c 2 = 17.5, and c 3 = 240.978 • C. Td and T are in • C, and the resulting VPD is in kPa.

Statistical Analysis
The mean air temperature, VPD, and precipitation data were analyzed at the annual (T y , VPD y , and P y ) and seasonal scales. We used the Sardinian hydrologic year that begins on October 1 (when runoff usually resumes) and ends on September 30 of the following year (when the dry season usually ends) [21]. We distinguished the seasons as (i) 'extended winter period' (December through March, T ew , VPD ew , and P ew ); (ii) 'shortened spring period' (April and May, T ss , VPD ss , and P ss ); (iii) 'extended summer period' (June through September, T es , VPD es , and P es ); and (iv) the 'shortened fall period' (October and November, T sf , VPD sf , and P sf ). Note that the use of the extended seasons is common in climate-change studies [21,23,66].
Trends were estimated by using the Mann-Kendall non-parametric test [67], while the slope of the linear trend was estimated with the Theil-Sen method [68,69], both of which are widely used in hydrologic applications [23,[70][71][72]. To reduce the impact of autocorrelation on the statistical significance of the Mann-Kendall test and on the Theil-Sen method, different prewhitening methods are proposed by the literature [73][74][75][76]. The algorithm, proposed by Collaud et al. (2020), was used here.
The Pettitt test, a method widely used in hydrology and climate analysis [23,77,78], was used to detect the presence of a change point in the time series of the average annual precipitation over Sardinia. The spatial distribution of the MAP values observed in the raingauge stations over the Sardinia region was estimated through an ordinary kriging model (linear semi-variogram, 250 m resolution, and 16 neighbors). Finally, the Pearson correlation coefficient, ρ, was used for estimating the linear correlation between two datasets.

Future Projections
For future climate scenarios of precipitation and air temperature, we considered the HadGEM2-AO climate model [79] predictions from the Coupled Model Intercomparison Project Phase 5 (CMIP5), and a model coming from the following generation of HadGEM configuration, the HadGEM3-GC31-LL [80,81], from the Phase 6 (CMIP6). Montaldo and Oren (2018) showed that the HadGEM2-AO model is the one that best represents the past variations in precipitation and air temperature in Sardinia. Future climate scenarios up to 2090 were used to predict the effect of future climate change on tree cover in the forests. Following Montaldo and Oren (2018), we considered two main representative concentration pathways (RCPs) of the HadGEM2-AO model (atmospheric horizontal resolution of 1.875 • × 1.25 • ), i.e., RCP 45 and RCP 85. In addition, two Shared Socioeconomic Pathways (SSPs) of the HadGEM3-GC31-LL model (atmospheric horizontal resolution of 1.875 • × 1.25 • ) were considered, i.e., SSP 245 and SSP 585). The multivariate bias-correction algorithm [82], which transfers all statistical characteristics of the land-observed meteorological data to the corresponding multivariate distribution of variables simulated by the GCMs, was applied for the downscaling and bias correcting of both HadGEM2-AO and HadGEM3-GC31-LL for Sardinian forests.

Long-Term Trends in Precipitation
Rainfall is not abundant, and on the coast, the MAP typically ranges from 415 mm (Cagliari station at 7 m asl) to 934 mm (Baunei C.ra station at 480 m asl) mm per year ( Figure 2a) and follows the Mediterranean pattern, so that rainfall is concentrated in the autumn (mean value of 180 mm) and winter (mean value of 305 mm), gradually decreases during the spring (mean value of 94 mm), and hits a low in the summer (mean value of 74 mm). In inland areas, rainfall locally exceeds 700 mm per year in hilly areas and 1000 mm in the mountains (Figure 2a). Considering the full rainfall dataset (1922-2019), the mean of the island's average annual rainfall is 729.08 mm (SD of 143.8 mm), and it varied between 469.5 mm (year 1999) and 1062 mm (year 1933). Nine of the top ten driest years that took place in this century occurred after 1980. Pettitt's test detects a changing point in the average annual precipitation over the island in 1980. A total of 51% of the observing stations showed a significant negative trend in annual precipitation, but the spatial variability is high (Figure 3a). Limiting the rainfall dataset to the stations within close proximity of the forests analyzed, it emerges that the annual precipitation decreases significantly in almost all forests, except for F3, F4, F11, and F17, and that the β Py < −3 mm per decade in the rainfall stations is close to F5 (−3.6 mm/y at Is Cannoneris Caserma), F7 (−3.72 mm/y at Monti Mannu Caserma), F8 (−3.1 mm/y at Marrubiu, and −3.31 at Baradili), and F12 and F13 (−3.12 mm/y at Ala dei Sardi). Seasonally, precipitation decreases were largest in the winter (Figure 3b Even in the extended winter average precipitation, we can observe an important spatial variability, with β Pew < −25 mm per decades in the rainfall stations within close proximity of F5 (−3.03 mm/y at Is Cannoneris Caserma), F7 (−2.70 mm/y at Monti Mannu Caserma), F8 (−2.43 mm/y at Marrubiu), and F13 (−2.41 mm/y at Ala dei Sardi). Variations in the Theil-Sen slope estimator of precipitation during the short spring (β Pss ) and the extended summer (β Pes ) across the island are in the range ± 0.65 mm/y, and ±1.00 mm/y for the short fall (β Psf ).
of the forests analyzed, it emerges that the annual precipitation decreases significantly in almost all forests, except for F3, F4, F11, and F17, and that the βPy < −3 mm per decade in the rainfall stations is close to F5 (−3.6 mm/y at Is Cannoneris Caserma), F7 (−3.72 mm/y at Monti Mannu Caserma), F8 (−3.1 mm/y at Marrubiu, and −3.31 at Baradili), and F12 and F13 (−3.12 mm/y at Ala dei Sardi). Seasonally, precipitation decreases were largest in the winter (Figure 3b), while little and no significant changes are observed for the other seasons. The island's average winter precipitation (Pew) was 344.47 mm (standard deviation, SD, of 98.18 mm), and it varied between 127.66 mm (year 1999) and 546.14 mm (year 1927). Even in the extended winter average precipitation, we can observe an important spatial variability, with βPew < −25 mm per decades in the rainfall stations within close proximity of F5 (−3.03 mm/y at Is Cannoneris Caserma), F7 (−2.70 mm/y at Monti Mannu Caserma), F8 (−2.43 mm/y at Marrubiu), and F13 (−2.41 mm/y at Ala dei Sardi). Variations in the Theil-Sen slope estimator of precipitation during the short spring (βPss) and the extended summer (βPes) across the island are in the range ± 0.65 mm/y, and ± 1.00 mm/y for the short fall (βPsf).     Figure  2c), it emerges that MAP is decreasing drastically in many parts of the island. In the period 1922-1979, areas with a MAP value of less than 700 mm occupied 26% of the island's sur-  Figure 2c), it emerges that MAP is decreasing drastically in many parts of the island. In the period 1922-1979, areas with a MAP value of less than 700 mm occupied 26% of the island's surface. This value rises to 62% in the period 1980-2019 and even includes some areas such as those in the southern part of the island and the Campidano Plain, with a MAP lower than 450 mm. At the 17 selected forests, we observed that MAP is decreasing everywhere by an average of −122.86 mm. F13 and F14 are the forests that recorded the highest (−243 mm) and the lowest (−84 mm) reduction values, respectively.

Trends in Tree Cover
A total of 155,518 MODIS pixels were sampled across the study area for monitoring Sardinian forests' tree cover (Figure 1), corresponding to about 28% of the entire island.  Table 3 provides the average Mann-Kendall τ and the Theil-Sen slope β values of the tree-canopy coverage, non-tree vegetation, and bare soil within each of the 17 forests.
Broad-leaved forests with significant negative trends in the Mann-Kendall τ values for tree-canopy coverage are clustered mainly in the forested areas covered by this study and particularly in F5, F8, F10, and F15. Ten of the thirteen forests that contain broad-leaved vegetation show a Theil-Sen slope β ranging between −0.30% TC per year and +0.12% TC per year. Forest F5, which is the second largest forest on the island, is the one with the lowest value of Mann-Kendall τ for the tree cover canopy (τTC), equal to −0.3 (Theil-Sen slope β = −0.38%TC per year). In this forest, as in F15 the decrease in τTC contrasts with an increase in the τ of the vegetation of less than 5 m (τNTV BLF ), while, for example, in F8 and F9, there was a growth in the value of bare soil (τBS). The Mann-Kendall τ of the tree-cover canopy of forests F3 and F4, which are in high elevation mountains in the central east of Sardinia, is, instead, not significant and moderately positive. The spatial distribution of trends in bushy sclerophyllous vegetation (Figure 4d-f) highlights an overall increase in treecover coverage, and a relative reduction in non-tree vegetation and non-vegetated areas. Bushy sclerophyllous vegetation in the analyzed forests shows a sensitive negative trend in τTC only in F3, and, on the contrary, 15 sites show a positive trend in τTC BSV , with βTC BSV ranging from +0.01 to +0.5% per year. shows a sensitive negative trend in τTC only in F3, and, on the contrary, 15 sites show a positive trend in τTCBSV, with βTCBSV ranging from +0.01 to +0.5% per year.  The relationship between tree canopy coverage and MAP 1922-2019 for BLF and BSV was investigated year by year. As an example, Figure 5 shows the relationship of MAP 1922-2019 versus tree-canopy coverage for the year 2008 for both BLF and BSV. We evaluated the breakpoint (MAP cp , identified by the regression analysis in Figure 5) at which the maximum tree cover is attained. The change in tree-canopy coverage as a function of MAP has been determined by using a 99th quantile piecewise linear regression with the identification of the breakpoint (the rainfall at which the maximum tree cover is attained), and of the maximum upper bound on tree cover (TC max ). We showed the scatter plot for the year 2008 in Figure  , with a SD of 3.53%, for bushy-sclerophyllous pixels ( Figure 6). There is a significant negative trend of TC max for broad-leaved vegetation (τ = −0.47, p-value < 0.05), while for bushy sclerophyllous vegetation, the trend is not significant (τ = −0.02 and p-value = 0.92). Table 3. The mean Mann-Kendall τ and the Theil-Sen slope β of the tree canopy coverage (TC), non-tree vegetation (NTV), and non-vegetated (bare) coverage (BS) for all the broad-leaved forests (BLFs) and bushy sclerophyllous vegetation (BSV) observed pixels falling within each of the 17 sites analyzed. The relationship between tree canopy coverage and MAP1922-2019 for BLF and BSV w investigated year by year. As an example, Figure 5 shows the relationship of MAP1922-2 versus tree-canopy coverage for the year 2008 for both BLF and BSV. We evaluated th breakpoint (MAPcp, identified by the regression analysis in Figure 5) at which the max mum tree cover is attained. The change in tree-canopy coverage as a function of MAP h been determined by using a 99th quantile piecewise linear regression with the identific tion of the breakpoint (the rainfall at which the maximum tree cover is attained), and the maximum upper bound on tree cover (TCmax). We showed the scatter plot for the ye 2008 in Figure 5

Long-Term Trends in Air Temperature and Vapor-Pressure Deficit
The average yearly air temperatures (Ty) in Sardinia have increased significantly by nearly 0.38 • C per decade (τ = +0.33 with p-value < 0.05) since 1990. Warming has been occurring significantly during all seasons but has been greatest during the short spring (τ = +0.33 with p-value < 0.05) and extended summer seasons (τ = +0.40 with p-value < 0.05), with a magnitude of increase of 0.48 and 0.51 • C per decade, respectively. The spatial variability of the Mann-Kendall τ values of the annual air temperature (τTy) on the island is displayed in Figure 7a. It shows an increasing trend over the entire island, with τ values ranging from 0.34 in the hilly areas of the island up to more than 0.5 in the northern and eastern part of the island and particularly close to the coastline. The Theil-Sen slope estimator (Figure 7b) reveals that Sardinia is warming on average by 0.4 • C per decade, with a gradient of increase from west (on average 0.035 • C/y) to east (on average 0.04 • C/y). Focusing on the 17 forests analyzed, all of them are subject to increasing annual and seasonal average temperatures. The maximum rises in average annual temperature are recorded in the forests located in the eastern side of the island (F9, F15, and F6, with a Theil-Sen slope estimator value of 0.52, 0.49, and 0.44 • C per decade, respectively), while by contrast, those located on the central western part of the island recorded the lowest rises (F2, F1, and F11, with a Theil-Sen slope estimator value of 0.52, 0.44, and 0.44 • C per decade, respectively). The Theil-Sen slopes (β) of the trend lines of air temperature during the short spring and extended summer are positive for all the forests, and the rise in T ss and T es is greater than 0.05 • C/y in 15 and 7 forests, respectively. Figure 8a,b show the spatial variability of the Theil-Sen slope estimator β on the island of the annual VPD (VPD y ) and in the 'extended summer season' (VPD es ). Overall, an increase of VPD was observed for most of the study area (Figure 8a). Increases in annual VPD have risen at a rate of approximately 0.003-0.007 kPa/y for 75% of the island, and this extrapolates to a total of 0.12 kPa for the 30-year period. The largest variations in annual VPD were in the western and northern portions of the island and the lower midsouth, while the lowest ones were in the inland hilly areas. During the extended summer season, the VPD increased over the largest areal extent of land, for which, 80% of the surface showed a Theil-Sen slope estimator in the range of 0.010-0.017 kPa/y. The most significant increases were in the mid-northern part of the island. VPD is increasing in all 17 of the forests analyzed, both annually and in every season. Increases in annual VPD have risen at a rate of approximately 0.005-0.009 kPa/y in seven forests, while for the remainder, the rate of growth is slightly lower. orded in the forests located in the eastern side of the island (F9, F15, and F6, with a Theil-Sen slope estimator value of 0.52, 0.49, and 0.44 °C per decade, respectively), while by contrast, those located on the central western part of the island recorded the lowest rises (F2, F1, and F11, with a Theil-Sen slope estimator value of 0.52, 0.44, and 0.44 °C per decade, respectively). The Theil-Sen slopes (β) of the trend lines of air temperature during the short spring and extended summer are positive for all the forests, and the rise in Tss and Tes is greater than 0.05 °C/y in 15 and 7 forests, respectively.

Evaluation of the Climate Impact on Forest Cover Changes
We found interesting relationships between the trends of the forest tree cover and some climate drivers, in particular, the MAP1990-2019 and the Theil-Sen slope estimator of the average annual VPD (Figure 9). Correlation coefficients (Table 4) whose magnitudes are between 0.5 and 0.7 indicate variables that can be considered moderately correlated,

Evaluation of the Climate Impact on Forest Cover Changes
We found interesting relationships between the trends of the forest tree cover and some climate drivers, in particular, the MAP 1990-2019 and the Theil-Sen slope estimator of the average annual VPD (Figure 9). Correlation coefficients (Table 4) whose magnitudes are between 0.5 and 0.7 indicate variables that can be considered moderately correlated, while correlation coefficients whose magnitudes are between 0.3 and 0.5 indicate variables that can be considered only slightly correlated. The Theil-Sen slope estimator of the percentage of tree cover of broad-leaved vegetation (βTC BLF ) was found to be positively correlated with MAP 1990-2019 (ρ = 0.71, p-value 0.004) and negatively with β VPDy (ρ = −0.61, p-value 0.02). Opposite to this is the behavior of the Theil-Sen slope estimator of the percentage of bare soil in the analyzed broad-leaved vegetation (βBS BLF ), which can be considered a measure of the forest mortality rate. βBS BLF shows a positive correlation with the annual change of annual VPD (β VPDy , ρ = 0.62, p-value 0.02) and a negative correlation with the mean annual precipitation (MAP 1990-2019 , ρ = −0.68, p-value 0.006). The correlations between the climatic variables analyzed and the variation in the percentage of tree cover and bare soil within the bushy sclerophyllous vegetation pixels of the forests analyzed are all low. Overall, therefore, the broad-leaved vegetation seems to be more affected by a variation in climate drivers than the bushy sclerophyllous vegetation. Figure 9 shows the scatterplots between the Theil-Sen slope estimators of the tree canopy (βTC BLF ) and of the bare soil (βBS BLF ) coverages and of the annual VPD(β VPDy ) and MAP 1990-2019 climate drivers. We fit (a) a linear function to describe the variability of βTC BLF to β VPDy (βTC BLF = −41.60 β VPDy + 0.092 (2), R 2 = 0.37; Figure 8a) Figure 8d). All of these four relationships are significant highlighting that the warmer conditions, expressed by both an increase of vapor-pressure deficit and a decrease of precipitation, are affecting tree cover in the Sardinian forests. In Figure 9a, most of the data are close to the βTC BLF (β VPDy ) linear regression line, except the F3 data, for which βTC BLF is moderately positive (but not significantly), and β VPDy is positive. Note that the F3 climate conditions are different from the other forests' conditions because F3 is the forest at the highest elevation ( Figure 1), with the lowest annual mean air temperature (~13 • C) and vapor-pressure deficit (~1.7 kPa) and the highest MAP (~900 mm). We used a sigmoidal function for βTC BLF (MAP 1990(MAP -2019 ) relationship for fitting the F3 and F4 βTC BLF positive values, clearly related to the particularly wet conditions (MAP > 750 mm).  Using the six future scenarios of the HadGEM2-AO of the CMIP5 and the six future scenarios of the HadGEM3-GC31-LL of the CMIP6 (Table 5), we estimated the future changes of tree canopy coverage (with respect to 2000-2019) in broad-leaved forests (Figure 10). All scenarios show a positive change in annual vapor-pressure deficit (Figure 10a), ranging from 4.6% in SC7 to 38.64% in SC12. A reduction in the percentage of tree cover is predicted by using the proposed βTCBLF (βVPDy) relationship for all the scenarios, ranging from 4.43% in SC7 to 38.52% in SC12 (peaking above −40% in the five forests, F3, F12, F4, F11, and F2; Figure 10c). Future changes of precipitation are predicted as both positive and negative for CMIP5 scenarios (SC1-SC6), and only negative for CMIP6 scenarios (SC7-SC12). SC6, S11, and S12 were the worst scenarios, with a reduction of 7.6%,12%, and 14.3%, which produced the maximum reduction of tree cover of 16.09%, 16.46%, and 17%, respectively (using the proposed βTCBLF (MAP1990-2019) relationship).  We used the estimated βTC BLF (β VPDy ) and βTC BLF (MAP 1990-2019 ) relationships of Figure 9 in the past for making a coarse prediction of the effects of future scenarios of precipitation and air temperature on forest tree covers.
Using the six future scenarios of the HadGEM2-AO of the CMIP5 and the six future scenarios of the HadGEM3-GC31-LL of the CMIP6 (Table 5), we estimated the future changes of tree canopy coverage (with respect to 2000-2019) in broad-leaved forests ( Figure 10). All scenarios show a positive change in annual vapor-pressure deficit (Figure 10a), ranging from 4.6% in SC7 to 38.64% in SC12. A reduction in the percentage of tree cover is predicted by using the proposed βTC BLF (β VPDy ) relationship for all the scenarios, ranging from 4.43% in SC7 to 38.52% in SC12 (peaking above −40% in the five forests, F3, F12, F4, F11, and F2; Figure 10c). Future changes of precipitation are predicted as both positive and negative for CMIP5 scenarios (SC1-SC6), and only negative for CMIP6 scenarios (SC7-SC12). SC6, S11, and S12 were the worst scenarios, with a reduction of 7.6%,12%, and 14.3%, which produced the maximum reduction of tree cover of 16.09%, 16.46%, and 17%, respectively (using the proposed βTC BLF (MAP 1990-2019 ) relationship).   Table 5 of (a) annual vaporpressure deficit (ΔVPDy) and (b) annual precipitation (ΔPy) in Sardinia. Variation in tree cover with respect to the historical value in the twelve scenarios due to the variation in VPDy alone (c) and in Py alone (d).

Discussion
Sardinian forests are subject to the effect of climate factors and are suffering, as demonstrated by the decrease of the canopy cover almost everywhere in recent decades. There is currently no tool to map the large-scale spatial variability and intensity of this phenomenon. This prompted us to investigate the possibility of using satellite data, and in fact, the analysis of the historical MODIS satellite data over Sardinia reveals the presence of hotspots of significant negative trends in the tree-cover percentage (Figure 4). These are mainly located in some of the forests analyzed by this study that are protected forests, not affected by fires during the period 2000-2020, and where the human impact is insignificant (no change in CLC between 2000 and 2018). The significant negative trends in the tree cover particularly affect the broad-leaved vegetation. A total of 86% of forests show a reduction in broad-leaved tree cover, and for seven forests, the Theil-Sen slope of the percentage of tree cover in broad-leaved vegetation is (βTCBLF) < 0.12%/y. The minimum value occurs in the forest of Monte Arcosu (F5), with a value of −0.30%/y, while the  Table 5 of (a) annual vaporpressure deficit (∆VPD y ) and (b) annual precipitation (∆P y ) in Sardinia. Variation in tree cover with respect to the historical value in the twelve scenarios due to the variation in VPDy alone (c) and in Py alone (d).

Discussion
Sardinian forests are subject to the effect of climate factors and are suffering, as demonstrated by the decrease of the canopy cover almost everywhere in recent decades. There is currently no tool to map the large-scale spatial variability and intensity of this phenomenon. This prompted us to investigate the possibility of using satellite data, and in fact, the analysis of the historical MODIS satellite data over Sardinia reveals the presence of hotspots of significant negative trends in the tree-cover percentage (Figure 4). These are mainly located in some of the forests analyzed by this study that are protected forests, not affected by fires during the period 2000-2020, and where the human impact is insignificant (no change in CLC between 2000 and 2018). The significant negative trends in the tree cover particularly affect the broad-leaved vegetation. A total of 86% of forests show a reduction in broad-leaved tree cover, and for seven forests, the Theil-Sen slope of the percentage of tree cover in broad-leaved vegetation is (βTC BLF ) < 0.12%/y. The minimum value occurs in the forest of Monte Arcosu (F5), with a value of −0.30%/y, while the maximum one occurs in F3, equal to 0.12%/y (Table 4). In some forests, such as F5 and F15, the decrease in βTC BLF is contrasted by an increase in shorter vegetation, and this is confirmed by the positive values of the Theil-Sen slope of the percentage of non-tree cover in broad-leaved vegetation (βNTV BLF ), while in other forests, such as F8 and F9, there is an increase in the percentage of bare soil (βBS BLF ). These results agree with the drought-triggered dieback observed in some of the analyzed forests in the last decade, visible also, for instance, from a comparison of aerial photography of the F5 forest in 2009 and 2019 ( Figure 2) and, in general, after the long drought of 2017 (Supplementary Figure S3). The behavior of the bushy sclerophyllous vegetation is the opposite: in 15 out of the 17 forests the Theil-Sen slope of the percentage of tree cover (βTC BSV ) is positive, and it reaches the maximum value in F12, F07, and F09. These results are in line with the findings of others, who reported a reduction in woody cover along all the Mediterranean Basin [83], or in specific regions such as Spain and Italy (see Gentilesca et al. (2017) [84] for a review) and in Portugal [85]. Definitely, a main reason for which the forest is suffering is the strong decrease of precipitation in the last century in Sardinia. Using the extended Sardinian rainfall database, we performed an accurate trend analysis of Sardinian rainfall that mainly highlighted a change point (1980) in rainfall after which there was an important reduction in the percentage of surfaces with an average value of rainfall of more than 700 mm, which decreases significantly from 74% in 1922-1979 to 38% in 1980-2019, and an increase from 1 to 8 (out of 17) of the number of forests that are below the same threshold. The estimated mean annual precipitation (MAP) threshold (Sankaran et al., 2008) that limits the tree cover in semi-arid forests is 700 mm. Hence, we estimated the first reason of the forest suffering, which is the drastic decrease of the rainfall in the last century in Sardinia. Starting from this new drier condition of the Sardinian forests, we analyzed the recent period, for which more data are hopefully available, and we investigated the climate factors that are affecting forests the most. The analysis of the 20 years of historical MODIS satellite data over Sardinia reveals the presence of hotspots of significant negative trends in tree-canopy coverage. We analyzed the historical precipitation pattern over Sardinia, highlighting a significant average negative trend in the average annual precipitation (the average Mann-Kendall τ over the island is −0.2, with a p-value ≤ 0.05, and the Theil-Sen slope (β Py ) is −13.81 mm per decade), which is due mainly to the negative trend of winter precipitation (average Mann-Kendall τ over the island is −0.21, with a p-value < 0.005). These alterations of the rainfall regime are a widespread phenomenon in the Mediterranean Basin [23,[30][31][32][33]35,37], due to climate-change effects [35]. In Sardinia, 61% of the rainfall stations showed a significant negative trend for the average annual precipitation, but the phenomenon has a marked variability over space (Figure 3a), with a decrease of the negative trend toward the mid-eastern coast (e.g., F3 and F4 forests), as anticipated by Montaldo and Sarigu (2017), who explained that this may be caused by the combined effect of the orography and of North Atlantic Oscillation (NAO) [86]. The reduction on annual total annual precipitation caused a reduction of MAP (Figure 2), and this has been recognized as one of the factors that most influence the spatial distributions of vegetation [87]. Sankaran et al. (2008) demonstrated that, when the MAP value falls below 700 mm, the upper bound of the percentage of tree cover may be limited mainly by water availability, while above this threshold, TC has little to no dependence on rainfall. The MAP analysis over Sardinia (Figure 2), performed by this study for before and after 1980 (change point in annual precipitation estimated with Petitt's method), highlighted an important reduction in the percentage of surfaces with an average value of rainfall of more than 700 mm, which decreases significantly from 74% in 1922-1979 to 38% in 1980-2019, and an increase from 1 to 8 (out of 17) of the number of forests that are below the same threshold (Table 3). More recently, in the period 1990-2019, the reduction in precipitation is even more marked, and in fact, in this period, 78% of the island records less than 700 mm of rainfall per year. A total of 16 out of 17 of the forests were affected by a rainfall reduction, peaking below −8.9% in F14, F13, and F8, and with F14 and F8 that dropped below the threshold of 600 mm of rainfall per year.
Having established that MAP plays a key role in determining the spatial distributions of vegetation, which is decreasing drastically in most of the island, and that this causes an increase in the number of historical forests recently located in water-limited areas (MAP < 700 mm), an attempt was made to assess the presence of any possible relationship with the spatiotemporal pattern of vegetation. The tree cover versus MAP relationships, which were found by other authors [40,42,56,88] in savanna-like ecosystems, are in line with those found in the Mediterranean context proposed by this study ( Figure 5). Indeed, our data indicate a strong dependence of tree cover on rainfall below the thresholds of 695.3 ± 51.65 mm and of 736.76 ± 48.3 mm for broad-leaved and bushy sclerophyllous vegetation, respectively, while above this, there is a higher probability of witnessing a woody plant encroachment. However, looking at the maximum tree cover reached by the vegetation over time (Figure 6), we noticed that the broad-leaved vegetation is showing a 3.4% per decade reduction in TC max , while no trend was found for bushy sclerophyllous vegetation. This means that something is reducing the ability of broad-leaved vegetation to return to its original state, and this has been recognized as an important indicator of forest degradation [89,90]. Although the analysis of the data ( Figure 5) suggested that precipitation is a key driver of tree cover stability, we investigated the further main climate factors that may explain the variability of the tree cover, such as global warming and the rising vapor-pressure deficit [91].
The average yearly air temperature (Ty) is increasing alarmingly in the entire island (with Mann-Kendall τ in the range 0.23-0.47) by about 0.39 • C per decade (see Figure 7), confirming climate-change-induced temperature-rise projections [92,93]. The annual mean temperature is growing in all forests; the smallest increases in βTy are found in the forests of the mountainous areas of the hinterland (F1 and F17) and along the west coast (F7 and F8), while forests located along the east coast record the greatest positive trends (F6, F9, and F15). Summer global warming is affecting the whole island ( Figure 7b) and, thus, all forests (>0.4 • C per decade), and this exacerbates the deleterious effects of leaf overheating and may affect a wide number of metabolic processes, including photosynthesis, respiration, water transport, and phenology [16]. Moreover, the increasing trend in the air temperature may aggravate the effects of drought, because temperature modulates plant water use through its effects on VPD [25,94]. A positive trend in the average annual vapor-pressure deficit has been found in the entire island (average β VPDy of about 0.05 kPa per decade), except for inland mountain areas, where the increase is moderate (i.e., F4 and F17 forests, Figure 8). The greatest increases in VPD y (β VPDy > 0.07 kPa/yr) affected the forests located in the proximity of the Campidano Plain (F14 and F8) and of the east coast of the island (F9). Plants respond to the rising vapor-pressure deficit regulating the gas exchange with the atmosphere through the closure of the stomata [94,95], and the stomatal sensitivity to the vapor-pressure deficit is highly variable across and within species [96]. The latter, such as Quercus ilex, which is the dominant species in the broad-leaved vegetation type (Figure 1), have a variable midday leaf water potential and keep their stomata open and a high photosynthetic rate for long periods, even in the presence of decreasing leaf water potential. This behavior might be beneficial when water is abundant, but under conditions of intense drought, it might endanger the plant [97]. The increase in the vapor-pressure deficit during the summer season is even greater in all forests (βVPD es = 0.007-0.015 kPa/y), making the phenomenon just described more acute and impacting plant transpiration and photosynthesis [95].
To sum up, the increasing vapor-pressure deficit substantially influences vegetation growth and forest mortality [91,98], and if high VPD conditions are combined with reduced precipitation, the risk of experiencing drought stress quickly increases [98]. Sardinia is experiencing a situation of this type. In the last three decades, on the one hand we are witnessing an increase in both average temperatures and the vapor-pressure deficit; and on the other, there is an important variation in MAP caused mainly by a drastic reduction in the winter average precipitation (see Figures 2 and 3), which is important for soil water recharge and is positively correlated with tree growth [20]. The combination of those climate drivers, as found by other authors [15,[99][100][101][102], may be the main cause of the shift from tall tree to shorter or no vegetation (Figure 4), and to the reduced ability of the vegetation to return to its original state (Figure 6), particularly for BLF vegetation. In closed-canopy forests, Bennett et al. (2015) [103] explain that taller trees face a more challenging hydraulic environment because they are exposed to higher solar radiation and leaf-to-air vapor-pressure deficits than those in the relatively buffered understory. This is what is most likely happening in F5, which is characterized by a MAP 1922-2019 of 697 mm, which dropped to 659 mm MAP 1990-2019 , and growing values of temperature (τT y = 0.33) and vapor-pressure deficit (τVPD y = 0.3) and where there is a decreasing trend in broadleaved tree canopy coverage (τTC BLF = −0.3) and a positive trend of broad-leaved non-tree vegetation coverage (τNTV BLF = +0.18). In F8, instead, the reduction in broad-leaved tree canopy coverage (τTC BLF = −0.19) is countered by an increase in bare soil (τBS BLF = +0.25). No significant trend was found in these forests for the percentage of tree cover in bushy sclerophyllous vegetation, demonstrating the greater resistance of Mediterranean maquis to climatic stresses, which is a clear step toward replacing forests with shrubland [15,104]. Areas with positive trend in the percentage of broad-leaved tree canopy coverage (F4 and F3) are mainly located in the inland mountain areas of the island, where water availability is not a problem, and the broad-leaved vegetation might benefit from the warmer weather.
The decreasing trend of broad-leaved tree canopy coverage (βTC BLF ) is positively related to MAP 1990-2019 (ρ = 0.71, p-value 0.004) and negatively related to the Theil-Sen slope of annual VPD (ρ = −0.61, p-value 0.02). On the contrary, the Theil-Sen slope of nonvegetated areas within broad-leaved plots (βBS BLF ) shows a moderate positive correlation with the Theil-Sen slope of annual VPD (ρ = −0.62, p-value 0.01) and a negative one with MAP 1990-2019 (ρ = −0.68, p-value 0.006). This means that, where we have a predominance of broad-leaved vegetation, as the average annual precipitation decreases or the annual VPD increases, the percentage of bare soil increases, while the tree cover decreases. Low or no correlation was found between bushy sclerophyllous vegetation, and the variables analyzed demonstrated low sensitivity at climate drivers. This is likely due to the fact that, drought-and heat-resistant species such as P. latifolia and Pistacia lentiscus, which are common in bushy sclerophyllous vegetation, tolerate the increase of temperature and dry conditions better than more mesic and less cold-sensitive species, such as Q. ilex [105,106].
The use of the proposed simplified βTC BLF (β VPDy ) and βTC BLF (MAP 1990-2019 ) relationships ( Figure 9) estimated from past data allowed us to coarsely predict the percentage of tree cover in the future climate scenarios of the GCMs and account for variations in average annual VPD and MAP 1990-2019 . Six future air-temperature and precipitation scenarios were derived from the HadGEM2-AO [79] predictions of the CMIP5, identified by Montaldo and Oren (2018) as the most suitable for Sardinia, and a further six scenarios from the following generation of HadGEM configuration of CMIP6, the HadGEM3-GC31-LL predictions. We considered three time periods of increasing distance in the future and two different atmospheric CO 2 concentrations (RCP 45 and RCP 85) for CMIP5, and two Shared Socioeconomic Pathways (SSP 245 and SSP 585) for CMIP6. All scenarios agree on the negative effect of the annual average VPD increase, and F2, F3, F4, F11, and F12 are the forests most affected. As expected, the higher the average annual VPD increase, the higher the tree-canopy-coverage decrease, with SC12 being the worst scenario. Future changes in precipitation were contrasting for CMIP5 scenarios; for example, SC6, in line with past precipitation trends, predicts a further drop in mean annual precipitation (−7.6%), which, based on the proposed model, would correspond to an average 2.1 times greater reduction in the tree cover (−16.09.5%). On the contrary, SC2 predicts a turnaround with an abrupt increase of precipitation (+21.5%) in the upcoming years, with a reduction in the number of forests in water-limited areas and an increase in the percentage of tree cover in almost all forests. Future changes in precipitation for CMIP6 scenarios agree on a precipitation reduction in the range of −3.4% (SC7) to −14.29% (S12). However, although the reduction in precipitation predicted in SC12 is almost double that predicted in SC6, the consequent average reduction in TC is comparable and stands at about −16%. Using the 'Space-for-time' substitution approach, a widely used technique to infer future trajectories of ecological systems from contemporary spatial patterns [107], we can predict a change of forest patterns in Sardinia, with more resistant tree species, such as those that compose the bushy sclerophyllous vegetation (P. latifolia, Pistacia lentiscus, Arbutus unedo, Mirtus, etc.). This already applied, for example, in a holm oak forest located in the Montserrat Mountain (Spain) [101], where, after drought-induced events of canopy die-off, there was a change in the composition and structure of the community, involving the increase in species, mostly small shrubs and herbaceous plants, as well as the arrival of non-dominant shrubs and short trees to the canopy. This approach, which follows the principles of adaptive forest management [84], can represent a valid alternative to the conservative management, based on the close-to-natural principles, that is currently required by forestry law in protected areas.

Conclusions
This study added evidence about the power of using satellite data for understanding the long-term tree-cover trend in forests and its relationships with climate drivers. We proposed an economic, user-friendly, and expeditious approach, which can be used by the forest management authorities in both developed and developing countries. The method allows users to map the current status of vegetation and historical trends and patterns and to identify future trajectories based on climate drivers. The results show that vegetation is responding to a warming and drier climate in different ways, as a function of species, mean annual precipitation (MAP), and vapor-pressure deficit (VPD) trends. Within broad-leaved vegetation (BLF) with a MAP lower than 700 mm, we are generally witnessing the decline of tall stands in favor of the increase of short stands, while under the same conditions, bushy sclerophyllous vegetation (BSV) showed an opposite behavior. On the contrary, within broad-leaved vegetation with a MAP higher than 700 mm, stands still have the potential to be stable even under the observed positive trends in temperature and VPD. Moreover, regardless of MAP, broad-leaved vegetation showed a drastic negative trend in TC max (τ = −0.61 and p-value < 0.05), demonstrating that it is not able to recover its original status after the climate disturbances of the last decades. We thus believe that Sardinian broad-leaved forests, as well as others Western Mediterranean oak forests, are declining [84,85,101,108]. In general, the most negative trends in tree-canopy coverage are observed in the broad-leaved vegetation of forests located in the southern-west part of the island, where there are important historical forests composed mainly of tall trees of Quercus ilex. This finding is in line with that of other authors who demonstrate that tall trees, and in particular holm oak, are at high risk of mortality during this changing climate [99,102,103,105]. The medium-and long-term forecasts, based on future climate scenarios and on the simplified models proposed, show that the situation is destined to worsen and deserves attention. It would be useful to demonstrate the compound effects of climate variables and the uncertainties of complex climate effects to vegetation; however, such an analysis is beyond the scope of this paper, which aims only to analyze historical precipitation data in Sardinia, show the actual state of the vegetation and the connection with historical climate data, and to provide simple models that may be easily used by managers and end-users of the forests. The observed overall tree-cover trends may help to identify necessary forest-management activities within the island, whether they are close-to-natural or adaptive. If management of the second type is chosen, the positive trends in tree-canopy coverage observed in most of the bushy sclerophyllous vegetation of the forests analyzed, rather than the absence of correlation with the climate drivers investigated, would lead us to encourage the use of species falling within this category of vegetation. Moreover, it is well-known that forest-cover dynamics have a substantial impact on local water resources [58], and given that Sardinia is experiencing a drastic reduction in runoff [23], future studies should consider the compound effects of those phenomena, even under climate-change scenarios.
Supplementary Materials: The following supporting information can be downloaded at: https://www. mdpi.com/article/10.3390/rs14194893/s1, Figure S1. The p-value of the Mann-Kendall τ (Figure 4) of the percentage of tree cover (a), non-tree vegetation (b) and bare soil (c) for all the observed pixels classified as broad-leaved forests (BLF) and of the percentage of tree cover (c), non-tree vegetation (d) and bare soil (e) for all the observed plots classified as Bushy sclerophyllous vegetation (BSV) over the 20 years period (2000-2019). The thin black outlines indicate the 17 forests analyzed. Figure S2

Acknowledgments:
The authors would like to thank NASA's Earth Science Data Systems (ESDS) program for providing MOD44B data, the European Centre for Medium-Range Weather Forecasts ECMWF for providing ERA5 2 m dewpoint temperature and 2 m temperature data, the climate modelling groups for producing and making available their model output, the Earth System Grid Federation (ESGF) for archiving the data and providing access and the multiple funding agencies who support CMIP5, CMIP6 and ESGF, and the Sardinian Regional Hydrographic Service of ARPAS, and particularly Michele Fiori, for providing rain and air-temperature data.

Conflicts of Interest:
The authors declare no conflict of interest.