SCIENCE ADVANCES | RESEARCH ARTICLE ENVIRONMENTAL STUDIES Copyright © 2020 The Authors, some Summer soil drying exacerbated by earlier spring rights reserved; exclusive licensee greening of northern vegetation American Association for the Advancement of Science. No claim to Xu Lian1, Shilong Piao1,2,3*, Laurent Z. X. Li4, Yue Li1, Chris Huntingford5, Philippe Ciais6, original U.S. Government Alessandro Cescatti7, Ivan A. Janssens8, Josep Peñuelas9,10, Wolfgang Buermann11,12, Works. Distributed Anping Chen13, Xiangyi Li1, Ranga B. Myneni14, Xuhui Wang1,6, Yilong Wang6, Yuting Yang15, under a Creative Zhenzhong Zeng1, Yongqiang Zhang16, Tim R. McVicar17,18 Commons Attribution NonCommercial Earlier vegetation greening under climate change raises evapotranspiration and thus lowers spring soil moisture, License 4.0 (CC BY-NC). yet the extent and magnitude of this water deficit persistence into the following summer remain elusive. We pro- vide observational evidence that increased foliage cover over the Northern Hemisphere, during 1982–2011, triggers an additional soil moisture deficit that is further carried over into summer. Climate model simulations independent- ly support this and attribute the driving process to be larger increases in evapotranspiration than in precipitation. This extra soil drying is projected to amplify the frequency and intensity of summer heatwaves. Most feedbacks operate locally, except for a notable teleconnection where extra moisture transpired over Europe is transported to central Siberia. Model results illustrate that this teleconnection offsets Siberian soil moisture losses from local Downloaded from http://advances.sciencemag.org/ on April 2, 2020 spring greening. Our results highlight that climate change adaptation planning must account for the extra summer water and heatwave stress inherited from warming-induced earlier greening. INTRODUCTION plex than just by modifying ET. For instance, more rapid spring soil A widespread shift in the phenology of northern vegetation under moisture losses raise moisture levels of the atmosphere, facilitate cloud global warming, toward earlier spring leaf emergence (i.e., earlier spring formation, and trigger extra precipitation (14–16). Any additional “greening”), has been unambiguously documented using indices of rainfall may compensate, at least partly, the greening-induced deficit satellite-monitored vegetation greenness (1–3) and ground observations of soil moisture. A further biophysical feedback is that soil moisture (4–6). Both modeling and observational studies have shown that the deficits result in less conversion of surface incident sunlight into water early onset of vegetation phenology causes an increase in ecosystem evaporation. Lower latent heat from the land may suppress the de- productivity through an extension of the active carbon uptake peri- velopment of clouds and subsequent precipitation (17, 18). Beyond od (7, 8). On the other hand, however, there is a growing concern local effects, these feedback mechanisms may also operate as telecon- that the earlier vegetation greening may induce soil moisture deficits nections between locations. Large-scale atmospheric circulation (14, 19) in late spring through increased evapotranspiration (ET), which may redistributes the extra recycled moisture across the Earth’s land sur- amplify summer droughts and have detrimental effects on summer face and contributes to additional precipitation in downwind regions. ecosystems (9–13). Direct observational evidence for this ET-centered The overall impact of these potential interseasonal feedbacks on summer legacy hypothesis is nonetheless very scarce. Furthermore, the influence water cycles, unfortunately, is still poorly characterized. of earlier spring greening on summer soil moisture can be more com- Quantifying this vegetation–soil moisture legacy effect must take full account of the complete loop of soil-vegetation-atmosphere ex- 1 changes of water and energy fluxes. Many feedback mechanisms, as Sino-French Institute for Earth System Science, College of Urban and Environmental Sciences, Peking University, Beijing 100871, China. 2Key Laboratory of Alpine Ecology, summarized in the above paragraph, are essential for accurate assess- Institute of Tibetan Plateau Research, Chinese Academy of Sciences, Beijing 100101, ments of global or regional soil moisture changes, yet are ignored by China. 3Center for Excellence in Tibetan Earth Science, Chinese Academy of Sciences, Beijing previous studies (9–12). In particular, without characterizing atmo- 100085, China. 4Laboratoire de Météorologie Dynamique, CNRS, Sorbonne Université, Ecole Normale Supérieure, Ecole Polytechnique, Paris, France. 5Centre for Ecology spheric feedbacks and including those that are interseasonal, it still and Hydrology, Wallingford, Oxfordshire OX10 8BB, UK. 6Laboratoire des Sciences remains an open question on the extent to which vegetation greening du Climat et de l’Environnement (LSCE), CEA CNRS UVSQ, 91191 Gif Sur Yvette, France. influences summer soil moisture. Here, we first analyze multiple strands 7 European Commission, Joint Research Centre, Directorate for Sustainable Resources, of direct observational data to characterize the spring-to-summer I-21027 Ispra (Varese), Italy. 8Department of Biology, University of Antwerp, Universite- itsplein 1, Wilrijk 2610, Belgium. 9CREAF, Cerdanyola del Valles, Barcelona, Catalonia vegetation–soil moisture feedbacks. We then perform factorial exper- 08193, Spain. 10CSIC, Global Ecology Unit CREAF-CSIC-UAB, Bellaterra, Barcelona, iments with a state-of-the-art coupled land-atmosphere climate model Catalonia 08193, Spain. 11Institute of the Environment and Sustainability, University to better understand the underlying mechanisms. of California, Los Angeles, Los Angeles, CA 90095, USA. 12Department of Geography, University of Augsburg, 86159 Augsburg, Germany. 13Department of Forestry and Natural Resources, Purdue University, West Lafayette, IN 46907, USA. 14Department of Earth and Environment, Boston University, Boston, MA 02215, USA. 15State Key RESULTS Laboratory of Hydro-Science and Engineering, Department of Hydraulic Engineer- ing, Tsinghua University, Beijing, China. 16Key Lab of Water Cycle and Related Observational evidence for summer drying induced by Land Surface Processes, Institute of Geographic Sciences and Natural Resources earlier spring greening Research, Chinese Academy of Sciences, Beijing 100101, China. 17CSIRO Land To examine spatiotemporal connections between recent spring pheno- and Water, Canberra, Australian Capital Territory, Australia. 18Australian Research logical changes (section S1 and fig. S1) and soil water content (SWC) Council Centre of Excellence for Climate Extremes, Canberra, Australian Capital Territory, Australia. dynamics in the following summer, we first analyze long-term satellite *Corresponding author. Email: slpiao@pku.edu.cn data of spring leaf area index (LAI; as a proxy for vegetation phenology Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 1 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE and growth; Materials and Methods) (20) and summer values of two spring LAI mode shows widespread positive anomalies (i.e., earlier observation-driven SWC-related datasets for the period of 1982–2011. greening) in the majority of the Eurasia (>80%) and eastern United The two SWC-related data include root-zone soil moisture from the States, while negative anomalies (i.e., delayed greening, also known as Global Land Evaporation Amsterdam Model (GLEAM; Materials and browning) cover vast areas of central and western North America Methods) (21) and the long-term reconstruction of terrestrial water (Fig. 1,C and E). The patterns of the leading LAI mode generally storage (TWS) anomalies based on Gravity Recovery and Climate resemble the trend pattern of satellite-observed spring LAI (fig. S1). As Experiment satellites (GRACE-REC) (22). Our study domain only a consequence, the leading summer SWC (or TWS) mode (Fig. 1, D and F) covers the middle to high latitudes of the Northern Hemisphere (25° to coupled to the spring LAI patterns demonstrates a direct response 90°N) because of their distinctive vegetation seasonality and particular- pattern of summer SWC to the emerging spring greening. The MCA- ly pronounced recent changes in springtime LAI phenology (fig. S1). based annual time series linked to the spring LAI and summer SWC We define spring across the domain as the period from March to May (or TWS) pattern are also highly correlated (r = 0.86 for GLEAM (MAM) and summer as the period from June to August (JJA). We SWC and r = 0.80 for GRACE-REC TWS, P < 0.05 for both; fig. S4). calculate partial correlations between spring LAI and summer SWC This reiterates the strong coupling between the two coupled patterns (or TWS), whereby covarying effects of summer LAI and climatic shown in Fig. 1 (C to F). In general, negative summer SWC anomalies variables (i.e., temperature, precipitation, and solar radiation) are (i.e., soil drying) are preceded by spring greening, whereas positive controlled for (Materials and Methods). As potential nonlocal pre- SWC anomalies (i.e., soil wetting) dominate the western North America cipitation feedbacks modulate the local relationship between spring with spring browning (Fig. 1, D and F versus C and E). The domi- LAI and summer SWC (14, 19), we also implement a lagged maximum nance of negative LAI-SWC coupling, as also seen by partial correla- covariance analysis (MCA; Materials and Methods, section S2) (23). tion patterns, implies an impact of earlier spring onset on local Downloaded from http://advances.sciencemag.org/ on April 2, 2020 This statistical approach isolates pairs of coupled spatial patterns that summer soil drying. The cropland-dominated central Europe has wetter explain a maximum fraction of the covariance between two space- summers during years with greener springs (Fig. 1, C to F), which is time data fields (in our case, spring LAI and summer SWC). Hence, in agreement with the positive correlation between spring greenness the MCA method has obvious advantages over routine correlation and summer soil moisture (Fig. 1, A and B). A notable exception analysis that only partially describes covariability between two fields occurs in the Siberia region, where summer SWC increases despite the (i.e., spatial correlation using geographical points only or temporal presence of a clear positive spring LAI anomaly and limited human correlation using time series only). interventions. This is also in sharp contrast to the negative partial Partial correlations between spring greenness and subsequent sum- correlations in central Siberia after removing covarying effects of pre- mer SWC show a widespread negative pattern across most northern cipitation (Fig. 1, A and B). Such difference implies that altered pre- lands (64% of land area for GLEAM SWC and 56% for GRACE-REC cipitation patterns driven by spring LAI changes cause extra moisture TWS) (Fig. 1, A and B), including significantly greening regions of to be brought into this area. Europe, eastern United States, East Asia, and Siberia (fig. S1). For temperate and semi-arid grasslands where LAI is coupled to shallow Simulating the spring-to-summer vegetation feedbacks moisture, this dominance of negative correlation between spring LAI Statistical analyses of observation-based data show that summertime and summer moisture is additionally supported by using micro- soil moisture anomalies negatively co-vary with the preceding spring wave-derived summer top-layer soil moisture (fig. S2). The wide- greening patterns. However, the detection of systematic legacy ef- spread negative correlation indicates that, when effects of interannual fects of spring LAI on summer soil moisture within data is partially variations of precipitation are factored out, earlier greening in spring affected by noise from weather variability, potential nonlinearity in generally results in local soil drying in subsequent summers. On the soil moisture responses, different rooting depth and moisture access other side, positive correlation is also apparent over regions including across ecosystems, and slowly varying drivers of both LAI and soil the North China Plain, central and southeastern Europe, and the U.S. Great moisture trends such as vapor pressure deficit and rising atmospheric Plains, especially for TWS that additionally captures storage varia- CO2. To confirm a causal link between spring greening and summer tions in canopy, snow, streams, lakes, and groundwater (Fig. 1, A soil drying, and to better understand the local and nonlocal mecha- and B). Many such regions are influenced by intensive agricultural nisms, we perform further analyses using the IPSL-CM4 coupled activities, which potentially decouples the relationship between spring land-atmosphere global climate model (GCM; Materials and Methods) LAI and summer SWC observed in natural ecosystems. By averag- (24). Coupled land-atmosphere models explicitly describe biophys- ing land surface undisturbed by cropland management, we find a ical interactions between the land surface and the atmosphere. This negative correlation between spring LAI and summer SWC that is includes local exchanges of water and energy fluxes and large-scale statistically significant (−0.44 for GLEAM SWC and −0.40 for GRACE- changes of atmospheric circulation and moisture convergence. Paired REC TWS, P < 0.05) (fig. S3). However, an increasing fraction of land time-evolving (i.e., transient) numerical experiments are run with the areas affected by croplands leads to a weakening of this negative cor- land surface component of IPSL-CM4, ORCHIDEE (Organizing Car- relation, which shifts robustly to being positive when croplands oc- bon and Hydrology In Dynamic Ecosystems; Materials and Methods), cupy over 50% of local land areas (fig. S3). forced with two contrasting spring LAI conditions. One set of sim- Figure 1 (C to F) presents the leading first MCA mode characterizing ulations (labeled as “​​LAI​MAM obs​ ​​”) is prescribed with annually varying spatial patterns of the coupling between spring LAI and summer LAI from satellite observations (20) during all seasons, including SWC. This mode, in particular, includes potential nonlocal impacts MAM. The second set of simulations (labeled as “​​LAI​MAM clim​ ​​”) is also of LAI-induced moisture transport anomalies on SWC at downwind prescribed with annually varying LAI during all seasons, except in regions. This mode pair explains 46 and 59% of the total squared MAM when it is fixed to the climatological mean LAI conditions covariance between the two variables for GLEAM SWC and GRACE- observed during 1982–2011 (Materials and Methods). Apart from REC TWS, respectively. For both SWC-related products, the first the spring LAI differences, the two model experiments are forced Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 2 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE Partial correlation coefficients (with GIMMS spring LAI) A GLEAM summer SWC B GRACE-REC summer TWS 0.6 0.6 0.4 0.4 0.2 0.2 0 0 −0.2 −0.2 −0.4 −0.4 −0.6 −0.6 1st MCA mode, GIMMS spring LAI versus GLEAM summer SWC, SFC = 46.3% C GIMMS spring LAI D GLEAM summer SWC Downloaded from http://advances.sciencemag.org/ on April 2, 2020 (m2 m −2) (m3 m −3) 0.04 0.04 0.03 0.03 0.02 0.02 0.01 0.01 0 0 −0.01 −0.01 −0.02 −0.02 −0.03 −0.03 −0.04 −0.04 1st MCA mode, GIMMS spring LAI versus GRACE-REC summer TWS, SFC = 58.9% E GIMMS spring LAI F GRACE-REC summer TWS (m2 m −2) (m3 m −3) 0.04 0.04 0.03 0.03 0.02 0.02 0.01 0.01 0 0 −0.01 −0.01 −0.02 −0.02 −0.03 −0.03 −0.04 −0.04 Fig. 1. Coupling between observed spring LAI and summer SWC during 1982–2011. Spatial pattern of partial correlation coefficients between Global Inventory Monitoring and Modeling Studies (GIMMS) spring LAI and (A) GLEAM summer SWC or (B) GRACE-REC summer TWS for 1982–2011. Black stipples indicate regions with a statistically significant correlation (P < 0.05). Heterogeneous regression maps of (C and E) GIMMS spring LAI and (D and F) GLEAM summer SWC (or GRACE-REC summer TWS), associated with the first MCA mode for 1982–2011. The squared fractional covariance (SFC) explained by the first MCA mode is 46.3 and 58.9% for the GLEAM and GRACE-REC datasets, respectively. identically with observed boundary conditions of sea-surface Model evaluation against global observational datasets supports temperature (SST), sea-ice fraction (SIF), and atmospheric CO2 the capacity of the IPSL-CM4 GCM to capture the annual and sea- concentration. This enables us to account for the vegetation bio- sonal features of key hydrological variables, including ET and pre- physical feedbacks stemming from interactions with changes in oceanic cipitation (fig. S5, A and B; further details in section S3). To avoid processes or CO2 radiative/physiological forcings. Multimember systematic biases from the diverging magnitude of modeled and paired ensembles are generated for each coupled model experiment observation-based SWC, SWC changes are instead compared as the by performing 30 repeated runs but with different initial conditions fraction of the present-day climatological mean (1982–2011) (25). At (Materials and Methods). the scale of the Northern Hemisphere (25° to 90°N), the IPSL-simulated Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 3 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE SWC shows a negative trend, whereas GLEAM SWC shows a weaker firm a strong coupling between spring LAI and summer SWC positive trend (fig. S5C). While this sign difference is possibly a anomalies (r = −0.51, P < 0.01) (Fig. 2A). Hence, on a year-to-year concern for global mean changes, we find zonally coherent patterns basis, the observed higher spring LAI values prescribed to our GCM cor- of the IPSL-simulated SWC trend compared to GLEAM SWC (fig. S5, respond to lower summer SWC estimates. We then examine whether D to I). The IPSL-simulated SWC broadly captures regions where there is a robust overall trend of decreasing SWC, driven by the soil wetting or drying occurs, at both annual and seasonal time scales significantly increasing spring LAI (greening). This effect is more (fig. S5, D to I). This regional consistency strengthens our confidence difficult to detect because the interannual trend of LAI over the past in examining both the observed zonally different fingerprint of spring 30 years is smaller than concurrent climate-driven variabilities. Aver- greening on summer soil moisture and the ability of our GCM to ex- age summer SWC across all the northern lands (25° to 90°N) tends plain such changes. to decrease at a rate of −0.11 ± 0.09 × 10−4 m3 m−3 year−1 (linear trend: The difference between the two factorial simulations with vari- mean ± 1 SD, P = 0.22) due to the spring LAI changes. This overall able and fixed spring LAI (​​LAI​MAM MAM obs​ ​​ − ​​LAIclim ​ ​ ​​) allows us to isolate summer SWC trend is negative, in line with our hypothesis that higher the fingerprint of observed spring LAI changes on summer soil green foliage cover is transpiring more in spring and leaving soils moisture (SWC) and other important hydrological variables. The drier in summer. The calculated trend is not statistically significant, symbol “” hereafter denotes the influence of spring LAI changes on partly due to the presence of contrasting effects on summer SWC in other target variables of interest. Our factorial model simulations con- regions with spring greening (i.e., positive LAI trends) versus spring A D Downloaded from http://advances.sciencemag.org/ on April 2, 2020 1 Spring LAI anomaly (m 2 m −2) −0.15 Summer SWC anomaly (%) 0.8 Spring LAI 0.6 Spring Browning Summer SWC Wetting 0.6 Fitted line for SWC −0.1 0.4 *** *** Summer ** 0.4 0.2 −0.05 0.2 Trend 0 0 0 n.s. −0.2 −0.2 ** n.s. n.s. Greening n.s. 0.05 Drying −0.4 ** −0.4 −0.6 0.1 −0.8 Corr: r = −0.55, P < 0.01 −0.6 Trend: −0.11 0.09 10 −4 m 3 m−3 yr −1, P = 0.22 −1 0.15 1985 1990 1995 2000 2005 2010 Spring LAI ET P Q SWC (0.01 year−1)(mm year−2)(mm year−2)(mm year−3)(10−4 m3 m−3 year−1) B E 1 Spring LAI anomaly (m 2 m −2) −0.15 Summer SWC anomaly (%) 0.8 Spring LAI 0.6 Spring Browning Summer SWC *** *** Wetting 0.6 Fitted line for SWC −0.1 0.4 *** Summer 0.4 0.2 −0.05 0.2 Trend 0 0 0 n.s. −0.2 −0.2 *** Greening 0.05 n.s. * ** Drying −0.4 *** −0.4 −0.6 0.1 −0.8 Corr: r = −0.7, P < 0.01 −0.6 Trend: −0.17 0.07 10 −4 m 3 m −3 yr −1, P < 0.05 −1 0.15 1985 1990 1995 2000 2005 2010 Spring LAI ET P Q SWC (0.01 year−1)(mm year−2)(mm year−2)(mm year−3)(10−4 m3 m−3 year−1) C F 2 −0.15 Spring LAI anomaly (m 2 m −2) Summer SWC anomaly (%) Spring LAI 0.6 Spring 1.5 Browning Summer SWC Wetting −0.1 0.4 Summer n.s. n.s. Fitted line for SWC 1 0.2 n.s. 0.5 −0.05 n.s. Trend 0 0 0 −0.5 −0.2 ** n.s. Greening 0.05 n.s. Drying n.s. −1 −0.4 n.s. n.s. 0.1 −1.5 Corr: r = −0.34, P = 0.07 −0.6 Trend: 0.15 0.19 10 −4 m 3 m −3 yr −1, P < 0.42 −2 0.15 1985 1990 1995 2000 2005 2010 Spring LAI ET P Q SWC (0.01 year−1)(mm year−2)(mm year−2)(mm year−3)(10−4 m3 m−3 year−1) Fig. 2. IPSL-simulated changes in summer soil moisture induced by spring LAI changes. (A to C) Interannual anomalies of the area-weighted average of spring LAI (green dotted lines) and summer SWC (blue solid lines) for (A) all northern latitudes (25° to 90°N), (B) regions with positive LAI trends (i.e., greening), and (C) regions with negative LAI trends (i.e., browning). The red lines indicate the least-squares linear regression of GCM-based SWC (straight lines) against time and the 95% confidence intervals (curves). Note that the right axes are reversed, so higher LAI values are toward the bottom of the plots. (D to F) Interannual trends in mean spring LAI and resultant changes in the spring and summer hydrological variables of ET, precipitation (P), runoff (Q), and SWC. The subplot maps at the bottom of each main panel display the corresponding averaged areas as gray [corresponding to (A) to (C)]. ***P < 0.01; **P < 0.05; *P < 0.1; n.s., P > 0.1. Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 4 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE browning (i.e., negative LAI trends). When instead considering only gions with earlier greening, toward regions with browning, within regions where spring LAI generally increases, a statistically significant the springtime. and much stronger summer drying trend in summer SWC (−0.17 ± Our factorial model simulations allow further investigation of the 0.07 × 10−4 m3 m−3 year−1, P < 0.05) is detected (Fig. 2B). On aver- soil moisture legacy effects into the summer following a spring with age, including areas with spring browning, this decreasing trend altered phenology and hydrology. We find that the spring soil mois- of summer SWC associated with spring LAI changes contributes an ture deficit caused by earlier greening is carried over into the summer additional 16% to the simulated decrease in summer SWC directly months (Fig. 2D, blue bars). The interseasonal soil moisture mem- caused by climate change during 1982–2011 (overall decrease, ory propagates this moisture deficit forward, and this also involves −0.68 × 10−4 m3 m−3 year−1 from the ​​LAI​MAM additional contribution of soil moisture–atmosphere biophysical feed- obs​ ​​ simulation). We also perform the MCA analysis on observed spring LAI- and backs. Given the short memory of days to weeks of the atmosphere, the atmospheric moisture anomaly from earlier greening is dissipat- IPSL-simulated summer SWC (from simulation ​​LAI​MAM obs​ ​​). The cou- ed mainly in the springtime through extra spring rainfall (Fig. 2D). pling patterns between spring LAI and summer SWC for the first The lower soil moisture levels in late spring due to a prolonged period mode (fig. S6, A and B) show high accordance with those purely of stronger decrease and insufficient rainfall resupply become less ac- based on observations (Fig. 1, C to F), confirming the ability of our cessible for plant root uptake and soil evaporation (27), particularly in GCM to simulate the interseasonal vegetation–soil moisture feedbacks. climatologically water-stressed regimes. Hence, summer ET decreases We further test whether insights based on our GCM-based factorial substantially across the entire study area (−0.10 ± 0.05 mm year−2), and simulations could be overly influenced by water-stressed regions at summer P (precipitation) decreases even more (−0.19 ± 0.08 mm year−2), low latitudes. In doing so, we calculate spring LAI and summer SWC leading to a slightly strengthened decreasing trend in summer SWC Downloaded from http://advances.sciencemag.org/ on April 2, 2020 anomalies weighted by both area and annual mean soil moisture (−0.11 ± 0.09 × 10−4 m3 m−3 year−1) (Fig. 2D). This summer drying values. Again, a significant negative correlation is found between spring signal is larger and statistically significant in areas with extra spring LAI and summer SWC (fig. S6, C and D). This spring LAI-induced greening (Fig. 2E). The significantly decreasing spring LAI (P < 0.05) drying trend of summer SWC (−0.09 ± 0.06 × 10−4 m3 m−3 year−1) in areas with spring browning do not, however, generate significant is consistent with that without assigning a lower weight to drier regions changes in any of the hydrological fluxes in summer (i.e., ET, P, run- (Fig. 2A), partly because earlier greening mainly occurs in regions off, and SWC; all P > 0.1) (blue bars in Fig. 2F). In these browning with favorable moisture conditions. Therefore, our GCM simulations regions, soil moisture anomalies in spring have little impact on sub- provide robust process-based evidence that the observed overall sequent summer water recycling (ET ≈ 0) and hence summer rain- earlier greening of northern land causes extra decreases in summer fall feedbacks. It is noteworthy that both spring greening and browning soil moisture. tend to decrease the amount of summer precipitation (Fig. 2, E Intra- and interseasonal soil moisture anomalies are determined by and F). This implies that the current spring LAI levels may have evolved the balance between changes in ET, precipitation, and runoff (schematic, to be an optimum that maximizes the precipitation retained on con- Fig. 3). To gain more insights into how spring vegetation changes tinents in the subsequent summer season. feedback on soil moisture, we further use the paired GCM experi- We further investigated the spatial pattern of the GCM-derived ments to examine changes in these major hydrological fluxes. High- soil moisture trend in SWC (Fig. 4; other hydrological fluxes in er spring LAI levels tend to accelerate water recycling from the land fig. S8). In general, the modeled summer SWC pattern and observed surface to the atmosphere through enhanced plant transpiration. As summer SWC anomalies associated with the first MCA mode (Fig. 1, transpiration dominates total terrestrial ET (average transpiration/ D and F) are highly similar, providing mechanistic support for our ET ratio of IPSL-CM4: 0.47 for spring and 0.76 for summer) (26), a observation-based findings. Besides, the summer SWC pattern higher spring transpiration directly leads to higher modeled spring closely resembles the spring SWC pattern (Fig. 4B versus Fig. 4A), ET (+0.28 ± 0.06 mm year−2, P < 0.01) (Fig. 2D). We also note that implying the dominant role of soil moisture memory in carrying this the overall trend of IPSL-simulated spring ET is not statistically sig- greening-induced water deficit into the summer. Areas with the high- nificant (P > 0.05), which is confirmed by two observation-based ET est rates of summer drying (>0.1 × 10−3 m3 m−3 year−1) occur main- products (fig. S5A). This is mainly because the significant ET increase ly in Europe, eastern Asia, western Asia, and eastern North America in regions with spring greening is compensated by the ET decrease (Fig. 4B), where pronounced advances of vegetation greening occur in regions with spring browning when integrated to the hemispheric during spring (fig. S1). Conversely, a summer wetting trend of soil scale (fig. S7). For this greening-enhanced ET, however, only a pro- moisture dominates regions with spring browning, such as central portion (71%) is subsequently recycled as additional land precipitation and northwestern North America (Fig. 4B versus fig. S1). One nota- during MAM (0.20 ± 0.09 mm year−2, P < 0.01) (Fig. 2D). During ble exception is again the central Siberian regions, where a strong MAM, this imbalance between increasing precipitation and ET due modeled soil wetting trend occurs in both spring and summer (>0.05 × to earlier phenology and LAI development results in a decreasing 10−3 m3 m−3 year−1) (Fig. 4), despite the presence of a clear spring trend in both soil moisture (−0.11 ± 0.09 × 10−4 m3 m−3 year−1) and greening tendency there that generally lowers summer soil moisture surface runoff (−0.08 ± 0.09 mm year−2). The larger increase in ET elsewhere (Fig. 4 and fig. S1). The GCM simulations show that this than in precipitation is much more pronounced when considering trend of regional wetting in Siberia is caused by different mechanisms only regions with spring greening, where it causes a strong signal of soil in the two seasons. In spring, the extra water vapor recycled into the drying (Fig. 2, B and E). However, in regions with spring browning, atmosphere in upstream regions (Europe) with strong greening is there is a stronger decrease in ET than in precipitation, thus causing carried to Siberia by westerly winds (fig. S9A). This imported extra soil wetting (Fig. 2, C and F). These contrasting responses, where ET precipitation moistens local soils and offsets the soil moisture deficit changes do not balance precipitation changes, suggest a geograph- originated from local spring greening (schematic, Fig. 3). During the ical redistribution of the extra atmospheric water vapor from re- following summer, when the prevailing westerly winds become weaker Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 5 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE Downloaded from http://advances.sciencemag.org/ on April 2, 2020 SWC Drying SWC = P − ET − Q Start Spring End Soil water soil water flux SWC Net inward accumulation Westerly winds Asian Monsoon Europe Siberia Eastern China Browning Greening Drying Wetting ET Precipitation Moisture transport Precipitation from local recycling Precipitation from moisture inflow Moisture outflow Fig. 3. Schematic of the effect of earlier greening on summer soil moisture. Earlier spring greening influences spring soil moisture by altering land-atmosphere water exchanges (via ET, P, and Q) and by the redistribution of atmospheric water vapor by atmospheric circulation. This spring soil moisture anomaly persists later into the following summer due to the carryover effects of soil moisture. The magnitude of this cross-seasonal vegetation feedback and the role of atmospheric circulation, how- ever, vary geographically. Three typical examples of the circulation-modulated vegetation feedback (Europe, Siberia, and eastern China) are displayed at the bottom of the schematic. (fig. S9B), the Siberian region operates as a more closed system with Implications for near-surface land climate in summer intensified local water recycling (fig. S8). Consequently, the imported Lower soil moisture can trigger the occurrence, or amplify the fre- extra soil moisture in springtime helps to sustain local positive quency and intensity of heatwaves, by altering the surface energy bal- feedbacks on summer precipitation and is ultimately retained in the ance (17, 28, 29). Our finding that earlier greening preconditions summer Siberian soils throughout the boreal summer. This additional sum- soils to be drier thus implies a shift in the distribution of summer mer wetting from the earlier greening of local and remote regions, surface air temperature toward even more hot extremes, adding to in total, accounts for over 90% of the overall summer soil wetting the warming-induced increase in extreme occurrence. To verify this trend derived from ​​LAI​MAM obs​ ​​, showing its much more impor­tant role mechanism with our paired GCM simulations, we calculate two ex- for Siberia than the ongoing climate change. treme hot temperature indices (30): (i) the summer maximum of daily Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 6 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE maximum temperature (TXx) and (ii) the number of summer days summer soil moisture. This, in turn, further aggravates the higher with maximum temperatures exceeding the local 90th percentile of frequency and intensity of hot summer extremes driven by rising a given base period, here as 1982–2011 (TX90p) (see Materials and greenhouse gas concentrations (31). Methods). We find that earlier spring greening increases both TXx and TX90p (Fig. 5 and fig. S10) across those land regions where greening depletes summer soil moisture more strongly (Fig. 4B). DISCUSSION Conversely, regions with additional summer soil wetting from Our combined use of satellite observations and model simulations sug- spring LAI changes exhibit fewer hot extremes in summer. On gests that a more intense water loss from land soils during spring, average, our paired GCM analysis shows that the overall soil drying due to earlier vegetation greening, lowers soil moisture in spring and causes an additional 23 and 26% increase (increasing rates of 0.8 days summer. This moisture depletion further reduces surface evapora- decade−1 for summer TX90p and 0.07°C decade−1 for summer tive cooling and thus increases the possibility of stronger summer TXx) in TX90p and TXx, respectively, which exacerbates the sum- heatwaves. Exacerbated soil water limitation, together with more ex- mer heatwaves caused directly by background global warming. This treme weather events in summer, may restrain vegetation growth finding implies the existence of a positive feedback: Human-induced and suppress ecosystem productivity, especially for climatologically global warming causes earlier spring greening and thus forces lower arid ecosystems that are sensitive to summer water availability (e.g., A Spring B Summer Downloaded from http://advances.sciencemag.org/ on April 2, 2020 (10 −4 m3 m−3 year−1) (10 −4 m3 m−3 year−1) 2 2 1.5 1.5 1 1 0.5 0.5 0 0 −0.5 −0.5 −1 −1 −1.5 −1.5 −2 −2 Fig. 4. Spatial patterns of IPSL-simulated trends in soil moisture induced by spring LAI changes. Spatial patterns of the linear trends in spring (left) and summer (right) SWC changes are shown. This change in soil moisture induced by earlier spring greening, SWC, is obtained as the difference between the two simulations with and without spring LAI changes. (i.e., ​​LAI​MAM MAM obs​ ​​–​​LAI​clim​ ​​). Black stipples indicate regions with a statistically significant linear trend (P < 0.05). 0.2 0.6 − 0.6 r = −0.27, P < 0.01 r = −0.41, P < 0.01 TX90p trend 0.15 0.3 0.3 TXx trend TXx or TX90p trend 0 0 0.1 −0.3 −0.3 −0.6 −0.6 0.05 −6 −4 −2 0 2 4 6 −6 −4 −2 0 2 4 6 SWC trend SWC trend 0 −0.05 TXx ( o C decade-1) −0.1 TX90p (days year -1) −0.15 1.6 −1.2 −0.8 −0.4 0.4–0 0–0.4 4–0.8 8–1.2 2–1.6 6–2.0 2 to − .6 to .2 to .8 to − 0. 0. 1. 1. − −1 −1 −0 SWC trend (10 m m year ) −4 3 −3 −1 Fig. 5. IPSL-simulated changes in summer extreme hot temperature indices related to SWC changes. The average linear trend of TXx and TX90p due to earlier greening for different summer SWC trends binned into 0.4 (×10−4 m3 m−3 year−1) intervals. The TXx and TX90p values are obtained as the difference between the two simulations with and without spring LAI changes (​​LAI​MAM MAM obs​ ​​ − ​​LAI​clim​ ​​). The insets are scatterplots of trends in TXx (left) and TX90p (right) versus trends in SWC for 1982–2011. Each colored dot in the scatterplots represents a grid cell in the pattern of trends in TXx or TX90p (fig. S10). Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 7 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE western United States). Our study has important implications for more into account in future assessment of hydrological impacts of large- accurate predictions of summer weather and climate features. For scale afforestation/reforestation. many parts of the world, climate models and measurements show that surface meteorological conditions are affected by both large- scale atmospheric circulation patterns and local soil moisture status in MATERIALS AND METHODS some regions with strong vegetation-atmosphere coupling (17, 32). Experimental design The results of this study highlight the necessity of also accounting We used satellite-derived measurements of LAI (as a proxy for veg- for spring plant phenology in weather and climate forecasting of sum- etation “greenness”) and reanalysis-based soil moisture to identify mer months, given their capacity to modulate summer soil moisture any significant relations between spring vegetation phenology and levels through the interseasonal biophysical feedbacks. summer soil moisture, using the MCA approach. We also used fac- One limitation of our analysis is the examination of diagnostics from torial experiments by the IPSL-CM Earth System Model to provide an ensemble based on a single climate model. This concern is offset by process understanding of observation-based greenness–to–soil mois- our findings of modeled spring-to-summer greenness–soil moisture ture linkages. linkages having consistency with results of parallel observation–based analysis (Fig. 1). However, for future studies, we encourage other mod- Observational datasets eling centers to perform similar factorial experiments to ours to ver- Vegetation greenness was characterized using a gridded biweekly LAI ify the lagged hydrological impacts of earlier vegetation phenology. dataset derived from the Global Inventory Monitoring and Model- This includes independent model evidence of our discovery that ing Studies (GIMMS) LAI3g at an 8-km spatial resolution and 15-day precipitation feedbacks can alleviate greening-induced soil drying for interval for the period 1982–2011 (20). This global LAI product was Downloaded from http://advances.sciencemag.org/ on April 2, 2020 some locations, involving the role of atmospheric teleconnections. generated from the GIMMS Advanced Very High Resolution Radi- Such assessment must rely on high-resolution GCMs that have ac- ometer (AVHRR) Normalized Difference Vegetation Index (NDVI) curate representations of large-scale atmospheric circulation pro- 3g data and Moderate Resolution Imaging Spectroradiometer (MODIS) cesses. Furthermore, root-zone soil moisture relevant to deep-rooted LAI data using an artificial neural network algorithm (20). plants is not directly measurable at large spatial scales. Hence, cur- The GIMMS LAI3g data have been demonstrated to have good rently, some uncertainties exist in soil moisture estimates from both quality through direct comparisons with field measurements and inter- observation-driven products and GCMs, in terms of the magnitude, comparisons with other satellite-derived LAI products (35). Monthly interannual variability, and trends. Future improvements of root-zone LAI maps were obtained by calculating the maximum value of the soil moisture may be fulfilled by integrating in situ soil moisture mea- two composites in each month at the original 8-km resolution and surements at different depths, microwave-based surface soil mois- were then aggregated to each 0.5° × 0.5° grid cell that was used for ture, and robust surface–to–root zone scaling approaches [e.g., (33)]. our model simulation during the study period. This will enable additional verification of our identified soil mois- Besides LAI, other observed boundary conditions prescribed for ture responses to the ongoing earlier greening and also aid enhanced IPSL-CM include SST and SIF maps and atmospheric CO2 concen- soil moisture parameterizations in the next generation of global climate trations, all covering the study period 1982–2011. Monthly SST and models (CGMs). While we are confident of the broad applicability SIF datasets were both obtained from the Atmospheric Model Inter- of our findings, such improved simulations of soil moisture will comparison Project (https://esgf-node.llnl.gov/news/list/input4mips/) allow refinement of magnitudes of our discovered feedbacks. with a 1° × 1° spatial resolution. Global annual atmospheric CO2 con- Human land-use changes may have strong imprints on the sea- centrations were derived from those used for transient simulations sonal linkage between natural plants and soil moisture. The phenology of the TRENDY project (“Trends and drivers of the regional scale of crops is strongly controlled by management actions, including fer- sources and sink of carbon dioxide”, http://dgvm.ceh.ac.uk/node/9). tilizer use, cropping calendars, multicropping, quick-growing hy- Two independent SWC-related datasets covering the period of brid cultivars, pest control, and, notably, irrigation. The latter can 1982–2011 were used for detecting coupled patterns of variability with effectively sustain soil moisture of irrigated croplands throughout the changes in spring LAI and for evaluating modeling results. The first growing season, causing different correlations between crop growth one is the satellite-derived root-zone soil moisture from GLEAM v3.2a and soil moisture. However, at present, many agricultural practices (21). The GLEAM dataset fully assimilates reanalysis-based tempera- are not explicitly represented in GCMs, which limits our capacity to ture and radiation and remotely sensed top-layer soil moisture, L-band simulate the hydrological consequences of cropland managements. vegetation optical depth (VOD), and snow water equivalents from As we illustrate how vegetation greening feedbacks can be strongly different satellite sensors (21) and also uses multisource precipita- modulated by human land use, this demonstrates the importance for tion forcings. The assimilation of surface soil moisture allows us to the next generation of GCMs to include dynamic crop growth and diagnose where irrigation in agricultural regions is perturbing the realistic irrigation schemes. This will enable more accurate simulation hydrological cycle and the vegetation–soil moisture feedbacks be- of vegetation–soil moisture interactions with human interventions. tween spring and summer. Meanwhile, the forced precipitation and In addition, the hydrological feedbacks of afforestation, reforesta- implied differences in hydrological balance provide data-based in- tion, and forest management practices over regions like China are dications of where atmospheric teleconnections are transporting ex- not assessed separately from those of natural vegetation greening. tra moisture. In addition, GLEAM incorporates VOD to calculate the Deliberately planted trees are expected to have similar variations in effects of water stress, heat stress, and phenological constraints on LAI as existing forests, causing similar alterations to the water cycle ET (21). In this way, it indirectly accounts for the influence of plant via changing transpiration and regional to large-scale precipitation phenological changes on soil moisture. VOD and LAI are intrinsically patterns. Nonetheless, the ability of plants to transpire water varies related, but this does not generate circularity as we test for the link with species composition and forest age (34), which should be taken between GLEAM summer SWC and spring LAI measurements from Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 8 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE optical satellite data. This is because VOD is controlled not only by with satellite-observed annually varying monthly LAI maps for plant structure (in part characterized by LAI) but also by plant hy- 1982–2011 (including MAM). The other experiment (abbreviated draulic strategies (i.e., stomatal and xylem regulation on water use) as ​​LAI​MAM clim​ ​​) was also prescribed with satellite-observed LAI maps (36). The seasonal development of L-band VOD and LAI is proven identically to the ​​LAI​MAM obs​ ​​experiment, except for spring months. to differ considerably for some biomes, such as woody savannas, of During MAM, the LAI maps for each year were replaced with average which water is mainly stored in stems or branches (36). The GLEAM monthly values over the entire study period (1982–2011). Spring SWC data have been extensively used in examining land-atmosphere greening in the real world may affect soil moisture by coupling to feedbacks and benchmarking climate model simulations [e.g., (16)]. changes in oceanic and atmospheric conditions; hence, both experi- Note that the GLEAM root-zone soil moisture, of which the maxi- ments were fed additionally with other boundary conditions of SST, mum depth ranges between 10 and 250 cm depending on local veg- SIF, and atmospheric CO2 concentration. With this numerical ex- etation types, was here treated as a substitute for total soil moisture. periment design, the fingerprint of earlier greening on summer hy- The other soil moisture–related data are GRACE-REC, a global drological variables can be calculated as the difference between the reconstruction of TWS anomalies based on the GRACE satellite retrievals (22). The past TWS anomalies were reconstructed using a sta- two experiments (i.e., ​​LAI​MAM MAM obs​ ​​–​​LAI​clim​ ​​). The original or MAM-­ tistical model forced by daily precipitation and temperature anoma- averaged LAI data–based time series (20) were used to prescribe lies, but trained with original TWS data from GRACE satellites for LAI for each GCM grid cell and for each PFT within a grid cell. 2002 to present (22). The TWS anomalies, expressed as the equiva- To reduce the effect of initial condition uncertainty due to internal lent water thickness (in centimeters), include not only soil moisture climate variability, we conducted an ensemble of 30 simulations for in the upper layers but also snowpack, canopy water, riverine and each experiment with different initial conditions. The initial condi- Downloaded from http://advances.sciencemag.org/ on April 2, 2020 lacustrine water, and groundwater. Removing these non–soil mois- tions were derived from the model outputs of an undisturbed 60-year ture components from the TWS anomalies is technically challeng- run (i.e., with average LAI for all seasons, average SST, SIF, and at- ing, so some caution is needed when interpreting the TWS-based mospheric CO2 concentrations): the first 30 years of the run for the results (Fig. 1, B, E, and F). spin-up of the simulation and then each of the last 30 years of the We also obtained a range of observation-based datasets to evalu- run for providing an alternative initial condition. All our model-based ate model performance for hydrological variables besides SWC (ET analyses were based on the 30-member ensemble mean of each tran- and P). No direct observations of global land ET have been reported, sient experiment. Only vegetated areas, defined as regions with mul- so we used two satellite-derived reconstructions of long-term land tiyear (1982–2011) average LAI ≥ 0.1, were considered to reduce the ET. The first, FLUXNET-MTE, is based on field eddy-covariance mea- spurious influence of bare and sparsely vegetated areas on the over- surements (37), and the second, WB-MTE, is based on information all trend in LAI. All quantities averaged over the Northern Hemi- for precipitation, river discharge, and changes in water mass balance sphere (25° to 90°N) were weighted by area, unless otherwise noted. (38). Both ET products used satellite measurements to scale up either site- or catchment-level ET observations to a 0.5° × 0.5° global grid Indices of extreme hot temperatures cell using a machine learning algorithm. Gridded data for global We used two metrics to measure modeled changes in summer ex- monthly precipitation were obtained from the Global Precipitation treme hot temperatures due to the earlier spring vegetation onset, Climatology Centre (39) version 6.0 and the Climate Research Unit following a previous study of extreme temperature detection (30). (CRU) TS3.22 datasets (40). One metric measures the frequency of extreme hot temperatures, ex- pressed as the number of summer days, with maximum daily tem- Model and experimental setup peratures exceeding the 90% percentile of a baseline period (TX90p). IPSL-CM (24) is a coupled land-atmosphere GCM, versions of which The percentiles were calculated for each calendar day over a 15-day continuously contribute to the Coupled Model Intercomparison moving window using the 30-year simulation forced by observed LAI Project initiatives, thus informing assessment reports of the Inter- and other boundary conditions (i.e., ​​LAI​MAM obs​ )​​ (1982–2011). The other governmental Panel on Climate Change. We used the IPSL-CM ver- metric measures the intensity of extreme hot temperature events, sion 4 GCM, which consists of the atmospheric component LMDZ expressed as the summer maximum of daily maximum temperature (the Laboratoire de Météorologie Dynamique general circulation model (TXx). Both indices were first calculated for both 30-ensemble mean with zooming capability) and the land surface component ORCHIDEE experiments (​​LAI​MAM MAM obs​ ​​ and ​​LAI​clim​ ​​), and changes in TX90p and TXx (41). The soil hydrology within ORCHIDEE consists of two mois- induced by earlier spring greening were then calculated as the dif- ture layers with a constant total depth of 2-m water content per unit ference between the indices derived from the two simulations with soil volume. Vegetation in this land surface model is described glob- and without spring LAI changes. ally with 13 plant functional types (PFTs), each using different gov- erning equations and taking different parameter values to represent Statistical analysis the onset and offset processes of the leafy season (41). Using the IPSL-CM We applied partial correlation analysis in this study to evaluate the coupled model, changes in hydrological fluxes due to vegetation relationship between spring LAI and summer soil moisture after statis- biophysical feedbacks were determined by replacing the modeled tically controlling for the covarying effects of summer LAI and cli- LAI in the ORCHIDEE land surface model with satellite observa- matic variables (including summer temperature, precipitation, and tions prescribed for each of the 13 PFTs. solar radiation). These climatic variables were all derived from the CRU To investigate how spring vegetation feedbacks modulate sum- TS3.22 datasets. We also used the MCA analysis to characterize the mer soil moisture patterns in the last 30 years, we conducted a pair coupling patterns between spring LAI and summer SWC. This ap- of transient GCM simulations prescribed with different spring LAI proach is also known as the singular value decomposition (SVD) (23), conditions. One experiment (abbreviated as ​​LAIMAM ​ bs​ ​​) was prescribed o which can be considered as a generalization of principal components Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 9 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE analysis. MCA has been widely used in environmental sciences to 9. S. Wolf, T. F. Keenan, J. B. Fisher, D. D. Baldocchi, A. R. Desai, A. D. Richardson, R. L. Scott, B. E. Law, M. E. Litvak, N. A. Brunsell, W. Peters, I. T. van der Laan-Luijkx, Warm spring detect teleconnection patterns in the ocean-atmosphere-land system reduced carbon cycle impact of the 2012 US summer drought. Proc. Natl. Acad. Sci. U.S.A. [e.g., (42)]. The MCA analysis looks for pairs of spatial patterns and 113, 5880–5885 (2016). associated time series (or called expansion coefficients) that efficiently 10. A. Angert, S. Biraud, C. Bonfils, C. C. Henning, W. Buermann, J. Pinzon, C. J. Tucker, I. Fung, explain the covariance between two space-time data fields. This was Drier summers cancel out the CO2 uptake enhancement induced by warmer springs. realized by performing SVD on the cross-covariance matrix between Proc. Natl. Acad. Sci. U.S.A. 102, 10823–10827 (2005). 11. W. Buermann, P. R. Bikash, M. Jung, D. H. Burn, M. Reichstein, Earlier springs decrease the two fields, yielding a set of singular vectors and singular values peak summer productivity in North American boreal forests. Environ. Res. Lett. 8, that are associated with pairs of spatial patterns. The time series for 024027 (2013). each pattern pair (fig. S4) were calculated by projecting one data field 12. W. Buermann, M. Forkel, M. O’Sullivan, S. Sitch, P. Friedlingstein, V. Haverd, A. K. Jain, onto the respective singular vector. The correlation coefficient be- E. Kato, M. Kautz, S. Lienert, D. Lombardozzi, J. E. M. S. Nabel, H. Tian, A. J. Wiltshire, D. Zhu, W. K. Smith, A. D. Richardson, Widespread seasonal compensation effects of spring tween the time series of the two fields indicates how strongly related warming on northern plant productivity. Nature 562, 110–114 (2018). the coupled patterns are. Heterogeneous regression maps (Fig. 1, C to F, 13. S. Sippel, M. Forkel, A. Rammig, K. Thonicke, M. Flach, M. Heimann, F. E. L. Otto, and fig. S6, A and B) were generated by regressing the gridded M. Reichstein, M. D. Mahecha, Contrasting and interacting changes in simulated spring values of one field to the time series of the pattern amplitude of the and summer carbon cycle extremes in European ecosystems. Environ. Res. Lett. 12, other. Detailed mathematical basis of the MCA algorithm is described 075006 (2017). 14. R. J. Van Der Ent, H. H. G. Savenije, B. Schaefli, S. C. Steele-Dunne, Origin and fate of in section S2. atmospheric moisture over continents. Water Resour. Res. 46, W09525 (2010). 15. Z. Zeng, S. Piao, L. Z. X. Li, T. Wang, P. Ciais, X. Lian, Y. Yang, J. Mao, X. Shi, R. B. Myneni, SUPPLEMENTARY MATERIALS Impact of Earth greening on the terrestrial water cycle. J. Climate 31, 2633–2650 (2018). Supplementary material for this article is available at http://advances.sciencemag.org/cgi/ 16. Y. Li, S. Piao, L. Z. X. Li, A. Chen, X. Wang, P. Ciais, L. Huang, X. Lian, S. Peng, Z. Zeng, Downloaded from http://advances.sciencemag.org/ on April 2, 2020 content/full/6/1/eaax0255/DC1 K. Wang, L. Zhou, Divergent hydrological response to large-scale afforestation and Supplementary Text vegetation greening in China. Sci. Adv. 4, eaar4182 (2018). Section S1. Earlier spring vegetation greening revealed by remote sensing 17. B. Quesada, R. Vautard, P. Yiou, M. Hirschi, S. I. Seneviratne, Asymmetric European Section S2. Mathematical basis of the MCA analysis summer heat predictability from wet and dry southern winters and springs. Nat. Clim. Chang. Section S3. Model evaluation using available observations 2, 736–741 (2012). Fig. S1. Spatial patterns of the trend in spring vegetation changes for 1982–2011. 18. A. Berg, K. Findell, B. Lintner, A. Giannini, S. I. Seneviratne, B. van den Hurk, R. Lorenz, Fig. S2. Grassland partial correlations between ESA-CCI summer SWC and GIMMS spring LAI A. Pitman, S. Hagemann, A. Meier, F. Cheruy, A. Ducharne, S. Malyshev, P. C. D. Milly, during 1982–2011. Land-atmosphere feedbacks amplify aridity increase over land under global warming. Fig. S3. Partial correlations between spring LAI and either summer SWC or summer TWS, as Nat. Clim. Chang. 6, 869–874 (2016). affected by agricultural extent. 19. A. L. S. Swann, I. Y. Fung, J. C. H. Chiang, Mid-latitude afforestation shifts general Fig. S4. Time series (or expansion coefficients) of the leading first MCA patterns. circulation and tropical precipitation. Proc. Natl. Acad. Sci. U.S.A. 109, 712–716 (2012). Fig. S5. Validation of IPSL-simulated hydrological variables. 20. Z. Zhu, J. Bi, Y. Pan, S. Ganguly, A. Anav, L. Xu, A. Samanta, S. Piao, R. R. Nemani, Fig. S6. Robustness test of the IPSL-simulated relationship between spring LAI and summer SWC. R. B. Myneni, Global data sets of vegetation leaf area index (LAI)3g and Fraction of Fig. S7. Comparison of trends of observation-based and IPSL-simulated ET. Photosynthetically Active Radiation (FPAR)3g derived from Global Inventory Modeling Fig. S8. Spatial patterns of IPSL-simulated trends in key hydrological variables induced by and Mapping Studies (GIMMS) Normalized Difference Vegetation Index (NDVI3g) for the spring LAI changes. period 1981 to 2011. Remote Sens. 5, 927–948 (2013). Fig. S9. IPSL-simulated average global prevailing wind at 10-m height. 21. B. Martens, D. G. Miralles, H. Lievens, R. van der Schalie, R. A. M. de Jeu, D. Fernández-Prieto, Fig. S10. Spatial patterns of IPSL-simulated trends in summer extreme hot temperature indices H. E. Beck, W. A. Dorigo, N. E. C. Verhoest, GLEAM v3: Satellite-based land evaporation induced by spring LAI changes. and root-zone soil moisture. Geosci. Model Dev. 10, 1903–1925 (2017). Reference (43) 22. V. Humphrey, L. Gudmundsson, S. I. Seneviratne, A global reconstruction of climate-driven subdecadal water storage variability. Geophys. Res. Lett. 44, 2300–2309 (2017). 23. C. S. Bretherton, C. Smith, J. M. Wallace, An intercomparison of methods for finding REFERENCES AND NOTES coupled patterns in climate data. J. Climate 5, 541–560 (1992). 1. J. Barichivich, K. R. Briffa, R. B. Myneni, T. J. Osborn, T. M. Melvin, P. Ciais, S. Piao, C. Tucker, 24. O. Marti, P. Braconnot, J.-L. Dufresne, J. Bellier, R. Benshila, S. Bony, P. Brockmann, Large-scale variations in the vegetation growing season and annual cycle of atmospheric P. Cadule, A. Caubel, F. Codron, N. de Noblet, S. Denvil, L. Fairhead, T. Fichefet, CO2 at high northern latitudes from 1950 to 2011. Glob. Chang. Biol. 19, 3167–3183 M.-A. Foujols, P. Friedlingstein, H. Goosse, J. Y. Grandpeix, E. Guilyardi, F. Hourdin, (2013). A. Idelkadi, M. Kageyama, G. Krinner, C. Lévy, G. Madec, J. Mignot, I. Musat, 2. S. Piao, J. Tan, A. Chen, Y. H. Fu, P. Ciais, Q. Liu, I. A. Janssens, S. Vicca, Z. Zeng, S.-J. Jeong, D. Swingedouw, C. Talandier, Key features of the IPSL ocean atmosphere model and its Y. Li, R. B. Myneni, S. Peng, M. Shen, J. Peñuelas, Leaf onset in the northern hemisphere sensitivity to atmospheric resolution. Climate Dynam. 34, 1–26 (2010). triggered by daytime temperature. Nat. Commun. 6, 6911 (2015). 25. A. Berg, J. Sheffield, P. C. D. Milly, Divergent surface and total soil moisture projections 3. Q. Liu, Y. H. Fu, Z. Zhu, Y. Liu, Z. Liu, M. Huang, I. A. Janssens, S. Piao, Delayed autumn under global warming. Geophys. Res. Lett. 44, 236–244 (2017). phenology in the Northern Hemisphere is related to change in both climate and spring 26. X. Lian, S. Piao, C. Huntingford, Y. Li, Z. Zeng, X. Wang, P. Ciais, T. R. McVicar, S. Peng, phenology. Glob. Chang. Biol. 22, 3702–3711 (2016). C. Ottlé, H. Yang, Y. Yang, Y. Zhang, T. Wang, Partitioning global land evapotranspiration 4. M. D. Schwartz, R. Ahas, A. Aasa, Onset of spring starting earlier across the Northern using CMIP5 models constrained by observations. Nat. Clim. Chang. 8, 640–646 (2018). Hemisphere. Glob. Chang. Biol. 12, 343–351 (2006). 27. M. Stéfanon, P. Drobinski, F. D’Andrea, N. De Noblet-Ducoudré, Effects of interactive 5. A. Menzel, T. H. Sparks, N. Estrella, E. Koch, A. Aasa, R. Ahas, K. Alm-Kübler, P. Bissolli, vegetation phenology on the 2003 summer heat waves. J. Geophys. Res. 117, D24103 (2012). O. Braslavská, A. Briede, F. M. Chmielewski, Z. Crepinsek, Y. Curnel, Å. Dahl, C. Defila, 28. M. Hirschi, S. I. Seneviratne, V. Alexandrov, F. Boberg, C. Boroneant, O. B. Christensen, A. Donnelly, Y. Filella, K. Jatczak, F. Måge, A. Mestre, Ø. Nordli, J. Peñuelas, P. Pirinen, H. Formayer, B. Orlowsky, P. Stepanek, Observational evidence for soil-moisture impact V. Remišová, H. Scheifinger, M. Striz, A. Susnik, A. J. H. Van Vliet, F.-E. Wielgolaski, S. Zach, on hot extremes in southeastern Europe. Nat. Geosci. 4, 17–21 (2011). A. Zust, European phenological response to climate change matches the warming 29. M. M. Vogel, R. Orth, F. Cheruy, S. Hagemann, R. Lorenz, B. J. J. M. Hurk, S. I. Seneviratne, pattern. Glob. Chang. Biol. 12, 1969–1976 (2006). Regional amplification of projected changes in extreme temperatures strongly controlled 6. Y. H. Fu, S. Piao, M. op de Beeck, N. Cong, H. Zhao, Y. Zhang, A. Menzel, I. A. Janssens, by soil moisture-temperature feedbacks. Geophys. Res. Lett. 44, 1511–1519 (2017). Recent spring phenology shifts in western Central Europe based on multiscale 30. S. E. Perkins, L. V. Alexander, On the measurement of heat waves. J. Climate 26, observations. Glob. Ecol. Biogeogr. 23, 1255–1263 (2014). 4500–4517 (2013). 7. C. D. Keeling, J. F. S. Chin, T. P. Whorf, Increased activity of northern vegetation inferred 31. N. Christidis, G. S. Jones, P. A. Stott, Dramatically increasing chance of extremely hot from atmospheric CO2 measurements. Nature 382, 146–149 (1996). summers since the 2003 European heatwave. Nat. Clim. Chang. 5, 46–50 (2014). 8. S. Piao, P. Friedlingstein, P. Ciais, N. Viovy, J. Demarty, Growing season extension and its 32. D. L. Swain, D. E. Horton, D. Singh, N. S. Diffenbaugh, Trends in atmospheric patterns impact on terrestrial carbon cycle in the Northern Hemisphere over the past 2 decades. conducive to seasonal precipitation and temperature extremes in California seasonal Global Biogeochem. Cycles 21, GB3018 (2007). precipitation and temperature extremes. Sci. Adv. 2, e1501344 (2016). Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 10 of 11 SCIENCE ADVANCES | RESEARCH ARTICLE 33. C. Albergel, C. Rüdiger, T. Pellarin, J. C. Calvet, N. Fritz, F. Froissard, D. Suquia, A. Petitpa, 42. A. Dai, Increasing drought under global warming in observations and models. B. Piguet, E. Martin, From near-surface to root-zone soil moisture using an exponential Nat. Clim. Chang. 3, 52–58 (2013). filter: An assessment of the method based on in-situ observations and model simulations. 43. S. Piao, Q. Liu, A. Chen, I. A. Janssens, Y. Fu, J. Dai, L. Liu, X. Lian, M. Shen, X. Zhu, Hydrol. Earth Syst. Sci. 12, 1323–1337 (2008). Plant phenology and global climate change: Current progresses and challenges. 34. S. T. Brantley, C. F. Miniat, P. V. Bolstad, Rainfall partitioning varies across a forest age Glob. Chang. Biol. 25, 1922–1940 (2019). chronosequence in the southern Appalachian Mountains. Ecohydrology 12, e2081 (2019). 35. H. E. Beck, T. R. McVicar, A. I. J. M. van Dijk, J. Schellekens, R. A. M. de Jeu, L. A. Bruijnzeel, Acknowledgments: We acknowledge the National Supercomputer Center in Tianjin, China Global evaluation of four AVHRR–NDVI data sets: Intercomparison and assessment (NSCC-JI, http://www.nscc-tj.cn/) and the National Computer Center IDRIS of CNRS in France for against Landsat imagery. Remote Sens. Environ. 115, 2547–2563 (2011). providing computing resources. Funding: This study was supported by the National Key Research 36. F. Tian, J.-P. Wigneron, P. Ciais, J. Chave, J. Ogée, J. Peñuelas, A. Ræbild, J.-C. Domec, and Development Program of China (grant 2018YFC1507704), the Second Tibetan Plateau X. Tong, M. Brandt, A. Mialon, N. Rodriguez-Fernandez, T. Tagesson, A. al-Yaari, Y. Kerr, Scientific Expedition and Research Program (Grant No. 2019QZKK0208), the National Natural C. Chen, R. B. Myneni, W. Zhang, J. Ardö, R. Fensholt, Coupling of ecosystem-scale plant Science Foundation of China (41561134016 and 41988101), the National Youth Top-notch Talent water storage and leaf phenology observed by satellite. Nat. Ecol. Evol. 2, 1428–1435 (2018). Support Program in China, and the 111 Project (B14001). C.H. acknowledges the Center for 37. M. Jung, M. Reichstein, P. Ciais, S. I. Seneviratne, J. Sheffield, M. L. Goulden, G. Bonan, Ecology and Hydrology (CEH) National Capability fund. P.C., I.A.J., and J.P. were funded by the A. Cescatti, J. Chen, R. de Jeu, A. J. Dolman, W. Eugster, D. Gerten, D. Gianelle, N. Gobron, European Research Council Synergy (grant SyG-2013-610028 IMBALANCE-P). Author J. Heinke, J. Kimball, B. E. Law, L. Montagnani, Q. Mu, B. Mueller, K. Oleson, D. Papale, contributions: S.P. designed the research. X.L. performed analysis. X.L. and S.P. drafted the paper. A. D. Richardson, O. Roupsard, S. Running, E. Tomelleri, N. Viovy, U. Weber, C. Williams, All authors contributed to the interpretation of the results and to the text. Competing interests: E. Wood, S. Zaehle, K. Zhang, Recent decline in the global land evapotranspiration trend The authors declare that they have no competing interests. Data and materials availability: All due to limited moisture supply. Nature 467, 951–954 (2010). methods needed to evaluate the conclusions in the paper are present in the paper and/or the 38. Z. Zeng, T. Wang, F. Zhou, P. Ciais, J. Mao, X. Shi, S. Piao, A worldwide analysis Supplementary Materials. All data related to this paper, along with code for interpretation of data of spatiotemporal changes in water balance-based evapotranspiration from 1982 and preparation of diagrams, are available on request from the corresponding author. to 2009. J. Geophys. Res. Atmos. 119, 1186–1202 (2014). 39. A. Becker, A. Meyer-Christoffer, B. Rudolf, K. Schamm, U. Schneider, M. Ziese, A description Submitted 15 February 2019 Downloaded from http://advances.sciencemag.org/ on April 2, 2020 of the global land-surface precipitation data products of the Global Precipitation Accepted 6 November 2019 Climatology Centre with sample applications including centennial (trend) analysis from Published 3 January 2020 1901–present. Earth Syst. Sci. Data 5, 71–99 (2013). 10.1126/sciadv.aax0255 40. I. Harris, P. D. Jones, T. J. Osborn, D. H. Lister, Updated high-resolution grids of monthly climatic observations—The CRU TS3.10 Dataset. Int. J. Climatol. 34, 623–642 (2014). Citation: X. Lian, S. Piao, L. Z. X. Li, Y. Li, C. Huntingford, P. Ciais, A. Cescatti, I. A. Janssens, 41. G. Krinner, N. Viovy, N. de Noblet-Ducoudré, J. Ogée, J. Polcher, P. Friedlingstein, P. Ciais, J. Peñuelas, W. Buermann, A. Chen, X. Li, R. B. Myneni, X. Wang, Y. Wang, Y. Yang, Z. Zeng, S. Sitch, I. C. Prentice, A dynamic global vegetation model for studies of the coupled Y. Zhang, T. R. McVicar, Summer soil drying exacerbated by earlier spring greening of northern atmosphere-biosphere system. Global Biogeochem. Cycles 19, GB1015 (2005). vegetation. Sci. Adv. 6, eaax0255 (2020). Lian et al., Sci. Adv. 2020; 6 : eaax0255 3 January 2020 11 of 11 Summer soil drying exacerbated by earlier spring greening of northern vegetation Xu Lian, Shilong Piao, Laurent Z. X. Li, Yue Li, Chris Huntingford, Philippe Ciais, Alessandro Cescatti, Ivan A. Janssens, Josep Peñuelas, Wolfgang Buermann, Anping Chen, Xiangyi Li, Ranga B. Myneni, Xuhui Wang, Yilong Wang, Yuting Yang, Zhenzhong Zeng, Yongqiang Zhang and Tim R. McVicar Sci Adv 6 (1), eaax0255. DOI: 10.1126/sciadv.aax0255 Downloaded from http://advances.sciencemag.org/ on April 2, 2020 ARTICLE TOOLS http://advances.sciencemag.org/content/6/1/eaax0255 SUPPLEMENTARY http://advances.sciencemag.org/content/suppl/2019/12/20/6.1.eaax0255.DC1 MATERIALS REFERENCES This article cites 43 articles, 5 of which you can access for free http://advances.sciencemag.org/content/6/1/eaax0255#BIBL PERMISSIONS http://www.sciencemag.org/help/reprints-and-permissions Use of this article is subject to the Terms of Service Science Advances (ISSN 2375-2548) is published by the American Association for the Advancement of Science, 1200 New York Avenue NW, Washington, DC 20005. The title Science Advances is a registered trademark of AAAS. Copyright © 2020 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution NonCommercial License 4.0 (CC BY-NC).