Amplified local cooling effect of forestation in warming Europe

Domain: earth_science
Published: 2025-09-25
Editorial outcome: minor_revision

Abstract
Forests exhibit local cooling or warming effects compared to adjacent openlands through biophysical processes. These temperature effects are predicted by earth system models to evolve in response to climate change. However, such temporal patterns remain unconstrained by observations and have not been detected in historical records. Here, by comparing the satellite observations of spatially nearby forests and openlands over the last two decades, we quantify temporal trends in local land surface temperature (LST) effects of forest change in Europe. During winter, the daytime warming effect of potential forestation weakens and reverses to cooling (−0.142 K/decade) with decreasing snow cover, as forests show less pronounced surface darkening trends than open lands. During summer, the daytime cooling effect intensifies (−0.188 K/decade) because forests remain more physiologically and hydrologically resilient to increasing soil dryness. These negative trends are broadly supported by state-of-the-art earth system models, though substantial inter-model variability persists. Given continued climate change, we emphasize the need to consider the dynamics of biophysical effects when comprehensive forest-related climate policies are formed.

Introduction
Forests provide numerous benefits to humans and the planet, such as producing food and energy, reducing soil erosion, and increasing water availability in downwind areas1–4. Increasing forest cover through restoration practices plays a crucial role in addressing global warming by regulating the hydrological cycle and enhancing cloud formation5–7. In addition to the global impacts, forests exhibit local cooling or warming effects due to their unique biophysical characteristics that regulate surface energy and water fluxes7,8. The biophysical effects of forests on local climate have garnered particular attention in recent years, as they can significantly intensify or offset the climate mitigation effects of carbon sequestration9,10.

Using numerical models or multi-scale observations, previous studies have thoroughly investigated the impacts of re/afforestation, deforestation, and forest degradation on local land surface temperature (LST) and their biophysical mechanisms11–13. At night, forests warm the land surface compared to openlands through two mechanisms. First, forest canopies absorb solar radiation and store heat during the daytime, subsequently release this energy at night. Second, the nighttime land surface is typically cooler than the near-surface air, forming a stable boundary layer where downward atmospheric heat fluxes warm the land surface. Forests, with their greater surface roughness, can intensify turbulent mixing and enhance this heat exchange process14,15. During the daytime, forests exhibit two competing effects compared to non-forest vegetation: (a) the non-radiative cooling effect driven by higher evapotranspiration (ET) rates (when water supply is sufficient)16, and (b) the radiative warming effect resulting from lower albedo17. The net temperature response depends on the balance between these opposing effects18, which is predominantly governed by the background climatic or environmental conditions (e.g., temperature, snow cover, and soil moisture)19. For instance, during boreal winter when ET is rendered inactive by low temperatures, the dense forest canopies obscure bright snow, which remains exposed in openlands. This snow masking effect creates a strong albedo contrast between forests and openlands and causes forests to absorb more solar radiation compared to surrounding openlands, resulting in the local warming effect4. In warm and wet tropics, forests demonstrate pronounced local cooling effects due to their enhanced ET rates relative to other vegetation types20. Limited soil moisture in arid regions or during drought events can constrain forests’ evaporative advantage, enabling radiative warming to dominate. Consequently, the magnitude of forests’ cooling effect at the seasonal or annual scale may depend on the duration of sufficient water availability21.

Given the profound influence of background climate on the biophysical processes, it is reasonable to infer that the regulation of forests on local temperature will shift in response to rising CO2 concentration, global warming, or varying water availability. In recent years, studies have revealed the role of rising temperature, changing aerosols or internal climate variability in the temporal dynamics of the biophysical temperature effects22–25. However, these findings are based on model simulations, which are subject to inaccuracies in the physical parameterization scheme and representation of land-atmosphere interaction26,27. Besides, prior simulations only show significant changes under high-emission or other ideal scenarios22–25. On the other hand, observational studies that assess the biophysical effects of forestation are generally performed in a relatively static manner13,18,28. These studies focus on spatial, seasonal or diurnal heterogeneity of the biophysical effect over a specific year or study period but ignore the potential long-term variations. Overall, direct observational evidence of temporal dynamics in local temperature effects of forests is still lacking in historical records. Evaluating these dynamic patterns and their driving mechanisms is useful for constraining model results and forming better forest-related climate policies.

This study aims to fill the knowledge gap regarding the temporal variations (especially trends) in the local LST effects of potential forestation and the underlying biophysical mechanisms. It should be noted that the evaluated LST effects are local effects and should be treated as the reference for climate adaptation rather than global climate mitigation. We focus on Europe for the following reasons. (1) Both radiative and non-radiative processes play important roles in the local LST effects of European forests, facilitating the exploration of different biophysical mechanisms underlying the long-term variation of the forests’ LST effects. (2) Thanks to the dense distribution of weather stations, high-quality gridded meteorological data in Europe can help to explore potential links between changes in biophysical processes and the background climate. (3) Europe is a hotspot of climate change and forestation projects (3 billion trees to be planted in the European Union by 2030 according to the European Green Deal29), making our assessment more relevant and instructive. Our study aims to answer the following three questions: (1) What are the temporal patterns of the observed local LST effects of European forests in the recent two decades? (2) What climate factors drive these temporal patterns? (3) What is the biophysical mechanism behind this climatic control? To address these questions, we first use observations from the Moderate Resolution Imaging Spectroradiometer (MODIS) and space-for-time method to estimate the local LST effects of forestation in Europe28. We then assess the temporal dynamics of the local LST effects and explore their relationships with different biophysical processes and potential climate drivers using statistical methods. To examine whether models can capture the observed temporal patterns, we evaluate the simulated trend of local LST effects of European forests in four state-of-the-art Earth system models against satellite observations. Our results could inform comprehensive assessments of climate effects induced by future forestation in Europe.

Results
Intensified local daytime cooling effects of forests
We first investigate the multi-year mean impacts of potential forestation on local LST in Europe winter and summer (ΔLSTw and ΔLSTs) (Supplementary Figs. 1 and 2). Results show that forests exhibit a slight winter warming effect (0.14±0.27 K, spatial mean±standard deviation) compared to spatially nearby openlands, primarily driven by nighttime warming (0.25±0.24 K), while the daytime signal remains negligible (0.01±0.46 K). In contrast, forests exhibit strong land surface cooling in summer (−0.75±0.41 K), as the daytime cooling (−2.11±0.98 K) is minimally offset by a much smaller nighttime warming (0.25±0.23 K). Our results confirm previous conclusions of the contrasting LST effects of mid-latitude forestation at seasonal and diurnal scales driven by different biophysical processes14,28.

We then focus on the temporal variations in ΔLSTw and ΔLSTs over the recent two decades. Our results show declining winter warming effects at the daily scale, with a trend of −0.052 K/decade (p=0.07) (Fig. 1a). This decrease in the daily mean ΔLSTw is dominated by a negative trend of the daytime ΔLSTw (−0.142 K/decade, p=0.04), while the nighttime ΔLSTw demonstrates only interannual variability with a negligible temporal trend (0.002 K/decade, p=0.97) (Fig. 1a). Spatially, 81% of grids show negative daytime ΔLSTw trends, with only 10% of grids mainly concentrated in southeastern Europe showing significant trends (p<0.05) (Fig. 1b).

During the summer nighttime, the warming effect shows an increasing trend of 0.040 K/decade (p=0.01), but the magnitude is only about one-fifth of the daytime ΔLSTs negative trend (−0.188 K/decade, p=0.01) (Fig. 1c). As a result, the daily mean cooling effects have been enhanced in the recent two decades (−0.083 K/decade, p=0.01) (Fig. 1c). Spatially, about 89% of grids show negative daytime ΔLSTs trends, and 34% of grids (mainly distributed in central Europe) show significant negative trends (p<0.05) (Fig. 1d). These results show that in both summer and winter, the daytime cooling effects of potential forestation have been significantly enhanced. Further analyses mainly focus on the biophysical mechanisms behind the negative trends in daytime temperature effects.

Impact of snow cover on daytime ΔLSTw dynamics
During winter, the sign and absolute value of ΔLSTw are generally dominated by the magnitude of the radiative process4,8,18. However, the negative trend in ΔLSTw could potentially arise from two mechanisms: (a) climate-driven alterations in net shortwave radiation absorption differences between forests and openlands (ΔNSRw), or (b) winter warming-induced enhancement in forests’ evaporative cooling advantage, represented by latent heat differences between the two land cover types (ΔLEw). To identify the dominant mechanism, we estimate the equivalent perturbations of the two processes on the surface energy balance (Supplementary Text 1).

Our analysis reveals that only ΔNSRw exhibits a statistically significant trend (−0.86 wm−2/decade, p=0.05), whereas the ΔLEw trend remains non-significant (−0.02 wm−2/decade, p=0.97) (Supplementary Fig. 3). Correlation analyses further support the greater importance of ΔNSRw affecting ΔLSTw dynamics, and confirm that this impact is mediated through modifications in the albedo difference (Δαw) between forests and openlands (Supplementary Figs. 4–6). Specifically, we show that forests are darker than openlands, as evidenced by negative Δαw values (−0.077, Supplementary Fig. 7). However, Δαw shows a positive trend (0.018 unitless/decade, p=0.04) (Fig. 2a), indicating that the albedo contrast between forests and openlands has been diminishing. Spatially, Δαw trends are positive in 90% areas, and about 15% of grids (concentrated in central and eastern Europe) are statistically significant (p<0.05) (Fig. 2b).

The revealed positive Δαw trend can be traced back to the stronger albedo decrease in openlands (multi-year mean: 0.236; trend: −0.036 unitless/decade, p=0.04) than in forests (multi-year mean: 0.159; trend: −0.017 unitless/decade, p=0.03) from 2004 to 2023 (Fig. 2a). This observed divergence in surface darkening trends may be associated with snow cover (SCw) reductions (Supplementary Fig. 8), and differences in the snow-albedo feedback intensity in different land covers30–34. We confirm that the surface darkening is more susceptible to SCw decrease in openlands than forests, as evidenced by the difference in spatial regression slopes between albedo trends and SCw trends (0.47 vs. 0.22, Fig. 2c). The possible mechanism is that snow tends to be masked by tree canopies and retained for longer durations in forest ecosystems, resulting in a lower sensitivity of forest albedo to snowpack. The divergent responses of albedo to SCw are also supported by the temporal regression results, which suggest that openland albedo is about two times more sensitive than forest albedo (Supplementary Fig. 9). Notably, the used SCw is a function of snow mass (represented by the snow water equivalent)35, which is more governed by the macroclimate and thus shows no difference between forests and nearby openlands. These results indicate that although forest albedo tends to decline due to SCw loss, the relative darkening effect of forests and the additional solar radiation absorbed by forests (compared to adjacent openlands) have diminished in the recent two decades.

The positive Δαw trend induced by snow dynamics further controls the temporal change of ΔLSTw. Specifically, in those grids with more pronounced SCw reductions, the Δαw trends is larger and the negative daytime ΔLSTw trends are also stronger (Fig. 2d). Meanwhile, the temporal relationships between SCw and ΔLSTw suggest that in those years with lower snow contents, the daytime cooling effect of forests is also stronger (slope=1.08±0.41), while the impact of SCw on the nighttime ΔLSTw is insignificant (slope=0.22±0.20) (Fig. 2e). Overall, our statistical evidence supports the control of SCw on daytime ΔLSTw dynamics through weakening the albedo difference between forests and openlands. Notably, the observed ΔLSTw trend should not be attributed to structural or physiological changes in forests, but rather to the more pronounced albedo reduction in openlands in a warming world. Furthermore, though not yet observed, rising temperatures may exceed a critical threshold that enhances cold-season ET in the future. This may strengthen the evaporative cooling effect of forests under future warming conditions, especially in coniferous forests with year-round transpiration activity.

Dominant role of increased soil dryness in daytime ΔLSTs dynamics
In contrast to the winter warming effect, European forests show significant summer daytime cooling compared to nearby openlands (Supplementary Fig. 2), driven by non-radiative mechanisms18. Over the past two decades, Europe has experienced a pronounced increase in summer air temperature (ATs), exacerbating atmospheric and soil dryness. Specifically, we find a significant upward trend in vapor pressure deficit (VPDs) and a concurrent decline in surface soil moisture (SSMs). Additionally, we show an increasing trend in summer downward solar radiation (DSRs) and decreasing trends in precipitation (Ps) and wind speeds (WSs) (Supplementary Fig. 10). We hypothesize that these climate changes may differentially influence turbulent fluxes in forests and openlands, thereby contributing to the observed temporal dynamics of ΔLSTs. Notably, the climate variables may themselves be influenced by the presence of forests (e.g., soil moisture retention and evaporative cooling of trees). However, it should be emphasized that the climatic drivers used in our analysis are the grid means at the coarse resolution, reflecting large-scale patterns rather than microclimates. Our analysis tends to capture the different responses of forests and openlands to background climate change.

Temporal correlation analyses using VPDs, ATs, DSRs, Ps and WSs as the potential drivers reveal strong relationships between daytime ΔLSTs and climatic variables except WSs (Fig. 3a, Supplementary Fig. 11). However, after accounting for covariates, the partial correlations between the daytime ΔLSTs and DSRs, Ps, and ATs are greatly reduced and even reverse in sign. In contrast, VPDs maintains a predominant negative correlation with daytime ΔLSTs, with the mean partial correlation coefficient of −0.42 (Fig. 3a). Notably, we can also find strong positive partial correlations between SSMs and ΔLSTs (Fig. 3b) when replacing VPDs with the SSMs, the other factor describing the dryness stress. This is due to the strong coupling between SSMs and VPDs at the annual scale (r=0.92, Supplementary Fig. 12). These findings raise the following question: whether variations in atmosphere or soil dryness are the primary factor governing changes in ΔLSTs over the past two decades?

To answer this question, we then use the binning regression method to further disentangle the respective contributions of VPDs and SSMs36. We conduct this analysis using the spatial samples of VPDs and SSMs trends with sufficient data for decoupling the two factors (Supplementary Fig. 13). We confirm the strong spatial correlation between SSMs and VPDs trends (r=−0.65) when all spatial samples are used. However, after stratifying either SSMs or VPDs trends into 10 quantile-based bins, correlations between the two variables show substantial attenuation (Fig. 3c), suggesting the two variables are decoupled. Then, we calculate the correlation between SSMs (or VPDs) trends and ΔLSTs trends in each VPDs (or SSMs) trend bin. We find the mean r between VPDs trends and ΔLSTs trends is only −0.11, whereas the mean r between SSMs trends and ΔLSTs trends reaches 0.40 (Fig. 3d). This result suggests that the change in SSMs is more likely to be the driver of ΔLSTs dynamics.

The critical role of SSMs is corroborated by non-linear statistical analysis using a Random Forest (RF) model (Supplementary Text 2). Variable importance analysis demonstrates that SSMs exerts a stronger influence than VPDs in predicting ΔLSTs dynamics. Contribution analysis also reveals that SSMs changes account for the largest portion (−0.152 K/decade, p<0.05) of the simulated negative ΔLSTs trend (−0.181 K/decade, p<0.05) (Supplementary Fig. 14), highlighting its dominant role in driving the observed enhanced forest cooling effect.

We also evaluate the sensitivity of summer latent heat (LEs) to VPDs in both forests and openlands using eddy covariance observations from European sites (Supplementary Text 3). The sensitivity analysis is based on the Penman-Monteith framework, which accounts for the influence of VPDs on surface resistance37. If VPDs is the main driver of enhanced forest cooling effect, we expect to observe greater positive LEs sensitivity in forests compared to openlands, in which case forests can produce stronger evaporative cooling effects in the context of rising VPDs. However, we find slightly higher positive sensitivity in openlands than forests, although the difference is not significant (Supplementary Table 1). This in-situ evidence, combined with our statistical results, suggests that VPDs may not explain the variation in daytime ΔLSTs, while the decline in SSMs appears to be the real driver.

We further conduct structural equation modeling (SEM) to elucidate the biophysical mechanisms underlying the impact of varying SSMs on daytime ΔLSTs dynamics (Fig. 4a). Our SEM incorporates two key biophysical mediators: the leaf area index difference (ΔLAIs) and the latent heat flux difference (ΔLEs) between nearby forests and openlands. The path effect (mean±standard error) is used to describe the strength and direction of the relationship. Following previous studies, we build SEMs at the grid scale and calculate the mean path effects of all models38,39. The analysis reveals three distinct pathways of SSMs affecting ΔLSTs:
(1) SSMs → ΔLAIs → ΔLSTs. SSMs exerts a significant negative effect on ΔLAIs (−0.41±0.15). Then, ΔLAIs demonstrates a negative effect on ΔLSTs (−0.19±0.09). This pathway suggests that reduced SSMs can amplify the greenness difference in forests and openlands, thereby enhancing forests’ cooling advantage through improved convective heat dissipation.
(2) SSMs → ΔLEs → ΔLSTs. SSMs shows a negative effect on ΔLEs (−0.25±0.11), followed by a strong negative influence of ΔLEs on ΔLSTs (−0.54±0.09). This pathway indicates that declining SSMs may impose less water limitation to forests compared to openlands, enabling forests to maintain stronger relative evaporative cooling through differential stomatal regulation (with greenness remaining constant).
(3) SSMs → ΔLAIs → ΔLEs → ΔLSTs. SSMs-induced change in ΔLAIs subsequently enhances ΔLEs (0.62±0.06), ultimately influencing ΔLSTs. This demonstrates that SSMs reduction amplifies evaporative cooling effects of forests by increasing the vegetation greenness contrast with openlands.

The pathway effects, representing the strength and direction of the mechanism, are quantified as the product of all constituent path coefficients (Fig. 4b). Results reveal that the second and third pathways (both mediated by ΔLEs) contribute equally to ΔLSTs variation (about 0.14±0.06), while the impact first pathway is weaker (0.08±0.05). These findings suggest SSMs can modulate ΔLSTs dynamics through multiple mechanisms, with forest evaporative cooling regulation emerging as the dominant process.

During summer, the vegetation greenness and ETrates are generally higher in forests than in adjacent openlands40. The negative effects of SSMs on ΔLAIs and ΔLEs revealed by the SEM analysis demonstrate that the decreased SSMs can amplify these greenness and ET differences between forests and openlands. Given the well-established constraint of soil moisture on vegetation growth41, we tend to attribute the negative effects of SSMs to forests’ stability and resilience to increasing soil dryness (i.e., lower sensitivity of forest greenness and ET to SSMs). Specifically, forests possess deeper root systems and enhanced soil moisture storage capacity, enabling them to sustain their cooling efficiency under climate change compared to openland vegetation.

Previous studies indicate that soil moisture depletion in Europe can exacerbate warming through land-atmosphere feedback mechanisms during heatwave events42,43. A comparative analysis of temperature responses across land cover types shows that openlands experience pronounced warming during long-lasting heatwaves, whereas forests exhibit moderate warming due to their conservative water-use strategies and sustained evaporative cooling44. Our results provide observational evidence that forests’ capacity to mitigate water deficit-induced warming persists not only during heatwave events but also across the interannual scales. However, it should be noted that our observational results represent the impact of soil moisture within a limited range under current climate regimes. If persistent warming or reduced precipitation drives soil moisture below a critical threshold, forests may lose their stability and the evaporative cooling capacity, thus failing to compensate for the warming effect induced by the dark canopy.

Evaluation of ΔLSTw and ΔLSTs trends in four CMIP6 models
To examine whether the observed trends of the biophysical LST effects of European forests can be reproduced by Earth System Models (ESMs), we use the output from the simulation of the historical period in the Coupled Model Intercomparison Project Phase 6 (CMIP6). The simulated effect of forests is isolated by a similar method for observational data, using differences between the forest and openland tiles (forest minus openland) within each model grid cell in four state-of-the-art ESMs that provide subgrid information (CESM2, CNRM-ESM2-1, GFDL-ESM4 and UKESM1-0-LL).

The multi-model mean suggests widespread declining trends of daily mean ΔLSTw (−0.068 K/decade, p<0.01) and ΔLSTs (−0.054 K/decade, p<0.01) over Europe during 1985–2014 (Fig. 5a, d). These trends align well with the observed ones (Fig. 1a, c), and the discrepancy in magnitudes can be partly attributed to differences in the period and spatial coverage of the observation and simulations (Supplementary Figs. 15 and 16). We should emphasize that the current models perform better in simulating negative trends in ΔLSTw than in ΔLSTs, as evidenced by the large intermodel differences in ΔLSTs trends simulated by four models. Specifically, the ΔLSTw trend values are similar across the four models (Fig. 5a). While for ΔLSTs, CESM2 simulates a significant declining trend of −0.133 K/decade, which dominates the multi-model mean signal. CNRM-ESM2-1, GFDL-ESM4 and UKESM1-0-LL simulate much weaker trends of −0.042, −0.018, and −0.027 K/decade, respectively (Fig. 5d). In addition, differences are found in the multi-year mean ΔLSTw or ΔLSTs from the four models (Fig. 5a, d), which can be attributed to the biases in surface energy partitioning and albedo responses to land cover changes45.

We further evaluate the model representation of biophysical processes behind the negative trends in ΔLSTw. Observational results indicate that the decreasing trend in the openland albedo is 2.10 times greater than that in the forest albedo, which is close to the multi-model mean of 1.94 (Fig. 5b). Moreover, we confirm that the observed slope (−5.51±1.64) between Δαw and ΔLSTw can be well reproduced in the multi-model mean (−5.57±1.39) (Fig. 5c). These findings demonstrate that the present model can produce reliable results in simulating the impact of snow dynamics on albedo and LST variations across vegetation types. However, considerable discrepancies exist across different models in these two relationships. For the summer results, we conduct partial correlation analysis between ΔLSTs and climate variables (Fig. 3a). The wind speed is excluded here considering the model data availability. Among all four models, only CNRM-ESM2-1 and UKESM1-0-LL capture the dominant role of SSMs on ΔLSTs dynamics (Fig. 5e). CESM, despite simulating the most pronounced negative trend in ΔLSTs, failed to reproduce the crucial role of SSMs. Overall, our comparison results suggest that while the multi-model ensemble can capture negative trends in ΔLSTw and ΔLSTs, discrepancies persist between observations and simulations (particularly in summer), along with intermodel differences in biophysical process representation. Our results are informative for improving current models for better simulations of local climate effects of European forests and their temporal dynamics in a warming world.

Methods
Satellite data
Satellite observations are adopted to estimate the potential biophysical effects of forestation. The satellite data used in this study are divided into two categories: MODIS products and other auxiliary land surface data. All satellite data are processed and output with a spatial resolution of 1 km.

LST data are obtained from the MYD11A1 dataset, which is retrieved by the split-window algorithm76. The MYD11 LST provides daily instantaneous observations at ~1:30 PM (daytime) and 1:30 AM (nighttime), corresponding to the daily maximum and minimum temperatures. Observations with LST errors greater than 2 K are filtered out according to the quality control flag. All high-quality observations are used to aggregate the summer (June, July, and August) and winter mean LST (December, January, and February) for both daytime and nighttime. We generate summer and winter LST from 2004 to 2023. It is worth noting that the winter LST for a given year is the mean value from January and February of the year and December of the previous year. For example, the 2023 winter LST is the mean value from December 2022 to February 2023. This aggregation principle is applied to all seasonal variables. We calculate the daily mean LST of summer and winter through the weighted combination of daytime and nighttime LST for both seasons80.

LE data are derived from the MOD16A2GF (MOD16A2 Gap-filled) dataset, with an original resolution of 8-day and 500 m. We select all observations retrieved by the main algorithm (Penman-Monteith model) and convert the LE unit to W/m². Then, all available records are spatially and temporally aggregated to 1 km seasonal means from 2004 to 2023. Leaf area index (LAI) data are obtained from the MCD15A3H dataset with a 4-day and 500 m resolution. All LAI values retrieved by the main Look-up-Table (LUT) method are aggregated to 1 km seasonal means from 2004 to 2023. The albedo data are from the MCD43A3 dataset, which contains daily 500 m black-sky and white-sky albedo for the shortwave band. We average the black-sky and white-sky albedo and then aggregate them to 1 km seasonal means from 2004 to 2023. Yearly land use/land cover data from 2004 to 2023 are derived from the MCD12Q1 dataset, with an original spatial resolution of 500 m. The International Geosphere-Biosphere Programme (IGBP) classification is used to identify forest and openland. The 500 m land cover maps are spatially aggregated to a 1 km scale based on the majority type.

Other auxiliary land surface data include forest loss data from the Global Forest Change (GFC) dataset81, impervious surface data from the Global Artificial Impervious Area (GAIA) dataset82, and digital elevation model (DEM) from the GMTED2010 dataset. The GFC dataset provides annual forest loss masks from 2000 to 2023 with 30 m spatial resolution based on time series analysis. All layers from 2004 to 2023 are composited, and a mask of forest loss during the study period is generated. The composited map is then spatially aggregated to 1 km resolution to match the MODIS data. The aggregated values range between 0 to 1, indicating the proportion of 30 m pixels where forest loss has occurred. GAIA data document the impervious surface mask from 1985 to 2018 with a spatial resolution of 30 m. Similar to the GFC datasets, the data from the last year (2018) are aggregated to a 1 km resolution impervious surface percentage map for further analysis. The DEM data of GMTED2010 has a spatial resolution of ~7.5 arc-seconds. We aggregate the original DEM to a 1 km spatial resolution for further analysis.

Reanalysis data
Meteorological variables are obtained from the monthly ERA5-Land reanalysis dataset83 with a spatial resolution of 0.1°. ERA5-Land dataset has been widely used for analysing the interaction between climate and terrestrial ecosystems. The used variables include 2m air temperature (AT), 2m dewpoint temperature (DT), precipitation, downward shortwave radiation, wind speed, surface soil moisture and snow cover. We aggregate these climatic variables into winter and summer means with a spatial resolution of 0.5°. Here, near-surface soil moisture data of layer 1 (0–7 cm depth) are used for analysis, as this part of the water can be utilized by herbaceous or woody vegetation, thereby containing the water limitation information for both forests and openlands. AT and DT are used to calculate the atmospheric vapor pressure deficit (VPD) using the following equations84:

VPD=e(AT) − e(DT) (1)

e(T)=6.112 × fw × exp(17.67 × T / (T+243.5)) (2)

fw=1.0007 + 3.46 × 10−6 × AP (3)

AP=1013.15 × ((T+273.16 − 0.0065 × Z) / (T+273.16))5.625 (4)

where AT and DT indicate the 2m air temperature and dewpoint temperature (Celsius), respectively; e(T) indicates the saturation water vapor (hPa) at a given temperature T (Celsius); and AP denotes the air pressure (hPa); Z indicates elevation (m).

Observed biophysical effects of forests on LST
This paper explores the temporal dynamics of forests’ biophysical processes under a changing climate. To ensure a sufficient sample size for robust results, we employ a space-for-time substitution strategy to estimate the a priori potential effect, rather than using the space-and-time approach for the a posteriori actual effect67,85. Specifically, we assume that all 1 km pixels within a 0.25° grid (comprising 30×30 pixels) share the same background climate12–14,28,86,87, and the potential effect of forests on the LST (ΔLST) can be estimated by the difference between the LST of forested areas (LSTf) and non-forest land cover (LSTo) within the grid:

ΔLST=LSTf − LSTo (5)

According to the MCD12Q1 land cover product, forest pixels should belong to evergreen needleleaf forests (ENFs), deciduous broadleaf forests (DBFs), evergreen broadleaf forests (EBFs), and mixed forests (MFs). Openland pixels include grasslands (GRA) and croplands (CRO). To minimize the potential impact of direct human activities, such as anthropogenic heat emissions and land use/land cover changes, on the ΔLST time series, we apply additional filtering criteria to the potential samples: (a) the land cover type remains unchanged during the study period; (b) the forest loss area should be less than 5% during the study period based on GFC data; and (c) the impervious surface percentage should be less than 5% according to GAIA data. For each 0.25° grid, we first determine the less dominant land cover type (forest or openland) based on the count of valid samples, then apply a ±100 m elevation threshold around its mean elevation to select comparable samples from the other type. We calculate ΔLST only if both the forest and openland sample sizes exceed 5 pixels. We also discard ΔLST calculations for grids where the mean elevation exceeds 500 m to avoid potential bias induced by mountainous terrain. The resultant mapping of the biophysical LST effect is spatially aggregated from the 0.25° scale to the 0.5° scale to match the meteorological data and mitigate the impact of potential outliers.

We calculate ΔLST for both winter and summer from 2004 to 2023. The criteria ensure that the forest and non-forest pixels remain consistent and stable during the study period, allowing us to attribute the temporal dynamics of ΔLST primarily to background climate change. By substituting the target variable, we also calculate the albedo, LE and LAI difference between forest and non-forest areas.

Simulated biophysical effects of forests on LST
To obtain the simulated biophysical effects of forests on LST over the historical period, we use the output from the historical simulation of the CMIP688. The historical simulation is forced by externally imposed and time-varying natural (e.g., solar variability and volcanic aerosols) and anthropogenic (e.g., greenhouse gases, aerosols, and land use/land cover changes) forcings. The last 30-year (i.e., 1985–2014) simulation output is used to calculate the biophysical effects of forests and their trends. Here, we choose the last 30 years period because it is close to our observational period (i.e., 2004–2023). Moreover, the 30-year time series can be considered long enough to minimize the influence of internal climate variabilities on the long-term trend.

A few models that participate in the historical simulation report the subgrid information following the data reporting protocol of the Land Use Model Intercomparison Project (LUMIP)89. Specifically, these models provide surface diagnostic variables outputs on four subgrid tiles, including primary and secondary land (psl, including trees, grasslands, barrens and vegetated wetlands), cropland (crp), pastureland (pst) and urban land, in each model grid cell. To create forest and openland comparisons within a grid cell for further analysis, we only retain the grid cells with the ratio of the tree cover fraction to the psl fraction exceeding 10%; this threshold is consistent with that used by the Food and Agriculture Organization to define forest cover90. For these selected grid cells, the psl tile can be considered as a forest tile. The crp and pst tiles are combined into an openland tile to keep consistency with the satellite observations. Since the forest and openland tiles in a grid cell are forced by the same atmospheric inputs (e.g., downward radiation, air temperature and precipitation), the background climate is the same between the two tiles within a grid cell. As such, differences in surface diagnostic variables (e.g., LST) between the two tiles are solely attributed to differences in land cover45,91. The local biophysical effects of forests on LST (ΔLST) in a grid cell can be expressed as follows45,91:

ΔLST=LSTpsl − LSTcrp=pst (6)

where LSTpsl denotes the forest LST value (psl tile), and LSTcrp=pst indicates the openland LST value (crp and pst combination). Specifically, if the outputs for both the crp and pst tiles are available, their arithmetic mean value is used; otherwise, the available one is used. LST can be replaced by other surface diagnostic variables (e.g., albedo) to quantify the biophysical effects of forests on corresponding variables.

It should be emphasized that the forest effects identified by this method should be interpreted as the result of a complete conversion from openlands to forests in a grid cell. In other words, the identified forest effects are insensitive to the area fraction of the forest or openland tile in a grid cell. Therefore, while historical land use and land cover changes are prevalent in the simulations, they do not influence the identified forest effects.

Four models (CESM2, CNRM-ESM2-1, GFDL-ESM4 and UKESM1-0-LL) are used for the subgrid analysis, as they provide the required output for the historical simulation. CESM2 and GFDL-ESM4 each have only one member that provides the subgrid output, and the available one (r10i1p1f1 for CESM2 and r1i1p1f1 for GFDL-ESM4) is used. CNRM-ESM2-1 and UKESM1-0-LL each have more than one member, and the first one (r1i1p1f2 for CNRM-ESM2-1 and r2i1p1f2 for UKESM1-0-LL) is used for the sake of keeping consistency with other models as well as saving computational cost.

Statistical analyses
We use the Theil-Sen estimator and Mann-Kendall test to analyse temporal trends92 (scipy and pymannkendall packages in Python). The Theil-Sen estimator is a non-parametric method for robust estimation of the slope by using the median value of a range of possible slopes. The Mann–Kendall test, which is also nonparametric, statistically assesses whether a temporal trend is monotonically upward or downward.

Correlation analyses are employed to examine the potential influence of given drivers on the target variables (numpy package in Python). These analyses are based on either spatial or temporal samples. Temporal correlation is conducted at both pixel and regional scales. For pixel-scale correlation, we calculate the spatial mean and standard deviation of correlation coefficients to represent regional results. For regional-scale correlation, we first calculate the mean value for the study area and then determine the correlation coefficient. In addition, partial temporal correlation analysis is used to explore the relative impacts of multiple climate drivers on ΔLSTs at the pixel scale (pingouin package in Python). The partial correlation can determine whether a climate variable exerts an independent effect on the target variable beyond its covariation with other climate drivers. For instance, the partial correlation can isolate the association between ΔLSTs and SSMs while accounting for the influence of other covariates (e.g., temperature and precipitation).

We also conduct spatial correlation analyses using trend values to examine whether long-term changes in the independent variables can induce the trend in the target variable. These analyses are conducted using either complete spatial samples or binned subsets of the data. The latter, referred to as the binning regression method, is used to disentangle the contribution of covarying VPDs and SSMs on ΔLSTs. Specifically, we first stratify either SSMs or VPDs trends into 10 percentile bins (0–10th, 10th–20th, …, 90th–100th) and then assess their correlation (Fig. 3c). Subsequently, we evaluate the correlation between SSMs (or VPDs) trends and ΔLSTs trends within each corresponding VPDs (or SSMs) trend bin (Fig. 3d). This binning approach can effectively decouple the two factors, allowing us to identify the real dominant driver by comparing the resulting correlation coefficient.

We use SEM to reveal the biophysical mechanism underlying the ΔLSTs dynamics (semopy package in Python). The SEM is designed based on previous studies93,94. Specifically, we assume that climate variables influence the temporal dynamics of ΔLSTs through three pathways: (P1) climate change affects ΔLSTs via ΔLAIs; (P2) climate change affects ΔLSTs via ΔLEs; and (P3) climate change affects ΔLAIs, which further alters ΔLEs and then ΔLSTs. For P1 and P3, the impact of climate variables on ΔLAIs arises from the diverse sensitivities of vegetation to climate change. Specifically, climate changes can differently promote or inhibit forest and non-forest vegetation growth, which further affects ΔLSTs through convective heat dissipation (P1) and the evaporative cooling effect (P3). Climate variables also directly modulate ΔLEs even when the vegetation greenness remains constant. This is due to inherent differences between different vegetation types44, resulting in diverse evaporative responses to climate change and consequently affecting ΔLSTs (P2). It should be noted that all drivers used in SEM (SSMs, ATs, DSRs, WSs and Ps) represent grid-level mean values, which primarily reflect background climatic and environmental conditions. These variables are unlikely to be influenced by internal variations within the grids (ΔLEs and ΔLAIs).

In practice, SEM is performed at the grid level. We use all the data within a 3×3 spatial moving window for each 0.5° grid to ensure a sufficient sample size for establishing a robust model. Notably, we first normalize the input data and ΔLSTs by subtracting the 20-year means for each grid to exclude the impact of spatial variability. SEM is conducted in those windows where the sample size is larger than 120. The data are normalized using the Z-score method before modeling. Following the previous study38, the goodness-of-fit index (GFI) is employed to evaluate the applicability and effectiveness of each model. We retain the models with GFI > 0.9 and calculate the mean path values along with the standard errors. The contributions of pathways P1, P2, and P3 are measured by their path effects, which are calculated by multiplying the contained path values. Uncertainties are estimated based on error propagation (error_propagation package in Python).

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "This manuscript is well-structured and well-written. The authors define ΔLST as the difference between forest LST and\nnearby open land LST. Their analysis reveals that ΔLST has intensified in recent decades, which they attribute to increasing\nVPD, favoring enhanced evaporative cooling in forests. I feel this manuscript is publishable if the authors can address the\nfollowing questions:\nLine 14: VPD increase will trigger stomatal closure, thus reducing ET cooling. How to reconcile this?\nLines 101-110: The observed forest cooling - could this be influenced by elevation effects? In the Methods section (Line\n577), is the assumption that the background climate remains uniform within a 0.25° × 0.25° area valid? I noticed that the\nauthors set a threshold of a mean elevation of 500 m to mitigate this bias. However, given that the typical lapse rate is\napproximately 6.5°C per 1,000 meters, this effect may be substantial enough to alter the direction of the temperature\ndifference (compare to the numbers you provided in Lines 101 to 110). In low-elevation hilly regions, forests may often be\nfound at higher elevations, while open land, if resulting from deforestation, tends to be at lower elevations due to easier\naccessibility. This elevation difference may contribute to the observed forest cooling effect.\nLine 161: In Fig 2c, the relationship is binned by the trend of snow cover. Can the authors reproduce the result if binned by\nthe trend of albedo?\nLine 165 – 168: not sure if the snow mass between forest and open land also can be affected by the elevation difference?\nLine 201 – 203: It is interesting to see a negative correlation between VPD and ΔLST_s. This suggests that when VPD is\nhigh, forests exhibit greater energy dissipation advantages compared to open lands. As the authors mentioned, VPD is\nhighly positively correlated to air temperature. So, this reminds me of an earlier heatwave study “Teuling, A., Seneviratne, S.,\nStöckli, R. et al. Contrasting response of European forest and grassland energy exchange to heatwaves. Nature Geosci 3,\n722–727 (2010). https://doi.org/10.1038/ngeo950”. Teuling et al find “initially, surface heating is twice as high over forest\nthan over grassland. Over grass, heating is suppressed by increased evaporation in response to increased solar radiation\nand temperature. Ultimately, however, this process accelerates soil moisture depletion and induces a critical shift in the\nregional climate system that leads to increased heating. … We conclude that the conservative water use of forest contributes\nto increased temperatures in the short term, but mitigates the impact of the most extreme heat and/or long-lasting events.” My\nquestion is whether the analysis in this Nature Communications manuscript represents short-term effects or longer-term\neffects?\nLine 272 - 275: The logic in this section feels somewhat unclear. The authors claim that higher VPD enhances ΔLAI_s and\nsuggest multiple similar cause-and-effect relationships with other variables. How can the causal inferences be reliably\nestablished? It seems that the reverse statement is also valid.\nIn the flowchart of Fig. 4, I notice a missing of surface resistance in the VPD-LE relationship. In the Penman-Monteith\nscheme, increasing VPD can have competing effects on ET: on one hand, increasing VPD enhances ET, while on the other,\nit also raises surface resistance, which suppresses ET. The results suggest that the VPD-driven enhancement of ET is more\ndominant. Could the authors clarify why this occurs when comparing ET between forests and open lands? I feel VPD is not\nthe driver of such findings, but the biome type - forests/open lands and their associated biophysical and physiological\nresponses - are the real drivers. Need to clarify this.\nSee Novick, K., Ficklin, D., Stoy, P. et al. The increasing importance of atmospheric demand for ecosystem water and\ncarbon fluxes. Nature Clim Change 6, 1023–1027 (2016). https://doi.org/10.1038/nclimate3114", "reviewer_2": "I find this manuscript has a lot of interesting data and discussion. Moreover, for me personally, it was highly motivating to\nread this discussion. I did find, especially at the beginning of the manuscript, that the work was not very well situated in the\nprincipal literature. All four of the first references to appear in the manuscript do not seem well suited to set up this general\ndiscussion. I have provided a number of suggestions in the attached submission pdf and am hopeful that some of these\nsuggestions will prove useful...\nI also found that the introduction was occasionally vague and would have benefitted from a somewhat more nuanced\ndiscussion... I have tried to point out places where I think this would be helpful.\nFinally, the discussion of radiative and non-radiative effects, as well as the discussion of winter and summertime effects, to\nmy mind, were not adequately introduced and motivated in the introductory parts of the paper. When these concepts are\nsuddenly announced, (the discussion about winter and summer time effects first comes up in the results section of the\npaper?), they are unexpected and seemingly unmotivated in the general context of the paper... This can be improved upon...\nI also find that the Authors do not adequately distinguish between global and local cooling... Both terms come up. But most\nof the time it sounds like the Authors are talking about local cooling... And the balance in terms of global cooling is never\nclearly estimated...? It is important not to conflate these terms, since local cooling is not necessarily equivalent to\nglobal cooling...?\nFurther, I would have like to see more discussion of the role of the availability of water and surface moisture. The cooling\npotential of forests strongly depends on the availability of water. As soon as this surface moisture disappears, the loss of\nlatent heat means the warming will suddenly take over... But there is little attention to this throughout the paper...\nFinally, I wonder why LEw never appears anywhere, but the effects of LEs are estimated... Clearly LEw will be significantly\nlower. But conifers can evapotranspire almost year round... So, to some degree, LEw should also be available and should\nhave an impact on surface temperatures... ???\nMore to the main point, however, the general finding that wintertime effects exhibit increased cooling due to the loss of snow\ncover, more or less independently of tree cover (the effect is stronger over open fields than over forested land), is an\ninteresting and valuable finding. However, it remains unclear why this phenomenon would be true in more forested areas...\nIn other words, why would the land surface become lighter under more forested areas, thus leading to wintertime\ncooling...??? I do not think this is adequately explained...?\nAll-in-all however, I did find the manuscript to be very interesting and compelling... I have included a number of comments in\nthe attached submissions pdf."}, "author_response": "Response to Reviewers’ Comments\nWe greatly appreciate the opportunity to revise our manuscript and thank all anonymous reviewers\nfor their constructive comments. The manuscript has been revised and improved accordingly. We\nhope that the revision will make it more acceptable for publication. Below are the point-by-point\nresponses to the comments, along with the revision of the manuscript (typed in red and Italic). The\nchanges are also marked in red in the manuscript file.\nReviewer #1 (Remarks to the Author):\nThis manuscript is well-structured and well-written. The authors define ΔLST as the difference\nbetween forest LST and nearby open land LST. Their analysis reveals that ΔLST has intensified\nin recent decades, which they attribute to increasing VPD, favoring enhanced evaporative\ncooling in forests. I feel this manuscript is publishable if the authors can address the following\nquestions:\nResponse: We appreciate the positive comments by the reviewer. We have carefully considered all\nthe comments and suggestions and made corresponding point-by-point responses. To address your\ncomment regarding the potential influence of elevation bias between forests and openlands, we have\nrefined our results by setting additional elevation restrictions when selecting pixels for comparison.\nAlso inspired by your comment regarding the short-term and long-term effects, we added new\nexperiments to explored the potential impact of soil moisture in summer (SM). Accordingly, we\ns\nupdated the results and refined the research findings. These refinements have significantly\nstrengthened the robustness of our findings, as detailed in our point-by-point responses below.\nLine 14: VPD increase will trigger stomatal closure, thus reducing ET cooling. How to reconcile\nthis?\nResponse: We appreciate your valuable comment about the impact of VPD. Based on your comment\nbelow, we have performed additional analyses to elucidate the mechanisms of enhanced forests\ncooling effects in summer (ΔLST). Our results demonstrates that while summer VPD (VPD)\ns s\nexhibits a strong correlation with ΔLST, increased atmospheric dryness does not appear to be the\ns\nactual factor driving the observed ΔLST patterns. Instead, the concurrent summer SM (SM)\ns s\ndepletion (tight coupling with VPD), should be the primary driver of the trends in ΔLST.\ns s\nDetailed\nanalyses supporting this conclusion are presented in our responses to your subsequent comments.\nTherefore, we have revised the relevant text in the abstract to more accurately reflect this\nmechanism (Line 12):\nIn the summertime, the vegetation physiological response to decreased soil moisture boosts the\ncooling effect of forests.\nLines 101-110: The observed forest cooling - could this be influenced by elevation effects? In\nthe Methods section (Line 577), is the assumption that the background climate remains\nuniform within a 0.25° × 0.25° area valid? I noticed that the authors set a threshold of a\nmean elevation of 500 m to mitigate this bias. However, given that the typical lapse rate is\napproximately 6.5°C per 1,000 meters, this effect may be substantial enough to alter the\ndirection of the temperature difference (compare to the numbers you provided in Lines 101\nto 110). In low-elevation hilly regions, forests may often be found at higher elevations, while\nopen land, if resulting from deforestation, tends to be at lower elevations due to easier\naccessibility. This elevation difference may contribute to the observed forest cooling effect.\nResponse: Thanks for your insightful comments regarding the impact of elevation. We fully agree\nthat elevation differences could have a potential impact on the evaluated LST effect of forests. We\nchecked our results and found a mean elevation difference of ~150 m between forest and grassland\nsamples under our original 500 m mean elevation constrain. This may lead to an overestimation of\nthe forest cooling effects, as you noted.\nTo mitigate this, we implemented additional screening criteria. First, we identified the less\ndominant land cover type (either forest or openland) within each 0.25° × 0.25° grid, then applied\na ±100 m elevation threshold around its mean elevation to select comparable samples from the\nother land cover type. Only when both forest and openland samples are larger than 5, the calculation\nof ΔLST was performed. In the revised manuscript, we revised corresponding methodological\ndescriptions (Line 618):\nFor each 0.25° grid, we first determine the less dominant land cover type (forest or openland)\nbased on the count of valid samples, then apply a ±100 m elevation threshold around its mean\nelevation to select comparable samples from the other type.\nThis refined approach reduced mean elevation differences to just ~ 20 m (vs. original ~150 m).\nBased on the refined method, we have updated the all the related numbers and figures. Given the\nextensive revisions made, we do not list all modifications in the response letter. Please refer to the\nrevised manuscript for detailed changes. However, we confirm that the original key finding,\nnegative trends in ΔLST and ΔLST , persisted (Fig. 1).\ns w\nFig 1. Temporal dynamics in winter and summer land surface temperature effects of European\nforests (ΔLST and ΔLST) in recent two decades.\nw s\nWe also appreciate your comment about the background climate assumption. While perfect\nclimate uniformity cannot be assumed at any spatial scale, our 0.25° grid setting was based on\nseveral considerations.\nFirst, previous research has shown that the decorrelation length of surface climate, the spatial\nscale at which these variables (e.g., temperature, precipitation) become relatively independent,\ntypically exceeds 30 km (Kim and Kosro, 2013; Larsen et al., 2013). This conclusion supports our\nassumption that sampled locations within each 0.25° (~30 km) grid cell share similar climatic\nconditions.\nSecond, the 0.25° window size is comparable with that in previous studies exploring the\nbiophysical effect of vegetation changes (Table R1). We have added these references in Line 604.\nSeveral studies also confirmed that the estimated LST effect is robust to the window size selection\n(Li et al., 2015; Zhang et al., 2024; Zhao et al., 2024).\nTable R1. Window size setting for the same background climate assumption in previous studies\nReference Window size Data\n(Lee et al., 2011) 28 km (mean) Site observation\n(Li et al., 2015) 50 × 28 km Remote sensing observation\n(Alkama and Cescatti, 2016) 50 km Remote sensing observation\n(Duveiller et al., 2018) 25 km Remote sensing observation\n(Chen et al., 2018) 21.6 km (mean) Site observation\n(Barnes et al., 2024) 30 km (maximum) Site observation\nThird, we also explored the background climate similarity using 43 paired and synchronous\nFLUXNET sites observations. The spatial distances of the paired sites range from 0 to 150 km, and\nthe elevation differences are less than 200 m. Here, we compared the air temperature (AT),\nprecipitation (P) and downward shortwave radiation (DSR) at different scales. For each variable,\nwe compared the annual means (Fig. R1) and the monthly ranges (the difference between the\nmonthly maximum and minimum values (Fig. R2). We also calculated the correlation coefficient\nbetween observations at the daily scale (Fig. R3).\nFig. R1. Comparison of annual temperature, precipitation and radiation means from spatial paired\nsite observations\nFig. R2. Comparison of monthly temperature, precipitation and radiation ranges from spatial\ntemperature, precipitation and radiation ranges from spatial\npairing site observations.\nFig. R3. Correlation of daily temperature precipitation and radiation between spatial paired site\nobservations\nFig. R1 shows that the paired sites within 30 km share similar annual mean AT, P and DSR,\nand the scatters are close to the 1:1 line. The average root mean square error (RMSE) for paired\nsites with distances less than 30 km is significantly smaller than for samples at distances greater\nthan 30 km, especially for precipitation (105 mm/y vs. 229 mm/y). Similar results are also found in\nmonthly ranges, which describe the seasonal variations of the climate variables (Fig. R2). For the\ndaily variation, the correlations between the paired site observations decrease with increasing\ndistance. Especially for precipitation, the mean correlation for paired sites with distance less than\n30 km can reach 0.81, while for those paired sites with distance larger than 30 km, the mean\ncorrelation is only 0.57. Overall, these evidences support the validity of the similar background\nclimate at the 0.25° scale.\nReference:\nAlkama, R. & Cescatti, A. Biophysical climate impacts of recent changes in global forest cover.\nScience (80-. ). 351, 600–604 (2016).\nBarnes, M. L. et al. A Century of Reforestation Reduced Anthropogenic Warming in the Eastern\nUnited States. Earth’s Futur. 12, (2024).\nChen, L., Dirmeyer, P. A., Guo, Z. & Schultz, N. M. Pairing FLUXNET sites to validate model\nrepresentations of land-use/land-cover change. Hydrol. Earth Syst. Sci. 22, 111–125 (2018).\nDuveiller, G., Hooker, J. & Cescatti, A. The mark of vegetation change on Earth’s surface energy\nbalance. Nat. Commun. 9, 64–75 (2018).\nKim, S. Y. & Kosro, P. M. Observations of near-inertial surface currents off Oregon: Decorrelation\ntime and length scales. J. Geophys. Res. Ocean. 118, 3723–3736 (2013).\nLarsen, M. A. D., Thejll, P., Christensen, J. H., Refsgaard, J. C. & Jensen, K. H. On the role of\ndomain size and resolution in the simulations with the HIRHAM region climate model. Clim.\nDyn. 40, 2903–2918 (2013).\nLee, X. et al. Observed increase in local cooling effect of deforestation at higher latitudes. Nature\n479, 384–387 (2011).\nLi, Y. et al. Local cooling and warming effects of forests based on satellite observations. Nat.\nCommun. 6, 1–10 (2015).\nZhang, Y. et al. Asymmetric impacts of forest gain and loss on tropical land surface temperature.\nNat. Geosci. 13, 823–831 (2024).\nZhao, J. et al. Forest fire size amplifies postfire land surface warming. Nature 633, 828–834 (2024).\nLine 161: In Fig 2c, the relationship is binned by the trend of snow cover. Can the authors\nreproduce the result if binned by the trend of albedo?\nResponse: Thank you for raising this concern. Fig. 2c (also Fig. R4a) reveals that surface darkening\n(α trends) is more sensitive to snow cover decline (SC trends) in openlands than in forests, as\nw w\nreflected by the steeper spatial regression slope for openlands.\nTo further validate these findings, we conducted additional analyses using α trend-based data\nw\nbinning as you suggested. The results presented in Fig. R4b confirm the robustness of the initial\ncorrelation under this alternative binning approach. We further verified that substituting SC with\nw\nΔα also yields a robust temporal correlation with ΔLST (Fig. R4c and R4d).\nw w\nFig. R4. Comparison of the results binned by SC and α . (a) Spatial relationships between SC\nw w w\ntrends and α trends binned by SC trends (same as Fig. 2c in the manuscript). (b) Binned by α\nw w w\ntrends. (c) Temporal relationships between SC and ΔLST (same as Fig. 2e in the manuscript).\nw w\n(b) Temporal relationships between Δα and ΔLST .\nw w\nLine 165 – 168: not sure if the snow mass between forest and open land also can be affected\nby the elevation difference?\nResponse: Thank you for raising the concern about snow mass. We agree that forests typically occur\nat higher elevations where lower temperatures may enhance snow retention compared to openlands.\nTo mitigate this impact, we have added additional constraint of elevation differences between\nsamples (please see the response to your comments about Lines 101-110). Based on the refined\nmethod, we have updated the all the related numbers and figures.\nWe found that the multi-year mean daytime ΔLST shifted from -0.20 K to 0.01 K after\nw\naccounting for elevation differences, while the negative temporal trend in ΔLST persists. These\nw\nresults suggest elevation induced snow differences may influence absolute value of ΔLST , but do\nw\nnot affect the key findings of the observed temporal trends.\nLine 201 – 203: It is interesting to see a negative correlation between VPD and ΔLST_s. This\nsuggests that when VPD is high, forests exhibit greater energy dissipation advantages\ncompared to open lands. As the authors mentioned, VPD is highly positively correlated to air\ntemperature. So, this reminds me of an earlier heatwave study “Teuling, A., Seneviratne, S.,\nStöckli, R. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci 3, 722–727 (2010). https://doi.org/10.1038/ngeo950”.\nTeuling et al find “initially, surface heating is twice as high over forest than over grassland.\nOver grass, heating is suppressed by increased evaporation in response to increased solar\nradiation and temperature. Ultimately, however, this process accelerates soil moisture\ndepletion and induces a critical shift in the regional climate system that leads to increased\nheating. … We conclude that the conservative water use of forest contributes to increased\ntemperatures in the short term, but mitigates the impact of the most extreme heat and/or\nlong-lasting events.” My question is whether the analysis in this Nature Communications\nmanuscript represents short-term effects or longer-term effects?\nResponse: Thank you for raising this concern. We have carefully read the reference (Teuling et al.,\n2010), and find the revealed mechanism very insightful. Teuling et al. used in-situ observations to\ndemonstrate that during a heatwave event, rising temperature (also rising VPD) and SM depletion\ndominate the turbulent fluxes of forests and grasslands at distinct stages: temperature (VPD)\ndominates initially, with grassland ET showing a stronger positive response; while SM becomes\ndominant later as water supply limits grassland ET cooling, and forests’ conservative water-use\nstrategy leads to lower temperature. Here, the temporal scale in Teuling et al.’s study focuses on\nheatwave events, while our manuscript evaluates interannual variations and trends. Therefore, we\nconsider our analysis to represent the long-term effect.\nIn the original version of our manuscript, we supposed that different ET responses to rising\ntemperature (also rising VPD) dominate the observed enhanced forest cooling. Specifically,\nincreased VPD may enhance forest ET more than grassland ET, thereby strengthening forest cooling.\nThis hypothesis appears contradictory to the findings of Teuling et al. at the first stage, which may\nbe due to the differences in the temporal scale.\nHowever, through reviewing Teuling et al.’s study, we recognized that SM can also influence\nthe cooling effect of forests, and this aspect was not sufficiently addressed in our original manuscript.\nWe also note that VPD and SM are strongly coupled at the annual scale. Thus, we conducted\nadditional analyses to identify the actual driver of the observed trends in ΔLST.\ns\nWe used binning\nregression method (Liu et al., 2020) to disentangle the impacts of VPD and SM.\ns\nThe results\ndemonstrated that SM shows stronger correlation with daytime ΔLSTₛ than VPD (Line 218):\ns s\nTemporal correlation analyses using VPD, AT, DSR, P and WS as the potential drivers\ns s s s s\nreveal strong relationships between daytime ΔLST and climatic variables except WS (Fig. 3a).\ns s\nHowever, after accounting for covariates, the partial correlations between the daytime ΔLST and\ns\nDSR, P, and AT are greatly reduced and even reverse in sign. In contrast, VPD maintain a\ns s s s\npredominant negative correlation with daytime ΔLST, with the mean partial correlation coefficient\ns\nof -0.42 (Fig. 3a). Notably, we can also find strong positive partial correlations between SM and\ns\nΔLST (Fig. 3b) when replacing VPD by the SM, the other factor describing the dryness stress.\ns s s\nThis is due to the strong coupling between SM and VPD at the annual scale (r = 0.92,\ns s\nSupplementary Fig. 10). These findings raise the following question: whether variations in\natmosphere or soil dryness are the primary factor governing changes in ΔLST over the past two\ns\ndecades?\nTo answer this question, we then use the binning regression method to further disentangle the\nrespective contributions of VPD and SM (Liu et al., 2020). We conduct this analysis using the\ns s\nspatial samples of VPD\ns\ntrend and SM\ns\ntrends with sufficient data for decoupling the two factors.\nWe confirm the strong spatial correlation between SM trend and VPD trends (r = -0.65) when all\ns s\nspatial samples are used. However, after stratifying either SM or VPD trends into 10 bins\ns s\naccording to percentiles, correlations between the two variables show substantial attenuation (Fig.\n3c), suggesting the two variables are decoupled. Then, we calculate the correlation between SM\ns\n(or VPD) trends and ΔLST trends in each VPD (or SM) trend bin. We find the mean r between\ns s s s\nVPD trends and ΔLST trends is only -0.11, whereas the mean r between SM trends and ΔLST\ns s s s\ntrends reach 0.40 (Fig. 3d). This binning regression analysis suggests that change SM is more likely\ns\nto be the driver of ΔLST dynamic.\ns\nFig 3. Dominant role of soil moisture in the temporal dynamics in the summer daytime cooling effect\n(ΔLST) of forests. (a) Regional mean correlation and partial correlation coefficients between five\ns\nclimate variables (vapor pressure deficit, downward solar radiation, precipitation, air temperature\nand wind speed) and the daytime ΔLST. The error bars indicate the standard deviations. (b) Similar\ns\nwith (a), but replace vapor pressure deficit by the soil moisture. (c) Correlation coefficients between\nsummer vapor pressure deficit (VPD) trends and summer soil moisture (SM) trends within either\ns s\nSM or VPD trends bins. (d) Correlation coefficients between VPD (or SM) trends and ΔLST\ns s s s s\ntrends within SM (or VPD) trends bins.\ns s\nFig. S10. Temporal correlation coefficient (r) between summer vapor pressure deficit (VPD) and\ns\nsurface soil moisture (SM). (a) Spatial pattern of r between VPD and SM. (b) Density of the r\ns s s\nvalues.\nWe also build a non-linear random forest model using both SM and VPD as the input. The\ns s\nresults also support that SM plays the more important role in predicting temporal dynamics of\ns\nΔLST (Line 255):\ns\nThe critical role of SM is corroborated by non-linear statistical analysis using a Random\ns\nForest (RF) model (Supplementary Text 1). Variable importance analysis demonstrates that SM\ns\nexerts a stronger influence than VPD in predicting ΔLST dynamics. Contribution analysis also\ns s\nreveals that SM changes account for the largest portion (0.152 K/decade, p < 0.05) of the simulated\ns\nnegative ΔLST trend (0.181 K/decade, p < 0.05), highlighting its dominant role in driving the\ns\nobserved enhanced forest cooling effect (Supplementary Fig. 11).\nSupplementary Text 1 Modelling ΔLST based on random forest\ns\nWe apply a non-linear machine learning method, Random Forest (RF), to evaluate the\nimportance and contribution of multiple potential climate drivers for predicting ΔLST (sklearn\ns\npackage in Python). Specifically, we use all spatial samples over 20 years to build the dataset. It is\nworth mentioning that we normalize the input climate variables and ΔLST by subtracting the 20-\ns\nyear means at the grid scale, since the model focus on the temporal dynamic of ΔLST rather than\ns\nthe spatial variation. We excluded air temperature from the model inputs due to its high correlation\nwith other climate variables.\nThe dataset is divided into training (80% of the samples) and test sets (20% of the samples),\nwith model parameters confirmed based on validation results from the test set (number of trees =\n100; maximum depth of the tree = 6). All the samples are then used to fit the best model. The relative\nimportance of each variable is indicated by the Gini coefficient. Moreover, we use the Shapley\nAdditive Explanations (SHAP) method to quantify the individual contribution of each variable in\neach sample (shap package in Python). SHAP is an approach to explain machine learning models\nbased on cooperative game theory, where SHAP values quantify the marginal contributions of\npredictors. We calculate the mean absolute SHAP values for each input variable as importance\nmetrics and explore the non-linear impact of climate variables on ΔLST by plotting their SHAP\ns\nvalues against input values.\nBased on the RF model, we set several scenarios to separate the relative contributions of\nclimate variables to the ΔLST trend. Specifically, we perform four experimental simulations: (S1)\ns\nvarying SM only; (S2) varying SM and VPD, while fixing DSR P and WS; (S3) varying SM,\ns s s s s s s\nVPD and DSR, while fixing P and WS; (S4) varying SM, VPD, DSR and P, while fixing WS;\ns s s s s s s s s\nand (S5) varying SM, VPD, DSR, P, and WS. The Theil-Sen slope and Mann-Kendall P-value of\ns s s s s\nsimulations S1, S2−S1, S3−S2, S4−S3, S5−S4 are used to evaluate the individual contributions of\nthe changing climate to the ΔLST trend.\ns\nThe model shows good accuracy on the test dataset (R2=0.56, RMSE=0.32 K), supporting our\nfurther analysis. The Gini importance and mean absolute SHAP values both suggest that SM is the\ns\nmost important driver, followed by VPD, P, DSR and WS (Supplementary Fig 11a). The marginal\ns s s s\ncontributions quantified by the SHAP value further reveal different impacts of SM and VPD\ns s\n(Supplementary Fig. 11b). SM exhibits the pronounce positive effect, with stronger impact when\ns\nSM anomaly is negative. While VPD only show slight negative impact under certain condition.\ns s\nOn the basis of the RF model, we estimate the contributions of climate variables to the long-\nterm trend of daytime ΔLST (Supplementary Fig. 11c). We first compared the reconstructed trends\ns\nwith all forcings and the observed trends, showing good consistency (observed: 0.188 K/decade,\np<0.05; simulated: 0.181 K/decade, p<0.05). By setting different scenarios, results suggest that\ndecreasing SM contributes the most to the negative trend in daytime ΔLST (-0.152 K/decade,\ns s\np<0.05), and increasing VPD shows lower contribution (-0.049 K/decade, p<0.05). The isolated\ns\ncontribution of DSR is also statistically significant, but the value is relatively small (0.012 K/decade,\ns\np<0.05). The contribution of WS and P are not significant.\ns s\nFig. S11. Impact of SM on the temporal dynamics in ΔLST based on the random forest model. (a)\ns s\nThe importance of climate variables is indicated by the absolute SHAP values and Gini coefficients.\n(b) Variation of SHAP values with changing SM or VPD anomaly. (c) Observed and modelled\ns s\n\nΔLST trends, and the isolated contributions from climate variables.\ns\nIn summary, we conclude that the SM depletion (similar with the second stage identified by\nTeuling et al.) is more likely to explain the observed intensification of summer forest cooling over\nEurope in recent decades.\nReference:\nLiu, L. et al. Soil moisture dominates dryness stress on ecosystem production globally. Nat.\nCommun. 11, 1–9 (2020).\nTeuling, A. J. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci. 3, 722–727 (2010).\nLine 272 - 275: The logic in this section feels somewhat unclear. The authors claim that higher\nVPD enhances ΔLAI_s and suggest multiple similar cause-and-effect relationships with other\nvariables. How can the causal inferences be reliably established? It seems that the reverse\nstatement is also valid.\nResponse: Thank you raising this concern about the causal relationship. We fully acknowledge that\nestablishing robust causal inferences is critical and fundamental in structural equation modeling\n(SEM) analysis. Our model design is based on established understanding of the biophysical\nprocesses of forests that influence local temperature and supported by previous research findings.\nFirst, we design the model by assuming that climate variables (temperature, precipitation…)\nexert different impacts on the vegetation growth and evapotranspiration (i.e. ΔLAI and ΔLE)\ns s\n(Zheng et al., 2022). Such differences in LAI or LE between forests and openlands are unlikely to\nalter background climate patterns. Additionally, the adopted space-for-time approach assumes that\nthe background climate remains constant within the window. This implies that changes in LAI and\nLE are solely driven by climatic factors and do not exert feedback effects on the atmosphere. Second,\nwe assume that ΔLAI and ΔLE can impact ΔLST, which has been proved by the previous research\ns s s\n(Prevedello et al., 2019).\nWe have clarified the logic of SEM in the revised manuscript (Lines 715) and included\nadditional references to support these assumptions:\nWe use SEM to reveal the biophysical mechanism underlying the ΔLST dynamics (semopy\ns\npackage in Python). The SEM is designed based on previous studies (Prevedello et al., 2019; Zheng\net al., 2022). Specifically, we assume that climate variables influence the temporal dynamics of\nΔLST through three pathways: (P1) climate change affects ΔLST via ΔLAI;\ns s s\n(P2) climate change\naffects ΔLST via ΔLE; and (P3) climate change affects ΔLAI, which further alters ΔLE and then\ns s s s\nΔLST. For P1 and P3, the impact of climate variables on ΔLAI arises from the diverse sensitivities\ns s\nof different vegetation types to climate changes. Specifically, climate changes can differently\npromote or inhibit forest and non-forest vegetation growth, which further affects ΔLST through\ns\nconvective heat dissipation (P1) and the evaporative cooling effect (P3). Climate variables also\ndirectly modulate ΔLE even the vegetation greenness remain constant. This is due to inherent\ns\ndifferences between different vegetation types (Teuling et al., 2010), resulting in diverse evaporative\nresponses to climate change and, consequently, affecting ΔLST (P2).\ns\nReference:\nPrevedello, J. A., Winck, G. R., Weber, M. M., Nichols, E. & Sinervo, B. Impacts of forestation and\ndeforestation on local temperature across the globe. PLoS One 14, 1–18 (2019).\nTeuling, A. J. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci. 3, 722–727 (2010).\nZheng, H. et al. Effects of Vegetation Changes and Multiple Environmental Factors on\nEvapotranspiration Across China Over the Past 34 Years. Earth’s Futur. 10, 1–16 (2022).\nIn the flowchart of Fig. 4, I notice a missing of surface resistance in the VPD-LE relationship.\nIn the Penman-Monteith scheme, increasing VPD can have competing effects on ET: on one\nhand, increasing VPD enhances ET, while on the other, it also raises surface resistance, which\nsuppresses ET. The results suggest that the VPD-driven enhancement of ET is more dominant.\nCould the authors clarify why this occurs when comparing ET between forests and open lands?\nI feel VPD is not the driver of such findings, but the biome type - forests/open lands and their\nassociated biophysical and physiological responses - are the real drivers. Need to clarify this.\nSee Novick, K., Ficklin, D., Stoy, P. et al. The increasing importance of atmospheric demand\nfor ecosystem water and carbon fluxes. Nature Clim Change 6, 1023 – 1027\n(2016). https://doi.org/10.1038/nclimate3114\nResponse: We sincerely appreciate the valuable comments. We have carefully read the referred\npaper and we fully agree with your perspective of the opposing effects of VPD on LE (or ET).\nHowever, it is difficult for our data-based analysis to explicitly incorporate the impact of surface\nresistance (r), due to the lack of reliable r data.\ns s\nTo address this issue, we performed additional analyses using in-situ observation and the\nPenman-Monteith framework, which can account for the impact of VPD on r.\ns\nThe results reveal\nthat when neglecting the r effect, VPD demonstrates a more pronounced positive effect on LE in\ns\nforests. However, when considering r responses to VPD, the sensitivity differences between forests\ns\nand openlands become statistically insignificant. This result implies that the negative effect resulting\nfrom the r response to rising VPD offsets the direct positive effect of VPD on forest\ns\nevapotranspiration.\nThese findings challenge our initial hypothesis that VPD exerts stronger stimulation on\nevaporative cooling in forests, highlighting the critical role of surface resistance dynamics in\nmediating these relationships. We have incorporated this supplementary analysis in the revised\nmanuscript (Line 262):\nWe also evaluate the sensitivity of summer latent heat (LE) to VPD in both forests and\ns s\nopenlands using eddy covariance observations from European sites (Supplementary Text 2). The\nsensitivity analysis is based on the Penman-Monteith framework, which accounts for the influence\nof VPD on surface resistance (Novick et al., 2016). If VPD is the main driver of enhanced forest\ns s\ncooling effect, we expect to observe greater positive LE sensitivity in forests compared to openlands,\ns\nin which case forests can produce stronger evaporative cooling effects in the context of rising VPD.\ns\nHowever, we find openlands show slightly higher positive sensitivity than forests, although the\ndifference is not significant. The direct regression analysis also supports the stronger LE responses\ns\nto VPD increases in openlands (Supplementary Table. 1). This in-situ observational evidence,\ns\ncombined with our statistical results (binning regression and random forest regression), suggests\nthat VPD may not explain the enhanced cooling effect of European forests, while the decline in SM\ns s\nappears to be the real dominant driver.\nSupplementary Text 2. Sensitivity of summer latent heat to increasing vapor pressure deficit\nElevated vapor pressure deficit (VPD) exerts dual opposing effects on latent heat (LE) flux:\nwhile it directly enhances LE through increased evaporative demand, it simultaneously elevates\nsurface resistance (rₛ), thereby suppressing LE (Novick et al., 2016). To investigate how VPD\ndifferentially influences LE in forests and openlands and consequently affects the cooling effect of\nforests, we analyze European eddy covariance observations from the FLUXNET2015 dataset and\nconduct sensitivity analyses.\nThe sensitivity of LE to VPD is assessed based on the Penman-Monteith framework, accounting\nfor both the aforementioned opposing effects (Lin et al., 2018):\n𝐿𝐸 (cid:3404)\n(cid:3057)(cid:4666)(cid:3019)(cid:3041)(cid:2879)(cid:3008)(cid:4667)(cid:2878)(cid:3096)(cid:3004)(cid:3291) (cid:4666)(cid:3023)(cid:3017)(cid:3005)/(cid:3045)(cid:3276) (cid:4667)\n(1)\n(cid:3057)(cid:2878)(cid:3082)(cid:4666)(cid:2869)(cid:2878)(cid:3045)(cid:3294)/(cid:3045)(cid:3276) (cid:4667)\n(cid:3008)(cid:3017)(cid:3017)\n𝐺 (cid:3404)𝐺 (cid:3397)𝐺 (2)\n(cid:3046) (cid:2868) (cid:2869) (cid:3023)(cid:3017)(cid:3005)(cid:3288)\nHere, Rn (W/m2) is the net radiation; G (W/m2) is the soil heat flux; Δ (kPa/K) is the slope of the\nsaturated vapor pressure curve; ρ is the air density (kg/m3); C (J/kg/K) is the specific heat of air;\np\nγ is the psychrometric constant (kPa/K); r (s/m) is the aerodynamic resistance, estimated using\na\nvegetation height; G (surface conductance) is the inverse of r; GPP is the gross primary\ns s\nproductivity. The model parameters (G , G , and m) is estimated via optimization (SciPy package\n0 1\nin Python) using daily-scale observations.\nWe analysis those observation in summer and we select sites with reliable LE accuracy (Nash-\nSutcliffe Efficiency > 0.5), resulting in 19 forested and 12 openland (grassland and cropland) sites.\nWe then introduce VPD perturbations into the observation data to estimate summer LE (LE)\ns s\nsensitivity to VPD (\n(cid:3031)(cid:3013)(cid:3006)(cid:3294)\n, LE response to a 1 hPa VPD increase). For comparison, we also\ns s s\n(cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294)\ncompute\n(cid:3031)(cid:3013)(cid:3006)(cid:3294)\n(i) without the impact of VPD on rₛ and (ii) via direct linear regression between\ns\n(cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294)\nVPD and LE (Supplementary Table 1).\ns s\nIf VPD is the primary driver of enhanced forest cooling effect, forests should exhibit a higher\ns\n(cid:3031)(cid:3013)(cid:3006)(cid:3294)\nthan openlands, implying stronger evaporative cooling under rising VPD. However, results\ns\n(cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294)\nindicate that this expectation holds only when neglecting the inhibitory effect of VPD on LE (via\ns s\nrₛ). When considering both opposing effects of VPD on LE, openlands show slightly higher mean\ns s\nsensitivity compared to forests. In addition, the direct regression analysis supports the stronger LE\ns\nresponses to VPD increases in openlands. These findings suggest that VPD may not explain the\ns s\nenhanced cooling effect observed in European forests.\nTable S1. Sensitivity of summer latent heat (LEs) flux to rising vapor pressure deficit (VPDs) based on FLUXNET\nobservations in Europe.\nRegression slope\n(cid:3031)(cid:3013)(cid:3006)(cid:3294) (w·m-2/hPa)\n(cid:3031)\n(cid:3031)\n(cid:3023)\n(cid:3013)\n(cid:3017)\n(cid:3006)\n(cid:3005)\n(cid:3294)\n(cid:3294)\n(w·m-2/hPa) (cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294) between VPDs and LEs\nignoring the impact of VPD on rs (w·m-2/hPa)\nForest (n=19) 2.08 ± 1.85 5.92 ± 2.95 3.87 ± 2.87\nOpenland (n=19) 2.48 ± 1.71 4.87 ± 1.32 6.44 ± 2.67\nThe values indicate the mean ± standard deviation across sites.\nCombining this sensitivity analysis result and the statistical evidence (please see the response\nto your comment above), we recognized that the decreasing SM appears to be the primary driver\ns\nof the observed daytime ΔLST trends. Thus, we updated the Fig. 4 and the related descriptions\ns\n(Line 276). The biophysical and physiological mechanisms corresponding to each path are also\nannotated to enhance conceptual clarity (Fig. 4a).\nWe further conduct structural equation modeling (SEM) to elucidate the biophysical\nmechanisms underlying the impact of varying SM on daytime ΔLST dynamics (Fig. 4a). Our SEM\ns s\nincorporates two key biophysical mediators: the leaf area index difference (ΔLAI) and latent heat\ns\nflux difference (ΔLE) between nearby forests and openlands. The path effect (mean ± standard error)\ns\nis used to describe strength and direction of the relationship. Following previous studies, we build\nSEMs at the local scale and calculate the mean path effects of all models (Shen et al., 2024; Zhong\net al., 2023). The analysis reveals three distinct pathways of SM affecting ΔLST.\ns s\n(1) SM→ΔLAI→ΔLST. SM exerts a significant negative effect on ΔLAI (-0.41 ± 0.15). Then,\ns s s s s\nΔLAI demonstrate a negative effect on ΔLST (-0.19 ± 0.09). This pathway suggests that reduced\ns s\nSM can amplify the greenness difference in forests and openlands, thereby enhancing forests’\ns\ncooling advantage through improved convective heat dissipation.\n(2) SM\ns\n→ΔLE\ns\n→ΔLST\ns\n. SM\ns\nhave a negative effect on ΔLE\ns\n(-0.25 ± 0.11), followed by a strong\nnegative influence of ΔLE on ΔLST (-0.54 ± 0.09). This pathway indicates that declining SM\ns s s\nimposes less water limitation to forests compared to openlands, enabling forests to maintain\nstronger relative evaporative cooling through differential stomatal regulation (with greenness\nremaining constant).\n(3) SM→ΔLAI→ΔLE→ΔLST. SM-induced changes in ΔLAI subsequently enhance ΔLE\ns s s s s s s\n(0.62 ± 0.06), ultimately influencing ΔLST. This demonstrates that SM reduction amplifies\ns s\nevaporative cooling effects of forests by increasing the vegetation greenness contrast with openlands.\nThe pathway effects, representing the strength and direction of the mechanism, are quantified\nas the product of all constituent path coefficients (Fig. 4b). Results reveal that the second and third\npathways (both mediated by ΔLE) contribute equally to ΔLST variation (about 0.14 ± 0.06), while\ns s\nthe impact first pathway is weaker (0.08 ± 0.05). These findings suggest SM can modulate ΔLST\ns s\ndynamics through multiple mechanisms, with forest evaporative cooling regulation emerging as the\ndominant process. Previous studies have demonstrated that reduced soil moisture in Europe can\namplify warming during heatwaves through land-atmosphere interactions (Dirmeyer et al., 2021;\nFischer et al., 2007). Our findings confirm forests exhibit greater stability and resistance to\nincreasing soil dryness compared to grasslands or croplands, which may mitigate land-atmosphere\ncoupling feedback and consequently enhance the relative cooling effect of forests (Teuling et al.,\n2010). This stability of forests is reflected in both more moderate evaporative limitation and\nsustained vegetation greenness under increasing soil dryness. Our results provide observational\nevidence suggesting that the capacity of forest ecosystems to mitigate water deficit-induced warming\npersists not only during heatwave events but also across interannual timescales.\nFig 4. Mechanisms underlying the effects of climate variables on the daytime LST effects of\nEuropean forests (ΔLST) in summer. (a) Structure equation models describing the biophysical\ns\nrelationships between climate variables and ΔLST. The numbers denote the path value (mean ±\ns\nstandard error). The arrow color indicates the sign of the path (red for positive and blue for\nnegative), and the thickness indicates the magnitude of the corresponding path value. (b) Bar plot\nof the pathway effects of climate variables affecting ΔLST. The error bars indicate the uncertainty\ns\nbased on error propagation.\nAccordingly, we also revised the discussion about the summer mechanism (Line 434).\nUsing a range of statistical methods, we find that declining SM is the predominant factor\ns\ncontributing to both the interannual variation and the negative trend in ΔLST. Path analysis shows\ns\nthat SM affects ΔLST dynamics through two mechanisms. First, we show ΔLAI is modulated by\ns s s\nSM, implying that as the soil becomes drier, the vegetation structure difference between forests and\ns\nopenlands becomes more evident. The higher ΔLAI could boost the cooling efficiency of forests\ns\nthrough stronger direct heat dissipation and evaporative cooling. The impact of SM on vegetation\ns\nhas been well-documented: increasing dryness pressure could trigger stomatal closure, inhibit\nphotosynthesis, and increase vegetation mortality (Li et al., 2022; Stocker et al., 2018). Our results\nconfirm a stronger negative effect of SM on non-forest vegetation than on forests (Hoek van Dijke\ns\net al., 2023). The potential reason is that forest ecosystems with high species richness and deep\nroots are more stable and resilient to dryness (Isbell et al., 2015; Tumber-Dávila et al., 2022).\nSecond, even under constant vegetation conditions, less soil water supply directly amplifies ΔLE\ns\nand enhances the cooling effect of forests. These findings are consistent with prior researches\n(Fu et al., 2022a, 2022b), indicating that soil moisture imposes stronger limitations on ET in\ngrasslands and croplands relative to forest ecosystems. This differential response is evidenced by\nboth the greater sensitivity of and the higher critical soil moisture thresholds in non-forest\necosystems. This phenomenon may be attributed to the capacity of forests to maintain ET through\ndeep soil water extraction, thereby exhibiting greater resilience to warming-induced surface soil\nmoisture deficits compared to openlands. Overall, the revealed impact of SM on ΔLST essentially\ns s\nreflects the diverse vegetation physiological responses to climate changes.\nReference:\nDirmeyer, P. A., Balsamo, G., Blyth, E. M., Morrison, R. & Cooper, H. M. Land‐Atmosphere\nInteractions Exacerbated the Drought and Heatwave Over Northern Europe During Summer\n2018. AGU Adv. 2, 1–16 (2021).\nFischer, E. M., Seneviratne, S. I., Lüthi, D. & Schär, C. Contribution of land-atmosphere coupling\nto recent European summer heat waves. Geophys. Res. Lett. 34, 1–6 (2007).\nFu, Z. et al. Critical soil moisture thresholds of plant water stress in terrestrial ecosystems. Sci. Adv.\n8, 1–12 (2022).\nFu, Z. et al. Uncovering the critical soil moisture thresholds of plant water stress for European\necosystems. Glob. Chang. Biol. 28, 2111–2123 (2022).\nHoek van Dijke, A. J. et al. Comparing forest and grassland drought responses inferred from eddy\ncovariance and Earth observation. Agric. For. Meteorol. 341, (2023).\nIsbell, F. et al. Biodiversity increases the resistance of ecosystem productivity to climate extremes.\nNature 526, 574–577 (2015).\nLi, W. et al. Widespread increasing vegetation sensitivity to soil moisture. Nat. Commun. 13, 3959\n(2022).\nLin, C. et al. Diel ecosystem conductance response to vapor pressure deficit is suboptimal and\nindependent of soil moisture. Agric. For. Meteorol. 250–251, 24–34 (2018).\nNovick, K. A. et al. The increasing importance of atmospheric demand for ecosystem water and\ncarbon fluxes. Nat. Clim. Chang. 6, 1023–1027 (2016).\nShen, P. et al. Biodiversity buffers the response of spring leaf unfolding to climate warming. Nat.\nClim. Chang. 14, 863–868 (2024).\nStocker, B. D. et al. Quantifying soil moisture impacts on light use efficiency across biomes. New\nPhytol. 218, 1430–1449 (2018).\nTeuling, A. J. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci. 3, 722–727 (2010).\nTumber-Dávila, S. J., Schenk, H. J., Du, E. & Jackson, R. B. Plant sizes and shapes above and\nbelowground and their interactions with climate. New Phytol. 235, 1032–1056 (2022).\nZhong, Z. et al. Disentangling the effects of vapor pressure deficit on northern terrestrial vegetation\nproductivity. Sci. Adv. 9, (2023).\nReviewer #1 (Remarks on code availability):\nI feel that this might be holding the authors to an excessively high standard. If Nature is asking\nwhether I can strictly replicate results using their code, I believe it would take a significant\namount of time, and I don't have the capacity to verify it. However, I feel is that it's unrealistic\nto expect the authors to cover every detail perfectly and provide a fully automated, one-click\nrunnable code.\nAdditionally, I believe that the provided code does not fully include all the components/scripts\nnecessary to reproduce all the results of this paper (e.g., the ESM analysis, MODIS data\nprocessing). Again, I feel this is holding the authors to an excessively high standard. They are\nnot software engineers.\nResponse: Thank you for the feedback concerning code availability. We fully recognize the\nimportance of transparent research practices. In response, we have reuploaded the code, including\nall processed datasets necessary to execute the analysis scripts\n(https://zenodo.org/records/15625044). These materials enable researchers to fully reproduce all\nstatistical analyses and figures presented in our manuscript.\nReviewer #2 (Remarks to the Author):\nI find this manuscript has a lot of interesting data and discussion. Moreover, for me personally,\nit was highly motivating to read this discussion. I did find, especially at the beginning of the\nmanuscript, that the work was not very well situated in the principal literature. All four of the\nfirst references to appear in the manuscript do not seem well suited to set up this general\ndiscussion. I have provided a number of suggestions in the attached submission pdf and am\nhopeful that some of these suggestions will prove useful...\nResponse: We appreciate your positive comments. We have revised our manuscript carefully and\nthe detailed corrections are listed below. Specifically, we have thoroughly reorganized the\nintroduction section by incorporating more appropriate foundational literature as suggested in your\ncomments. According to your concern regarding the water availability, we conducted additional\nexperiments to explore the impact of soil moisture, and added related discussions. We hope the\nrevised paper will be more suitable for publication.\nI also found that the introduction was occasionally vague and would have benefitted from a\nsomewhat more nuanced discussion... I have tried to point out places where I think this would\nbe helpful.\nResponse: Thanks for the careful reading and pointing out the confusing statements. Your detailed\nsuggestions have helped us significantly improve the manuscript. We provide point-by-point\nresponses to each of your specific comments about the introduction in the PDF file, along with the\ncorresponding revisions we have implemented (please see the response below).\nFinally, the discussion of radiative and non-radiative effects, as well as the discussion of winter\nand summertime effects, to my mind, were not adequately introduced and motivated in the\nintroductory parts of the paper. When these concepts are suddenly announced, (the\ndiscussion about winter and summer time effects first comes up in the results section of the\npaper?), they are unexpected and seemingly unmotivated in the general context of the paper...\nThis can be improved upon...\nResponse: We appreciate this constructive suggestion. We have incorporated background on\nradiative/non-radiative mechanisms and seasonal (winter/summer) effects in the introduction,\nensuring these concepts are properly motivated before appearing in later sections (Line 39).\nDuring the daytime, forests exhibit two competing effects compared to non-forest vegetation:\n(a) the non-radiative cooling effect driven by higher evapotranspiration (ET) rates (Mildrexler et\nal., 2011), and (b) the radiative warming effect resulting from lower albedo (R.A.Betts, 2000).\nThe net temperature response depends on the balance between the opposing radiative and non-radiative\neffects (Bright et al., 2017), and their relative dominance is largely governed by background climatic\nconditions (Pitman et al., 2011). For instance, forests in mid-latitude regions (e.g., North America,\nEurope and China) demonstrate seasonally contrasting LST effects (Li et al., 2015; Tang et al.,\n2018). During cold winters, bright snow covers short vegetation but is masked by dense forest\ncanopies. As a result, forests absorb more solar radiation and show a local warming effect compared\nto surrounding openlands. Conversely, during the peak growing season (warm summer), forests\nexert cooling effects under sufficient soil moisture conditions, owing to their higher ET rates\ncompared to other vegetation types. However, in arid regions, limited soil moisture can constrain\nforest evaporative cooling despite abundant radiation and optimal temperature, causing radiative\nwarming to dominate (Rohatyn et al., 2023).\nReference:\nBright, R. M. et al. Local temperature response to land cover and management change driven by\nnon-radiative processes. Nat. Clim. Chang. 7, 296–302 (2017).\nLi, Y. et al. Local cooling and warming effects of forests based on satellite observations. Nat.\nCommun. 6, 1–10 (2015).\nMildrexler, D. J., Zhao, M. & Running, S. W. A global comparison between station air temperatures\nand MODIS land surface temperatures reveals the cooling role of forests. J. Geophys. Res.\nBiogeosciences 116, 1–15 (2011).\nPitman, A. J. et al. Importance of background climate in determining impact of land-cover change\non regional climate. Nat. Clim. Chang. 1, 472–475 (2011).\nRohatyn, S., Rotenberg, E., Tatarinov, F., Carmel, Y. & Yakir, D. Large variations in afforestation-\nrelated climate cooling and warming effects across short distances. Commun. Earth Environ.\n4, 1–10 (2023).\nR.A.Betts. Offset of the potential carbon sink from boreal forestation by decreases in surface albedo.\nNature 408, 187–190 (2000).\nTang, B., Zhao, X. & Zhao, W. Local effects of forests on temperatures across Europe. Remote Sens.\n10, 1–24 (2018).\nI also find that the Authors do not adequately distinguish between global and local cooling...\nBoth terms come up. But most of the time it sounds like the Authors are talking about local\ncooling... And the balance in terms of global cooling is never clearly estimated...? It is\nimportant not to conflate these terms, since local cooling is not necessarily equivalent to\nglobal cooling...?\nResponse: We appreciate this important criticism. We fully agree with that global and local effects\nare inherently different. As you mentioned, our manuscript focusses on the local effect. We have\nadded explicit discussion distinguishing between local and global cooling effects of forests in the\nrevised manuscript, with particular attention to their differing mechanisms (Line 477):\nThere are several potential caveats or issues when interpreting our results. First, the adopted\n“space-for-time” method provides priori estimates of local temperature effects of forests (Wang et\nal., 2023), which is different from the global effects. At the planetary scale, forests can exert more\ncomplicated effects on climate through multiple processes such as changes in cloud formation and\nlarge-scale atmospheric circulations (Ellison et al., 2024; Xu et al., 2022). Thus, our results should\nbe treated as the reference for local climate adaptation rather than global climate mitigation.\nWe also distinguished between the global and local impacts in the background statement (Line\n25):\nThrough conservation, proper management, and restoration practices, forest ecosystems could\nplay a crucial role in addressing global warming by regulating the hydrological cycle and\nenhancing cloud formation (Cerasoli et al., 2021; Ellison et al., 2024; Luo et al., 2024). In addition\nto the global impact, forests exhibit significant local cooling or warming effects through altering\nsurface energy and water fluxes (Ellison et al., 2024; Lawrence et al., 2022).\nReference:\nCerasoli, S., Yin, J. & Porporato, A. Cloud cooling effects of afforestation and reforestation at\nmidlatitudes. Proc. Natl. Acad. Sci. U. S. A. 118, 1–7 (2021).\nEllison, D., Pokorný, J. & Wild, M. Even cooler insights: On the power of forests to (water the Earth\nand) cool the planet. Glob. Chang. Biol. 30, 1–20 (2024).\nLawrence, D., Coe, M., Walker, W., Verchot, L. & Vandecar, K. The Unseen Effects of Deforestation:\nBiophysical Effects on Climate. Front. For. Glob. Chang. 5, 1–13 (2022).\nLuo, H., Quaas, J. & Han, Y. Decreased cloud cover partially offsets the cooling effects of surface\nalbedo change due to deforestation. Nat. Commun. 15, 7345 (2024).\nWang, H., Yue, C. & Luyssaert, S. Reconciling different approaches to quantifying land surface\ntemperature impacts of afforestation using satellite observations. Biogeosciences 20, 75–92\n(2023).\nXu, R. et al. Contrasting impacts of forests on cloud cover based on satellite observations. Nat.\nCommun. 13, 670 (2022).\nFurther, I would have like to see more discussion of the role of the availability of water and\nsurface moisture. The cooling potential of forests strongly depends on the availability of water.\nAs soon as this surface moisture disappears, the loss of latent heat means the warming will\nsuddenly take over... But there is little attention to this throughout the paper...\nResponse: Thank you for raising this essential issue. We full acknowledge that both atmospheric\nwater demand (VPD) and soil water supply (soil moisture, SM) can have important effects on the\nevaporative cooling of vegetation. Given the strong annual-scale correlation between VPD and SM,\ntheir respective roles may have been conflated in the original manuscript. Thus, we conducted\nadditional analyses to identify the actual driver of the observed trends in ΔLST. Specifically, we\ns\nused binning regression method (Liu et al., 2020) to decouple summer VPD (VPD) and summer\ns\nSM (SM). The results demonstrated that it is SM, rather than VPD, dominates with the long-term\ns s\ntrend of daytime ΔLSTₛ, which largely addresses this comment (Line 218):\nTemporal correlation analyses using VPD, AT, DSR, P and WS as the potential drivers\ns s s s s\nreveal strong relationships between daytime ΔLST and climatic variables except WS (Fig. 3a).\ns s\nHowever, after accounting for covariates, the partial correlations between the daytime ΔLST and\ns\nDSR, P, and AT are greatly reduced and even reverse in sign. In contrast, VPD maintain a\ns s s s\npredominant negative correlation with daytime ΔLST, with the mean partial correlation coefficient\ns\nof -0.42 (Fig. 3a). Notably, we can also find strong positive partial correlations between SM and\ns\nΔLST (Fig. 3b) when replacing VPD by the SM, the other factor describing the dryness stress.\ns s s\nThis is due to the strong coupling between SM and VPD at the annual scale (r = 0.92,\ns s\nSupplementary Fig. 10). These findings raise the following question: whether variations in\natmosphere or soil dryness are the primary factor governing changes in ΔLST over the past two\ns\ndecades?\nTo answer this question, we then use the binning regression method to further disentangle the\nrespective contributions of VPD and SM (Liu et al., 2020). We conduct this analysis using the\ns s\nspatial samples of VPD\ns\ntrend and SM\ns\ntrends with sufficient data for decoupling the two factors.\nWe confirm the strong spatial correlation between SM trend and VPD trends (r = -0.65) when all\ns s\nspatial samples are used. However, after stratifying either SM or VPD trends into 10 bins\ns s\naccording to percentiles, correlations between the two variables show substantial attenuation (Fig.\n3c), suggesting the two variables are decoupled. Then, we calculate the correlation between SM\ns\n(or VPD) trends and ΔLST trends in each VPD (or SM) trend bin. We find the mean r between\ns s s s\nVPD trends and ΔLST trends is only -0.11, whereas the mean r between SM trends and ΔLST\ns s s s\ntrends reach 0.40 (Fig. 3d). This binning regression analysis suggests that change SM is more likely\ns\nto be the driver of ΔLST dynamic.\ns\nFig 3. Dominant role of soil moisture in the temporal dynamics in the summer daytime cooling effect\n(ΔLST) of forests. (a) Regional mean correlation and partial correlation coefficients between five\ns\nclimate variables (vapor pressure deficit, downward solar radiation, precipitation, air temperature\nand wind speed) and the daytime ΔLST. The error bars indicate the standard deviations. (b) Similar\ns\nwith (a), but replace vapor pressure deficit by the soil moisture. (c) Correlation coefficients between\nsummer vapor pressure deficit (VPD) trends and summer soil moisture (SM) trends within either\ns s\nSM or VPD trends bins. (d) Correlation coefficients between VPD (or SM) trends and ΔLST\ns s s s s\ntrends within SM (or VPD) trends bins.\ns s\nFig. S10. Temporal correlation coefficient (r) between summer vapor pressure deficit (VPD) and\ns\nsurface soil moisture (SM). (a) Spatial pattern of r between VPD and SM. (b) Density of the r\ns s s\nvalues.\nWe then build a non-linear random forest model using both SM and VPD as the input. The\ns s\nresults also support that SM plays the more important role in predicting temporal dynamics of\ns\nΔLST (Line 255):\ns\nThe critical role of SM is corroborated by non-linear statistical analysis using a Random\ns\nForest (RF) model (Supplementary Text 1). Variable importance analysis demonstrates that SM\ns\nexerts a stronger influence than VPD in predicting ΔLST dynamics. Contribution analysis also\ns s\nreveals that SM changes account for the largest portion (0.152 K/decade, p < 0.05) of the simulated\ns\nnegative ΔLST trend (0.181 K/decade, p < 0.05), highlighting its dominant role in driving the\ns\nobserved enhanced forest cooling effect (Supplementary Fig. 11).\nSupplementary Text 1 Modelling ΔLST based on random forest\ns\nWe apply a non-linear machine learning method, Random Forest (RF), to evaluate the\nimportance and contribution of multiple potential climate drivers for predicting ΔLST (sklearn\ns\npackage in Python). Specifically, we use all spatial samples over 20 years to build the dataset. It is\nworth mentioning that we normalize the input climate variables and ΔLST by subtracting the 20-\ns\nyear means at the grid scale, since the model focus on the temporal dynamic of ΔLST rather than\ns\nthe spatial variation. We excluded air temperature from the model inputs due to its high correlation\nwith other climate variables.\nThe dataset is divided into training (80% of the samples) and test sets (20% of the samples),\nwith model parameters confirmed based on validation results from the test set (number of trees =\n100; maximum depth of the tree = 6). All the samples are then used to fit the best model. The relative\nimportance of each variable is indicated by the Gini coefficient. Moreover, we use the Shapley\nAdditive Explanations (SHAP) method to quantify the individual contribution of each variable in\neach sample (shap package in Python). SHAP is an approach to explain machine learning models\nbased on cooperative game theory, where SHAP values quantify the marginal contributions of\npredictors. We calculate the mean absolute SHAP values for each input variable as importance\nmetrics and explore the non-linear impact of climate variables on ΔLST by plotting their SHAP\ns\nvalues against input values.\nBased on the RF model, we set several scenarios to separate the relative contributions of\nclimate variables to the ΔLST trend. Specifically, we perform four experimental simulations: (S1)\ns\nvarying SM only; (S2) varying SM and VPD, while fixing DSR P and WS; (S3) varying SM,\ns s s s s s s\nVPD and DSR, while fixing P and WS; (S4) varying SM, VPD, DSR and P, while fixing WS;\ns s s s s s s s s\nand (S5) varying SM, VPD, DSR, P, and WS. The Theil-Sen slope and Mann-Kendall P-value of\ns s s s s\nsimulations S1, S2−S1, S3−S2, S4−S3, S5−S4 are used to evaluate the individual contributions of\nthe changing climate to the ΔLST trend.\ns\nThe model shows good accuracy on the test dataset (R2=0.56, RMSE=0.32 K), supporting our\nfurther analysis. The Gini importance and mean absolute SHAP values both suggest that SM is the\ns\nmost important driver, followed by VPD, P, DSR and WS (Supplementary Fig 11a). The marginal\ns s s s\ncontributions quantified by the SHAP value further reveal different impacts of SM and VPD\ns s\n(Supplementary Fig. 11b). SM exhibits the pronounce positive effect, with stronger impact when\ns\nSM anomaly is negative. While VPD only show slight negative impact under certain condition.\ns s\nOn the basis of the RF model, we estimate the contributions of climate variables to the long-\nterm trend of daytime ΔLST (Supplementary Fig. 11c). We first compared the reconstructed trends\ns\nwith all forcings and the observed trends, showing good consistency (observed: 0.188 K/decade,\np<0.05; simulated: 0.181 K/decade, p<0.05). By setting different scenarios, results suggest that\ndecreasing SM contributes the most to the negative trend in daytime ΔLST (-0.152 K/decade,\ns s\np<0.05), and increasing VPD shows lower contribution (-0.049 K/decade, p<0.05). The isolated\ns\ncontribution of DSR is also statistically significant, but the value is relatively small (0.012 K/decade,\ns\np<0.05). The contribution of WS and P are not significant.\ns s\nFig. S11. Impact of SM on the temporal dynamics in ΔLST based on the random forest model. (a)\ns s\nThe importance of climate variables is indicated by the absolute SHAP values and Gini coefficients.\n(b) Variation of SHAP values with changing SM or VPD anomaly. (c) Observed and modelled\nsn\nΔLST trends, and the isolated contributions from climate variables.\ns\nWe also performed additional analyses using in-situ observation and the Penman-Monteith\nframework to evaluate the sensitivity of LE to VPD. If VPD is the main driver of enhanced forest\ns s s\ncooling effect, we expect to observe greater positive LE sensitivity in forests compared to openlands.\ns\nHowever, the results do not support this hypothesis (Line 262).\nWe also evaluate the sensitivity of summer latent heat (LE) to VPD in both forests and\ns s\nopenlands using eddy covariance observations from European sites (Supplementary Text 2). The\nsensitivity analysis is based on the Penman-Monteith framework, which accounts for the influence\nof VPD on surface resistance (Novick et al., 2016). If VPD is the main driver of enhanced forest\ns s\ncooling effect, we expect to observe greater positive LE sensitivity in forests compared to openlands,\ns\nin which case forests can produce stronger evaporative cooling effects in the context of rising VPD.\ns\nHowever, we find openlands show slightly higher positive sensitivity than forests, although the\ndifference is not significant. The direct regression analysis also supports the stronger LE responses\ns\nto VPD increases in openlands (Supplementary Table. 1). This in-situ observational evidence,\ns\ncombined with our statistical results (binning regression and random forest regression), suggests\nthat VPD may not explain the enhanced cooling effect of European forests, while the decline in SM\ns s\nappears to be the real dominant driver.\nSupplementary Text 2. Sensitivity of summer latent heat to increasing vapor pressure deficit\nElevated vapor pressure deficit (VPD) exerts dual opposing effects on latent heat (LE) flux:\nwhile it directly enhances LE through increased evaporative demand, it simultaneously elevates\nsurface resistance (rₛ), thereby suppressing LE (Novick et al., 2016). To investigate how VPD\ndifferentially influences LE in forests and openlands and consequently affects the cooling effect of\nforests, we analyze European eddy covariance observations from the FLUXNET2015 dataset and\nconduct sensitivity analyses.\nThe sensitivity of LE to VPD is assessed based on the Penman-Monteith framework, accounting\nfor both the aforementioned opposing effects (Lin et al., 2018):\n𝐿𝐸 (cid:3404)\n(cid:3057)(cid:4666)(cid:3019)(cid:3041)(cid:2879)(cid:3008)(cid:4667)(cid:2878)(cid:3096)(cid:3004)(cid:3291) (cid:4666)(cid:3023)(cid:3017)(cid:3005)/(cid:3045)(cid:3276) (cid:4667)\n(1)\n(cid:3057)(cid:2878)(cid:3082)(cid:4666)(cid:2869)(cid:2878)(cid:3045)(cid:3294)/(cid:3045)(cid:3276) (cid:4667)\n(cid:3008)(cid:3017)(cid:3017)\n𝐺 (cid:3404)𝐺 (cid:3397)𝐺 (2)\n(cid:3046) (cid:2868) (cid:2869) (cid:3023)(cid:3017)(cid:3005)(cid:3288)\nHere, Rn (W/m2) is the net radiation; G (W/m2) is the soil heat flux; Δ (kPa/K) is the slope of the\nsaturated vapor pressure curve; ρ is the air density (kg/m3); C (J/kg/K) is the specific heat of air;\np\nγ is the psychrometric constant (kPa/K); r (s/m) is the aerodynamic resistance, estimated using\na\nvegetation height; G (surface conductance) is the inverse of r; GPP is the gross primary\ns s\nproductivity. The model parameters (G , G , and m) is estimated via optimization (SciPy package\n0 1\nin Python) using daily-scale observations.\nWe analysis those observation in summer and we select sites with reliable LE accuracy (Nash-\nSutcliffe Efficiency > 0.5), resulting in 19 forested and 12 openland (grassland and cropland) sites.\nWe then introduce VPD perturbations into the observation data to estimate summer LE (LE)\ns s\nsensitivity to VPD (\n(cid:3031)(cid:3013)(cid:3006)(cid:3294)\n, LE response to a 1 hPa VPD increase). For comparison, we also\ns s s\n(cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294)\ncompute\n(cid:3031)(cid:3013)(cid:3006)(cid:3294)\n(i) without the impact of VPD on rₛ and (ii) via direct linear regression between\ns\n(cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294)\nVPD and LE (Supplementary Table 1).\ns s\nIf VPD is the primary driver of enhanced forest cooling effect, forests should exhibit a higher\ns\n(cid:3031)(cid:3013)(cid:3006)(cid:3294)\nthan openlands, implying stronger evaporative cooling under rising VPD. However, results\ns\n(cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294)\nindicate that this expectation holds only when neglecting the inhibitory effect of VPD on LE (via\ns s\nrₛ). When considering both opposing effects of VPD on LE, openlands show slightly higher mean\ns s\nsensitivity compared to forests. In addition, the direct regression analysis supports the stronger LE\ns\nresponses to VPD increases in openlands. These findings suggest that VPD may not explain the\ns s\nenhanced cooling effect observed in European forests.\nTable S1. Sensitivity of summer latent heat (LEs) flux to rising vapor pressure deficit (VPDs) based on FLUXNET\nobservations in Europe.\nRegression slope\n(cid:3031)(cid:3013)(cid:3006)(cid:3294) (w·m-2/hPa)\n(cid:3031)\n(cid:3031)\n(cid:3023)\n(cid:3013)\n(cid:3017)\n(cid:3006)\n(cid:3005)\n(cid:3294)\n(cid:3294)\n(w·m-2/hPa) (cid:3031)(cid:3023)(cid:3017)(cid:3005)(cid:3294) between VPDs and LEs\nignoring the impact of VPD on rs (w·m-2/hPa)\nForest (n=19) 2.08 ± 1.85 5.92 ± 2.95 3.87 ± 2.87\nOpenland (n=19) 2.48 ± 1.71 4.87 ± 1.32 6.44 ± 2.67\nThe values indicate the mean ± standard deviation across sites.\nCombined the sensitivity analysis result and the statistical evidence, we recognized that the\ndecreasing SM appears to be the primary driver of the observed daytime ΔLST trends. Thus, we\ns s\nupdated the Fig. 4 and the related mechanism descriptions accordingly (Line 276).\nWe further conduct structural equation modeling (SEM) to elucidate the biophysical\nmechanisms underlying the impact of varying SM on daytime ΔLST dynamics (Fig. 4a). Our SEM\ns s\nincorporates two key biophysical mediators: the leaf area index difference (ΔLAI) and latent heat\ns\nflux difference (ΔLE) between nearby forests and openlands. The path effect (mean ± standard error)\ns\nis used to describe strength and direction of the relationship. Following previous studies, we build\nSEMs at the local scale and calculate the mean path effects of all models (Shen et al., 2024; Zhong\net al., 2023). The analysis reveals three distinct pathways of SM affecting ΔLST.\ns s\n(1) SM→ΔLAI→ΔLST. SM exerts a significant negative effect on ΔLAI (-0.41 ± 0.15). Then,\ns s s s s\nΔLAI demonstrate a negative effect on ΔLST (-0.19 ± 0.09). This pathway suggests that reduced\ns s\nSM can amplify the greenness difference in forests and openlands, thereby enhancing forests’\ns\ncooling advantage through improved convective heat dissipation.\n(2) SM\ns\n→ΔLE\ns\n→ΔLST\ns\n. SM\ns\nhave a negative effect on ΔLE\ns\n(-0.25 ± 0.11), followed by a strong\nnegative influence of ΔLE on ΔLST (-0.54 ± 0.09). This pathway indicates that declining SM\ns s s\nimposes less water limitation to forests compared to openlands, enabling forests to maintain\nstronger relative evaporative cooling through differential stomatal regulation (with greenness\nremaining constant).\n(3) SM→ΔLAI→ΔLE→ΔLST. SM-induced changes in ΔLAI subsequently enhance ΔLE\ns s s s s s s\n(0.62 ± 0.06), ultimately influencing ΔLST. This demonstrates that SM reduction amplifies\ns s\nevaporative cooling effects of forests by increasing the vegetation greenness contrast with openlands.\nThe pathway effects, representing the strength and direction of the mechanism, are quantified\nas the product of all constituent path coefficients (Fig. 4b). Results reveal that the second and third\npathways (both mediated by ΔLE) contribute equally to ΔLST variation (about 0.14 ± 0.06), while\ns s\nthe impact first pathway is weaker (0.08 ± 0.05). These findings suggest SM can modulate ΔLST\ns s\ndynamics through multiple mechanisms, with forest evaporative cooling regulation emerging as the\ndominant process. Previous studies have demonstrated that reduced soil moisture in Europe can\namplify warming during heatwaves through land-atmosphere interactions (Dirmeyer et al., 2021;\nFischer et al., 2007). Our findings confirm forests exhibit greater stability and resistance to\nincreasing soil dryness compared to grasslands or croplands, which may mitigate land-atmosphere\ncoupling feedback and consequently enhance the relative cooling effect of forests (Teuling et al.,\n2010). This stability of forests is reflected in both more moderate evaporative limitation and\nsustained vegetation greenness under increasing soil dryness. Our results provide observational\nevidence suggesting that the capacity of forest ecosystems to mitigate water deficit-induced warming\npersists not only during heatwave events but also across interannual timescales.\nFig 4. Mechanisms underlying the effects of climate variables on the daytime LST effects of\nEuropean forests (ΔLST) in summer. (a) Structure equation models describing the biophysical\ns\nrelationships between climate variables and ΔLST. The numbers denote the path value (mean ±\ns\nstandard error). The arrow color indicates the sign of the path (red for positive and blue for\nnegative), and the thickness indicates the magnitude of the corresponding path value. (b) Bar plot\nof the pathway effects of climate variables affecting ΔLST. The error bars indicate the uncertainty\ns\nbased on error propagation.\nWe also revised the discussion about the summer mechanism (Line 434).\nUsing a range of statistical methods, we find that declining SM is the predominant factor\ns\ncontributing to both the interannual variation and the negative trend in ΔLST. Path analysis shows\ns\nthat SM affects ΔLST dynamics through two mechanisms. First, we show ΔLAI is modulated by\ns s s\nSM, implying that as the soil becomes drier, the vegetation structure difference between forests and\ns\nopenlands becomes more evident. The higher ΔLAI could boost the cooling efficiency of forests\ns\nthrough stronger direct heat dissipation and evaporative cooling. The impact of SM on vegetation\ns\nhas been well-documented: increasing dryness pressure could trigger stomatal closure, inhibit\nphotosynthesis, and increase vegetation mortality (Li et al., 2022; Stocker et al., 2018). Our results\nconfirm a stronger negative effect of SM on non-forest vegetation than on forests (Hoek van Dijke\ns\net al., 2023). The potential reason is that forest ecosystems with high species richness and deep\nroots are more stable and resilient to dryness (Isbell et al., 2015; Tumber-Dávila et al., 2022).\nSecond, even under constant vegetation conditions, less soil water supply directly amplifies ΔLE\ns\nand enhances the relative cooling effect of forests. These findings are consistent with prior\nresearches (Fu et al., 2022a, 2022b), indicating that soil moisture imposes stronger limitations on\nET in grasslands and croplands relative to forest ecosystems. This differential response is evidenced\nby both the greater sensitivity of and the higher critical soil moisture thresholds in non-forest\necosystems. This phenomenon may be attributed to the capacity of forests to maintain ET through\ndeep soil water extraction, thereby exhibiting greater resilience to warming-induced surface soil\nmoisture deficits compared to openlands. Overall, the revealed impact of SM on ΔLST essentially\ns s\nreflects the diverse vegetation physiological responses to climate changes.\nReference:\nDirmeyer, P. A., Balsamo, G., Blyth, E. M., Morrison, R. & Cooper, H. M. Land‐Atmosphere\nInteractions Exacerbated the Drought and Heatwave Over Northern Europe During Summer\n2018. AGU Adv. 2, 1–16 (2021).\nFischer, E. M., Seneviratne, S. I., Lüthi, D. & Schär, C. Contribution of land-atmosphere coupling\nto recent European summer heat waves. Geophys. Res. Lett. 34, 1–6 (2007).\nFu, Z. et al. Critical soil moisture thresholds of plant water stress in terrestrial ecosystems. Sci. Adv.\n8, 1–12 (2022).\nFu, Z. et al. Uncovering the critical soil moisture thresholds of plant water stress for European\necosystems. Glob. Chang. Biol. 28, 2111–2123 (2022).\nHoek van Dijke, A. J. et al. Comparing forest and grassland drought responses inferred from eddy\ncovariance and Earth observation. Agric. For. Meteorol. 341, (2023).\nIsbell, F. et al. Biodiversity increases the resistance of ecosystem productivity to climate extremes.\nNature 526, 574–577 (2015).\nLi, W. et al. Widespread increasing vegetation sensitivity to soil moisture. Nat. Commun. 13, 3959\n(2022).\nLin, C. et al. Diel ecosystem conductance response to vapor pressure deficit is suboptimal and\nindependent of soil moisture. Agric. For. Meteorol. 250–251, 24–34 (2018).\nLiu, L. et al. Soil moisture dominates dryness stress on ecosystem production globally. Nat.\nCommun. 11, 1–9 (2020).\nNovick, K. A. et al. The increasing importance of atmospheric demand for ecosystem water and\ncarbon fluxes. Nat. Clim. Chang. 6, 1023–1027 (2016).\nShen, P. et al. Biodiversity buffers the response of spring leaf unfolding to climate warming. Nat.\nClim. Chang. 14, 863–868 (2024).\nStocker, B. D. et al. Quantifying soil moisture impacts on light use efficiency across biomes. New\nPhytol. 218, 1430–1449 (2018).\nTeuling, A. J. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci. 3, 722–727 (2010).\nTumber-Dávila, S. J., Schenk, H. J., Du, E. & Jackson, R. B. Plant sizes and shapes above and\nbelowground and their interactions with climate. New Phytol. 235, 1032–1056 (2022).\nZhong, Z. et al. Disentangling the effects of vapor pressure deficit on northern terrestrial vegetation\nproductivity. Sci. Adv. 9, (2023).\nFinally, I wonder why LEw never appears anywhere, but the effects of LEs are estimated...\nClearly LEw will be significantly lower. But conifers can evapotranspire almost year round... So,\nto some degree, LEw should also be available and should have an impact on surface\ntemperatures... ???\nResponse: Thank you for your insightful comment. We fully agree that evergreen forests can sustain\ntranspiration throughout the year. To investigate this potential effect, we calculated the winter latent\nheat flux (LE) difference between forests and openlands (ΔLE ), and then examined the relationship\nw\nbetween the winter land surface temperature effect of forests (ΔLST ) and ΔLE . Results indicate\nw w\nthat the correlation between ΔLE and ΔLST is statistically insignificant (Line 200). This finding\nw w\nsupports previous studies reporting that the albedo effect dominates the local temperature effect of\nboreal forests in the cold season (Alkama and Cescatti, 2016; Lawrence et al., 2022).\nSimilarly, we find negligible correlation between ΔLST and non-radiative processes\nw\nrepresented by winter latent heat flux differences (ΔLE\nw\n) (Supplementary Fig. 7).\nFig. S7. Temporal correlation coefficient (r) between LE and ΔLST . (a to c) Spatial pattern of r\nw w\nbetween LE and daily mean, daytime and nighttime ΔLST . (d) Density of the r values.\nw w\nReference:\nAlkama, R. & Cescatti, A. Biophysical climate impacts of recent changes in global forest cover.\nScience (80-. ). 351, 600–604 (2016).\nLawrence, D., Coe, M., Walker, W., Verchot, L. & Vandecar, K. The Unseen Effects of Deforestation:\nBiophysical Effects on Climate. Front. For. Glob. Chang. 5, 1–13 (2022).\nMore to the main point, however, the general finding that wintertime effects exhibit increased\ncooling due to the loss of snow cover, more or less independently of tree cover (the effect is\nstronger over open fields than over forested land), is an interesting and valuable finding.\nHowever, it remains unclear why this phenomenon would be true in more forested areas... In\nother words, why would the land surface become lighter under more forested areas, thus\nleading to wintertime cooling...??? I do not think this is adequately explained...?\nResponse: We sincerely appreciate the opportunity to clarify this important point. It is important to\nnote that forests have become “relatively” brighter compared to openlands during the past two\ndecades of snow decline. This relatively brighter trend is evidenced by the increasing Δα trend\nw\n(brown line in Fig. 2a), despite both ecosystems exhibiting absolute darkening trends. The relative\nbrightening stems from forests (green line) darkening at a slower rate than openlands (yellow line).\nThe mechanism is that snow tends to be masked by tree canopies in forest ecosystems, resulting in\nthe greater sensitivity of openland α to snow changes compared to forest α .\nw w\nFig 2a. Interannual variations of winter albedo effect of European forests (Δα ) and the individual\nw\nwinter albedo (α ) for forests and openlands.\nw\nTo prevent misinterpretation, we have modified the result description (Line 166):\nThese results indicate that although forest albedo tends to decline due to snow cover loss, the\nrelative darkening effect of forests, and the additional solar radiation absorbed by forests (compared\nto adjacent openlands) have diminished in the recent two decades.\nAll-in-all however, I did find the manuscript to be very interesting and compelling... I have\nincluded a number of comments in the attached submissions pdf.\nResponse: Thank you for your positive feedback. Below are the detailed responses to the comments\nas provided in the PDF document.\nReviewer #2’s comments in the PDF documents:\nLine 17: incomplete sentence...\nResponse: Thanks for pointing this out. We have revised the sentence and now it reads (Line 14):\nThe negative trends in LST effects of European forests are broadly supported by four state-of-\nthe-art earth system models, though substantial inter-model variability persists in the multi-year\nmeans.\nLine 24: There are many much better references for supporting this point... For example,\nneither reference noted here talks about advantages of forests for water availability in\ndownwind areas...???\nA few options to consider:\nCool Insights https://bit.ly/40zlLv2\nOn the Forest Cover- Water Yield Debate https://bit.ly/40t3fUM\nFrom Myth to Concept and Beyond https://bit.ly/4hAenGp\nResponse: Thank you for pointing this out. We have replaced the original reference with those you\nrecommended (Line 23).\nForests provide numerous benefits to humans and the planet, such as producing food and\nenergy, reducing soil erosion, and increasing water availability in downwind areas (Bonan, 2008;\nBorrelli et al., 2017; Ellison et al., 2017, 2012).\nReference:\nEllison, D., Futter, M. N. & Bishop, K. On the forest cover-water yield debate: From demand- to\nsupply-side thinking. Global Change Biology vol. 18 806–820 (2012).\nBorrelli, P. et al. An assessment of the global impact of 21st century land use change on soil erosion.\nNat. Commun. 8, (2017).\nEllison, D. et al. Trees, forests and water: Cool insights for a hot world. Glob. Environ. Chang. 43,\n51–61 (2017).\nBonan, G. B. Forests and climate change: Forcings, feedbacks, and the climate benefits of forests.\nScience (80-. ). 320, 1444–1449 (2008).\nLine 26: Again, the cited references here do not really discuss and defend impacts global\nclimate cooling... But there are several appropriate references that do...\nEven Cooler Insights https://bit.ly/46PTpzi\nhttps://doi.org/10.1073/pnas.2026241118\nhttps://www.nature.com/articles/s41467-024-51783-y\nResponse: Thank you for your valuable suggestion. We have replaced the original reference with\nthose you recommended (Line 25).\nThrough conservation, proper management, and restoration practices, forest ecosystems could\nplay a crucial role in addressing global warming by regulating the hydrological cycle and\nenhancing cloud formation (Cerasoli et al., 2021; Ellison et al., 2024; Luo et al., 2024).\nReference:\nCerasoli, S., Yin, J. & Porporato, A. Cloud cooling effects of afforestation and reforestation at\nmidlatitudes. Proc. Natl. Acad. Sci. U. S. A. 118, 1–7 (2021).\nLuo, H., Quaas, J. & Han, Y. Decreased cloud cover partially offsets the cooling effects of surface\nalbedo change due to deforestation. Nat. Commun. 15, 7345 (2024).\nEllison, D., Pokorný, J. & Wild, M. Even cooler insights: On the power of forests to (water the Earth\nand) cool the planet. Glob. Chang. Biol. 30, 1–20 (2024).\nLine 27: While certainly the role of carbon is \"overemphasized\" in this literature, there is\nin fact a far more extensive literature... Much of it, for example, is referenced and detailed in\nthis paper:\nEven Cooler Insights https://bit.ly/46PTpzi\nI see that a small part of this literature does appear in the following sentences, which\nsuggests that the sentence itself here is primarily misleading...\nResponse: Thank you for point this out. For clarity, we have removed the detailed statements about\nthe carbon sink effects of forests. Now the background part reads (Line 23):\nForests provide numerous benefits to humans and the planet, such as producing food and\nenergy, reducing soil erosion, and increasing water availability in downwind areas. Through\nconservation, proper management, and restoration practices, forest ecosystems could play a crucial\nrole in addressing global warming by regulating the hydrological cycle and enhancing cloud\nformation. In addition to the global impacts, forests exhibit local cooling or warming effects due to\ntheir unique biophysical characteristics that regulate surface energy and water fluxes. The\nbiophysical effects of forests on local climate have garnered particular attention in recent years, as\nit can significantly intensify or offset the climate mitigation effects of carbon sequestration.\nLine 45: Vague..., in what ways...? I would have put the emphasis on the availability of water\nand how that might change with some climate driven processes and temperature\nchange...???\nAgain, see: Even Cooler Insights https://bit.ly/46PTpzi\nResponse: Thanks for the comment. We have revised the discussions about the background climate\nand added sentence about the impact of water availability (Line 42):\nThe net temperature response depends on the balance between the opposing radiative and\nnon-radiative effects (Bright et al., 2017), and their relative dominance is largely governed by\nbackground climatic conditions (Pitman et al., 2011). For instance, forests in mid-latitude regions\n(e.g., North America, Europe and China) demonstrate seasonally contrasting LST effects (Li et al.,\n2015; Tang et al., 2018). During cold winters, bright snow covers short vegetation but is masked\nby dense forest canopies. As a result, forests absorb more solar radiation and show a local\nwarming effect compared to surrounding openlands. Conversely, during the peak growing season\n(warm summer), forests exert cooling effects under sufficient soil moisture conditions, owing to\ntheir higher ET rates compared to other vegetation types. However, in arid regions, limited soil\nmoisture can constrain forest evaporative cooling despite abundant radiation and optimal\ntemperature, causing radiative warming to dominate (Rohatyn et al., 2023).\nReference:\nBright, R. M. et al. Local temperature response to land cover and management change driven by\nnon-radiative processes. Nat. Clim. Chang. 7, 296–302 (2017).\nPitman, A. J. et al. Importance of background climate in determining impact of land-cover change\non regional climate. Nat. Clim. Chang. 1, 472–475 (2011).\nLi, Y. et al. Local cooling and warming effects of forests based on satellite observations. Nat.\nCommun. 6, 1–10 (2015).\nTang, B., Zhao, X. & Zhao, W. Local effects of forests on temperatures across Europe. Remote Sens.\n10, 1–24 (2018).\nRohatyn, S., Rotenberg, E., Tatarinov, F., Carmel, Y. & Yakir, D. Large variations in afforestation-\nrelated climate cooling and warming effects across short distances. Commun. Earth Environ.\n4, 1–10 (2023).\nLine 56: This statement is quite broad...??? Platitude...??? Can you be more precise, more\nnuanced...?\nResponse: Thank for the comment. We have revised the sentence and now it reads (Line 55):\nGiven the profound influence of background climate on forest biophysical processes, it is\nreasonable to infer that the regulation of forests on local temperature will shift in response to rising\nCO concentration, global warming, or varying water availability.\ndiscussions of radiative and non-radiative impacts more clearly delineated and categorized\nfrom the beginning...\nResponse: Thanks for the suggestion. We have clarified the radiative and non-radiative impacts in\nthe introduction section (Line 39):\nDuring the daytime, forests exhibit two competing effects compared to non-forest vegetation: (a)\nthe non-radiative cooling effect driven by higher evapotranspiration (ET) rates (Mildrexler et al.,\n2011), and (b) the radiative warming effect resulting from lower albedo (R.A.Betts, 2000).\nReference:\nMildrexler, D. J., Zhao, M. & Running, S. W. A global comparison between station air temperatures and\nMODIS land surface temperatures reveals the cooling role of forests. J. Geophys. Res.\nBiogeosciences 116, 1–15 (2011).\nR.A.Betts. Offset of the potential carbon sink from boreal forestation by decreases in surface albedo.\nNature 408, 187–190 (2000).\nLine 104: Up to now, nighttime effects have NOT been previously mentioned... However,\nyou note that such effects have been found previously, and that there is a literature... Can\nyou motivate these effects...? And shouldn't they also be discussed in to introductory\nsections above...? The drivers of these effects need to be carefully founded and\ndescribed...???\nResponse: Thanks for pointing this out. We have added description about the nighttime LST effect\nof forests in the introduction section (Line 36):\nAt night, forests warm the land surface through two mechanisms: (1) the release of daytime\nheat storage during nighttime, and (2) reduced aerodynamic resistance that enhances turbulent\nmixing, facilitating downward heat transfer from warmer air aloft (Lee et al., 2011; Schultz et al.,\n2017).\nReference:\nLee, X. et al. Observed increase in local cooling effect of deforestation at higher latitudes. Nature\n479, 384–387 (2011).\nSchultz, N. M., Lawrence, P. J. & Lee, X. Global satellite data highlights the diurnal asymmetry of\nthe surface temperature response to deforestation. J. Geophys. Res. Biogeosciences 122, 903–\n917 (2017).\nLine 150: Could this relationship also be confirmed via other methods, such as NDVI...? And\nmore generally, why would forested areas become lighter in color over time, if snow cover is\ndeclining...? This seems counter-intuitive and contradictory...???\nResponse: Sorry for the misleading. Our results demonstrate that while both forested and open areas\nexperience surface darkening due to snow cover decline, the darkening effect is more pronounced\nin openlands. This differential response means forests have become \"relatively brighter\" compared\nto adjacent openlands over time. We have clarified this point in the revised manuscript (Line 166):\nThese results indicate that although forest albedo tends to decline due to snow cover loss, the\nrelative darkening effect of forests, and the additional solar radiation absorbed by forests (compared\nto adjacent openlands) have diminished in the recent two decades.\nLine 205: VPD is of course also tied to the availability of water on the land surface... If less\nwater can evaporate, VPD will rise...\nResponse: Thanks for the professional comment. Since VPD and availability of water (soil moisture)\nare strongly coupled at the annual scale, we have added new experiments disentangle the\ncontributions of increasing VPD and decreasing soil moisture in the revised manuscript. Please see\nthe response to your major comment about soil moisture above.\nLine 226: Clearly some of these variables, in particular VPD and T, exhibit significant\ncolinearity... How is this handled... And how do you know which of these two variables is the\ndominant factor...???\nResponse: Thank you for raising this methodological consideration. We fully acknowledge the\ncollinearity between climate variables, particularly VPD and T. To address this inherent covariation\namong meteorological drivers, we employed partial correlation analysis, which is a robust approach\nthat isolates the unique relationship between the driver and the target variable (e.g., VPD and ΔLST)\ns\nwhile statistically controlling for confounding effects of covariates (e.g., T, precipitation,\nradiation…). By removing the influence of covariance, the residual VPD-ΔLST correlation reflects\ns\nmoisture-driven effects independent of thermal or radiation contributions, as well as any non-linear\nor interactive relationships not captured by simple linear covariation with temperature.\nUsing partial correlation, we reveal a critical distinction (Fig. 3a-b): while direct correlations\n(light bars) showed significant relationships between ΔLST and all climate variables (except wind\ns\nspeed), partial correlations (dark bars) demonstrated that VPD or soil moisture (SM) maintain the\nstrongest independent associations with ΔLSTs even after accounting for other drivers. It is\nimportant to note that partial correlation cannot fully disentangle the effects of SM and VPD when\nboth two variables are included in one model simultaneously, as they exhibit near-perfect\ncollinearity at the annual scale (r = 0.92). That’s why we add other experiments to disentangle their\ncontributions.\nFig 3a and 3b. Regional mean correlation and partial correlation coefficients between climate\nvariables and the daytime ΔLSTs.\nWe have enhanced the methods section to provide a detailed explanation of our partial temporal\ncorrelation analysis approach (Line 696):\nIn addition, partial temporal correlation analysis is used to explore the relative impacts of multiple\nclimate drivers on ΔLST at the pixel scale (pingouin package in Python). The partial correlation can\ns\ndetermine whether a climate variable can exert an independent effect on target variable beyond its\ncovariation with other climate drivers. For instance, the partial correlation can isolate the association\nbetween ΔLST and SM while accounting for the influence of other covariates (e.g., temperature and\ns s\nprecipitation) on SM.\ns\nLine 266: Are there any estimates of LEw...? Evergreen forests will evapotranspire almost\nyear round...??? So what does this relationship look like...?\nResponse: Thank you for your insightful comment. We examined the relationship between ΔLST\nw\nand ΔLE (Fig. R5). Results indicate that the correlation between ΔLE and ΔLST is statistically\nw w w\ninsignificant (Fig. R5a). In contrast, Δα exhibits a significant correlation with ΔLST (Fig. R5b).\nw w\nThese findings suggest that temporal variations in radiative processes, rather than\nevapotranspiration-related mechanisms, are more likely to drive the observed changes in ΔLST .\nw\nFig. R5. Temporal relationship between ΔLST and ΔLE or Δα .\nw w w\nIn addition, we observed an unexpected negative ΔLE , indicating that forest latent heat flux\nw\nwas slightly lower than that of openlands (multi-year mean: -1.3 w/m2). This finding is aligned with\nprevious study based on land surface model simulation and remote sensing data (Breil et al., 2021;\nZan et al., 2024). To further assess this result, we compared winter LE observations from European\nforest and openland (grassland and cropland) sites using FLUXNET2015 dataset (Fig. R6). The\nresults revealed that, on average, openland LE exceeded forest LE, though the difference was not\nstatistically significant. This result ensures the reliability of our remote sensing-based results.\nFig. R6. Winter latent heat flux (LE) at European forest and openland sites based on FLUXNET\nobservations. Each data point represents the site mean. Red asterisk indicates the muti-site mean.\nIn summary, our results indicate that while wintertime evapotranspiration persists in forests,\nits magnitude does not differ significantly from openlands. The temporal variation in ΔLST appears\nw\nto be driven primarily by radiative processes rather than LE-related mechanisms.\nReference:\nBreil, M., Davin, E. L. & Rechid, D. What determines the sign of the evapotranspiration response\nto afforestation in European summer? Biogeosciences 18, 1499–1510 (2021).\nZan, B. et al. Spatiotemporal inequality in land water availability amplified by global tree\nrestoration. Nat. Water (2024) doi:10.1038/s44221-024-00296-5.\nLine 278: The availability of water, however, is also key here... If the water disappears, then\nthese effects will flip... So, how does one integrate the effects of water availability into these\nmeasures...???\nResponse: Thanks for pointing this issue. Please see our response to your comment above, in which\nthe role of water availability (soil moisture) is explored and discussed in detail.\nLine 391: This could be better expressed... One could just as easily substitute southwards\nhere, since you do not clearly specify distinctions between cooling and warming effects...???\nResponse: Thanks for your suggestion about the expression. We have modified the sentence and\nnow it reads (Line 413):\nThis result implies that forests in mid- to high-latitude regions could exhibit similar reduced\nwarming effects or even reversed net cooling in a warming world. Such changes may push the\nlatitudinal boundary between forest warming and cooling effects poleward, primarily due to\ndiminished snow-related albedo feedbacks."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have made significant improvements, conducted substantial additional analyses, and even revised their\nconclusions regarding the role of VPD. I am satisfied with most of these efforts and have only one remaining comment.\n\"In the summertime, the vegetation physiological response to decreased soil moisture boosts the cooling effect of forests.\"\nThis statement may still be confusing or even misleading. It implies that forests have a stronger cooling effect than\ngrasslands solely due to their physiological response to dry soils.\nIn reality, while physiological processes such as stomatal regulation do influence surface energy balance, the observed\ncooling effect in forests may also stem from hydrological advantages, including deeper soil water storage and more\nconsistent precipitation recharge in forested areas.", "reviewer_2": "I find that a lot of interesting points come up in this analysis. I am, however, a bit concerned that the paper could still benefit\nfrom some significant editing and rewriting. I raise a number of questions in the attached submission pdf. I guess that,\nprimarily, I am unconvinced by the forest albedo logic used herein and would be tempted to argue in a different direction on\nmany points. Please see the related comments in the text. Further, I do think that the paper seems inclined to somehow lay\nthe blame for these changes on forests, when actually, most of the important changes seem to be occurring elsewhere (in\nopen fields that have reduced amounts of forest cover, and not in heavily forested areas) and are primarily driven by climate\nwarming (and not any change in forests themselves)....\nI do think that with another round of editing, this paper can be significantly improved... I am uncertain, however, how one\nshould resolve the problem of over-emphasising forest albedo effects. This would ultimately require some different\nmeasures. Perhaps the Authors can respond to this point...? What do they think is the best procedure here...? Personally, I\nwould be inclined to take this in a somewhat different direction. Moreover, while the forest albedo literature considers albedo\neffects, it entirely neglects cloud related effects... Since these are presumably also tied to forests, in part, what does this say\nabout the methodology employed in the analysis...???"}, "author_response": "Response to Reviewers’ Comments\nWe greatly appreciate the anonymous reviewers for the constructive comments. The manuscript has\nbeen revised and improved accordingly. We hope that the revision will make it more acceptable for\npublication. Below are the point-by-point responses to the comments, along with the revision of the\nmanuscript (typed in red and Italic). The changes are also marked in red in the manuscript file.\nReviewer #1 (Remarks to the Author):\nThe authors have made significant improvements, conducted substantial additional analyses,\nand even revised their conclusions regarding the role of VPD. I am satisfied with most of these\nefforts and have only one remaining comment.\nResponse: We appreciate the constructive and valuable comments by the reviewer during the review\nprocess.\n\"In the summertime, the vegetation physiological response to decreased soil moisture boosts\nthe cooling effect of forests.\"\nThis statement may still be confusing or even misleading. It implies that forests have a\nstronger cooling effect than grasslands solely due to their physiological response to dry soils.\nIn reality, while physiological processes such as stomatal regulation do influence surface\nenergy balance, the observed cooling effect in forests may also stem from hydrological\nadvantages, including deeper soil water storage and more consistent precipitation recharge\nin forested areas.\nResponse: Thanks for pointing this out. We have revised the sentence to accurately describe this\nmechanism (Line 13):\nDuring summer, the daytime cooling effect intensifies (-0.188 K/decade) because forests remain\nmore physiologically and hydrologically resilient to increasing soil dryness.\nReviewer #2 (Remarks to the Author):\nI find that a lot of interesting points come up in this analysis. I am, however, a bit concerned\nthat the paper could still benefit from some significant editing and rewriting. I raise a number\nof questions in the attached submission pdf. I guess that, primarily, I am unconvinced by the\nforest albedo logic used herein and would be tempted to argue in a different direction on\nmany points. Please see the related comments in the text.\nResponse: We sincerely appreciate your positive and constructive feedback. We understand your\nconcerns regarding our interpretation of the winter land surface temperature effect of conversion\nfrom openlands to forests (ΔLST ). In response to your feedback, we have conducted additional\nw\nexperiments to quantitatively compare the contributions of albedo- and ET-related processes to the\nΔLST trend. We identify the albedo-related process, rather than the ET-related process, as the main\nw\ndriver of ΔLST change. Please see the response to your later comments on line 115.\nw\nFurther, I do think that the paper seems inclined to somehow lay the blame for these changes\non forests, when actually, most of the important changes seem to be occurring elsewhere (in\nopen fields that have reduced amounts of forest cover, and not in heavily forested areas) and\nare primarily driven by climate warming (and not any change in forests themselves)....\nResponse: Thank you for pointing this out. Our study analyzes differences between forests and\nadjacent openlands (Δx) to derive a priori estimates of potential forestation effects (openland-to-\nforest conversion). We fully agree with your view that the change in Δx should not be exclusively\nattributed to forests. The Δx variations may stem from either forest modifications or alterations in\nthe reference openlands. To ensure clarity, we employ the term “forestation” to emphasize that our\nanalysis specifically examines forest-openland comparisons. We have carefully reviewed the\nmanuscript and modified ambiguous descriptions accordingly.\nFor instance, we have revised the title (Line 1):\nAmplified Local Cooling Effect of Forestation in Warming Europe\nWe also revised the sentence in the abstract (Line 7):\nBy comparing the satellite observations of spatially nearby forests and openlands over the last\ntwo decades, we quantify the temporal trends in local land surface temperature (LST) effects of\npotential forestation in Europe.\nWe added a sentence to highlight that the change is more significant in openlands in winter\nresults (Line 189):\nNotably, the observed ΔLST trend should not be attributed to structural or physiological\nw\nchanges in forests, but rather to the more pronounced albedo reduction in openlands in a warming\nworld.\nFurthermore, we added a discussion to clarify that the change in ΔLST may originate from\nboth forests or openlands (Line 505):\nFirst, the observed negative trends in ΔLST and ΔLST reflect variations in LST difference\nw s\nwhen comparing forests and openlands, rather than the inherent capacity of existing forests to\nreverse warming trends. These patterns may originate from either modifications in forest\ncharacteristics or changes in the adjacent openland used as reference. Fundamentally, the negative\ntrends emerge because forests exhibit a weaker temperature response than openlands in a warming\nclimate, and this divergence stems from distinct biophysical and physiological properties between\nforest and openland vegetation.\nWe hope these revisions will enhance the readability.\nI do think that with another round of editing, this paper can be significantly improved... I am\nuncertain, however, how one should resolve the problem of over-emphasising forest albedo\neffects. This would ultimately require some different measures. Perhaps the Authors can\nrespond to this point...? What do they think is the best procedure here...? Personally, I would\nbe inclined to take this in a somewhat different direction.\nResponse: We sincerely appreciate your insightful suggestion regarding the potential mechanism\nof the ΔLST variations. In response to your valuable feedback, we have added experiments to\nw\ncompare the albedo- and ET-related effect based on the energy balance equation. This analysis is to\nquantitatively assess the relative importance of these two processes in ΔLST changes. We also added\na discussion about the ET-related mechanism. Please see the detailed response to your later\ncomments on line 115.\nMoreover, while the forest albedo literature considers albedo effects, it entirely neglects cloud\nrelated effects... Since these are presumably also tied to forests, in part, what does this say\nabout the methodology employed in the analysis...???\nResponse: We appreciate your insightful critique regarding cloud-related effects. We fully agree\nthat cloud feedback is an important pathway of forest affecting the climate system. As you rightly\npointed out, the space-for-time substitution method employed in our study primarily captures forests’\ndirect modulation of surface albedo rather than planetary albedo. We would like to note that several\nstudies have qualitatively assessed forests’ impact on cloud formation using similar methodologies\ncombined with satellite observations (Duveiller et al., 2021; Xu et al., 2022). However,\nquantitatively evaluating the equivalent radiative forcing of forest-cloud feedback remains\nchallenging.\nIf the remote sensing-based results are properly interpreted, we believe that our methodological\napproach and findings are scientifically valid and useful. First, while the method focuses on the\nlocal effects, it provides accurate assessments that can serve as valuable benchmarks for improving\nclimate models and reanalysis datasets regarding the land-atmosphere interaction, which is a\nfundamental prerequisite for understanding forests’ influence on cloud formation. Previous\nmodeling studies examining forest cover impacts on cloud-radiative effects have shown\nconsiderable uncertainty in their results (Cerasoli et al., 2021; Luo et al., 2024). Second, the used\nmethod captures the direct impacts of forest on local LST, which are effective across spatial scales.\nIn contrast, forest-cloud feedback operates at the regional scale and depends on the intensity of\nforest change. For instance, even forestation/deforestation in a small area can immediately affect\nLST, while cloud impacts may be masked by the atmospheric advection. Thus, our assessments\nprovide insights for climate adaptation strategies at the local scale. These findings complement\nrather than contradict other regional or larger-scale modeling studies. We have added a discussion\nabout this issue (Line 504):\nSecond, the adopted “space-for-time” method provides a priori estimates of the local\ntemperature effect. The influence of forests on regional cloud formation and subsequent downward\nradiation represents a significant yet overlooked mechanism in current assessments. This omission\nmay lead to a systematic underestimation of forests’ net cooling effect. Meanwhile, forests can exert\nmore complicated effects on climate through multiple processes such as changes in large-scale\natmospheric circulations. Coupled models can effectively estimate these effects, notwithstanding the\nconsiderable uncertainties. Thus, our results should be treated as the reference for local climate\nadaptation rather than global climate mitigation.\nReference:\nCerasoli, S., Yin, J. & Porporato, A. Cloud cooling effects of afforestation and reforestation at\nmidlatitudes. Proc. Natl. Acad. Sci. U. S. A. 118, 1–7 (2021).\nDuveiller, G. et al. Revealing the widespread potential of forests to increase low level cloud cover.\nNat. Commun. 12, 1–15 (2021).\nLuo, H., Quaas, J. & Han, Y. Decreased cloud cover partially offsets the cooling effects of surface\nalbedo change due to deforestation. Nat. Commun. 15, 7345 (2024).\nXu, R. et al. Contrasting impacts of forests on cloud cover based on satellite observations. Nat.\nCommun. 13, 670 (2022).\nLine 11: Since this observation runs in the face of quite a large amount of literature which\nsuggests the opposite relationship, perhaps this is worth highlighting, or perhaps this should\nbe stated in a more nuanced fashion....? See for example some of the literature highlighted in\nEllison et al., 2024, perhaps, in particular, Davin et al.,\nhttp://journals.ametsoc.org/doi/10.1175/2009JCLI3102.1 or, Lawrence et al., 2022\nResponse: Thank you for pointing this out. Our findings reveal a winter daytime warming effect in\nforests compared to openlands initially. However, this daytime warming gradually diminishes and\neventually reverses, while the nighttime warming effect persists. For clarity, we have revised the\nrelevant sentence (Line 7):\nBy comparing the satellite observations of spatially nearby forests and openlands over the last\ntwo decades, we quantify the temporal trends in local land surface temperature (LST) effects of\npotential forestation in Europe. During winter, the daytime warming effect of forestation diminishes\nand reverses to cooling (-0.142 K/decade) due to weaker surface darkening trends in forests than\nopenlands, induced by decreasing snow cover.\nLine 26: Could perhaps highlight relative importance of enhancing/increasing forest cover...,\nsince only change in forest cover is likely to help redress global warming...\nResponse: Thank you for the suggestion. We have revised the sentence and now it reads (Line 26):\nIncreasing forest cover through restoration practices could play a crucial role in addressing\nglobal warming by regulating the hydrological cycle and enhancing cloud formation.\nLine 39: ???...\nResponse: Thank you for the comment. We have expanded the description of this mechanism (Line\n37):\nAt night, forests warm the land surface compared to open lands through two mechanisms. First,\nforest canopies absorb solar radiation and store heat during the daytime and subsequently release\nthis energy at night. Second, the nighttime land surface is typically cooler than the near-surface air,\nforming a stable boundary layer where downward atmospheric heat fluxes warm the land surface.\nForests, with their greater surface roughness, can intensify turbulent mixing, enhance this heat\nexchange process, leading to higher nighttime LST observed in forested areas.\nLine 41: I think this is not entirely correct... I agree this is the way we commonly think about\nthis relationship. But, consider the following: (b) only or primarily occurs in the absence of (a).\nThe net temperature effect is a function of the availability of soil moisture., not of the net\nradiative balance. ...?\nWarming occurs when soil moisture is depleted and is no longer able to mitigate albedo\neffects, otherwise the surface is kept from warming through ET production (energy absorption,\nlatent heat) and the transfer of heat energy into the lower atmosphere..\nSo warming/cooling is not measured by the net radiative effects, but rather of the relative\nduration of time during which soil moisture is available and not...??? (I am sure this can be\nbetter expressed...).\nThis is especially true at the micro/local LST scale at which you measure climate impacts here,\nbut ignores larger scale global climate impacts....\nResponse: Thanks for your insightful comment. We fully recognize the fundamental role of soil\nmoisture (SM) in governing vegetation water availability and its consequent influence on\nevapotranspiration-mediated cooling effects. However, we hold the view that this SM-related\nmechanism remains consistent with the established paradigm that “the net temperature response is\ndetermined by the balance between competing radiative and non-radiative processes”, which is\nbased on the surface energy balance (SEB) framework (Bright et al., 2017; Duveiller et al., 2018;\nRigden and Li, 2017). In the SEB equation, LST governs upward longwave radiation (LW↑), which\ncan be expressed as:\n𝐿𝑊 ↑(cid:3404)𝜀𝜎𝐿𝑆𝑇(cid:2872) (cid:3404)𝑆𝑊 ↓(cid:4666)1(cid:3398)𝛼(cid:4667)(cid:3397)𝐿𝑊 ↓(cid:3398)𝐿𝐸(cid:3398)𝐻(cid:3398)𝐺 (R1)\nwhere ε represents emissivity, σ denotes the Stefan-Boltzmann constant, SW↓ indicates downward\nshortwave radiation, α is albedo, LW↓ represents downward longwave radiation, and LE, H, and G\ncorrespond to latent heat, sensible heat, and ground heat flux, respectively.\nAssuming limited impact of forest cover on background downward radiation and neglecting\nthe emissivity differences between forests and open lands, the LST effect (ΔLST, forest minus open\nlands) can be expressed in differential form (with Δ prefix representing forest minus open land\ndifferences):\n(cid:2869)\n𝛥𝐿𝑆𝑇 (cid:3404) (cid:4666)(cid:3398)𝛥𝛼(cid:3401)𝑆𝑊 ↓(cid:3398)𝛥𝐿𝐸(cid:3398)𝛥(cid:4666)𝐻(cid:3397)𝐺(cid:4667)(cid:4667) (R2)\n(cid:2872)(cid:3097)(cid:3084)(cid:3013)(cid:3020)(cid:3021)(cid:3119)\nEquation R2 decomposes ΔLST into three components on the right side. The radiative effect is\ncaptured by the first, while non-radiative effects are represented by the remaining two terms, with -\nΔLE being the dominant component.\nThe control of SM can be well integrated into this framework. Specifically, SM is a critical\nenvironmental factor modulating transpiration rates of vegetation, and thus serves as a key driver of\nΔLE (i.e. 𝛥𝐿𝐸 (cid:3404)𝑓(cid:4666)𝑆𝑀(cid:4667) ). Through affecting ΔLE, SM can further influence both the sign and\nmagnitude of ΔLST. Specifically, forests maintain higher transpiration rates than non-forest\nvegetation (positive ΔLE) across a broad spectrum of soil moisture, resulting in a significant cooling\neffect of forests (negative ΔLST). However, if soil moisture falls below the wilting point of trees,\nthe forests’ evaporative advantage should be constrained, allowing albedo-driven radiative effects\nto dominate ΔLST.\nTherefore, we believe that the impact of SM does not conflict with the mechanism of the\nbalance between radiative and non-radiative effects. We emphasize that SM’s control operates\nthrough its influence on the radiative/non-radiative balance by regulating ΔLE. This creates a\nunified framework where background climatic or environmental conditions determine the feasibility\nof evaporative cooling, and the energy balance further quantifies the ΔLST expression.\nBased on your comment, we have revised the paragraph to highlight the role of background\nclimate and environmental conditions, including SM (Line 44):\nDuring the daytime, forests commonly exhibit two competing effects compared to non-forest\nvegetation: (a) the non-radiative cooling effect driven by higher evapotranspiration (ET) rates\n(when water supply is sufficient), and (b) the radiative warming effect resulting from lower albedo.\nThe net temperature response depends on the balance between these opposing effects, which is\npredominantly governed by the background climatic or environmental conditions (e.g., temperature,\nsnow cover, and water availability). For instance, during boreal winter when ET is rendered inactive\nby low temperatures, the dense forest canopies obscure bright snow, which remains exposed in\nopenlands. This snow masking effect creates a strong albedo contrast between forests and openlands\nand causes forests to absorb more solar radiation compared to surrounding openlands, resulting in\nthe local warming effect. In warm and wet tropics, forests demonstrate pronounced local cooling\neffects due to their enhanced ET rates relative to other vegetation types. Limited soil moisture in\narid regions or during drought events can constrain forests’ evaporative advantage, enabling\nradiative warming to dominate. Consequently, the magnitude of forests’ cooling effect at seasonal\nor annual scale may depends on the duration of sufficient water availability.\nRegarding the scale-dependency of climate effects, we acknowledge that our remote sensing-\nbased approach has limitations in capturing forest-climate effect at larger scales. We have explicitly\naddressed this constraint in our discussion and emphasized the importance of proper interpretation\nof our findings (Line 517):\nMeanwhile, forests can exert more complicated effects on climate through multiple processes\nsuch as changes in large-scale atmospheric circulations. Coupled models can effectively estimate\nthese effects, notwithstanding the considerable uncertainties. Thus, our results should be treated as\nthe reference for local climate adaptation rather than global climate mitigation.\nReference:\nBright, R. M. et al. Local temperature response to land cover and management change driven by\nnon-radiative processes. Nat. Clim. Chang. 7, 296–302 (2017).\nDuveiller, G., Hooker, J. & Cescatti, A. The mark of vegetation change on Earth’s surface energy\nbalance. Nat. Commun. 9, 64–75 (2018).\nRigden, A. J. & Li, D. Attribution of surface temperature anomalies induced by land use and land\ncover changes. Geophys. Res. Lett. 44, 6814–6822 (2017).\nLine 49: At least as long as ET is rendered inactive by lower temperatures.... This could be\nmade more explicit... Note that coniferous trees evapotranspire for much longer periods of\ntime on the annual scale... ???\nResponse: Thanks for the professional comment. We have revised the sentence and now it reads\n(Line 51):\nFor instance, during boreal winter when ET is rendered inactive by low temperatures, the dense\nforest canopies obscure bright snow, which remains exposed in openlands. This snow masking effect\ncreates a strong albedo contrast between forests and openlands and causes forests to absorb more\nsolar radiation compared to surrounding openlands, resulting in the local warming effect.\nLine 115: I tend to think that these estimates should be paired with indications of the changes\nin wintertime ET behavior. One might expect, with coniferous species at least, that warming\nwill increase the potential for wintertime ET production, and therefore wintertime cooling\neffects from forests...?\nResponse: Thanks for the insightful comment. We agree with your opinion that in the warmer winter,\nthe evaporative cooling effect of forests (compared to openlands) may be intensified and lead to the\nobserved ΔLST trends. In response to your suggestion, we have conducted additional experiments\nw\nexamining this process (Supplementary text 1). Our results confirm that albedo-driven variations,\nrather than evapotranspiration-related processes, are primarily responsible for the observed trend in\nΔLST (Line 147).\nw\nDuring winter, the sign and absolute value of ΔLST are generally dominated by the magnitude\nw\nof the radiative process. However, the negative trend in ΔLST could potentially arise from two\nw\nmechanisms: (a) climate-driven alterations in net shortwave radiation absorption differences\nbetween forests and openlands (ΔNSR ), or (b) winter warming-induced enhancement in forests’\nw\nevaporative cooling advantage, represented by latent heat differences between the two land cover\ntypes (ΔLE ). To identify the dominant mechanism, we estimate the equivalent perturbations of the\nw\ntwo processes on the surface energy balance (Supplementary Text 1). Our analysis reveals that only\nΔNSR exhibits a statistically significant trend (-0.86 w·m-2/decade, p = 0.05), whereas the ΔLE\nw w\ntrend remains non-significant (-0.02 w·m-2/decade, p = 0.97) (Supplementary Fig. 3). Correlation\nanalyses of interannual variability further support the greater importance of radiative processes,\nand confirm that this effect is mediated through modifications in the albedo difference (Δα ) between\nw\nforests and openlands (Supplementary Figs. 4 to 6).\nSupplementary Text 1. Identifying the dominant process of ΔLST change\nw\nWe employ the surface energy balance (SEB) framework to identify the primary processes\ndriving variations in ΔLST . This framework has been widely adopted for decomposing the effects\nw\nof vegetation changes on local LST (Bright et al., 2017; Duveiller et al., 2018; Rigden and Li, 2017).\nIn the SEB equation, LST governs upward longwave radiation (LW↑), which can be expressed as:\n𝐿𝑊 ↑(cid:3404)𝜀𝜎𝐿𝑆𝑇(cid:2872) (cid:3404)𝑆𝑊 ↓(cid:4666)1(cid:3398)𝛼(cid:4667)(cid:3397)𝐿𝑊 ↓(cid:3398)𝐿𝐸(cid:3398)𝐻(cid:3398)𝐺 (1)\nwhere ε represents emissivity, σ denotes the Stefan-Boltzmann constant, SW↓ indicates downward\nshortwave radiation, α is albedo, LW↓ represents downward longwave radiation, and LE, H, and G\ncorrespond to latent heat, sensible heat, and ground heat flux, respectively.\nAssuming minimal forest cover influence on background downward radiation and negligible\nemissivity differences between forests and openlands, the LST effect (ΔLST, forest minus openlands)\ncan be expressed in differential form (with Δ prefix representing forest minus open land differences):\n(cid:2869)\n𝛥𝐿𝑆𝑇 (cid:3404) (cid:4666)(cid:3398)𝛥𝛼(cid:3401)𝑆𝑊 ↓(cid:3398)𝛥𝐿𝐸(cid:3398)𝛥(cid:4666)𝐻(cid:3397)𝐺(cid:4667)(cid:4667) (2)\n(cid:2872)(cid:3097)(cid:3084)(cid:3013)(cid:3020)(cid:3021)(cid:3119)\nEquation 2 decomposes ΔLST into three components on the right side. The radiative effect is\ncaptured by the first, while non-radiative effects are represented by the remaining two terms, with -\nΔLE being the dominant component.\nWe calculate winter differences in net solar radiation (ΔNSR = -Δα ·SW↓) and latent heat\nw w\n(ΔLE ) between spatially adjacent forests and openlands for each year, analogous to ΔLST . As\nw w\nshown in Supplementary Fig. 3a, positive ΔNSR values indicate a warming effect from potential\nw\nforestation. Surprisingly, we consistently observe negative ΔLE values (positive -ΔLE) across all\nyears, suggesting slightly lower latent heat flux in forests compared to openlands during winter.\nThis finding aligns with previous modeling and remote sensing studies (Breil et al., 2021; Zan et al.,\n2024), indicating that forests exhibit no evaporative cooling advantage over openlands in winter.\nTemporal analysis reveals a statistically insignificant trend in ΔLE (p = 0.97), while ΔNSR\nw w\nshowed a significant negative trend (p = 0.05) (Supplementary Fig. 3a). Similar results are obtained\nwhen comparing evergreen needleleaf forests (which maintain winter evapotranspiration) with\nopenlands (Supplementary Fig. 3b). These findings strongly suggest that the observed trend in\nΔLST is driven primarily by decreasing shortwave radiation absorption differences between forests\nw\nand openlands, rather than by changes in ΔLE .\nw\nFurthermore, we compute annual correlation coefficients between ΔLST (daily, daytime, and\nw\nnighttime) and both ΔNSR and ΔLE . A positive correlation emerged between daytime ΔLST and\nw w w\nΔNSR (r = 0.42; Supplementary Fig. 4), while correlations between ΔLST and ΔLE are\nw w w\nstatistically insignificant (Supplementary Fig. 5). These interannual correlation patterns\ncorroborate our trend-based conclusions.\nNotably, ΔNSR represents the product of both background downward solar radiation (DSR )\nw w\nand albedo differences (Δα ) between forests and openlands. Additional analysis reveals weak\nw\ncorrelation between DSR and ΔLST (Supplementary Fig. 6). These results suggest that Δα , rather\nw w w\nthan DSR , should serve as the principal driver of ΔLST variation.\nw w\nSupplementary Fig. 3. Temporal variations of winter net solar radiation difference (ΔNSR ) and\nw\nlatent heat difference (ΔLE ) between forests and openlands in Europe. (a) All forest samples are\nw\nused. (b) Only evergreen needleleaf forest samples are used.\nSupplementary Fig. 4. Temporal correlation coefficient (r) between ΔNSR and ΔLST . (a to c)\nw w\nSpatial pattern of r between ΔNSR and daily mean, daytime and nighttime ΔLST . (d) Density of\nw w\nthe r values.\nSupplementary Fig. 5. Temporal correlation coefficient (r) between ΔLE and ΔLST . (a to c)\nw w\nSpatial pattern of r between ΔLE and daily mean, daytime and nighttime ΔLST . (d) Density of the\nw w\nr values.\nSupplementary Fig. 6. Temporal correlation coefficient (r) between DSR and ΔLST . (a to c)\nw w\nSpatial pattern of r between DSR and daily mean, daytime and nighttime ΔLST . (d) Density of the\nw w\nr values.\nMoreover, we concur with your perspective that warmer winters could theoretically enhance\nthe evaporative cooling effect of forests relative to openlands, potentially contributing to ΔLST\nw\nvariations. Our results did not detect such signals in the current climate regime. This absence may\nbe attributed to the nonlinear response of LE to temperature increases. Specifically, current winter\nwarming trends in Europe remain insufficient to significantly stimulate transpiration in forests. In\nthe future, continued warming may eventually cross a critical threshold where temperature increases\nbecome sufficient to activate wintertime evapotranspiration processes. Under such a scenario, ET-\nrelated mechanisms could potentially emerge as a more dominant factor in governing winter land\nsurface temperature differences. We have added a discussion about this mechanism in the discussion\n(Line 199):\nFurthermore, though not yet observed, rising temperatures may exceed a critical threshold\nthat enhances cold-season ET in the future. This may strengthen the evaporative cooling effect of\nforests under future warming conditions, especially in coniferous forests with year-round\ntranspiration activity.\nReference:\nBreil, M., Davin, E. L. & Rechid, D. What determines the sign of the evapotranspiration response\nto afforestation in European summer? Biogeosciences 18, 1499–1510 (2021).\nBright, R. M. et al. Local temperature response to land cover and management change driven by\nnon-radiative processes. Nat. Clim. Chang. 7, 296–302 (2017).\nDuveiller, G., Hooker, J. & Cescatti, A. The mark of vegetation change on Earth’s surface energy\nbalance. Nat. Commun. 9, 64–75 (2018).\nRigden, A. J. & Li, D. Attribution of surface temperature anomalies induced by land use and land\ncover changes. Geophys. Res. Lett. 44, 6814–6822 (2017).\nZan, B. et al. Spatiotemporal inequality in land water availability amplified by global tree restoration.\nNat. Water (2024) doi:10.1038/s44221-024-00296-5.\nMoreover, since these effects are small overall, one wonders how much they really matter in\nthe larger picture... Should we worry more about summertime warming (and the related tree\ncooling potential), or should we be more concerned about wintertime changes overall...???\nWhich of these two phenomena has the larger, more important impact on climate change\nphenomena...?\nResponse: We appreciate your insightful comments. We acknowledge that the observed ΔLSTtrend\nis limited compared to the magnitude of climate change, as we discussed in the manuscript (Line\n513):\nThird, the enhanced forest local cooling effect implies that planting trees in proper areas\nremains a promising local solution against the risk of warming, especially in highly populated\nregions. Nonetheless, tree restoration is not a panacea for climate change. Although our results\nhighlight the negative biophysical feedback to the warming trend in forests, this feedback is unable\nto offset the warming trend driven by rising atmospheric CO .\nmagnitude of ΔLST and ΔLST changes are equivalent to about 2~6% and 6~9% climate warming\nw s\ntrend in winter and summer, respectively (Supplementary Text 4). We recognize that summer effects\nmay remain more directly relevant for human comfort and safety in a warming world. However, the\nwinter effects operate through distinct mechanisms requiring separate consideration, and carry\nimportant ecological implications that should not be overlooked.\nFirst, the revealed relationship between snow and ΔLST could provide constraints for climate\nw\nmodels predicting vegetation-snow interactions. Second, there is an ongoing debate regarding the\nsuitability of afforestation for mid- to high latitudes due to the undesired warming effect of\nafforestation in the cold season. Our result suggests the potential for expanded forestation feasibility\nin such regions. Specifically, the weakening winter warming signal suggests we may need to re-\nevaluate traditional limitations on forestation in mid- or high-latitudes, as projected snow cover\nreductions under climate change could further diminish albedo-related positive forcing and\npotentially expand suitable areas for climate mitigation through forestation (as we discussed in Line\n449). The observed changes in ΔLST provide valuable insights into how forestation strategies\nw\nmight need to adapt to changing snow regimes under global warming. Overall, the understanding\nof seasonal forest-climate interactions contributes to more geographically targeted and feasible\nforestation approaches.\nLine 145: From my perspective, as noted also above, it is not obvious that this should be true.\nResponse: Thanks for the comment. We have added additional results to identify the mechanism.\nPlease see the detailed response to your comments on line 115.\nLine 149: Empirically, these findings correspond well with what I would expect...\nResponse: We appreciate the reviewers’ recognition of our findings. We refined this sentence to\nidentify the direction of change in Δα (Line 162):\nw\nSpecifically, we show that forests are darker than openlands, as evidenced by the negative Δα\nw\nvalue (-0.077, Supplementary Fig. 7). However, Δα shows a positive trend (0.018 unitless/decade,\nw\np=0.04) (Fig. 2a), indicating that the albedo contrast between forests and openlands has been\ndiminishing.\nLine 155: I think one still has to demonstrate that these differences are not true in the longer\nterm...? Has the direction of change actually changed, or only become more pronounced....?\nResponse: We appreciate your comments regarding the long-term dynamics of snow-albedo effects.\nOur satellite-based analysis demonstrates a consistent pattern over the past two decades (2004-2023):\nboth forests and openlands have exhibited significant albedo declines, with a less pronounced\nreduction in forests (-0.017 unitless/decade, p<0.05) compared to openlands (-0.036 unitless/decade,\np<0.05). This divergence results in a significant positive trend in Δα (forest minus openland),\nw\nreducing the albedo contrast between these land cover types.\nWe understand your potential concerns about the long-term persistence of these snow-albedo\neffects. While our MODIS-based analysis is constrained by the temporal limitations of the satellite\nrecord, our findings align well with established literature. Previous studies have shown stronger\nsnow-albedo feedback mechanisms in openlands compared to vegetated surfaces (Alessandri et al.,\n2020; Yu et al., 2022), supporting our observed patterns. This difference arises primarily from the\ngreater snow-masking effect of forest canopies, which diminishes snow’s contribution to surface\nalbedo relative to openlands. We have incorporated these supporting references and refined our\ndescription accordingly (Line 168):\nThe revealed positive Δα trend can be traced back to the stronger albedo decrease in\nw\nopenlands (multi-year mean: 0.236; trend: -0.036 unitless/decade, p=0.04) than in forests (multi-\nyear mean: 0.159; trend: -0.017 unitless/decade, p=0.03) from 2004 to 2023 (Fig. 2a). This\nobserved divergence in surface darkening trends may be associated with snow cover (SC )\nw\nreductions (Supplementary Fig. 8), and differences in the snow-albedo feedback intensity in different\nland covers (Alessandri et al., 2020; Atlaskina et al., 2015; Li et al., 2017; Manninen et al., 2019;\nYu et al., 2022).\nReference:\nAlessandri, A., Catalano, F., De Felice, M., van den Hurk, B. & Balsamo, G. Varying snow and\nvegetation signatures of surface albedo feedback on the Northern Hemisphere land warming.\nEnviron. Res. Lett. 16, (2020).\nAtlaskina, K., Berninger, F. & de Leeuw, G. Satellite observations of changes in snow-covered land\nsurface albedo during spring in the Northern Hemisphere. Cryosph. 9, 1879–1893 (2015).\nLi, Q., Ma, M., Wu, X. & Yang, H. Snow Cover and Vegetation-Induced Decrease in Global Albedo\nFrom 2002 to 2016. J. Geophys. Res. Atmos. 124–138 (2017) doi:10.1002/2017JD027010.\nManninen, T. et al. Monitoring changes in forestry and seasonal snow using surface albedo during\n1982-2016 as an indicator. Biogeosciences 16, 223–240 (2019).\nYu, L., Leng, G. & Python, A. Attribution of the spatial heterogeneity of Arctic surface albedo\nfeedback to the dynamics of vegetation, snow and soil properties and their interactions. Environ.\nRes. Lett. 17, (2022).\nLine 201: this is interesting and could be highlighted more... But it is important to weigh the\nwarming effect against this negligible correlation...???\nResponse: Thank you for raising this important point. We suppose that our correlation analysis may\nbe insufficient to prove the conclusion of the limited role of ΔLE in the observed ΔLST trends.\nw w\nTo address this, we conducted additional experiments that clearly show albedo-driven processes,\nrather than LE-related mechanisms, dominate the winter ΔLST variations. Please see the detailed\nw\nresponse to your comments on line 115.\nLine 212: And this lack of moisture should also correlate strongly with the increasing\nimportance and impact of albedo effects, since there is no, or there is a reduced, latent heat\neffect to minimize them....???\nResponse: We appreciate your insightful view. As you correctly pointed out, decreased water\navailability can constrain vegetation’s evaporative cooling capacity and amplify the relative\nimportance of albedo effects. However, we would like to clarify that our study examines the trend\nand impacts of surface SM (0-7 cm, based on the definition of ERA5 data), which has been widely\nused in previous observation- and model-based climate change studies (Hsu and Dirmeyer, 2023;\nLiu et al., 2023; Yao et al., 2025). We have revised the abbreviation to SSM (summer surface soil\ns\nmoisture) in the manuscript for clarity.\nIn forest ecosystems, SSM represents only a fraction of the total available water, as trees can\ns\naccess deeper water reserves through their root systems. In contrast, herbaceous plants in openlands\nprimarily rely on SM in the upper layers. Therefore, we suppose that the reductions in SSM would\ns\nmore severely limit evaporative cooling in openland vegetation than in forests.\nOur subsequent findings support this hypothesis, demonstrating that SSM negatively affects\ns\nthe latent heat and leaf area index differences (ΔLE and ΔLAI) between forests and openlands (Fig.\ns s\n4), which means decreasing SSM can have a stronger limitation on short vegetation, resulting in\ns\nthe amplification of the difference between forests and openlands. It should be noted that we did not\ninclude deep layer SM in the analysis, mainly due to its high uncertainty.\nThis differential response has been documented in a previous study (Teuling et al., 2010),\nwhich investigated heatwave-induced temperature variations across vegetation types, indicating the\nmore significant warming in grassland than forests due to the rapid SM loss. Our work extends these\nfindings by validating this mechanism over the annual scale.\nFurthermore, we acknowledge the potential nonlinearity of the SM’s impacts, as you\nmentioned in the comment. Specifically, under extreme water scarcity (perhaps beyond a critical\nthreshold), even forests may lose their evaporative cooling advantage, and the albedo warming\neffect may become dominant. However, currently available observational data have not yet captured\nthese nonlinear effects. We have added a discussion about this mechanism (Line 349):\nHowever, it should be noted that our observational results represent the impact of soil moisture\nwithin a limited range under current climate regimes. If persistent warming or reduced precipitation\ndrives soil moisture below a critical threshold, forests may lose their stability and the evaporative\ncooling capacity, thus failing to compensate for the warming effect induced by the dark canopy.\nReference:\nHsu, H. & Dirmeyer, P. A. Soil moisture-evaporation coupling shifts into new gears under increasing\nCO2. Nat. Commun. 14, (2023).\nLiu, Y., Yang, Y. & Song, J. Variations in Global Soil Moisture During the Past Decades: Climate or\nHuman Causes? Water Resour. Res. 59, (2023).\nTeuling, A. J. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci. 3, 722–727 (2010).\nYao, L. et al. Emergent constraints on global soil moisture projections under climate change.\nCommun. Earth Environ. 6, 1–8 (2025).\nLine 237: we would of course expect this decoupling as SM tends to zero...?\nResponse: Thank you for the comment. We agree with your view regarding the decoupling between\nsummer vapor pressure deficit (VPD) and surface soil moisture (SSM) as soil moisture approaches\ns s\nzero. However, our decoupling methodology is not specifically designed to address this particular\ncondition.\nWe should first highlight the need for a decoupling process because VPD and SSM exhibit a\ns s\nstrong correlation, and both are significantly correlated with ΔLST. This high collinearity makes it\ns\nchallenging to determine which factor primarily drives ΔLST variability based on temporal\ns\ncorrelation analysis (as suggested by Fig. 3a and 3b, not shown here).\nTo address this, we drew upon a previous study (Liu et al., 2020), and adopted a binning\nregression method. Here, our analysis is based on spatial samples of VPD trends and SSM trends\ns s\nwith sufficient data to decouple the two factors. For clarity, we have added a schematic illustration\nof this method in the supplementary materials (Supplementary Fig. 13a).\nIn this framework, we first stratified VPD trends into 10 quantile-based groups to constrain\ns\ntheir variation. Within each subset, we then examined the relationship between SSM trends and\ns\nΔLST trends using identical samples (Supplementary Fig. 13c), thereby isolating the independent\ns\neffect of SSM by controlling for VPD variability. We subsequently reversed this process, grouping\ns s\ndata by SSM trends quantiles to evaluate the VPD-ΔLST relationship (Supplementary Fig. 13b\ns s s\nand 3d).\nThe results reveal a clear pattern: SSMs maintain robust correlations with ΔLST when\ns\naccounting for VPD variation, and the VPD-ΔLST relationships became statistically insignificant\ns s s\nwhen controlling for SSM (Fig. 3d, not shown here). This analysis provides evidence that SSM\ns s\nvariation is the primary driver of ΔLST trends, resolving the initial ambiguity caused by the\ns\ncollinearity between VPD and SSM.\ns s\nSupplementary Fig. 13. Schematic illustration of the binning regression method for isolating the\neffects of summer vapor pressure deficit (VPD) and surface soil moisture (SSM) on the land surface\ns s\ntemperature effect of potential forestation (ΔLST). (a) Stratification of VPD trends into 10 quantile-\ns s\nbased groups. (b) Stratification of SSM trends into 10 quantile-based groups. (c) SSM-ΔLST\ns s s\nrelationship within the given constrained VPD trend group. (d) VPD-ΔLST relationship within the\ns s s\ngiven constrained SSM trend group.\ns\nReference:\nLiu, L. et al. Soil moisture dominates dryness stress on ecosystem production globally. Nat.\nCommun. 11, 1–9 (2020).\nLine 310: What is missing for me, here, is the pathway effect of LAI on SM...? Of the two\nvariables, forest effects in my view are primary, while SM effects are dependent upon the\npresence/absence of forest....??? How do you address this...?\nResponse: Thank you for raising this concern about the causal relationship. We agree with your\nview that there may be a bidirectional relationship between SM and LAI. We should first emphasize\nthe definition of the variables used to build the model. The summer surface SM (SSM) is the grid-\ns\nlevel average that includes all land cover types (forest, openland, etc.). Thus, SSM primarily reflects\ns\nbackground environmental conditions driven by climate change, rather than vegetation-specific\ncharacteristics. On the other hand, ΔLAI specifically captures the difference in LAI between forest\ns\nand openland samples within the grid.\nOur modeling approach is based on the logic that climate-driven increase or decline in SSM\ns\nmay affect forest and openland vegetation differently. For instance, in dry years, forests with their\ndeeper root systems may maintain higher LAI compared to more sensitive openland vegetation. This\nwould increase ΔLAI.\ns\nWe recognize your point about potential feedback mechanisms. However, we believe that in\nour framework, the influence of ΔLAI on grid-scale SSM is likely minimal. This is because ΔLAI\ns s s\nrepresents relative differences in vegetation within stable pixels (where land cover type has not\nchanged). This internal difference may not affect the overall SM of the grid. We have added a\ndescription to clarify this causality issue (Line 764):\nIt should be noted that all drivers used in SEM (SSM, AT, DSR, WS and P) represent grid-\ns s s s s\nlevel mean values, which primarily reflect background climatic and environmental conditions.\nThese variables are unlikely to be influenced by internal variations within the grids (ΔLE and\ns\nΔLST).\ns\nLine 319: I did not see this tested in any way in the modeled pathways...? See related comment\nabove...???\nResponse: Thank you for identifying this logical issue. Our structural equation modeling (SEM)\nresults demonstrate that SSM negatively affect both the LE and LAI differences (ΔLE and ΔLAI)\ns s s\nbetween forests and openlands. This negative relationship indicates that declining SSM can amplify\ns\nthe ET and greenness difference between these land cover types. Building upon our findings and\nexisting literature about soil moisture deficit impacts on warming, we conclude that the observed\nnegative effect likely stems from forests’ greater stability and resilience to soil dryness. We have\nrevised this section to improve clarity and better illustrate this mechanism (Line 332):\nDuring summer, the vegetation greenness and ET rates are generally higher in forests than in\nadjacent openlands (Tang et al., 2018). The negative effects of SSM on ΔLAI and ΔLE revealed by\ns s s\nthe SEM analysis demonstrate that the decreased SSM can amplify these greenness and ET\ns\ndifferences between forests and openlands. Given the well-established constraint of soil moisture on\nvegetation growth (Kannenberg et al., 2024), we tend to attribute the negative effects of SSM to\ns\nforests’ stability and resilience to increasing soil dryness (i.e., lower sensitivity of forest greenness\nand ET to SSMs), enabling forests to better maintain their cooling capacity under climate change\ncompared to openland vegetation.\nPrevious studies indicate that soil moisture depletion in Europe can exacerbate warming\nthrough land-atmosphere feedback mechanisms during heatwave events (Dirmeyer et al., 2021;\nFischer et al., 2007). A comparative analysis of temperature responses across land cover types\nshows that openlands experience pronounced warming during long-lasting heatwaves, whereas\nforests exhibit moderate warming due to their conservative water-use strategies and sustained\nevaporative cooling (Teuling et al., 2010). Our results provide observational evidence that forests’\ncapacity to mitigate water deficit-induced warming persists not only during heatwave events but\nalso across the interannual scales.\nReference:\nDirmeyer, P. A., Balsamo, G., Blyth, E. M., Morrison, R. & Cooper, H. M. Land-Atmosphere\nInteractions Exacerbated the Drought and Heatwave Over Northern Europe During Summer\n2018. AGU Adv. 2, 1–16 (2021).\nFischer, E. M., Seneviratne, S. I., Lüthi, D. & Schär, C. Contribution of land-atmosphere coupling\nto recent European summer heat waves. Geophys. Res. Lett. 34, 1–6 (2007).\nKannenberg, S. A., Anderegg, W. R. L., Barnes, M. L., Dannenberg, M. P. & Knapp, A. K. Dominant\nrole of soil moisture in mediating carbon and water fluxes in dryland ecosystems. Nat. Geosci.\n(2024) doi:10.1038/s41561-023-01351-8.\nTang, B., Zhao, X. & Zhao, W. Local effects of forests on temperatures across Europe. Remote Sens.\n10, 1–24 (2018).\nTeuling, A. J. et al. Contrasting response of European forest and grassland energy exchange to\nheatwaves. Nature Geosci. 3, 722–727 (2010).\nLine 321: Of course, at some point, one might expect this \"stability\" to slip as warming\nsurpasses the ability of a forest to store adequate soil moisture and as ET begins to exceed\nthe supply of water (due to long-term rising temperatures and declining rainfall)...???\nResponse: Thanks for raising this important issue. We have added a discussion about this\nmechanism (Line 349):\nHowever, it should be noted that our observational results represent the impact of soil moisture\nwithin a limited range under current climate regimes. If persistent warming or reduced precipitation\ndrives soil moisture below a critical threshold, forests may lose their stability and the evaporative\ncooling capacity, thus failing to compensate for the warming effect induced by the dark canopy.\nLine 334: ???\nResponse: Thanks for the careful reading. We have corrected the typo.\nLine 401: ??? Why not say that the radiative forcing of European forests is virtually\nunchanged...??? The changes you highlight occur elsewhere, not in the forested areas...???\nI do not think it makes sense to confuse forests with this phenomenon...??? And it certainly\ndoesn't make sense to blame forests (or suggest they are in any way responsible) for this\noutcome...???\nResponse: Thank you for the constructive criticism. We apologize for the confusion caused by the\nwording. We have revised the relevant sentence to improve clarity (Line 435):\nThus, the diminished SC could result in a more pronounced reduction in openland albedo but\nw\nhave a limited impact on forest albedo, thereby weakening the relative radiative warming effect of\nforests compared to openlands.\nLine 416: Same comment as above... And this statement should include a stronger sense of\nthe directionality implied in the effects...?\nResponse: Thank you for your feedback. After consideration, we have removed this statement from\nthe revised version for two reasons: (1) large uncertainties remain concerning both the local climatic\neffects of mid-latitude forestation and the precise transitional latitude between warming and cooling\neffects; and (2) the original wording’s potential for misinterpretation. We note that this particular\npoint was not central to our overall argument, and its removal does not affect the conclusions.\nResponse to Reviewer #2’s Comments\nWe greatly appreciate the anonymous reviewer #2 for the constructive comments. The manuscript\nhas been revised and improved accordingly. Below are the point-by-point responses to the\ncomments, along with the revision of the manuscript (typed in red and Italic). The changes are also\nmarked in red in the manuscript file.\nReviewer #2 (Remarks to the Author):\nI am impressed by the lengths to which the Authors have gone to address the many\ncomments raised. Though ideally, I would like to continue this discussion, I do\nnot wish to burden the Authors unnecessarily. I do think that this manuscript\nalready provides a lot of interesting and highly relevant analysis. I did still\nmake a number of comments in the attached submission pdf and would encourage the\nAuthors to consider these. Most of the comments are relatively minor in character\nand can easily be addressed (or ignored where appropriate). Thus, I leave it up\nto the Authors whether they adopt these suggestions or not.\nI do think the manuscript could probably benefit from one last round of editing\nand perhaps some minor shortening. I do hope the provided comments can prove\nuseful here.\nResponse: We sincerely appreciate your positive feedback. We have considered your comments\nraised in the PDF file and made corresponding revisions. Please see the following responses for\ndetails.\nLine 9: And/or the consequences of more extensive forest loss... I think the\nphrasing matters here, in particular because the conventional view is that more\nforest cover means more warming, due to the albedo effects... But you are saying\nthe opposite here... So I would try to emphasize this more prominently... You\nare not really criticizing forestation, but rather its opposite...\nResponse: Thanks for the advice. We have revised the wording to avoid ambiguity (Line 43):\nBy comparing the satellite observations of spatially nearby forests and openlands over the last\ntwo decades, we quantify temporal trends in local land surface temperature (LST) effects of forest\nchange in Europe.\nLine 12: Some minor editing here and there would still be helpful. Here, I would\nmove this phrase to just after \"the daytime warming effect of forestation induced\nby...\nResponse: Thanks for the comment. We have revised the sentence as you suggested (Line 45):\nDuring winter, the daytime warming effect of potential forestation weakens and reverses to\ncooling (-0.142 K/decade) with decreasing snow cover, as forests show less pronounced surface\ndarkening trends than openlands.\nLine 144: Rephrase slightly...? In contrast, forests exhibit strong land surface\ncooling in summer (), at the strong daytime cooling effect () is minimally offset\nby a much smaller nighttime warming effect ()... \nResponse: Thanks for the advice. We have revised the sentence as you suggested (Line 148):\nIn contrast, forests exhibit strong land surface cooling in summer (-0.75 ± 0.41 K), as the\ndaytime cooling (-2.11 ± 0.98 K) is minimally offset by a much smaller nighttime warming (0.25\n± 0.23 K).\nLine 120: Our/The results...\nResponse: Revised as you suggested.\nLine 125: with only...?\nResponse: Revised as you suggested.\nLine 156: new paragraph...?\nResponse: Revised as you suggested.\nLine 171: It may also be worth noting that snow cover under forests disappears\nmore slowly/is more easily retained for longer periods of time than snow cover\nin open areas...\nResponse: Thanks for the suggestion. We have added related descriptions (Line 205):\nThe possible mechanism is that snow tends to be masked by tree canopies and retained for\nlonger durations in forest ecosystems, resulting in a lower sensitivity of forest albedo to snowpack.\nLine 175: suggested, evidenced by, supported by...\nResponse: Revised as you suggested.\nLine 177: Or mention the point raised above here... Or further condense...\nResponse: We have added related descriptions in Line 205. Please see the above response for the\nrevision.\nLine 255: I tend to think this analysis should also still consider the differences\nacross openlands and forested lands...? I would expect to see significant\ndifferences here...?\nResponse: Thank you for raising this concern. We admit that soil moisture may have difference in\nforested area and openland. We have highlight this point in Line 336:\nSpecifically, forests possess deeper root systems and enhanced soil moisture storage capacity,\nenabling them to sustain their cooling efficiency under climate change compared to openland\nvegetation.\nLine 307: good... I see this point is covered here... Perhaps above one could\nmention that it will be discussed further below... Above, it just feels like\nthis point is missing...\nResponse: Thanks for the advice. We have mentioned this point above as you suggested (Line 245):\nNotably, the climate variables may themselves be influenced by the presence of forests (e.g.,\nsoil moisture retention and evaporative cooling of trees). However, it should be emphasized that the\nclimate drivers used in our analysis are the grid means at the coarse resolution, reflecting large-\nscale patterns rather than microclimates. Our analysis tends to capture the distinct responses of\nforests and openlands to background climate change.\nLine 336: There are presumably two things going on here: 1) forests have deeper\nroot systems and can dig deeper for water... and 2) forested areas presumably\nstore more soil moisture and thus are not as quickly or easily affected by soil\nmoisture decline...\nDoes one see higher soil moisture measurements, in general, across the forested\nareas...?\nYou also mention below that forests may be better at conserving water (e.g.,\nshutting down stomata)... So this would be a third mechanism...\nResponse: We fully agree with the point you mentioned. We have highlighted the two mechanisms\nyou mentioned here (Line 336):\nSpecifically, forests possess deeper root systems and enhanced soil moisture storage capacity,\nenabling them to sustain their cooling efficiency under climate change compared to openland\nvegetation.\nLine 504: This paragraph needs to be broken up into smaller pieces... Way too\nlong...\nResponse: The paragraph has been sperated into 3 parts as you suggested.\nLine 520: I appreciate this point... It might be wise to find ways to introduce\nthis restriction earlier on in the discussion, and perhaps even the abstract...\nJust to say that the analysis is focused strictly on local level effects...\nResponse: Thanks for the advice. We have highlighted the restriction at the beginning of the\nmanuscript (Line 114):\nIt should be noted that the evaluated LST effects are local effects and should be treated as the\nreference for local climate adaptation rather than global climate mitigation.\nLine 548: Is it forestation that exhibits negative trends, or effects on open\nlands...???\nResponse: Thanks for pointing this out. We have revised the sentence for clarity (Line 534):\nOverall, our findings provide solid satellite evidence that the cooling effects of forest\nconversion from openlands in Europe have been amplified in the last two decades."}]

Decision letter