Manuscript 1 Title 2  Full title: Summer soil drying exacerbated by earlier spring greening of northern 3 vegetation 4  Short title: Inter-seasonal vegetation-soil moisture feedbacks 5 6 Authors 7 Xu Lian1, Shilong Piao1,2,3,*, Laurent Z. X. Li4, Yue Li1, Chris Huntingford5, Philippe 8 Ciais6, Alessandro Cescatti7, Ivan A. Janssens8, Josep Peñuelas9,10, Wolfgang 9 Buermann11,12, Anping Chen13, Xiangyi Li1, Ranga B. Myneni14, Xuhui Wang1,6, Yilong 10 Wang6, Yuting Yang15, Zhenzhong Zeng1, Yongqiang Zhang16, Tim R. McVicar17,18 11 12 Affiliations 1 13 Sino-French Institute for Earth System Science, College of Urban and Environmental 14 Sciences, Peking University, Beijing 100871, China. 2 15 Key Laboratory of Alpine Ecology, Institute of Tibetan Plateau Research, Chinese 16 Academy of Sciences, Beijing 100101, China. 3 17 Center for Excellence in Tibetan Earth Science, Chinese Academy of Sciences, Beijing 18 100085, China. 4 19 Laboratoire de Météorologie Dynamique, CNRS, Sorbonne Université, Ecole Normale 20 Supérieure, Ecole Polytechnique, Paris, France. 5 21 Centre for Ecology and Hydrology, Wallingford, Oxfordshire, OX10 8BB, U.K. 6 22 Laboratoire des Sciences du Climat et de l’Environnement (LSCE), CEA CNRS UVSQ, 23 91191 Gif Sur Yvette, France. 7 24 European Commission, Joint Research Centre, Directorate for Sustainable Resources, I- Science Advances Manuscript Page 1 of 34 25 21027 Ispra (VA), Italy. 8 26 Department of Biology, University of Antwerp, Universiteitsplein 1, Wilrijk 2610, 27 Belgium. 9 28 CREAF, Cerdanyola del Valles, Barcelona, Catalonia 08193, Spain. 10 29 CSIC, Global Ecology Unit CREAF- CSIC-UAB, Bellaterra, Barcelona, Catalonia 30 08193, Spain. 11 31 Department of Geography, University of Augsburg, Augsburg, Germany. 12 32 Institute of the Environment and Sustainability, University of California, Los Angeles, 33 Los Angeles, CA, USA. 13 34 Department of Forestry and Natural Resources, Purdue University, West Lafayette, IN 35 46907, USA. 14 36 Department of Earth and Environment, Boston University, Boston, Massachusetts 37 02215, USA. 15 38 State Key Laboratory of Hydro-Science and Engineering, Department of Hydraulic 39 Engineering, Tsinghua University, Beijing, China. 16 40 Key Lab of Water Cycle and Related Land Surface Processes, Institute of Geographic 41 Sciences and Natural Resources Research, Chinese Academy of Sciences, Beijing 42 100101, China. 17 43 CSIRO Land and Water, Canberra, Australian Capital Territory, Australia. 18 44 Australian Research Council Centre of Excellence for Climate Extremes, Canberra, ACT, 45 Australia. 46 *Corresponding author. Email: slpiao@pku.edu.cn 47 48 Other information on this paper: Science Advances Manuscript Page 2 of 34 49 Keywords: soil moisture, earlier spring greening, vegetation biophysical feedbacks, 50 terrestrial water cycle, summer heatwaves 51 One-sentence summary: Extra water loss from earlier spring greening dries soils locally 52 in summer, but can moisten soils in areas downwind. 53 Abbreviations: AMIP, Atmospheric Model Intercomparison Project; AVHRR, Advanced 54 Very High Resolution Radiometer; CMIP, Coupled Model Intercomparison Project; 55 CRU, Climate Research Unit; GCM, global climate model; GIMMS, Global Inventory 56 Monitoring and Modeling Studies; GLEAM, Global Land Evaporation Amsterdam 57 Model; ET, evapotranspiration; GPCC, Global Precipitation Climatology Centre; IPCC, 58 Intergovernmental Panel on Climate Change; GRACE-REC, reconstructed TWS 59 anomalies from Gravity Recovery and Climate Experiment satellites; IPSL-CM4, Institut 60 Pierre Simon Laplace Climate Model version 4; JJA, June, July and August; LAI, Leaf M AM M AM 61 Area Index; L A I obs , numerical experiment forced with observed MAM LAI; L A I c lim , 62 numerical experiment forced with climatological MAM LAI; LMDZ, the Laboratoire de 63 Météorologie Dynamique general circulation model with Zooming capability; MAM, 64 March, April and May; MCA, maximum covariance analysis; MODIS, Moderate 65 Resolution Imaging Spectroradiometer; FLUXNET-MTE, the ET product based on eddy- 66 covariance measurements and the model tree ensemble algorithm; WB-MTE, the ET 67 product based on water mass balance and the model tree ensemble algorithm; NDVI, 68 Normalized Difference Vegetation Index; ORCHIDEE, Organising Carbon and 69 Hydrology In Dynamic Ecosystems; P, precipitation; PFTs, plant functional types; Q, 70 surface runoff; SIF, sea-ice fraction; SST, sea-surface temperature; SVD, singular value 71 decomposition; SWC, soil water content; TRENDY, trends and drivers of the regional 72 scale sources and sinks of carbon dioxide; TWS, terrestrial water storage; TXx, the Science Advances Manuscript Page 3 of 34 73 summer maximum of daily maximum temperature; TX90p, the number of summer days 74 with maximum temperatures exceeding the local 90th percentile of a given base period. 75 Science Advances Manuscript Page 4 of 34 76 Abstract 77 Earlier vegetation greening under climate change raises evapotranspiration and 78 thus lowers spring soil moisture, yet the extent and magnitude of this water deficit 79 persistence into the following summer remain elusive. We provide observational 80 evidence that increased foliage cover over the Northern Hemisphere, during 1982- 81 2011, triggers an additional soil moisture deficit that is further carried over into 82 summer. Climate model simulations independently support this, and attribute the 83 driving process to be larger increases in evapotranspiration than in precipitation. 84 This extra soil drying is projected to amplify the frequency and intensity of summer 85 heatwaves. Most feedbacks operate locally, except for a significant teleconnection 86 where extra moisture transpired over Europe is transported to central Siberia. 87 Model results illustrate that this teleconnection offsets Siberian soil moisture losses 88 from local spring greening. Our results highlight that climate change adaptation 89 planning must account for the extra summer water and heatwave stress inherited 90 from warming-induced earlier greening. 91 Science Advances Manuscript Page 5 of 34 92 MAIN TEXT 93 Introduction 94 A widespread shift in the phenology of northern vegetation under global warming, 95 towards earlier spring leaf emergence (i.e. earlier spring “greening”), has been 96 unambiguously documented using indices of satellite-monitored vegetation greenness (1– 97 3) and ground observations (4–6). Both modeling and observational studies have shown 98 that the early onset of vegetation phenology causes an increase in ecosystem productivity, 99 through an extension of the active carbon uptake period (7, 8). On the other hand, 100 however, there is a growing concern that the earlier vegetation greening may induce soil 101 moisture deficits in late spring through increased evapotranspiration (ET), which may 102 amplify summer droughts and have detrimental effects on summer ecosystems (9–13). 103 Direct observational evidence for this ET-centered legacy hypothesis is nonetheless very 104 scarce. Furthermore, the influence of earlier spring greening on summer soil moisture can 105 be more complex than just via modifying ET. For instance, more rapid spring soil 106 moisture losses raise moisture levels of the atmosphere, facilitates cloud formation, and 107 triggers extra precipitation (14–16). Any additional rainfall may compensate, at least 108 partly, the greening-induced deficit of soil moisture. A further biophysical feedback is 109 that soil moisture deficits result in less conversion of surface incident sunlight into water 110 evaporation. Lower latent heat from the land may suppress the development of clouds and 111 subsequent precipitation (17, 18). Beyond local effects, these feedback mechanisms may 112 also operate as teleconnections between locations. Large-scale atmospheric circulation 113 (14, 19) redistributes the extra recycled moisture across the Earth’s land surface and 114 contributes to additional precipitation in downwind regions. The overall impact of these 115 potential inter-seasonal feedbacks on summer water cycles, unfortunately, is still poorly 116 characterized. Science Advances Manuscript Page 6 of 34 117 118 Quantifying this vegetation-soil moisture legacy effect must take full account of the 119 complete loop of soil-vegetation-atmosphere exchanges of water and energy fluxes. 120 Many feedback mechanisms, as summarized in the above paragraph, are essential for 121 accurate assessments of global or regional soil moisture changes, yet are ignored by 122 previous studies (9–12). In particular, without characterising atmospheric feedbacks and 123 including those that are inter-seasonal, it still remains an open question on the extent to 124 which vegetation greening influences summer soil moisture. Here, we first analyse 125 multiple strands of direct observational data to characterise the spring-to-summer 126 vegetation-soil moisture feedbacks. We then perform factorial experiments with a state- 127 of-the-art coupled land-atmosphere climate model to better understand the underlying 128 mechanisms. 129 130 Results 131 Observational Evidence for Summer Drying Induced by Earlier Spring Greening 132 To examine spatiotemporal connections between recent spring phenological changes 133 (section S1 and fig. S1) and soil water content (SWC) dynamics in the following summer, 134 we first analyse long-term satellite data of spring Leaf Area Index (20) (LAI, as a proxy 135 for vegetation phenology and growth, Methods), and summer values of two observation- 136 driven SWC-related data sets, for the period of 1982–2011. The two SWC-related data 137 include root-zone soil moisture from the Global Land Evaporation Amsterdam Model 138 (21) (GLEAM; Methods) and the long-term reconstruction of terrestrial water storage 139 (TWS) anomalies based on Gravity Recovery and Climate Experiment satellites (22) 140 (GRACE-REC). Our study domain only covers the middle to high latitudes of the 141 Northern Hemisphere (25–90°N), because of their distinctive vegetation seasonality and Science Advances Manuscript Page 7 of 34 142 particularly pronounced recent changes in springtime LAI phenology (fig. S1). We define 143 spring across the domain as the period from March to May (MAM), and summer as the 144 period from June to August (JJA). We calculate partial correlations between spring LAI 145 and summer SWC (or TWS), whereby covarying effects of summer LAI and climatic 146 variables (i.e. temperature, precipitation and solar radiation) are controlled for (Methods). 147 As potential non-local precipitation feedbacks modulate the local relationship between 148 spring LAI and summer SWC (14, 19), we also implement a lagged maximum covariance 149 analysis (23) (MCA; Methods, section S2). This statistical approach isolates pairs of 150 coupled spatial patterns that explain a maximum fraction of the covariance between two 151 space-time data fields (in our case, spring LAI and summer SWC). Hence, the MCA 152 method has obvious advantages over routine correlation analysis that only partially 153 describes covariability between two fields (i.e. spatial correlation using geographical 154 points only, or temporal correlation using time series only). 155 156 Partial correlations between spring greenness and subsequent summer SWC show a 157 widespread negative pattern across most northern lands (64% of land area for GLEAM 158 SWC and 56% for GRACE-REC TWS) (Fig. 1, A and B), including significantly 159 greening regions of Europe, eastern U.S., east Asia and Siberia (fig. S1). For temperate 160 and semi-arid grasslands where LAI is coupled to shallow moisture, this dominance of 161 negative correlation between spring LAI and summer moisture is additionally supported 162 by using microwave-derived summer top-layer soil moisture (fig. S2). The widespread 163 negative correlation indicates that, when effects of inter-annual variations of precipitation 164 are factored out, earlier greening in spring generally results in local soil drying in 165 subsequent summers. On the other side, positive correlation is also apparent over regions 166 including the North China Plain, central and southeastern Europe and the U.S. Great Science Advances Manuscript Page 8 of 34 167 Plains, especially for TWS that additionally captures storage variations in canopy, snow, 168 streams, lakes and groundwater (Fig. 1, A and B). Many such regions are influenced by 169 intensive agricultural activities, which potentially decouples the relationship between 170 spring LAI and summer SWC observed in natural ecosystems. By averaging land surface 171 undisturbed by cropland management, we find a negative correlation between spring LAI 172 and summer SWC that is statistically significant (-0.44 for GLEAM SWC and -0.40 for 173 GRACE-REC TWS, p < 0.05) (fig. S3). However, an increasing fraction of land areas 174 affected by croplands leads to a weakening of this negative correlation, which shifts 175 robustly to being positive when croplands occupy over 50% of local land areas (fig. S3). 176 177 Fig. 1, C-F presents the leading first MCA mode characterising spatial patterns of the 178 coupling between spring LAI and summer SWC. This mode in particular includes 179 potential non-local impacts of LAI-induced moisture transport anomalies on SWC at 180 downwind regions. This mode pair explains 46% and 59% of the total squared covariance 181 between the two variables, for GLEAM SWC and GRACE-REC TWS respectively. For 182 both SWC-related products, the first spring LAI mode shows widespread positive 183 anomalies (i.e. earlier greening) in the majority of the Eurasia (> 80%) and eastern U.S., 184 while negative anomalies (i.e. delayed greening, or browning) covers vast areas of central 185 and western North America (Fig. 1, C and E). The patterns of the leading LAI mode 186 generally resemble the trend pattern of satellite-observed spring LAI (fig. S1). As a 187 consequence, the leading summer SWC (or TWS) mode (Fig. 1, D and F) coupled to the 188 spring LAI patterns demonstrates a direct response pattern of summer SWC to the 189 emerging spring greening. The MCA-based annual time series linked to the spring LAI 190 and summer SWC (or TWS) pattern are also highly correlated (r = 0.86 for GLEAM 191 SWC, r = 0.80 for GRACE-REC TWS, p < 0.05 for both; fig. S4). This reiterates the Science Advances Manuscript Page 9 of 34 192 strong coupling between the two coupled patterns shown in Fig. 1, C-F. Generally, 193 negative summer SWC anomalies (i.e. soil drying) are preceded by spring greening, 194 whereas positive SWC anomalies (i.e. soil wetting) dominate the western North America 195 with spring browning (Fig. 1, D and F versus Fig. 1, C and E). The dominance of 196 negative LAI-SWC coupling, as also seen by partial correlation patterns, implies an 197 impact of earlier spring onset on local summer soil drying. The cropland-dominated 198 central Europe have wetter summers during years with greener springs (Fig. 1, C-F), 199 which is in agreement with the positive correlation between spring greenness and summer 200 soil moisture (Fig. 1, A and B). A notable exception occurs in the Siberia region, where 201 summer SWC increases despite the presence of a clear positive spring LAI anomaly and 202 limited human interventions. This is also in sharp contrast to the negative partial 203 correlations in central Siberia after removing covarying effects of precipitation (Fig. 1, A 204 and B). Such difference implies that altered precipitation patterns driven by spring LAI 205 changes cause extra moisture to be brought into this area. 206 207 Simulating the Spring-to-summer Vegetation Feedbacks 208 Statistical analyses of observation-based data show that summertime soil moisture 209 anomalies negatively co-vary with the preceding spring greening patterns. However, the 210 detection of systematic legacy effects of spring LAI on summer soil moisture within data 211 is partially impacted by noise from weather variability, potential non-linearity in soil 212 moisture responses, different rooting depth and moisture access across ecosystems, and 213 slowly varying drivers of both LAI and soil moisture trends such as vapour pressure 214 deficit and rising atmospheric CO2. To confirm a causal link between spring greening and 215 summer soil drying, and to better understand the local and non-local mechanisms, we 216 perform further analyses using the IPSL-CM4 coupled land-atmosphere global climate Science Advances Manuscript Page 10 of 34 217 model (24) (GCM) (Methods). Coupled land-atmosphere models describe explicitly 218 biophysical interactions between the land surface and the atmosphere. This includes local 219 exchanges of water and energy fluxes, and large-scale changes of atmospheric circulation 220 and moisture convergence. Paired time-evolving (i.e. transient) numerical experiments 221 are run with the land surface component of IPSL-CM4, ORCHIDEE (Methods), forced 222 with two contrasting spring LAI conditions. One set of simulations (labelled as “ L A I oMbAs M 223 ”) is prescribed with annually varying LAI from satellite observations (20) during all 224 seasons, including MAM. The second set of simulations (labelled as “ L A I cMlimA M ”) is also 225 prescribed with annually varying LAI during all seasons, except in MAM when it is fixed 226 to the climatological mean LAI conditions observed during 1982–2011 (Methods). Apart 227 from the spring LAI differences, the two model experiments are forced identically with 228 observed boundary conditions of sea-surface temperature, sea-ice fraction and 229 atmospheric CO2 concentration. This enables us to account for the vegetation biophysical 230 feedbacks stemming from interactions with changes in oceanic processes or CO2 231 radiative/physiological forcings. Multi-member paired ensembles are generated for each 232 coupled-model experiment by performing 30 repeated runs but with different initial 233 conditions (Methods). 234 235 Model evaluation against global observational datasets supports the capacity of the IPSL- 236 CM4 GCM to capture the annual and seasonal features of key hydrological variables, 237 including ET and precipitation (fig. S5, A and B, further details in section S3). To avoid 238 systematic biases from the diverging magnitude of modelled and observation-based 239 SWC, SWC changes are instead compared as the fraction of the present-day 240 climatological mean (1982–2011) (25). At the scale of the Northern Hemisphere (25– 241 90°N), the IPSL-simulated SWC shows a negative trend whereas GLEAM SWC shows a Science Advances Manuscript Page 11 of 34 242 weaker positive trend (fig. S5C). Whilst this sign difference is possibly a concern for 243 global mean changes, we find zonally coherent patterns of the IPSL-simulated SWC trend 244 compared to GLEAM SWC (fig. S5, D-I). The IPSL-simulated SWC broadly captures 245 regions where soil wetting or drying occurs, at both annual and seasonal timescales (fig. 246 S, D-I). This regional consistency strengthens our confidence in examining both the 247 observed zonally different fingerprint of spring greening on summer soil moisture, and 248 the ability of our GCM to explain such changes. 249 250 The difference between the two factorial simulations with variable and fixed spring LAI ( – L A I c lim M AM M AM 251 L A I obs ) allows us to isolate the fingerprint of observed spring LAI changes on 252 summer soil moisture (SWC), as well as other important hydrological variables. The 253 symbol “” hereafter denotes the influence of spring LAI changes on other target 254 variables of interest. Our factorial model simulations confirm a strong coupling between 255 spring LAI and summer SWC anomalies (r = –0.51, p < 0.01) (Fig. 2A). Hence, on a 256 year-to-year basis, the observed higher spring LAI values prescribed to our GCM 257 correspond to lower summer SWC estimates. We then examine if there is a robust overall 258 trend of decreasing SWC, driven by the significantly increasing spring LAI (greening). 259 This effect is more difficult to detect because the inter-annual trend of LAI over the past 260 30 years is smaller than concurrent climate-driven variabilities. Average summer SWC 261 across all the northern lands (25–90°N) tends to decrease at a rate of –0.11±0.09×10-4 m3 262 m-3 yr-1 (linear trend: mean ± one standard derivation, p = 0.22) due to the spring LAI 263 changes. This overall summer SWC trend is negative, in line with our hypothesis that 264 higher green foliage cover is transpiring more in spring and leaving soils drier in summer. 265 The calculated trend is, though, not statistically significant, partly due to the presence of 266 contrasting effects on summer SWC in regions with spring greening (i.e. positive LAI Science Advances Manuscript Page 12 of 34 267 trends) versus spring browning (i.e. negative LAI trends). When instead considering only 268 regions where spring LAI generally increases, a statistically significant and much 269 stronger summer drying trend in summer SWC (–0.17±0.07×10-4 m3 m-3 yr-1, p < 0.05) 270 is detected (Fig. 2B). On average, including areas with spring browning, this decreasing 271 trend of summer SWC associated with spring LAI changes contributes an additional 16% 272 to the simulated decrease in summer SWC directly caused by climate change during 273 1982–2011 (overall decrease: –0.68×10-4 m3 m-3 yr-1 from the M AM L A I obs simulation). 274 275 We also perform the MCA analysis on observed spring LAI and IPSL-simulated summer M AM 276 SWC (from simulation L A I obs ). The coupling patterns between spring LAI and summer 277 SWC for the first mode (fig. S6, A and B) show high accordance with those purely based 278 on observations (Fig. 1, C-F), confirming the ability of our GCM to simulate the inter- 279 seasonal vegetation-soil moisture feedbacks. We further test whether insights based on 280 our GCM-based factorial simulations could be overly influenced by water-stressed 281 regions at low-latitudes. In doing so, we calculate spring LAI and summer SWC 282 anomalies weighted by both area and annual mean soil moisture values. Again, a 283 significant negative correlation is found between spring LAI and summer SWC (fig. S6, 284 C and D). This spring LAI-induced drying trend of summer SWC (–0.09±0.06×10-4 m3 285 m-3 yr-1) is consistent with that without assigning a lower weight to drier regions (Fig. 286 2A), partly because earlier greening mainly occurs in regions with favorable moisture 287 conditions. Therefore, our GCM simulations provide robust process-based evidence that 288 the observed overall earlier greening of northern land causes extra decreases in summer 289 soil moisture. 290 Science Advances Manuscript Page 13 of 34 291 Intra- and inter-seasonal soil moisture anomalies are determined by the balance between 292 changes in ET, precipitation and runoff (Schematic, Fig. 3). To gain more insights into 293 how spring vegetation changes feedback on soil moisture, we further use the paired GCM 294 experiments to examine changes in these major hydrological fluxes. Higher spring LAI 295 levels tend to accelerate water recycling from the land surface to the atmosphere through 296 enhanced plant transpiration. As transpiration dominates total terrestrial ET (26) (average 297 transpiration/ET ratio of IPSL-CM4: 0.47 for spring and 0.76 for summer), a higher 298 spring transpiration directly leads to higher modelled spring ET (+0.28±0.06 mm yr-2, p < 299 0.01) (Fig. 2D). We also note that the overall trend of IPSL-simulated spring ET is not 300 statistically significant (p > 0.05), which is confirmed by two observation-based ET 301 products (fig. S5A). This is mainly because the significant ET increase in regions with 302 spring greening is compensated by ET decrease in regions with spring browning, when 303 integrated to the hemispheric scale (fig. S7). For this greening-enhanced ET, however, 304 only a proportion (71%) is subsequently recycled as additional land precipitation during 305 MAM (0.20±0.09 mm yr-2, p < 0.01) (Fig. 2D). During MAM, this imbalance between 306 increasing precipitation and ET due to earlier phenology and LAI development results in 307 a decreasing trend in both soil moisture (–0.11±0.09×10-4 m3 m-3 yr-1) and surface runoff 308 (–0.08±0.09 mm yr-2). The larger increase in ET than in precipitation is much more 309 pronounced when considering only regions with spring greening, where it causes a strong 310 signal of soil drying (Fig. 2, B and E). However, in regions with spring browning, there is 311 a stronger decrease in ET than in precipitation, thus causing soil wetting (Fig. 2, C and F). 312 These contrasting responses, where ET changes do not balance precipitation changes, 313 suggest a geographical redistribution of the extra atmospheric water vapour from regions 314 with earlier greening, towards regions with browning, within the springtime. 315 Science Advances Manuscript Page 14 of 34 316 Our factorial model simulations allow further investigation of the soil moisture legacy 317 effects into the summer following a spring with altered phenology and hydrology. We 318 find that the spring soil moisture deficit caused by earlier greening is indeed carried over 319 into the summer months (Fig. 2D, blue bars). The inter-seasonal soil moisture memory 320 propagates this moisture deficit forward, and this also involves additional contribution of 321 soil moisture-atmosphere biophysical feedbacks. Given the short memory of days to 322 weeks of the atmosphere, the atmospheric moisture anomaly from earlier greening is 323 dissipated mainly in the springtime, through extra spring rainfall (Fig. 2D). The lower soil 324 moisture levels in late spring due to a prolonged period of stronger decrease and 325 insufficient rainfall re-supply, become less accessible for plant root uptake and soil 326 evaporation (27), particularly in climatologically water-stressed regimes. Hence, summer 327 ET decreases substantially across the entire study area (–0.10±0.05 mm yr-2), and 328 summer P (precipitation) decreases even more (–0.19±0.08 mm yr-2), leading to a 329 slightly strengthened decreasing trend in summer SWC (–0.11±0.09×10-4 m3 m-3 yr-1) 330 (Fig. 2D). This summer drying signal is larger and statistically significant in areas with 331 extra spring greening (Fig. 2E). The significantly decreasing spring LAI (p < 0.05) in 332 areas with spring browning do not, however, generate significant changes in any of the 333 hydrological fluxes in summer (i.e. ET, P, runoff and SWC; all p > 0.1) (blue bars in Fig. 334 2F). In these browning regions, soil moisture anomalies in spring have little impact on 335 subsequent summer water recycling (ET » 0) and hence summer rainfall feedbacks. It is 336 noteworthy that both spring greening and browning tend to decrease the amount of 337 summer precipitation (Fig. 2, E and F). This implies that the current spring LAI levels 338 may have evolved to be an optimum that maximizes the precipitation retained on 339 continents in the subsequent summer season. 340 Science Advances Manuscript Page 15 of 34 341 We further investigated the spatial pattern of the GCM-derived soil moisture trend in 342 SWC (Fig. 4; other hydrological fluxes in fig. S8). Generally, the modelled summer 343 SWC pattern and observed summer SWC anomalies associated with the first MCA 344 mode (Fig. 1, D and F) are highly similar, providing mechanistic support for our 345 observation-based findings. Besides, the summer SWC pattern closely resembles the 346 spring SWC pattern (Fig. 4B versus Fig. 4A), implying the dominant role of soil 347 moisture memory in carrying this greening-induced water deficit into the summer. Areas 348 with the highest rates of summer drying (> 0.1×10-3 m3 m-3 yr-1) occur mainly in Europe, 349 eastern Asia, western Asia and eastern North America (Fig. 4B), where pronounced 350 advances of vegetation greening occur during spring (fig. S1). Conversely, a summer 351 wetting trend of soil moisture dominates regions with spring browning, such as central 352 and northwestern North America (Fig. 4B versus fig. S1). One notable exception is again 353 the central Siberian regions where a strong modelled soil wetting trend occurs in both 354 spring and summer (> 0.05×10-3 m3 m-3 yr-1) (Fig. 4), despite the presence of a clear 355 spring greening tendency there that generally lowers summer soil moisture elsewhere 356 (Fig. 4 and fig. S1). The GCM simulations show that this trend of regional wetting in 357 Siberia is caused by different mechanisms in the two seasons. In spring, the extra water 358 vapour recycled into the atmosphere in upstream regions (Europe) with strong greening, 359 is carried to Siberia by westerly winds (fig. S9A). This imported extra precipitation 360 moistens local soils and offsets the soil moisture deficit originated from local spring 361 greening (Schematic Fig. 3). During the following summer, when the prevailing westerly 362 winds become weaker (fig. S9B), the Siberian region operates as a more closed system 363 with intensified local water recycling (fig. S8). Consequently, the imported extra soil 364 moisture in springtime helps to sustain local positive feedbacks on summer precipitation 365 and is ultimately retained in the Siberian soils throughout the boreal summer. This Science Advances Manuscript Page 16 of 34 366 additional summer wetting from the earlier greening of local and remote regions, in total, M AM 367 accounts for over 90% of the overall summer soil wetting trend derived from L A I obs , 368 showing its much more important role for Siberia than the ongoing climate change. 369 370 Implications for Near-surface Land Climate in Summer 371 Lower soil moisture can trigger the occurrence, or amplify the frequency and intensity of 372 heatwaves, by altering the surface energy balance (17, 28, 29). Our finding that earlier 373 greening pre-conditions summer soils to be drier thus implies a shift in the distribution of 374 summer surface air temperature towards even more hot extremes, adding to the warming- 375 induced increase in extreme occurrence. To verify this mechanism with our paired GCM 376 simulations, we calculate two extreme hot temperature indices (30): i) the summer 377 maximum of daily maximum temperature (TXx) and ii) the number of summer days with 378 maximum temperatures exceeding the local 90th percentile of a given base period, here as 379 1982–2011 (TX90p) (see Methods). We find that earlier spring greening increases both 380 TXx and TX90p (Fig. 5 and fig. S10) across the land regions where it depletes summer 381 soil moisture more strongly (Fig. 4B). Conversely, regions with additional summer soil 382 wetting from spring LAI changes exhibit fewer hot extremes in summer. On average, our 383 paired GCM analysis shows that the overall soil drying causes an additional 23% and 384 26% increase (increasing rates of 0.8 days decade-1 for summer TX90p and 0.07 °C 385 decade-1 for summer TXx) in TX90p and TXx, respectively, which exacerbates the 386 summer heatwaves caused directly by background global warming. This finding implies 387 the existence of a positive feedback: human-induced global warming causes earlier spring 388 greening and thus forces lower summer soil moisture. This, in turn, aggravates further the 389 higher frequency and intensity of hot summer extremes driven by rising greenhouse gas 390 concentrations (31). Science Advances Manuscript Page 17 of 34 391 392 Discussion 393 Our combined use of satellite observations and model simulations suggest that a more 394 intense water loss from land soils during spring, due to earlier vegetation greening, 395 lowers soil moisture in spring and summer. This moisture depletion further reduces 396 surface evaporative cooling and thus increases the possibility of stronger summer 397 heatwaves. Exacerbated soil water limitation, together with more extreme weather events 398 in summer, may restrain vegetation growth and suppress ecosystem productivity, 399 especially for climatologically arid ecosystems that are sensitive to summer water 400 availability (e.g. western U.S.). Our study has important implications for more accurate 401 predictions of summer weather and climate features. For many parts of the world, climate 402 models and measurements show that surface meteorological conditions are affected by 403 both large-scale atmospheric circulation patterns, as well as by local soil moisture status 404 in some regions with strong vegetation-atmosphere coupling (17, 32). The results of this 405 study highlight the necessity of also accounting for spring plant phenology in weather and 406 climate forecasting of summer months, given their capacity to modulate summer soil 407 moisture levels through the inter-seasonal biophysical feedbacks. 408 409 One limitation of our analysis is the examination of diagnostics from an ensemble based 410 on a single climate model. This concern is offset by our findings of modelled spring-to- 411 summer greenness-soil moisture linkages having consistency with results of parallel 412 observation-based analysis (Fig. 1). However, for future studies, we encourage other 413 modelling centres to perform similar factorial experiments to ours, to verify the lagged 414 hydrological impacts of earlier vegetation phenology. This includes independent model 415 evidence of our discovery that precipitation feedbacks can alleviate greening-induced soil Science Advances Manuscript Page 18 of 34 416 drying for some locations, involving the role of atmospheric teleconnections. Such 417 assessment must rely on high-resolution GCMs that have accurate representations of 418 large-scale atmospheric circulation processes. Furthermore, root-zone soil moisture 419 relevant for deep-rooted plants are not directly measurable at large spatial scales. Hence, 420 currently some uncertainties exist in soil moisture estimates from both observation-driven 421 products and GCMs, in terms of the magnitude, inter-annual variability and trends. Future 422 improvements of root-zone soil moisture may be fulfilled by integrating in situ soil 423 moisture measurements at different depths, microwave-based surface soil moisture and 424 robust surface-to-root-zone scaling approaches (e.g. ref. 33). This will enable additional 425 verification of our identified soil moisture responses to the ongoing earlier greening, and 426 also aid enhanced soil moisture parameterizations in the next generation of CGMs. Whilst 427 we are confident of the broad applicability of our findings, such improved simulations of 428 soil moisture will allow refinement of magnitudes of our discovered feedbacks. 429 430 Human land-use changes may have strong imprints on the seasonal linkage between 431 natural plants and soil moisture. The phenology of crops is strongly controlled by 432 management actions, including fertilizer use, cropping calendars, multi-cropping, quick- 433 growing hybrid cultivars, pest control and notably, irrigation. The latter can effectively 434 sustain soil moisture of irrigated croplands throughout the growing season, causing 435 different correlations between crop growth and soil moisture. However, at present, many 436 agricultural practices are not explicitly represented in GCMs, which limits our capacity to 437 simulate the hydrological consequences of cropland managements. As we illustrate how 438 vegetation greening feedbacks can be strongly modulated by human land-use, this 439 demonstrates the importance for the next generation of GCMs to include dynamic crop 440 growth and realistic irrigation schemes. This will enable more accurate simulation of Science Advances Manuscript Page 19 of 34 441 vegetation-soil moisture interactions with human interventions. In addition, the 442 hydrological feedbacks of afforestation, reforestation and forest management practices 443 over regions like China are not assessed separately from those of natural vegetation 444 greening. Deliberately planted trees are expected to have similar variations in LAI as 445 existing forests, causing similar alterations to the water cycle via changing transpiration 446 and regional to large-scale precipitation patterns. Nonetheless, the ability of plants to 447 transpire water varies with species composition and forest age (34), which should be 448 taken into account in future assessment of hydrological impacts of large-scale 449 afforestation/reforestation. 450 451 Materials and Methods 452 Experimental Design 453 We utilised satellite-derived measurements of LAI (as a proxy for vegetation 454 “greenness”) and reanalysis-based soil moisture to identify any significant relations 455 between spring vegetation phenology and summer soil moisture, using the Maximum 456 Covariance Analysis (MCA) approach. We also used factorial experiments by the IPSL- 457 CM Earth System Model to provide process understanding of observation-based 458 greenness-to-soil moisture linkages. 459 460 Observational Data Sets 461 Vegetation greenness was characterized using a gridded bi-weekly LAI data set derived 462 from the Global Inventory Monitoring and Modeling Studies (GIMMS) LAI3g at an 8- 463 km spatial resolution and 15-day interval for the period 1982–2011 (20). This global LAI 464 product was generated from the GIMMS Advanced Very High Resolution Radiometer 465 (AVHRR) Normalized Difference Vegetation Index (NDVI) 3g data and Moderate Science Advances Manuscript Page 20 of 34 466 Resolution Imaging Spectroradiometer (MODIS) LAI data using an artificial neural 467 network algorithm (20). The GIMMS LAI3g data have been demonstrated to have good 468 quality through direct comparisons with field measurements and inter-comparisons with 469 other satellite-derived LAI products (35). Monthly LAI maps were obtained by 470 calculating the maximum value of the two composites in each month at the original 8-km 471 resolution, and were then aggregated to each 0.5°×0.5° grid cell that was used for our 472 model simulation during the study period. 473 474 Besides LAI, other observed boundary conditions prescribed for IPSL-CM include SST 475 and SIC maps and atmospheric CO2 concentrations, all covering the study period 1982– 476 2011. Monthly SST and SIC data sets were both obtained from the Atmospheric Model 477 Intercomparison Project (AMIP; https://esgf-node.llnl.gov/news/list/input4mips/) with a 478 1°×1° spatial resolution. Global annual atmospheric CO2 concentrations were derived 479 from those used for transient simulations of the TRENDY project (“Trends and drivers of 480 the regional scale sources and sink of carbon dioxide”, http://dgvm.ceh.ac.uk/node/9). 481 482 Two independent SWC-related data sets covering the period of 1982–2011 were used for 483 detecting coupled patterns of variability with changes in spring LAI and for evaluating 484 modeling results. The first one is the satellite-derived zoot-zone soil moisture from the 485 Global Land Evaporation Amsterdam Model (21) (GLEAM v3.2a). The GLEAM data set 486 fully assimilates reanalysis-based temperature and radiation, and remotely-sensed top- 487 layer soil moisture, L-band vegetation optical depth (VOD) and snow water equivalents 488 from different satellite sensors (21), and also uses multi-source precipitation forcings. 489 The assimilation of surface soil moisture allows us to diagnose where irrigation in 490 agricultural regions is perturbing the hydrological cycle and the vegetation-soil moisture Science Advances Manuscript Page 21 of 34 491 feedbacks between spring and summer. Meanwhile, the forced precipitation and implied 492 differences in hydrological balance provide data-based indications of where atmospheric 493 teleconnections are transporting extra moisture. In addition, GLEAM incorporates VOD 494 to calculate the effects of water stress, heat stress and phenological constraints on ET 495 (21). In this way, it indirectly accounts for the influence of plant phenological changes on 496 soil moisture. VOD and LAI are intrinsically related, but this does not generate circularity 497 as we test for the link between GLEAM summer SWC and spring LAI measurements 498 from optical satellite data. This is because VOD is controlled not only by plant structure 499 (in part characterized by LAI), but additionally by plant hydraulic strategies (i.e. stomatal 500 and xylem regulation on water use) (36). The seasonal development of L-band VOD and 501 LAI is proven to differ considerably for some biomes, such as woody savannas of which 502 water is mainly stored in stems or branches (36). The GLEAM SWC data have been 503 extensively used in examining land-atmosphere feedbacks and benchmarking climate 504 model simulations (e.g. ref. 16). Note that the GLEAM root-zone soil moisture, of which 505 the maximum depth ranges between 10 and 250 cm depending on local vegetation types, 506 was here treated as a substitute for total soil moisture. 507 508 The other soil moisture-related data is GRACE-REC, a global reconstruction of terrestrial 509 water storage (TWS) anomalies based on the GRACE satellite retrievals (22). The past 510 TWS anomalies are reconstructed using a statistical model forced by daily precipitation 511 and temperature anomalies, but trained with original TWS data from GRACE satellites 512 for 2002-present (22). The TWS anomalies, expressed as the equivalent water thickness 513 (units: cm), include not only soil moisture in the upper layers, but also snowpack, canopy 514 water, riverine and lacustrine water and groundwater. Removing these non-soil moisture Science Advances Manuscript Page 22 of 34 515 components from the TWS anomalies is technically challenging, so some caution is 516 needed when interpreting the TWS-based results (Fig. 1, B, E and F). 517 518 We also obtained a range of observation-based data sets to evaluate model performance 519 for hydrological variables besides SWC (ET and P). No direct observations of global land 520 ET have been reported, so we used two satellite-derived reconstructions of long-term land 521 ET. The first, FLUXNET-MTE, is based on field eddy-covariance measurements (37), 522 and the second, WB-MTE, is based on information for precipitation, river discharge and 523 changes in water mass balance (38). Both ET products use satellite measurements to scale 524 up either site- or catchment-level ET observations to a 0.5°×0.5° global grid cell using a 525 machine learning algorithm. Gridded data for global monthly precipitation were obtained 526 from the Global Precipitation Climatology Centre (39) (GPCC) version 6.0 and the 527 Climate Research Unit (CRU) TS3.22 data sets (40). 528 529 Model and Experimental Setup 530 IPSL-CM (24) is a coupled land-atmosphere global climate model (GCM), versions of 531 which continuously contribute to the Coupled Model Intercomparison Project (CMIP) 532 initiatives, thus informing assessment reports of the Intergovernmental Panel on Climate 533 Change (IPCC). We employed the IPSL-CM version 4 GCM, which consists of the 534 atmospheric component LMDZ (the Laboratoire de Météorologie Dynamique general 535 circulation model with Zooming capability) and the land surface component ORCHIDEE 536 (Organizing Carbon and Hydrology In Dynamic Ecosystems (41)). The soil hydrology 537 within ORCHIDEE consists of two moisture layers with a constant total depth of 2 m 538 water content per unit soil volume. Vegetation in this land surface model is described 539 globally with 13 plant functional types (PFTs), each using different governing equations Science Advances Manuscript Page 23 of 34 540 and taking different parameter values to represent the onset and offset processes of the 541 leafy season (41). Using the IPSL-CM coupled model, changes in hydrological fluxes due 542 to vegetation biophysical feedbacks were determined by replacing the modeled LAI in the 543 ORCHIDEE land surface model with satellite observations prescribed for each of the 13 544 PFTs. 545 546 To investigate how spring vegetation feedbacks modulate summer soil moisture patterns 547 in the last 30 years, we conducted a pair of transient GCM simulations prescribed with M AM 548 different spring LAI conditions. One experiment (abbreviated as L A I obs ) was prescribed 549 with satellite-observed annually varying monthly LAI maps for 1982–2011 (including M AM 550 MAM). The other experiment (abbreviated as L A I c lim ) was also prescribed with satellite- M AM 551 observed LAI maps identically to the L A I obs experiment, except for spring months. 552 During MAM, the LAI maps for each year were replaced with average monthly values 553 over the entire study period (1982–2011). Spring greening in the real world may affect 554 soil moisture via a coupling to changes in oceanic and atmospheric conditions, hence both 555 experiments were fed additionally with other boundary conditions of sea-surface 556 temperature (SST), sea-ice fraction (SIF) and atmospheric CO2 concentration. With this 557 numerical experiment design, the fingerprint of earlier greening on summer hydrological 558 variables can be calculated as the difference between the two experiments (i.e. M AM L A I obs – M AM 559 L A I c lim ). The original or MAM-averaged LAI data-based time series (20) were used to 560 prescribe LAI for each GCM grid cell and for each PFT within a grid cell. To reduce the 561 effect of initial conditions uncertainty due to internal climate variability, we conducted an 562 ensemble of 30 simulations for each experiment with different initial conditions. The 563 initial conditions were derived from the model outputs of an undisturbed 60-year run (i.e. 564 with average LAI for all seasons, average SST, SIF and atmospheric CO2 concentrations): Science Advances Manuscript Page 24 of 34 565 the first 30 years of the run for the spin-up of the simulation and then each of the last 30 566 years of the run for providing an alternative initial condition. All our model-based 567 analyses were based on the 30-member ensemble mean of each transient experiment. 568 Only vegetated areas, defined as regions with multi-year (1982–2011) average LAI  0.1, 569 were considered to reduce the spurious influence of bare and sparsely vegetated areas on 570 the overall trend in LAI. All quantities averaged over the Northern Hemisphere (25– 571 90°N) are weighted by area, unless otherwise noted. 572 573 Indices of Extreme Hot Temperatures 574 We used two metrics to measure modelled changes in summer extreme hot temperatures 575 due to the earlier spring vegetation onset, following a previous study of extreme- 576 temperature detection (30). One metric measures the frequency of extreme hot 577 temperatures, expressed as the number of summer days with maximum daily 578 temperatures exceeding the 90% percentile of a baseline period (TX90p). The percentiles 579 were calculated for each calendar day over a 15-day moving window using the 30-year M AM 580 simulation forced by observed LAI and other boundary conditions (i.e. L A I obs ) (1982– 581 2011). The other metric measures the intensity of extreme hot temperature events, as 582 expressed as the summer maximum of daily maximum temperature (TXx). Both indices 583 were first calculated for both 30-ensemble mean experiments ( L A I oMbAs M and M AM L A I c lim ), and 584 changes in TX90p and TXx induced by earlier spring greening were then calculated as 585 the difference between the indices derived from the two simulations with and without 586 spring LAI changes. 587 588 Statistical Analysis Science Advances Manuscript Page 25 of 34 589 We applied partial correlation analysis in this study to evaluate the relationship between 590 spring LAI and summer soil moisture, after statistically controlling for the covarying 591 effects of summer LAI and climatic variables (including summer temperature, 592 precipitation and solar radiation). These climatic variables are all derived from the CRU 593 TS3.22 data sets. We also used the MCA analysis to characterize the coupling patterns 594 between spring LAI and summer SWC. This approach is also known as the singular value 595 decomposition (SVD) (23), which can be considered as a generalization of principle 596 component analysis. MCA has been widely used in environmental sciences to detect 597 teleconnection patterns in the ocean-atmosphere-land system (e.g. ref. 42). The MCA 598 analysis looks for pairs of spatial patterns and associated time series (or called expansion 599 coefficients) that efficiently explain the covariance between two space-time data fields. 600 This is realized by performing SVD on the cross-covariance matrix between the two 601 fields, yielding a set of singular vectors and singular values that are associated with pairs 602 of spatial patterns. The time series for each pattern pair (fig. S4) are calculated by 603 projecting one data field onto the respective singular vector. The correlation coefficient 604 between the time series of the two fields indicates how strongly related the coupled 605 patterns are. Heterogeneous regression maps (Fig. 1, C-F and fig. S6, A and B) are 606 generated by regressing the gridded values of one field to the time series of the pattern 607 amplitude of other. Detailed mathematical basis of the MCA algorithm is described in 608 section S2. 609 610 611 Supplementary Materials 612 Supplementary Text 613 Section. S1. Earlier spring vegetation greening revealed by remote sensing Science Advances Manuscript Page 26 of 34 614 Section. S2. Mathematical basis of the MCA analysis 615 Section. S3. Model evaluation using available observations 616 Figures S1-S10 617 Fig. S1. Spatial patterns of the trend in spring vegetation changes for 1982-2011. 618 Fig. S2. Grassland partial correlations between ESA-CCI summer SWC and GIMMS 619 spring LAI during 1982-2011. 620 Fig. S3. Partial correlations between spring LAI and either summer SWC or summer 621 TWS, as affected by agricultural extent. 622 Fig. S4. Time series (or expansion coefficients) of the leading first MCA patterns. 623 Fig. S5. Validation of IPSL-simulated hydrological variables. 624 Fig. S6. Robustness test of the IPSL-simulated relationship between spring LAI and 625 summer SWC. 626 Fig. S7. Comparison of trends of observation-based and IPSL-simulated ET. 627 Fig. S8. Spatial patterns of IPSL-simulated trends in key hydrological variables induced 628 by spring LAI changes. 629 Fig. S9. IPSL-simulated average global prevailing wind at 10-m height. 630 Fig. S10. Spatial patterns of IPSL-simulated trends in summer extreme hot temperature 631 indices induced by spring LAI changes. 632 633 634 References and Notes 635 1. J. Barichivich et al., Large-scale variations in the vegetation growing season and annual cycle 636 of atmospheric CO2 at high northern latitudes from 1950 to 2011. Glob. Chang. Biol. 19, 3167– 637 3183 (2013). 638 2. S. Piao et al., Leaf onset in the northern hemisphere triggered by daytime temperature. Nat. Science Advances Manuscript Page 27 of 34 639 Commun. 6, 6911 (2015). 640 3. Q. Liu et al., Delayed autumn phenology in the Northern Hemisphere is related to change in 641 both climate and spring phenology. Glob. Chang. Biol. 22, 3702–3711 (2016). 642 4. M. D. Schwartz, R. Ahas, A. Aasa, Onset of spring starting earlier across the Northern 643 Hemisphere. Glob. Chang. Biol. 12, 343–351 (2006). 644 5. A. Menzel et al., European phenological response to climate change matches the warming 645 pattern. Glob. Chang. Biol. 12, 1969–1976 (2006). 646 6. Y. H. Fu et al., Recent spring phenology shifts in western Central Europe based on multiscale 647 observations. Glob. Ecol. Biogeogr. 23, 1255–1263 (2014). 648 7. C. D. Keeling, J. F. S. Chin, T. P. Whorf, Increased activity of northern vegetation inferred from 649 atmospheric CO2 measurements. Nature. 382, 146 (1996). 650 8. S. Piao, P. Friedlingstein, P. Ciais, N. Viovy, J. Demarty, Growing season extension and its 651 impact on terrestrial carbon cycle in the Northern Hemisphere over the past 2 decades. Glob. 652 Biogeochem. Cycles. 21, 1–11 (2007). 653 9. S. Wolf et al., Warm spring reduced carbon cycle impact of the 2012 US summer drought. Proc. 654 Natl. Acad. Sci. U. S. A. 113, 5880–5885 (2016). 655 10. A. Angert et al., Drier summers cancel out the CO2 uptake enhancement induced by warmer 656 springs. Proc. Natl. Acad. Sci. U. S. A. 102, 10823–10827 (2005). 657 11. W. Buermann, P. R. Bikash, M. Jung, D. H. Burn, M. Reichstein, Earlier springs decrease peak 658 summer productivity in North American boreal forests. Environ. Res. Lett. 8, 24027 (2013). 659 12. W. Buermann et al., Widespread seasonal compensation effects of spring warming on northern 660 plant productivity. Nature. 562, 110–114 (2018). 661 13. S. Sippel et al., Contrasting and interacting changes in simulated spring and summer carbon 662 cycle extremes in European ecosystems. Environ. Res. Lett. 12, 75006 (2017). 663 14. R. J. Van Der Ent, H. H. G. Savenije, B. Schaefli, S. C. Steele-Dunne, Origin and fate of Science Advances Manuscript Page 28 of 34 664 atmospheric moisture over continents. Water Resour. Res. 46, W09525 (2010). 665 15. Z. Zeng et al., Impact of Earth greening on the terrestrial water cycle. J. Clim. 31, 2633–2650 666 (2018). 667 16. Y. Li et al., Divergent hydrological response to large-scale afforestation and vegetation 668 greening in China. Sci. Adv. 4, eaar4182 (2018). 669 17. B. Quesada, R. Vautard, P. Yiou, M. Hirschi, S. I. Seneviratne, Asymmetric European summer 670 heat predictability from wet and dry southern winters and springs. Nat. Clim. Chang. 2, 736– 671 741 (2012). 672 18. A. Berg et al., Land-atmosphere feedbacks amplify aridity increase over land under global 673 warming. Nat. Clim. Chang. 6, 869–874 (2016). 674 19. A. L. S. Swann, I. Y. Fung, J. C. H. Chiang, Mid-latitude afforestation shifts general circulation 675 and tropical precipitation. Proc. Natl. Acad. Sci. U. S. A. 109, 712–716 (2012). 676 20. Z. Zhu et al., Global data sets of vegetation leaf area index (LAI)3g and Fraction of 677 Photosynthetically Active Radiation (FPAR)3g derived from Global Inventory Modeling and 678 Mapping Studies (GIMMS) Normalized Difference Vegetation Index (NDVI3g) for the period 679 1981 to 2011. Remote Sens. 5, 927–948 (2013). 680 21. B. Martens et al., GLEAM v3: Satellite-based land evaporation and root-zone soil moisture. 681 Geosci. Model Dev. 10, 1903–1925 (2017). 682 22. V. Humphrey, L. Gudmundsson, S. I. Seneviratne, A global reconstruction of climate-driven 683 subdecadal water storage variability. Geophys. Res. Lett. 44, 2300–2309 (2017). 684 23. C. S. Bretherton, C. Smith, J. M. Wallace, An intercomparison of methods for finding coupled 685 patterns in climate data. J. Clim. 5, 541–560 (1992). 686 24. O. Marti et al., Key features of the IPSL ocean atmosphere model and its sensitivity to 687 atmospheric resolution. Clim. Dyn. 34, 1–26 (2010). 688 25. A. Berg, J. Sheffield, P. C. D. Milly, Divergent surface and total soil moisture projections under Science Advances Manuscript Page 29 of 34 689 global warming. Geophys. Res. Lett. 44, 236–244 (2017). 690 26. X. Lian et al., Partitioning global land evapotranspiration using CMIP5 models constrained by 691 observations. Nat. Clim. Chang. 8, 640–646 (2018). 692 27. M. Stéfanon, P. Drobinski, F. D’Andrea, N. De Noblet-Ducoudré, Effects of interactive 693 vegetation phenology on the 2003 summer heat waves. J. Geophys. Res. Atmos. 117, D24103 694 (2012). 695 28. M. Hirschi et al., Observational evidence for soil-moisture impact on hot extremes in 696 southeastern Europe. Nat. Geosci. 4, 17–21 (2011). 697 29. M. M. Vogel et al., Regional amplification of projected changes in extreme temperatures 698 strongly controlled by soil moisture‐temperature feedbacks. Geophys. Res. Lett. 44, 1511–1519 699 (2017). 700 30. S. E. Perkins, L. V. Alexander, On the measurement of heat waves. J. Clim. 26, 4500–4517 701 (2013). 702 31. N. Christidis, G. S. Jones, P. A. Stott, Dramatically increasing chance of extremely hot 703 summers since the 2003 European heatwave. Nat. Clim. Chang. 5, 46–50 (2014). 704 32. D. L. Swain et al., Trends in atmospheric patterns conducive to seasonal precipitation and 705 temperature extremes in California seasonal precipitation and temperature extremes. Sci. Adv. 706 2, e1501344 (2016). 707 33. Albergel, C. et al., From near-surface to root-zone soil moisture using an exponential filter: an 708 assessment of the method based on in-situ observations and model simulations. Hydrol. Earth 709 Syst. Sci. 12, 1323–1337 (2008). 710 34. S. T. Brantley, C. F. Miniat, P. V. Bolstad, Rainfall partitioning varies across a forest age 711 chronosequence in the southern Appalachian Mountains. Ecohydrology, 12, e2081 (2019). 712 35. H. E. Beck et al., Global evaluation of four AVHRR–NDVI data sets: Intercomparison and 713 assessment against Landsat imagery. Remote Sens. Environ. 115, 2547–2563 (2011). Science Advances Manuscript Page 30 of 34 714 36. F. Tian et al., Coupling of ecosystem-scale plant water storage and leaf phenology observed by 715 satellite. Nat. Ecol. Evol. 2, 1428–1435 (2018). 716 37. M. Jung et al., Recent decline in the global land evapotranspiration trend due to limited 717 moisture supply. Nature. 467, 951–954 (2010). 718 38. Z. Zeng et al., A worldwide analysis of spatiotemporal changes in water balance-based 719 evapotranspiration from 1982 to 2009. J. Geophys. Res. Atmos. 119, 1186–1202 (2014). 720 39. A. Becker et al., A description of the global land-surface precipitation data products of the 721 Global Precipitation Climatology Centre with sample applications including centennial (trend) 722 analysis from 1901–present. Earth Syst. Sci. Data. 5, 71–99 (2013). 723 40. I. Harris, P. D. Jones, T. J. Osborn, D. H. Lister, Updated high-resolution grids of monthly 724 climatic observations – the CRU TS3.10 Dataset. Int. J. Climatol. 34, 623–642 (2014). 725 41. G. Krinner et al., A dynamic global vegetation model for studies of the coupled atmosphere- 726 biosphere system. Glob. Biogeochem. Cycles. 19, 1–33 (2005). 727 42. A. Dai, Increasing drought under global warming in observations and models. Nat. Clim. 728 Chang. 3, 52–58 (2013). 729 43. S. Piao et al., Plant phenology and global climate change: Current progresses and challenges. 730 Glob. Chang. Biol. 25, 1922–1940 (2019). 731 732 733 Acknowledgments 734 General: We acknowledge the National Supercomputer Center in Tianjin, China (NSCC- 735 JI, http://www.nscc-tj.cn/), and the National Computer Center IDRIS of CNRS in France 736 for providing computing resources. Funding: This study was supported by the National 737 Key Research and Development Program of China (Grant-2018YFC1507704), the 738 Second Tibetan Plateau Scientific Expedition and Research Program (Grant No. 739 2019QZKK0208), National Natural Science Foundation of China (41561134016), Science Advances Manuscript Page 31 of 34 740 National Youth Top-notch Talent Support Program in China and the 111 Project 741 (B14001). C. H. acknowledges the Center for Ecology and Hydrology (CEH) National 742 Capability fund. P. C., I. J. and J. P. were funded by the European Research Council 743 Synergy grant SyG-2013-610028 IMBALANCE-P. Author contributions: S. P. 744 designed the research; X. L. performed analysis; X. L. and S. P. drafted the paper; and all 745 authors contributed to the interpretation of the results and to the text. Competing 746 interests: The authors declare that they have no competing interests. Data and materials 747 availability: All methods needed to evaluate the conclusions in the paper are present in 748 the paper and/or the Supplementary Materials. All data related to this paper, along with 749 code for interpretation of data and preparation of diagrams is available on request from 750 the corresponding author. 751 752 753 Figure captions 754 755 Fig. 1. Coupling between observed spring LAI and summer SWC during 1982-2011. 756 Spatial pattern of partial correlation coefficients between GIMMS spring LAI and 757 (A) GLEAM summer SWC, or (B) GRACE-REC summer TWS for 1982-2011. 758 Black stipples indicate regions with a statistically significant correlation (p < 759 0.05). Heterogeneous regression maps of (C, E) GIMMS spring LAI and (D, F) 760 GLEAM summer SWC (or GRACE-REC summer TWS), associated with the first 761 MCA mode for 1982-2011. The squared fractional covariance (SFC) explained by 762 the first MCA mode is 46.3% and 58.9% for the GLEAM and GRACE-REC data 763 sets, respectively. 764 Science Advances Manuscript Page 32 of 34 765 Fig. 2. IPSL-simulated changes in summer soil moisture induced by spring LAI 766 changes. (A–C) Inter-annual anomalies of the area-weighted average of spring 767 LAI (green dotted lines) and summer SWC (blue solid lines) for (A) all northern 768 latitudes (25–90°N), (B) regions with positive LAI trends (i.e. greening) and (C) 769 regions with negative LAI trends (i.e. browning). The red lines indicate the least- 770 squares linear regression of GCM-based SWC (straight lines) against time and 771 the 95% confidence intervals (curves). Note that the right axes are reversed, so 772 higher LAI values are towards the bottom of the plots. (D–F) Inter-annual trends 773 in mean spring LAI and resultant changes in the spring and summer hydrological 774 variables of ET, precipitation (P), runoff (Q) and SWC. The subplot maps at the 775 bottom of each main panel display the corresponding averaged areas as gray (the 776 same as panels A–C). 777 778 Fig. 3. Schematic of the effect of earlier greening on summer soil moisture. Earlier 779 spring greening influences spring soil moisture by altering land-atmosphere water 780 exchanges (via ET, P and Q) and by the redistribution of atmospheric water 781 vapour by atmospheric circulation. This spring soil moisture anomaly persists 782 later into the following summer due to the carry-over effects of soil moisture. The 783 magnitude of this cross-seasonal vegetation feedback and the role of atmospheric 784 circulation, however, vary geographically. Three typical examples of the 785 circulation-modulated vegetation feedback (Europe, Siberia and eastern China) 786 are displayed at the bottom of the schematic. 787 788 Fig. 4. Spatial patterns of IPSL-simulated trends in soil moisture induced by spring 789 LAI changes. Shown are spatial patterns of the linear trends in spring (left) and Science Advances Manuscript Page 33 of 34 790 summer (right) SWC changes. This change in soil moisture induced by earlier 791 spring greening, SWC, is obtained as the difference between the two simulations 792 with and without spring LAI changes. (i.e. M AM L A I obs – L A I cMlimA M ). Black stipples 793 indicate regions with a statistically significant linear trend (p < 0.05). 794 795 Fig. 5. IPSL-simulated changes in summer extreme hot temperature indices related 796 to SWC changes. The average linear trend of (A) TXx and (B) TX90p due to 797 earlier greening for different summer SWC trends binned into 0.4 (  10-4 m3 m-3 798 yr-1) intervals. The TXx and TX90p values are obtained as the difference 799 between the two simulations with and without spring LAI changes ( L A I oMbAs M – 800 M AM L A I c lim ). The insets are scatterplots of trends in TXx (left) and TX90p (right) 801 versus trends in SWC for 1982–2011. Each colored dot in the scatterplots 802 represents a grid cell in the pattern of trends in TXx or TX90p (fig. S10). Science Advances Manuscript Page 34 of 34