Data<\/p>\n
TWS data come from the GRACE and GRACE-FO mass concentration solutions (RL06) developed by the Center for Space Research at the University of Texas25<\/a>. The dataset represents monthly total water mass variations (expressed as equivalent water height) at a gridded resolution of 0.25\u00b0. We omit data before August 2002 (due to sensor calibration at the outset of the GRACE mission) and for water-year 2017 (due to incomplete data during the mission gap between GRACE and GRACE-FO). We use GRACE TWS because it is an observational data product that holistically measures land water, capturing all surface and underground stocks.<\/p>\n To account for observational uncertainty, we use three daily gridded precipitation data products: the Global Precipitation Climatology Project (GPCP) Daily Precipitation Analysis v.1.3 (ref.\u200921<\/a>), the Global Precipitation Climatology Center (GPCC) Daily Analysis v.2022 (ref.\u200923<\/a>) and the National Oceanic and Atmospheric Administration Climate Prediction Center (CPC) Unified Gauge-Based Analysis22<\/a>. GPCC is a station-based product (1\u00b0 resolution) with coverage over 1982\u20132020, whereas CPC (0.5\u00b0) and GPCP (1\u00b0) integrate station and satellite observations over 1979\u20132022 and 1997\u20132022, respectively. We repeated our analysis for the more recent GPCP v.3.3 and find consistent main effects56<\/a> (Supplementary Fig. 4<\/a>), but use the more complete earlier version for our main analysis (GPCP v.3.3 contains missing daily values at high latitudes after 2020).<\/p>\n We use two net all-sky surface shortwave radiation datasets: the Global Energy and Water Exchanges Surface Radiation Budget (GEWEX-SRB, limited here to 2002\u20132017)52<\/a> daily dataset and the NASA Clouds and the Earth\u2019s Radiant Energy System Energy Balanced and Filled monthly dataset (NASA-EBAF, limited here to 2002\u20132022)53<\/a>, each at 1\u00b0 resolution. We use monthly evapotranspiration from the GLEAM v.3.8a, at 0.25\u00b0 resolution for the overlapping GRACE period (2002\u20132020; ref.\u200954<\/a>). Main river basin boundaries are from the HydroBASINS v.1.0 dataset, which is derived from Shuttle Radar Topography Mission elevation data collected in 2000 (ref.\u200951<\/a>). All data URLs are enumerated in Supplementary Table 1<\/a>.<\/p>\n Daily climate model data are from the land-hist experiment of the Coupled Model Intercomparison Project, Phase 6, Land-Surface, Snow and Soil Moisture Intercomparison Project (CMIP6-LS3MIP)55<\/a>. In the land-hist experiment, the land-surface components of climate models are forced by atmospheric forcings from the Global Soil Wetness Project phase 3, whose precipitation data are derived from the NOAA-CIRES-DOE Twentieth Century Reanalysis (20CR). We include results from the two models that report all required variables: MIROC6 with land component MATSIRO and VISIT-e57<\/a>, and CNRM-CM6-1 with land component ISBA-CTRIP58<\/a>. We use the mrtws (total water storage), tas (mean surface air temperature) and pr (precipitation rate) variables. No participating models report these required variables plus evaporation, precluding extending our climate model analysis towards mechanisms.<\/p>\n Data aggregation and precipitation Gini coefficient<\/p>\n All data are temporally aggregated to the water-year scale, using October to September for the Northern Hemisphere water-year and July to June in the Southern Hemisphere. On this water-year scale, snow and ice storage changes are negligible. We further limit our spatial domain to regions without permanent snow and ice accumulation. We compute daily mean temperature data using the average of daily maximum and minimum temperatures from CPC. We also compute climatological mean precipitation as the average of the annual total over the full precipitation data record.<\/p>\n We linearly detrend GRACE data to isolate interannual variability, which is largely governed by hydroclimate variability, from long-term trends that include anthropogenic water use and management26<\/a>,31<\/a>. The rationale for this focus on interannual variability is to establish whether more concentrated daily precipitation distributions lead to cumulative water balance changes, averaged across the water-year, absent the influence of coincident economic, demographic and mean warming trends.<\/p>\n Correspondingly, and because GRACE only estimates TWS relative anomalies, we detrend and demean all other hydroclimate variables (except climatological precipitation, which is time-invariant). Finally, we interpolate all datasets to a common 0.5\u00b0 grid as a compromise across the resolutions of the datasets. We repeated our analyses at alternative coarser common resolutions (1\u00b0, 2\u00b0 and 3\u00b0) and found consistent results (Supplementary Fig. 5<\/a>). We present long-term linear trends in terms of changes per 20\u2009years to ease interpretation.<\/p>\n To quantify the annual concentration of daily precipitation, we apply the Gini coefficient, a common economic measure of income inequality59<\/a>, to the three daily precipitation datasets. The Gini coefficient is based on the Lorenz curve, which plots the cumulative share of the total quantity (that is, cumulative share of annual precipitation) against the cumulative share of the population (that is, days per year), sorted from least to greatest (Fig. 1d<\/a>). The Gini coefficient is defined as the area between the Lorenz curve and the 1:1 line of equality. For more concentrated distributions, such as all of the annual precipitation falling in 1\u2009day, the Lorenz curve deviates further from the line of equality, giving a Gini coefficient closer to 1. If all days are equally rainy, the Gini yields a value of 0. We use a computationally efficient discrete formulation to calculate the annual Gini coefficient of sorted daily precipitation (GP):<\/p>\n $${G}_{{\\rm{P}}}=\\frac{1}{n}\\,\\left(n+1-2\\frac{{\\Sigma }_{i}^{k}{\\Sigma }_{i}^{n}{P}_{i}}{{\\Sigma }_{i}^{n}{P}_{i}}\\right)$$<\/p>\n \n (1)\n <\/p>\n where n is the number of days per year and P is a sorted vector of precipitation values for a given water-year and location. The cumulative sum of P approximates the Lorenz curve (inner numerator double summation), and the normalized integral of that curve is estimated by summing the cumulative values (outer numerator double summation) and dividing by the total sum (denominator summation). This formulation is analytically identical to more common Gini formulations, but is computationally faster. As with the other variables, we linearly detrend and demean GP to isolate local interannual anomalies.<\/p>\n GP is\u2014by formulation\u2014sensitive to both incremental dry days and precipitation extremes. These joint sensitivities are mathematically coupled: in normalized terms, incremental dry days necessitate a greater concentration of annual precipitation into the remaining wet days, whereas incremental precipitation extremes leave less annual precipitation to be distributed across the remaining days of the year. This coupling maps conceptually with expectations from the physics of hydrologic intensification, in which daily precipitation responds more to warming than its annual total. We confirm this joint sensitivity of GP by regressing its annual anomalies on annual occurrence of extreme precipitation days (>99th percentile of nonzero precipitation days), dry days (<0.1\u2009mm), and total annual precipitation as well as on annual mean temperature.<\/p>\n Statistical models of TWS effect<\/p>\n For each precipitation dataset, the data aggregation and GP calculations yield a panel of annual anomalies in TWS, mean temperature (T), total precipitation (P) and GP, plus long-term (absolute) climatological mean precipitation (\\(\\bar{P}\\)). The three panels vary in sample size because of differing time availability of the precipitation datasets (Supplementary\u00a0Table 2<\/a>).<\/p>\n To assess the influence of GP on TWS for each panel, we apply linear fixed-effects models, a type of panel regression that enables controlling for time-invariant unobserved spatial heterogeneity using spatial fixed effects, as well as global shocks in any year using time fixed effects60<\/a>. Our approach minimizes confounding by spatial or temporal trends in environmental, climate or anthropogenic factors that may covary with GP, without explicit inclusion in the model. Spatial fixed effects remove all static differences across grid cells, such as baseline differences in soils, hydrogeology, vegetation and groundwater use or irrigation infrastructure. Water-year fixed effects absorb global shocks that affect all regions simultaneously (notably, El Ni\u00f1o events). By detrending all variables, we also exclude long-term groundwater depletion trends, ensuring that the model captures only interannual variability.<\/p>\n This identification strategy isolates the emergent interannual influence of precipitation concentration on TWS. The models are of the form<\/p>\n $$\\,{\\mathrm{TWS}}^{i,t}=\\theta {T}^{i,t}+\\pi {P}^{i,t}+\\gamma {G}_{{\\rm{P}}}^{i,t}+\\chi {G}_{{\\rm{P}}}^{i,t}{\\bar{P}}^{i}+{\\sigma }^{i}+{\\tau }^{t}+{\\varepsilon }^{i,t}$$<\/p>\n \n (2)\n <\/p>\n where \u03b8 and \u03c0 are coefficients of temperature and annual precipitation, \u03c3 and \u03c4 are spatial and time fixed effects, superscripts i and t denote grid cell and water-year, and \u03b5 is the residual. These terms are included as controls for the influence of other hydroclimate or unobserved variables. We estimate the effect of GP as an interaction with climatological mean precipitation (using linear coefficient \u03b3 and interaction coefficient \u03c7) to explicitly assess variation in the effect across the global aridity gradient. As such, the marginal effect of an anomaly \u0394GP on TWS is conditional on climatological precipitation and is given by<\/p>\n $$\\Delta \\mathrm{TWS}=(\\gamma +\\chi \\bar{P})\\Delta {G}_{{\\rm{P}}}$$<\/p>\n \n (3)\n <\/p>\n We interpret this marginal effect as the generalizable response of local interannual TWS variation to local anomalies of precipitation concentration, conditioned on local climatological wetness. We present the marginal effect of GP on TWS per a 10 percentage-point change in GP (that is, a GP anomaly of 0.1) as GP values are bounded between 0 and 1, and vary interannually by a smaller portion of that range. For each precipitation dataset, we estimate the coefficients in equation (2<\/a>) separately using absolute anomalies (to quantify the effect in hydrologically meaningful units) and standardized anomalies (to enable comparison of effect sizes across predictors). Coefficients are estimated using the lfe package in R. Global Pearson correlation among predictors (pooled across space and time) is under 0.5 (r2 under 0.25), indicating little scope for influence of collinearity on regression estimates (Extended Data Fig. 10<\/a>). To ease visualization of standardized coefficients, we present the conditional GP effect at the global population-weighted climatological mean precipitation value of 1,000\u2009mm.<\/p>\n As the marginal effect of GP on TWS is conditioned on climatological precipitation, its standard error depends on the variance of coefficients \u03c7 and \u03b3 as well as the climatological precipitation. We estimate this standard error (SE)\u00a0using the variance of linear combinations of estimators:<\/p>\n $$\\mathrm{SE}(\\gamma +\\chi \\bar{P})=\\sqrt{\\mathrm{var}(\\gamma )+{\\bar{P}}^{2}\\mathrm{var}(\\chi )+2\\bar{P}\\cdot \\mathrm{cov}(\\gamma ,\\chi )}$$<\/p>\n \n (4)\n <\/p>\n When estimating the variances and covariance of \u03c7 and \u03b3, we cluster standard errors at the main river basin scale to account for spatial and temporal autocorrelation of gridded data51<\/a>\u2014this step yields more conservative standard errors and avoids type-I errors due to the non-independence of observations. We also tested Conley spatial standard errors61<\/a> with distance cutoffs of 300\u2009km (roughly the effective resolution of GRACE) and 750\u2009km (the average length scale of global river basins in HydroSHEDS; Supplementary Fig. 6a\u2013c<\/a>), and Driscoll\u2013Kraay temporal standard errors with time cutoffs of 3\u2009years and 5\u2009years (Supplementary Fig. 6d\u2013f<\/a>). All approaches yield similar confidence intervals. At a given value of climatological precipitation, we consider the marginal effect of GP on TWS as statistically significant, for which its median estimate plus or minus two standard errors excludes zero.<\/p>\n We also estimate a reduced form of equation (2<\/a>) at the main river basin scale to check the robustness of the global effect at finer spatial scales51<\/a>, and to characterize its spatial variability. We omit the interaction between GP and climatological precipitation, as variation in climatological precipitation within basins is smaller than across the globe. For simplicity, we present these basin-scale TWS effects as the average across the three products, wherever they are significant (P\u2009<\u20090.05) and of consistent sign in at least two of the three precipitation data products. Basin-scale TWS effects are fairly consistent across the three data products (Extended Data Fig. 5<\/a>).<\/p>\n Although spatial fixed effects control for time-invariant spatial differences, they do not control for confounders that vary interannually in response to GP itself, such as irrigation responses to year-to-year variation in TWS change, which might deplete groundwater. To isolate the influence of irrigation, we add an interaction term, I, for the grid cell area equipped for irrigation to equation (2<\/a>) (using data from the UN Food and Agriculture AQUASTAT Global Information System on Water and Agriculture):<\/p>\n $${{\\rm{T}}{\\rm{W}}{\\rm{S}}}^{i,t}=\\theta {T}^{i,t}+\\pi {P}^{i,t}+\\gamma {G}_{{\\rm{P}}}^{i,t}+\\chi {G}_{{\\rm{P}}}^{i,t}{\\bar{P}}^{i}+\\varphi {G}_{{\\rm{P}}}^{i,t}{I}^{i}+\\delta {I}^{i}{G}_{{\\rm{P}}}^{i,t}{\\bar{P}}^{i}+{\\sigma }^{i}+{\\tau }^{t}+{\\varepsilon }^{i,t}$$<\/p>\n \n (5)\n <\/p>\n In this model, the influence of GP on TWS depends on both climatological mean precipitation and irrigation area, which we evaluate at varying levels. The approach enables our isolation of irrigation influence, despite 95% of grid cells in our analysis containing <1% irrigated area.<\/p>\n For the CMIP6-LS3MIP land-hist climate model data, we apply the identical data aggregation and statistical modelling in equations (2<\/a>\u20134<\/a>) to the simulated TWS, T, P and GP. To test the robustness of the estimated TWS effect to temporal and spatial variability in gauge coverage, we regress the basin-level TWS effects on the on basin mean areal gauge density (gauges per basin) and on the coefficient of variation of the number of gauges per basin over time (standard deviation of gauges per basin divided by the average), based on the GPCC gauge count data. The relationships show no significant dependence of the basin-level TWS effect on areal station density or the degree of year-to-year variation in gauge density (Extended Data Fig. 6b,c<\/a>). Empirical results using only GRACE (through water 2016) and excluding GRACE-FO (2018\u20132022) are consistent with those using the complete GRACE record62<\/a> (Supplementary Table 3<\/a>).<\/p>\n Statistical models of effect mechanism<\/p>\n To clarify the mechanism of the effect of GP on TWS, we develop two companion panel regressions similar to those used to assess the TWS effect, but using GLEAM evapotranspiration (E) and shortwave radiation (S) anomalies from GEWEX-SRB and NASA EBAF in place of TWS in equation (2<\/a>). To isolate the intensity-partitioning mechanism of the TWS effect, we add a regressor in the form of shortwave radiation anomalies in equation (2<\/a>) to explicitly control for the radiative effect of additional dry days:<\/p>\n $${\\mathrm{TWS}}^{i,t}={\\theta }^{{\\prime} }{T}^{i,t}+{\\pi }^{{\\prime} }{P}^{i,t}+{\\gamma }^{{\\prime} }{G}_{{\\rm{P}}}^{i,t}+{\\chi }^{{\\prime} }{G}_{{\\rm{P}}}^{i,t}{\\bar{P}}^{i}+\\zeta {S}^{i,t}+\\psi {S}^{i,t}{\\bar{P}}^{i}+{\\sigma }^{{\\prime} i}+{\\tau }^{{\\prime} t}+{\\varepsilon }^{{\\prime} i,t}$$<\/p>\n \n (6)\n <\/p>\n We model the effect of S on TWS, including an interaction with climatological precipitation to account for variable effects across gradients of mean aridity (and thus energy compared with moisture limitation). As such, \u03b6 is the TWS sensitivity to shortwave radiation anomalies and \u03c8 is the interaction effect with mean precipitation. The prime symbols indicate coefficients estimated with the S control, to distinguish from coefficients in equation (2<\/a>).<\/p>\n Introducing the S control removes the TWS variance attributable to shortwave anomalies, which are positively correlated with GP, from the estimated GP effect. We consider this positive S\u2013GP correlation as a facet of the hydrologic effects of precipitation concentration, so we omit this control when estimating the full TWS effect. However, we include the S control in equation (6<\/a>) to isolate the land-surface hydrologic partitioning effects of daily precipitation distribution on TWS, apart from the fact that more concentrated precipitation years are sunnier. We interpret this isolated effect as the influence of more intense precipitation on surface partitioning, such as infiltration-excess or saturation-excess ponding of precipitation on the surface.<\/p>\n For simplicity, we present this partitioning effect as an average of effects across the three precipitation datasets, with uncertainty characterized as the maximum and minimum of the full two standard error interval across the three precipitation and two radiation data products. This approach is conservative, as standard errors would be smaller if we first averaged the precipitation data, rather than averaging effects. However, results for individual datasets are consistent with those from the three-product average (Supplementary Fig. 1<\/a>). To estimate the relative fractional contribution of the partitioning effect to the full effect (fD), we divide the partitioning effect by the marginal effect estimated using equations (2<\/a>) and (3<\/a>):<\/p>\n $${f}_{{\\rm{D}}}=\\frac{{\\gamma }^{{\\prime} }+{\\chi }^{{\\prime} }\\bar{P}}{\\gamma +\\chi \\bar{P}}$$<\/p>\n \n (7)\n <\/p>\n Because it is simpler to visualize, we present this fractional distribution effect both for the average across data products and for each individual product (based on the data in Supplementary Fig. 1<\/a>).<\/p>\n Reanalysis models contain appreciable biases in the simulation of daily precipitation, especially for extreme values that contribute disproportionately to annual precipitation and strongly influence concentration63<\/a>. As such, we favour an empirical approach based on observational hydroclimate and TWS data. But understanding the mechanism of the observed effect requires assessing surface water budget terms that are poorly observed globally. We use GLEAM evapotranspiration data to establish the plausibility of evapotranspiration responses to GP. But as GLEAM assimilates precipitation data (and thus biases) from reanalysis, it probably underrepresents daily precipitation concentration effects. We similarly use the FLUXCOM CERES-GPCP flux-tower based energy and water flux dataset (2001\u20132014) to provide an in situ window on evaporative responses to precipitation concentration64<\/a>; however, as its observational constraint is limited to the approximately 100\u2009m radius surrounding the tower, these data are not comparable to GLEAM, estimated on a 0.25\u00b0 grid. To overcome these observational limitations on water-budget terms (evaporation and runoff), we extend a simple, process-based land\u2013atmosphere model with hydrology65<\/a>, driven by the daily precipitation datasets at a set of grid points sampling the global precipitation gradient.<\/p>\n Idealized land-surface model<\/p>\n We use an idealized process-based land-surface model to provide additional evidence for the negative TWS effect of GP and to clarify its mechanism (the underlying surface water balance changes and their drivers). Although more complex land-surface contain more complete process representations, the ultimate energy\u2013moisture balance rests on nested and coupled parameterizations, making a mechanistic diagnosis difficult, if not impossible.<\/p>\n Embracing recent community calls for robust and simple theories of land surface climate11<\/a>, we run land-surface model simulations that parsimoniously capture variations in surface precipitation partitioning and land\u2013atmosphere coupling. We use experimental design that enables direct comparison with the observed TWS effect. First, as the observed TWS effect is due to GP anomalies and because GRACE TWS measures anomalies, we focus on simulating the sensitivity of TWS to changes in GP, rather than mean TWS itself. Second, our statistical model includes space fixed effects to capture time-invariant differences in land cover and soil characteristics. Our idealized model mirrors this approach by using median parameters globally, rather than calibrating parameters locally. In short, rather than developing and validating a complete and realistic model of all of the processes driving TWS, we develop the minimal physics required to plausibly explain our observational results and assess the most important underlying processes.<\/p>\n Our process-based model parsimoniously represents the coupled surface partitioning of precipitation and energy\u2013moisture balance. It extends the one-dimensional simple energy\u2013moisture balance (SEMB) model, developed in ref.\u200965<\/a>, to represent precipitation partitioning at the land surface and simulate a surface water stock (Extended Data Figs. 9<\/a> and 10<\/a>). In the original model formulation, precipitation infiltrated instantaneously into the soil moisture stock (MS), irrespective of its intensity. On saturation, excess runoff was treated as an instantaneous water flux outside the model domain, precluding evaporation from surface water within the model. Our extended model (SEMB with hydrology, SEMB-H) includes intensity-dependent surface partitioning of precipitation towards infiltration compared with infiltration and saturation excess, the latter two being retained within the model using a surface water stock (ML), whose outputs include evaporation and runoff. Physically, these developments allow the model to more completely represent surface precipitation partitioning processes: when precipitation meets the land, water can pool at the surface if the infiltration capacity of the soil or maximum saturation have been reached, and this ponded surface water can run off or evaporate.<\/p>\n The surface energy balance is given in terms of time-evolving surface air temperature (T) as<\/p>\n $$C\\frac{{\\rm{d}}T}{{\\rm{d}}t}={F}_{\\mathrm{SW}}-D-L{E}_{{\\rm{S}}}={F}_{\\mathrm{SW}}-\\alpha (T-\\bar{{T}_{{\\rm{D}}}})-\\frac{L{\\rho }_{{\\rm{a}}}{M}_{{\\rm{S}}}}{{r}_{{\\rm{S}}}}[{q}_{\\mathrm{sat}}(T)-\\bar{q}]$$<\/p>\n \n (8)\n <\/p>\n where C is the effective heat capacity of the soil layer, FSW is the net surface shortwave radiation, D is the sum of net longwave radiation and the ground and sensible heat fluxes (that is, the \u2018dry\u2019 energetic response), ES is evaporation from soils and L is the enthalpy of vapourization of water. Following ref.\u200965<\/a>, D is parameterized as a function of the dewpoint depression (\\(T-\\bar{{T}_{{\\rm{D}}}}\\)) using a dry energetic parameter \u03b1, representing the temperature-damping efficiency of non-evaporative fluxes. We derive the value of \u03b1 from ref.\u200965<\/a>, estimated by regressing the sum of longwave, sensible heat and ground heat fluxes on the surface dewpoint depression temperature in ERA-Interim reanalysis. Evaporation is parameterized as a function of vapour pressure deficit (\\({q}_{\\mathrm{sat}}(T)-\\bar{q}\\)) by rS, which is the bulk soil surface resistance (\u03c1a is the mean density of air). Atmospheric humidity variables (\\(\\bar{q}\\) and \\(\\bar{{T}_{{\\rm{D}}}}\\)) are fixed to climatological mean values, as local variability in atmospheric evaporative demand is more strongly related to temperature than specific humidity variation. Evaporation is also a linear function MS (that is, soil wetness), expressed in terms of saturation state.<\/p>\n The change in soil moisture stock is precipitation less soil evaporation, infiltration and saturation-excess ponding rate of water, and soil drainage. We simulate it as<\/p>\n $${\\mu }_{{\\rm{S}}}\\frac{{\\rm{d}}{M}_{{\\rm{S}}}}{{\\rm{d}}t}=P-{E}_{{\\rm{S}}}-{Q}_{{\\rm{S}}}(P)-{Q}_{{\\rm{D}}}({M}_{{\\rm{S}}})$$<\/p>\n \n (9)\n <\/p>\n where P is daily precipitation, QS is ponded surface water (infiltration excess and surface excess), QD is soil drainage (also known as subsurface runoff in many land-surface models) and \u03bcS is a geometric parameter to convert between water equivalent heights and saturation state. QS is parameterized using the US Soil Conservation Service Curve Number scheme developed by the United States Department of Agriculture66<\/a>,67<\/a>,68<\/a>, which posits a quadratic increase in the surface ponding rate (saturation and infiltration excess) as a function of daily precipitation intensity and a parameter, S, that captures the maximal retention of soil water:<\/p>\n $${Q}_{{\\rm{S}}}(P)=\\frac{{(P-\\lambda S)}^{2}}{P-(1-\\lambda )S}$$<\/p>\n \n (10)\n <\/p>\n $$S=25.4\\,\\left(\\frac{\\mathrm{1,000}}{\\mathrm{CN}}-10\\right)$$<\/p>\n \n (11)\n <\/p>\n In this parametrization, \u03bb is a parameter representing the ratio between initial abstraction of precipitation and potential maximum soil water retention (S), which determines the minimum daily precipitation at which ponding occurs. Beyond this minimum, surface ponding rate increases according to the curve number (CN) parameter (a value between 0 and 100), which reflects the surface ponding rate sensitivity to precipitation intensity. We use CN values of 30 (dense forest cover), 50 (sparse woody vegetation on sandy loam soil) and 80 (grassland or impervious surfaces) in our model experiments67<\/a>. The soil drainage (or subsurface runoff term), QD, is parametrized in terms of unsaturated hydraulic conductivity:<\/p>\n $${Q}_{{\\rm{D}}}({M}_{{\\rm{S}}})={k}_{\\mathrm{sat}}\\,{{M}_{{\\rm{S}}}}^{c}$$<\/p>\n \n (12)\n <\/p>\n where ksat is the saturated hydraulic conductivity and c is parameter controlling the nonlinear response of soil water leakage to soil saturation states above the field capacity69<\/a>. Parameter values are for median sandy loam soils (Supplementary Table 3<\/a>).<\/p>\n A key innovation in the model is to assess the potential for precipitation to pool on the land surface and evaporate before it contributes to total runoff. We model ponded surface water and soil drainage as contributing to the surface water stock (ML), which is governed by<\/p>\n $${\\mu }_{{\\rm{L}}}\\frac{{\\rm{d}}{M}_{{\\rm{L}}}}{{\\rm{d}}t}=P-{E}_{{\\rm{L}}}+{Q}_{{\\rm{S}}}(P)+{Q}_{{\\rm{D}}}({M}_{{\\rm{S}}})-{Q}_{{\\rm{L}}}({M}_{{\\rm{L}}})$$<\/p>\n \n (13)\n <\/p>\n where EL is evaporation form the surface water stock, QL is runoff out of the surface water stock, and \u03bcL is a geometric parameter to convert between water equivalent heights and saturation state (that is, fractional \u2018fullness\u2019 of the reservoir). In our specification, QD, or soil drainage, contributes to surface water stocks, reflecting the dominant contribution of soil drainage to streamflow through subsurface flow70<\/a>,71<\/a>. EL is parametrized similarly to ES, but as a function of ML instead of MS and with a lower bulk surface resistance (rL). This lower resistance captures stomatal control on transpiration over land, and the fact that water evaporates more easily from a free surface (Supplementary Table 3<\/a>). As our model stipulates a minority areal fraction covered by the surface water stock (10% in the median parameters; Supplementary Table 3<\/a>), we assume that the temperature (and thus vapour pressure deficit) over surface water is dictated by the energy balance over soils. This simplification enables us to simulate EL without simulating horizontal transfers within the atmosphere, but results are similar using a simple bulk temperature-homogenization scheme (Supplementary Fig. 7<\/a>). QL is parametrized as a quadratic function of ML and a maximum fractional daily outflow rate parameter (\u03c9), such that QL\u2009=\u2009\u03c9ML2. In this formulation, runoff occurs only at the boundary of the model, enabling the model to track whether ponded surface water evaporates before leaving the grid cell. We use a global mean range of \u03c9 values based on the mean turnover rate of freshwater in lakes and reservoirs, estimated as the ratio of mean runoff (5\u2009\u00d7\u2009104\u2009km3\u2009yr\u22121; ref.\u200972<\/a>) to global total surface freshwater volume (2\u2009\u00d7\u2009105\u2009km3; ref.\u200973<\/a>), or around 0.1% per day. We also run the model with larger \u03c9 values (0.5\u20131% per day) to account for the volumetric dominance of large, slowly overturning lakes.<\/p>\n We initialize our model at a surface air temperature of 280\u2009K and with a saturation state of 0.5 for the soil and surface water stocks. In contrast to more complex land-surface models featuring deep soil water stocks, all prognostic variables are governed by fast first-order dynamics with characteristic e-folding times of the order of days, as defined by the geometry, parameters and forcings of the model. As such, we run our idealized model experiments without a multi-year spin-up common for complex land-surface models.<\/p>\n Land-surface model experimental design<\/p>\n To simulate the TWS effect using SEMB-H, we run the model using P and FSW forcings from CPC and GEWEX-SRB at 1,400 points randomly sampled across the global range of climatological precipitation. At each location, we run the model twice, using forcings for the highest and lowest GP year available in the observational CPC data. This approach captures variation in both daily P and S and how they, in turn, shape GP. To isolate the effect of GP from that of annual precipitation, we scale the daily P time series from each of these two years to ensure that they have equal totals (Extended Data Fig. 8<\/a>). This eliminates confounding TWS effects of low and high precipitation years, analogously to our control for annual precipitation in equation (2<\/a>). We estimate the simulated TWS effect (\u0394TWSsim) as the difference between maximum and minimum GP years in volume-weighted sums of ML and MS. To match the water-year scale of our empirical analysis, we take this difference as an average across the water-year (denoted by overbars):<\/p>\n $${\\Delta {\\rm{T}}{\\rm{W}}{\\rm{S}}}_{{\\rm{s}}{\\rm{i}}{\\rm{m}}}={\\mu }_{{\\rm{S}}}{a}_{{\\rm{S}}}(\\overline{{{M}_{{\\rm{S}}}|}_{{G}_{{\\rm{P}}}^{\\text{max}}}}-\\overline{{{M}_{{\\rm{S}}}|}_{{G}_{{\\rm{P}}}^{\\text{min}}}})+{\\mu }_{{\\rm{L}}}{a}_{{\\rm{L}}}(\\overline{{{M}_{{\\rm{L}}}|}_{{G}_{{\\rm{P}}}^{\\text{max}}}}-\\overline{{{M}_{{\\rm{L}}}|}_{{G}_{{\\rm{P}}}^{\\text{min}}}})$$<\/p>\n \n (14)\n <\/p>\n where aS and aL are the areal fractions of soil compared with surface water stocks. As these simulated TWS responses reflect a wide range of actual GP differences, we linearly scale \u0394TWSsim to a uniform 0.1\u0394GP as in our empirical analysis. Finally, we express the simulated TWS effect as a global statistical relationship in terms of climatological precipitation, for comparison to our empirical results:<\/p>\n $$\\Delta {{{\\rm{TWS}}}_{\\mathrm{sim}}}^{i}={\\beta }_{0}+{\\beta }_{1}{\\bar{P}}^{i}+{\\varepsilon }^{i}$$<\/p>\n \n (15)\n <\/p>\n where i denotes points in our sample. As we keep annual precipitation and surface parameters constant in our simulations, equation (15<\/a>) is analogous to equation (2<\/a>), such that \u03b20 and \u03b21 are the simulated equivalents of \u03b3 and \u03c7. We estimate uncertainty in equation (15<\/a>) analogously to in equation (4<\/a>), and compare the confidence intervals of the empirical and simulated TWS effects. We also track the mean differences in water-year totals of EL, ES, QS and QL between the maximum and minimum GP years to attribute the TWS effect to surface water budget changes. EL reflects evaporation from ponded surface water as well as from initial moisture stocks, as EL is initialized at 0.5 (that is, half full), and as such is not necessarily limited by QS.<\/p>\n The main free parameters in SEMB-H are the dry radiative coefficient (\u03b1), the CN, the fractional daily outflow rate from the surface water stock (\u03c9) and the geometric parameters governing the potential height and fractional area of the soil and surface water stocks (h and a). To test the robustness of our simulation results to values of these free parameters, we rerun the simulations using 243 combinatoric parameter sets at the 1,400 points (a total of about 365,000 simulations; Supplementary Table 3<\/a>), and compare the range of simulated TWS effects (equation (15<\/a>)) to the confidence intervals of the empirical effect. We further characterize the relative influence of free parameter values on the simulated TWS effect using Sobol\u2019s indices, a variance-based sensitivity analysis across the full parameter space74<\/a>. This analysis allows us to assess how uncertainty in each parameter, individually and jointly (through interactions), affects uncertainty in the simulated TWS effect.<\/p>\n Imputed TWS changes<\/p>\n Using the global coefficients from equation (2<\/a>) and local climatological precipitation across global land area, we estimate the TWS impact of GP trends (\\(\\Delta {\\mathrm{TWS}}_{\\mathrm{hist}}\\)) over 1980\u20132022 using<\/p>\n $$\\Delta {\\overline{{\\rm{T}}{\\rm{W}}{\\rm{S}}}}_{\\mathrm{hist}}^{i}=(\\gamma +\\chi {\\bar{P}}^{i})\\Delta {{G}_{{\\rm{P}}}}_{\\mathrm{hist}}^{i\\,}$$<\/p>\n \n (16)\n <\/p>\n In this formulation, TWS impacts depend on climatological precipitation and concentration trends at a given grid cell (superscripts i). This approach assumes that TWS effects of concentration, as estimated from detrended interannual variability, operate similarly over multi-decadal trends. We present the imputed historical TWS impact as an average across data products for simplicity.<\/p>\n To project future TWS impacts of concentration under continued global warming, we first use a simple thermodynamic model to project daily precipitation distribution changes under mean warming across global land45<\/a>. This model enables efficient and transparent estimates of future concentration globally, isolated from the influence of long-term annual precipitation changes and avoids consistent daily precipitation biases in general circulation models. The model intensifies observed daily precipitation distributions (I) at a fractional scaling rate of \u03c6 per degree of local warming (\u0394T):<\/p>\n $${I}^{i{\\prime} }={(1+\\varphi )}^{\\Delta {T}^{i}}{I}^{i}$$<\/p>\n \n (17)\n <\/p>\n For \u03c6, we use the Clausius\u2013Clapeyron saturation vapour pressure sensitivity to temperature of 7%\u2009\u00b0C\u22121, drawing from observational literature linking daily precipitation intensity to temperature8<\/a>,75<\/a>. Our goal is to briefly contextualize the TWS implications of median concentration expectations, so here we do not explore uncertainty in this parameter. For the initial observed intensity, I, we use daily precipitation distributions for the water-year with GP closest to the long-term median. Precipitation intensification is projected at an illustrative 1\u2009\u00b0C of global mean surface air temperature warming beyond the study period (or about 2\u2009\u00b0C above the pre-industrial limit), in line with current plausible mitigation scenarios for the late-twenty-first century. The local surface air temperature warming pattern at 1\u2009\u00b0C of global mean warming is derived from the average of 26 CMIP6 models76<\/a>, run under the SSP5-8.5 scenario.<\/p>\n Annual precipitation generally changes by less than 7%\u2009\u00b0C\u22121, implying an increase in dry days. Consistent with our empirical model formulation, we enforce constant annual precipitation to isolate the influence of concentration changes under warming. To achieve this, we offset the increased precipitation from intensification by converting the lowest-M precipitation days to dry days, where M is the number of days for which the total precipitation equals the intensification-induced annual increase. We then difference the GP of this final concentrated distribution from the GP of the initial, observed distribution to estimate GP changes under warming. These projected GP changes are of similar order to historical trends (1\u201310 basis points; Fig. 1g<\/a> and Extended Data Fig. 3<\/a>).<\/p>\n Finally, we use this GP-change field to project changes in TWS, as in equation (16<\/a>). Apart from mapping local TWS impacts, we also estimate global population exposure to standardized mean TWS changes using population data from the Intersectoral Impact Model Intercomparison Project77<\/a>. We first project the TWS impacts onto the population data to derive a frequency distribution of exposure, then quantify the proportion of the global population experiencing standardized anomalies exceeding 0.5 standard deviations\u2014an illustrative threshold indicative of abnormally dry conditions10<\/a>. We simply provide a median estimate for these population exposure estimates as a first-order characterization the global importance of these effects. Importantly, our approach assumes that spatial and temporal patterns of global irrigation remain constant in the future. Because intensively irrigated areas experience a stronger TWS effect from precipitation concentration (Extended Data Fig. 6a<\/a>), expansion or intensification of irrigated areas could lead to locally stronger future drying. Conversely, improvements in irrigation efficiency may reduce withdrawals per unit area in response to precipitation concentration.<\/p>\n","protected":false},"excerpt":{"rendered":"Data TWS data come from the GRACE and GRACE-FO mass concentration solutions (RL06) developed by the Center for…\n","protected":false},"author":2,"featured_media":44205,"comment_status":"","ping_status":"","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[2],"tags":[17,2754,12226,1553,12227,431],"class_list":{"0":"post-44204","1":"post","2":"type-post","3":"status-publish","4":"format-standard","5":"has-post-thumbnail","7":"category-canada","8":"tag-canada","9":"tag-climate-change","10":"tag-humanities-and-social-sciences","11":"tag-hydrology","12":"tag-multidisciplinary","13":"tag-science"},"_links":{"self":[{"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/posts\/44204","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/users\/2"}],"replies":[{"embeddable":true,"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/comments?post=44204"}],"version-history":[{"count":0,"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/posts\/44204\/revisions"}],"wp:featuredmedia":[{"embeddable":true,"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/media\/44205"}],"wp:attachment":[{"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/media?parent=44204"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/categories?post=44204"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/www.europesays.com\/canada\/wp-json\/wp\/v2\/tags?post=44204"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}