Enhanced spring warming in a Mediterranean mountain by atmospheric circulation

We analyzed trends of air temperature across the Cévennes National Park in Southern France, a mid-altitude coastal mountain experiencing a rapid spread of forests at the expense of rangelands and submitted to Mediterranean Sea influences and so, impacted by local and regional processes of climate change. Since 1980, April to June warming trend reached a maximum temperature increase of + 0.124 °C year−1 and uniform whatever the altitude. Minimum temperature increased by + 0.058 °C year−1 at 500 m altitude and + 0.089 °C year−1 at 1500 m. Concomitantly, forest cover is increasing by + 0.51% year−1. Using an intrinsic biophysical mechanism model, we demonstrated that, at monthly scale, the forest surface is 1.7–3.1 °C cooler than that of nearby grasslands. As a result, the decrease in albedo corresponding to the conversion from grasslands to dense forests, translates into a cooling of maximum air temperatures of 0.023 °C year−1 which contributes to slow down the warming rate enhancement. Spring warming trends co-varied with negative WeMO phases associated with a low in the Gulf of Cádiz and an anticyclone in Central Europe. An east to west pressure gradient increases atmospheric humidity leading to a strong water vapor feedback, enhancing the forcing of thermal long wave radiations and hence the rise in temperature.

their topographic features induce modifications to atmospheric flows and exchange processes which preclude theoretical frameworks developed for horizontally flat homogeneous terrain to be applied with confidence 6 . Third, the local near-surface atmospheric flows strongly influence the mesoscale circulation thus further complicating the use of common up-and downscaling approaches 22 . The resulting uncertainties limit our ability to project land-atmosphere interactions in mountain areas under changing future climate and land use 23 .
In addition to these local scale effects, regional atmospheric circulation patterns can be altered with global warming, in turn modifying mass flows at terrestrial/oceanic interfaces. Variations in sea-level pressure dipoles have been shown to influence climate northeast of the Iberian peninsula 24 or to reduce and enhance water-limited areas influences on Mediterranean coastal areas 25 . They are also key factors affecting local climates, independently of the general warming environment.
In this work we focus on the Cévennes National Park in Southern France as a study case, combining a midaltitude coastal mountainous situation in the Mediterranean hot spot of climate change, affected by land cover changes and at the interface between Atlantic and Mediterranean influences. As such, it is submitted to multifactorial impacts of local and regional scale processes of climate change. We specifically address the following questions: (1) In term of strength and seasonality, what are the ongoing elevation-dependent warming (EDW) rates for both minimal and maximal air temperature we observed over the whole study area? (2) Do land use changes contribute to accelerate or to slow down the enhancement of warming rates? (3) How do large-scale climate drivers indirectly influence strength and seasonality of EDW rates?

Results
Trends in monthly temperatures. For the period considered, i.e. from 1980 to 2015, we calculated for each month and for all stations, the least-square linear trend of monthly average maximum temperatures (t x ) and minimum temperatures (t n ). Figure 1 plots the trends obtained for the Mont Aigoual station. Only trends in April, May and June were significant. For example, for t x the June one reached 0.110 (95% confidence intervals 0.026, 0.029) °C year −1 or 1.1 °C per decade. For April, these trends reached 0.100 (0.025, 0.028) °C year −1 and 0.088 (0.027, 0.030) °C year −1 for t x and t n respectively. This pattern of change was observable for all stations. As a result, we aggregated months from April to June and analyzed the average temperatures over this period. For Mont Aigoual, t x and t n trends were then 0.103 (0.041, 0.044) °C year −1 (p < 0.001) and 0.076 (0.031, 0.040) °C year −1 (p < 0.001). Over the complete series of observations available for this station since 1901, Fig. 2 shows that the starting point of our analysis, i.e. 1980, appeared to be in the change point region detected by 26 .
These significant warming trends were not exclusive to the highest station. Larger values were observed, for example, at Nasbinal and Brenoux-Mende, with 0.100 (0.049, 0.037) and 0.109 (0.029, 0.024) °C year −1 for t x and 0.130 (0.044, 0.033) and 0.096 (0.026, 0.025) °C year −1 for t n respectively. These stations are located at altitudes of 1284 m and 1019 m respectively. One low elevation station did not have a significant trend for their t n : Arphy. Altitude played a highly significant role in t n trend (Fig. 3a), which increased linearly with a slope of 0.0032 ± 0.0008 °C year −1 per 100 m and a y-intercept of 0.0418 ± 0.0057 °C year −1 (R 2 = 0.36; p-value = 0.0003).
For t x, no significant trends were observed (Fig. 3b). We can consider this trend as spatially uniform and equal to 0.124 ± 0.007 °C year −1 , corresponding to a + 1.24 °C warming per decade. The effect of altitude was also tested by separating the stations according to their aspects. For t n and t x , the slopes and intercepts were not significantly different depending on whether the aspect was either north-facing or south-facing slopes.
The long monthly precipitation time series available, 1901-2014, for 13 of the 32 stations do not show significant trends for annual totals or seasonal ones from April to June, either for the whole period or only for the period ranging from 1980 to 2014. Graphs for the wettest station, Mont Aigoual (mean annual precipitation MAP = 1792 ± 462 mm) and the driest one, Générargues (MAP = 1217 ± 331 mm) are reported in Fig SI3. For these two stations we analyzed the co-variation between rain totals from April to June and the corresponding J a n F e b M a r A p r i l M a y J u n J u l A u  Figure 4a and b present the within-year patterns of change in white sky albedo WSA for the four target areas: coniferous and broadleaf forest, low shrubland and natural grassland. These patterns cover the twelve years (2000-2011) of 8-day mean ± SD. For the two forested areas, winter albedos were very close with 0.085 ± 0.005 for the coniferous canopy and 0.087 ± 0.006 for the broadleaved one. Then albedo increased to peak early May at 0.126 ± 0.007 and 0.136 ± 0.005 followed by slight summer declines and then finally went   . Linear warming trends observed on (a) mean minimal (blue circle ± 95% CI) and (b) maximal (red circle ± 95% CI) April to June temperatures (°C year −1 ) and elevation (m asl) for the 32 studied stations.  (Fig. 4a). The WSA temporal patterns for low shrublands and natural grasslands were affected by snow cover within the cold period. On average, snow began late November-early December for the two vegetation types and disappeared end of February in grasslands and later, end of March, in the low shrublands (Fig. 4b). During the vegetative period, albedos reached their maximum values at the beginning of July for natural grasslands (0.186 ± 0.011) and then mid-July for low shrublands (0.156 ± 0.005). Standard deviation was higher in natural grasslands because they were more sensitive to local variations of soil water limitations and erratic stormy rainfall. When pooling the three forest types in the forest land cover type, a significant linear relationship was observed between April to June WSA and percent forest cover (PFC) surrounding the 32 weather stations: WSA = − 7 × 10 -4 ± 7 × 10 -5 PFC + 0.186 ± 0.005 (R 2 = 0.75, p-value < 0.0001). This means that the albedo decreased by 7 × 10 -3 when the percentage of forest increased by 10%. The albedo level corresponding to the y-intercept was equal to the one previously observed on the target zone of pure natural grassland plotted in Fig. 4b. This local albedo around the stations, which is assumed to be slightly constant over the observation period 2000-2011, played a significant role in warming trends. For t x , this warming trend increased linearly with WSA with a slope of 0.385 ± 0.182 °C year −1 albedo unit −1 and a y-intercept of 0.063 ± 0.027 °C year −1 (R 2 = 0.15; p-value = 0.04) (Fig. 5). This means that transition from natural grasslands with an albedo of 0.185 to a forest with an albedo of 0.125 dampened the warming trends by 0.23 °C per decade. For t n , no significant relationship between warming trend and WSA was observed. The slope was lower than the one obtained for t x trend, with a value of 0.189 ± 0.217 °C year −1 albedo unit −1 (R 2 = 0.03; p-value = 0.39). The dampening effect of forests on the warming trend for minimum temperatures was less than for maximum temperatures. This result tends to demonstrate the likely effect of biophysical variables such as albedo on changes in the regional field of air temperatures. However, this site-specific result does not contribute in explaining the slightly decreasing, statistically non-significant, elevation-dependent warming of spring maximum temperatures.

IBPM simulations.
In our afforestation scenario, we applied the IBPM model. We examined the three components of change in surface temperature separately. In April, the applied albedos were respectively 0.180 ± 0.002 and 0.110 ± 0.006 for natural grasslands and forests. Those in July were 0.184 ± 0.002 and 0.124 ± 0.002. Over the 6-month period analyzed, changes in albedo were rather constant and weak, ranging from 0.059 ± 0.009 and 0.070 ± 0.003. As a consequence, changes in daily net shortwave radiation were limited ranging from + 24.6 ± 9.8 W m −2 in September to + 37.7 ± 10.1 W m −2 in July. The shortwave radiative forcing effect due to albedo change yielded a warming between + 0.48 ± 0.22 °C in April and + 1.06 ± 0.39 °C in July. The main effect was the result of changes in roughness. Figure SI4 plots the monthly field of wind velocity. Mean daily wind speed at 2-m height was the highest in April with a value of 2.64 ± 1.05 m s −1 and the lowest in August at 2.00 ± 0.75 m s −1 . The corresponding aerodynamic resistances for grasslands and forests were respectively   Monthly cooling induced by changes in the energy redistribution factor attributable to changes in surface roughness ranged from − 1.98 ± 0.75 °C in April to − 4.03 ± 1.04 °C in July. The role of energy redistribution associated with changes in Bowen ratio was limited. Figure SI5 represents the time courses of Bowen ratios of both ecosystems. During the well-watered month of April, Bowen ratio for natural grassland was 1.79 ± 0.41 and 1.38 ± 0.15 for forest. With the occurrence of the summer drought, this ratio increased for both ecosystems respectively to 7.48 ± 1.09 and 5.25 ± 0.70 in July and then 7.10 ± 1.66 and 5.34 ± 0.79 in August. The Bowen ratio cooling only reached − 0.17 ± 0.07 °C in April and − 0.15 ± 0.04 °C in July. Finally, by summing up the three terms of the surface temperature change, we obtained a decrease between − 1.67 ± 0.78 °C in April and − 3.12 ± 1.11 °C in July ( Fig. 6). We demonstrate here that the drivers of this unexpected spring warming are not local and surely not related to landscape changes, so we studied further atmospheric patterns known to have a significant impact on the Mediterranean climate.
Large-scale influences. By considering Mont Aigoual as our reference station, we analyzed the correlation between the monthly temperatures t x and t n from 1980 to 2015 with the corresponding values of three teleconnection indices, i.e. WEMO, NAO and MO. The time courses of the nonparametric Spearman's correlation have been reported on Table 2. For WEMO, Spearman r was highly significant (p < 0.001) in April and May for both t x and t n . The correlation declined in June but with p-values remaining nevertheless < 0.01. The second time-period for which we observed a significant statistical link between WEMO and temperatures was from September to November with a peak in September (p < 0.001) followed by a gradual decrease in the significance. Correlations during the other months were not significant. For NAO index, highly significant correlations were observed for A p r M a y J u n J u l A u g S e p Figure 6. Components of monthly changes in surface temperature (red bar ± SD) in the afforestation scenario: shortwave radiative forcing effect due to albedo change (brown bar ± SD), energy redistribution associated with changes both roughness (blue bar ± SD) and Bowen ratio (green bar ± SD). By averaging indices and temperatures for April to June, months with highly significant warming trends (see Fig. 1), we found, as expected, highly significant correlations of WEMO with both t n (r = − 0.618, p < 0.001) and t x (r = − 0.698, p < 0.001). For the same period, the correlations were not significant for either t n or t x both with NAO and MO. For extending such analysis to the whole-data set, which is the 32 station data, we performed a principal component analysis (PCA). The first component of the PCA accounted for 43.7% and 32.8% respectively of the total variance for t x and t n , and the second component only described 6.4% and 7.6%. The Spearman correlations between the three indices and the first eigenvectors were only significant for WEMO with r = − 0.54 (p < 0.001) and r = − 0.67 (p < 0.001) for t n et t x respectively. The link between climate indices and more specifically the WEMO and maximum and minimum temperatures in spring raises questions about the underlying mechanisms.

Discussion
We have highlighted ongoing warming trends that are largely higher than those observed at global scale for the same time step and for very close change points. For the three-month spring period from April to June, the warming trend reached for maximum air temperature 0.124 °C year −1 , i.e. 1.24 °C per decade. This warming was uniform whatever the altitude. For the same period, minimum air temperature increased by 0.058 °C year −1 at 500 m altitude and this warming rate reaches 0.089 °C year −1 at 1500 m. All other months displayed the weakest warming trends which were not statistically significant. For comparison, global rates over 1979-2010 ranged from 0.141 ± 0.015 to 0.175 ± 0.012 °C/decade depending on data sources 1 .
Our results agreed with previous reported warming rates in the north western Mediterranean basin. They appeared to be the highest across the world and contributed to place the Mediterranean region as a climate change hot-spot 27 . Over Spain, the warming rate is particularly significant for mean air temperature in spring and summer and close for both seasons to 0.6 °C/decade 12 . However, here exist regional variations. In Catalonia (Spain), Reference 28 reported only significant warming rates in spring with similar magnitude to our observations. Over Italy, Reference 29 observed a peak trend occurring in summer with mean air temperature increasing from ca. 0.65 °C/decade. While a consensus is emerging that fall and winter are not significantly affected by warming, the other two seasons are those with the highest elevated warming rates. However, the fact that either spring or summer or both may be affected suggests that underlying meso-scale climate processes are active. If we shift our focus to mid-altitude mountainous regions, the scarcity of observations does not allow us to generalize our results. However, sparse evidences do indicate an increase in the warming rates with altitude, such as that observed more often in high mountains 10 . In Spain, Reference 15 observed it in the Montseny natural park (Catalan Pre-Coastal Range) and 16 in the Sierra de Cazorla (Betic mountain range). Similarly, Reference 30 noted that the increase in air temperature is higher with altitude in the Abruzzo Region (central Italy).
Among the causes of enhanced warming of such magnitudes, the degree of landscape changes should be considered. Land Use and Land Cover (LULC) change is an important driver of environmental change, occurring at the same time as, and interacting with, climate change 31,32 . In the Cévennes National Park, we observed a rapid spreading of forests and woodlands at the expense of rangelands. Such changes have led to a linear increase in the percentage of forest cover from 38 to 68% between 1955 and 2015, at the expense of non-forested areas which declined from 62 to 25%. This landscape closure occurred in parallel with an increase in percentage of coniferous forests within the forest land cover type. In our study area, mechanisms of invasion by Pinus spp. have been detailed in 33 . These mechanisms have a rather broad nature 34 . Afforestation has already led to a significant decrease in water yield flowing from catchments whose headwaters are located within the CNP 35 and significant losses of biodiversity 36 . During the last centuries, Mediterranean mountain landscapes have experienced  www.nature.com/scientificreports/ a gradually decreasing intensity of human impacts due to declines of small-scale agriculture, rangelands, and forest use, particularly in low productivity areas 20,37,38 . The increasing abandonment of mountain and rural areas, generally triggered by the decline of livestock grazing induced a natural expansion of tree cover arising from secondary succession or closure of pre-existing woodlands. Natural reforestation is a heterogeneous and site-dependent process that is driven by topographic, climate, and socio-economic factors 33,39 . So, we identified common mechanisms in the CNP, in Spain 19,40 and in Italy 41 .
Recently, attention has focused on understanding the potential for LULC change to affect surface temperature by modifying biophysical properties, energy acquisition and partitioning 42 . In this work, we limited our analysis in assessing consequences of changing vegetation cover with the highest contrast: natural grasslands and dense coniferous forests. Surface temperature is a state variable reflecting land surface-atmosphere feedback driven by water and energy exchange, and it is linked to the roughness of vegetation cover and the water exchange processes, which have been synthesized here by the Bowen ratio. As a result of our simulations, forests with a lower albedo absorb more solar irradiation, leading to a warming component ranging from 0.5 to 1.1 °C from April to July when compared to nearby grassland. The enhanced turbulent fluxes of latent and sensible heat of forest have a combined cooling effect, from − 2.2 to − 4.2 °C that offsets the albedo warming effect. Finally, at monthly time scale, the surface of the forest is 1.7-3.1 °C cooler than the grassland, indicating a substantial cooling effect. Our results are close to observations obtained in comparing paired flux towers 32,42 or from remotely-sensed products 43 . However, assessing climate mitigation potential of forests requires an understanding of how LULC change affects not only surface temperature but also air temperature. To this point, extending effects of surface temperature to effects on air above the surface has been difficult 44 . For this reason, apparently puzzling results of warmer surfaces under colder air and vice versa have been reported 45 . This paradox extends to regional scale assessments 46 . The impact of a surface forcing is determined by the effective heat capacity of the atmosphere, which is defined by the depth of the planetary boundary layer. This can vary by an order of magnitude on contrasted time and spatial scales, and so this results in a much boosted temperature response in shallow boundary layers compared to a thicker one. Over mountainous terrain, multi-scale thermally driven flows play a key role in developing the planetary boundary layer and understanding of the underlying processes are still limited 47 . Therefore, it is difficult to infer consequences in terms of air temperature trends of the large LULC changes that drastically affect the CNP landscapes. However, the analysis we conducted in mapping the albedo surrounding the meteorological stations of our study area gives us a reasonable answer. The decrease in albedo from 0.185 to 0.125, i.e. the one corresponding to the change from natural grasslands to dense coniferous forests, translates into a cooling of the maximum temperatures of 0.26 °C/decade which contributes to slow down the warming rate enhancement. This result is completely in line with the estimates given in the synthesis of 48 .
As we demonstrated that the drivers of this unexpected spring warming are not local and surely not related to landscape changes, we studied three climate indices known to have significant impact on the Mediterranean climate: the North-Atlantic Oscillation (NAO), the Mediterranean Oscillation index (MO) and the Western Mediterranean Oscillation index (WeMO). Difficulties in the interpretation arise in using climate indices to the fact that another level of complexity is added, namely the link between climate indices and local climate. Several studies have shown that climate indices are linked with precipitation regimes and temperature trends 24,29 . European climates are particularly influenced by NAO, a north-south sea-level pressure dipolar pattern, with one of the centers located over Iceland and the other one approximately over the Azores Islands that controls the location and strength of storm tracks over the North-Atlantic and drives patterns in temperature and precipitation. Over the last decades, the phase of the NAO has been shifting from mostly negative to mostly positive index values 49 . As a consequence, this results in an intensified westerly flow during the positive phase, that brings warm maritime air to Europe leading to a warming of southern Europe. Reference 50 observed "the similarity of the departure patterns of temperature anomalies over these regions are strongly related to the persistent and exceptionally strong positive phase of the NAO index since the early 1980s". However, in our study its impact remains rather limited, particularly with respect to spring temperatures. Other teleconnections are known to also significantly influence the Mediterranean climate, particularly the WeMO index. A positive phase of the WeMO typically shows an anticyclone located in the Gulf of Cádiz and a low-pressure area in Genoa Gulf whereas a negative WeMO phase will show a low in the Gulf of Cádiz and an anticyclone in Central Europe. During the positive phase, the prevailing winds are typically west and northwest, originating from the North Atlantic area. In contrast, a negative WeMO phase is associated to airflow which has traveled over the Mediterranean Sea. This configuration is at the origin of the high spring warming that we have evidenced in the Cevennes National Park. WeMO has been retained as a proxy for local environmental conditions such as Sea Surface Temperature (SST) with negative values associated with high SST 51 .
The evolution of near-surface air temperature is mostly impacted by radiation, heat advection, and land-atmosphere interactions. Refernce 52 demonstrated that the larger part of year-to-year variation and longterm trend of daily temperature can be explained by the radiation. The effects of incoming shortwave and thermal longwave radiations have to be analyzed apart. A widespread reduction of incoming shortwave or global solar radiation between the 1950s and 1980s has been observed and since late 1980s a reversal in this trend has been recorded in many regions of the world 53,54 . Changes in transmissivity of the atmosphere due to change in anthropogenic aerosol emissions and cloudiness are considered to be the most likely causes 53,55,56 . For mainland Europe, daily brightening may reach 4.5 W m −2 per decade with a significant amplification in spring for the southern regions 57,58 . For the north western Mediterranean basin, this seasonality covariates with large-scale atmospheric circulation 59 . Brightening contributes in part to the rapid temperature increase since 1980 we observed over the CNP and that is considerably larger than the rise expected from greenhouse gas increases alone 60 .
The negative WeMO spring phases are associated with a strong east to west positive gradient in atmospheric humidity caused by evaporation as air travels across the Mediterranean Sea, leading to a strong water vapor feedback enhancing the forcing and temperature rise by about a factor of three 61 www.nature.com/scientificreports/ increases in response to increasing air specific humidity, i.e. the concentration of water vapor expressed as the mass of water vapor per unit mass of moist air. It is proportional to the partial pressure of water vapor and inversely proportional to the atmospheric pressure. This relationship is nonlinear with substantially higher sensitivities at low levels of air specific humidity, especially below 5 g kg -163 . These conditions occur more frequently during the cool spring season in the CNP area and might cause the observed elevation-dependent warmings 7 .
Here, we suspect that additional heat release from water vapor acts as the main amplifying mechanism for climate change. This influence should be carefully studied in the future. We have seen that air and sea interact on a wide range of scales, shaping climate and influencing weather. Coupled regional climate system models can help us to unravel these complex interactions. Recent simulations show that the atmospheric components correctly reproduce both large-scale and local features of the Mediterranean climate, although large temperature biases were observed. The sea component correctly reproduces the sea surface properties along with their interannual variation and circulations 64 . In spring, it is likely the air-sea interactions in the Mediterranean Sea that exert a significant influence on the regional biases 65 . Despite ongoing model improvements, it is suggested that higher spatial resolution could be needed for improving simulation of large-scale atmospheric circulations that impact the climate of the Cévennes National Park. The southern part corresponds to a Mediterranean-type climate turning gradually into a continental climate as it moves northward, leading to contrasted weather patterns. Its topography results in precipitation yearly totals ranging from 800 mm to almost 1800 mm. From 1901 to 2014, Mont Aigoual received an averaged amount of 1792 mm with 25% and 75% percentiles of 904 mm and 2117 mm. Annual average temperature over the wholearea ranged from 3 to 16 °C. Landscapes were shaped until the middle of the twentieth century by agriculture and pastoralism. Since 1950, the region has been subjected to a severe human depopulation, which has led to a severe agricultural decline and drastic changes in the landscape. Today the mosaic of vegetation units is more complex with coniferous forests (Pinus sylvestris, P. nigra), broadleaf deciduous (Fagus sylvatica, Castanea sativa) and evergreen (Quercus ilex) forests but also shrublands and grasslands.

Material and methods
Climate analysis. Trends in monthly averaged daily minimum and maximum temperature data collected at 32 meteorological surface stations (henceforth stations) were analyzed with a non-parametric statistics, the Theil-Sen estimator 66 (Sen 1968), more accurate than simple linear regression for skewed or heteroscedastic data. Data were downloaded from the MétéoFrance National Climatic Data Center (https:// publi theque. meteo. fr). Stations were selected to maximize completeness, minimize station changes, and optimize spatial coverage over Cévennes National Park (central protected plus adhesion zones) and its close surrounding. We didn't retain stations near or in urban areas to minimize the effects of urbanization on temperature.
Selected stations had at least 90% of data actually registered for the analysis period of 35 years, from January 1980 to December 2015 (see station list on Table SI1 and map on Fig. SI1). Each station was characterized by its elevation and aspect: north-or south-facing slope. Analyses were performed on time series of individual months for each station to determine if trends occurred during particular seasons. We choose the highest instrumented station, Mont Aigoual, as a reference site thanks to more than a century of monthly data. Additionally, we analyzed co-variations of minimum and maximum temperatures with monthly rainfall data from 13 stations among the 32. The CORINE Land Cover maps were produced in 2000 (CLC2000) and in 2006 (CLC2006) as part of a European Community program to generate digital land-use/land-cover maps covering the European continent (https:// land. coper nicus. eu/ pan-europ ean/ corine-land-cover). The maps were produced using satellite Landsat7 ETM + data and other ancillary data (Digital Elevation Model and aerial photos). The maps have 44 classes and a spatial detail comparable to that of the paper map or a 1/100,000.
Additionally, we used two other digital maps: the National Forest Inventory map for 1999 (NFI1999) was based on a regular 1 km × 1 km samples with detailed vegetation inventory on each grid point scaled up with the help of false color orthorectified aerial photos (scale 1/100,000) and the GISLR2006 built for the Occitanie region (where CNP belongs to) from supervised classification of satellite Landsat7 ETM + data (scale 1/100,000) available at OpenIG: https:// www. siglr. org. For further analysis we grouped together vegetation classes in three major Surface albedo. In order to characterize regional pattern of albedo changes, we used MODIS (Moderate Resolution Imaging Spectrometry) MCD43A3 albedo products (Collection 5) available at the geoportal https:// lpdaac. usgs. gov/ data_ access. Albedo products are derived from surface reflectances, corrected for atmospheric scattering, ozone absorption, and aerosols. They represent an eight-day 500-m composite combining reflectance data from Terra and Aqua satellites. As albedo depends not only on surface properties, but also on solar angle and atmospheric conditions, we considered in this work white-sky albedo to isolate the dependence of albedo on surface characteristics. White-sky albedo (WSA) or diffuse albedo, represents bi-hemispherical reflectance under isotropic skylight illumination; therefore, it is independent of sky conditions. WSA is particularly useful for characterizing biophysical variables of contrasting land covers at pixel scales. They were used primarily to evaluate differences in surface albedo between contrasted vegetation types over 12 to 16 pixels for target zones of pure vegetation class relatively unchanged from 2000 to 2011. We also evaluated the mean surface albedo of area of 36 pixels centered on the 32 weather stations concurrently with the percent cover of the 'forest' type.
Change in surface temperature with the intrinsic biophysical mechanism model. We limited this analysis in assessing consequences of changing vegetation cover from natural grasslands to dense coniferous forests on the energy balance and consequently, on surface temperature. For this purpose, we applied the intrinsic biophysical mechanism (IBPM) model proposed by 31 in prognostic mode (details are available in the supporting information). The change in land surface temperature depends on three terms: (1) the shortwave radiative forcing effect due to albedo change and the role of energy redistribution associated with changes in both (2) roughness and (3) Bowen ratio. We focused here only on monthly averaged daily temperature for the period ranging from April to September. The daily climate data came from the Lafage-INRA experimental station (lat. 3°05′E, long. 43°55′N, 804 m asl, see Table SI1): maximum and minimum temperatures, precipitation, global radiation and wind speed. Albedo variations are derived from MODIS WSA (see previous paragraph). The aerodynamic resistances (r a ) of the two vegetation covers (grassland and coniferous forest) are calculated monthly for a wind velocity which is the average velocity observed for each month of the selected time window (see Fig. SI3). The Bowen ratios β are obtained by simulating a simple soil water balance model based on the calculation of the Priestley-Taylor potential evapotranspiration (PET) and a relationship between the ratio of actual evapotranspiration (AET) on PET and the relative soil water content proposed by 68 . This model has been validated for natural grasslands on measurements of soil water storage (see 33  Large-scale influences: the role of teleconnection. The relationship between teleconnection patterns and air temperature was studied by means of the Pearson's correlation coefficient. The selected teleconnection patterns were the most prominent for Western Europe, the NAO, and two Mediterranean patterns: the MO and the WeMO. The North Atlantic Oscillation index (NAO) values are from the Gibraltar-Iceland dipole 69 . The Mediterranean Oscillation index (MO) values are defined by the dipole Gibraltar-Lod Airport (Israel) 70 . Lastly, the Western Mediterranean Oscillation index (WeMO) has recently been suggested for South Western Europe 71 using the dipole Cadiz (southern Spain)-Padua (central Italy), with the Mediterranean coastland as main area of influence. In their dipoles, the three indices share values of Gibraltar or its surrounding areas. The three indices were calculated using normalized sea level pressure. Daily MO and monthly NAO index were obtained from the Climate Research Unit, East Angelia University, UK (http:// www. cru. uea. ac. uk/). MO data were further aggregated on a monthly time step. The WeMO monthly data were made available by the Climatology Group, University of Barcelona (http:// www. ub. es/). The Pearson correlation coefficient has been used to quantify the relationships between large-scale circulations and both seasonal temperature time series and main eigenvectors of a Principal Component Analysis (PCA) processed on the whole-CNP t x and t n data set.