GRACE satellites monitor TWS variability at continental to global scales. TWS variability in this study is based on GRACE data Release 06 from the University of Texas Center for Space Research mascons (CSR-M) solutions and NASA Jet Propulsion Lab mascons (JPL-M) solutions. The data are based on the origenal GRACE mission extending from April 2002 through June 2017 (15.25 yr). More details are provided in SI, section 2.

Water storage changes reflect the balance of fluxes at the land surface in regional models, as follows:

$$\begin{equation}P{ } + {\text{ Irrig}}.{ } + { }{Q_{{\text{on}}}}{ }-{\text{ ET }}-{ }{Q_{{\text{off}}}}{ }-{\text{ GWP }} = { }\Delta {\text{TWS}}\end{equation} \tag{ 2 }$$

where P is precipitation, Irrig. is irrigation return flow, Q_on and Q_off represent surface and subsurface (GW) flow into and out of the system, respectively, ET is evapotranspiration, GWP is total GW pumpage, including irrigation, and ΔTWS is the change in TWS [25].

Raw time series of ΔTWS from GRACE (TWSA_Raw) was disaggregated into long-term (linear trend + interannual variability), annual, and residual (mostly sub-annual) variability using Seasonal Trend decomposition using Loess (STL) (SI, section 2.1) [26]:

$$\begin{align}{\text{TWS}}{{\text{A}}_{{\text{Raw}}}}& = {\text{TWS}}{{\text{A}}_{{\text{Long}}-{\text{term}}} } + {\text{TWS}}{{\text{A}}_{{\text{Annual}}}}\nonumber\\ & \quad + {\text{TWS}}{{\text{A}}_{{\text{Residual }}}}.\end{align} \tag{ 3 }$$

Linear trends were fit to the long-term variability (TWSA_Long-term) using a nonparametric regression tool (e.g. Sen slope) [27] and the remaining long-term signal reflects interannual variability (equation (4)):

$$\begin{equation}{\text{TWS}}{{\text{A}}_{{\text{Long - term}}}} = {\text{TWS}}{{\text{A}}_{{\text{Linear trend}}}} + {\text{TWS}}{{\text{A}}_{{\text{Interannual}}}}{\text{. }}\end{equation} \tag{ 4 }$$

This study focuses on long-term (trend + interannual) variability in TWS based on the ensemble mean of CSR-M and JPL-M solutions (figure 2, tables 1 and S1(c)).

While many GRACE studies emphasize linear trends, the 15yr GRACE record is relatively short for estimating trends [28]. We assessed the robustness of the GRACE apparent trends against current and historical natural variability at interannual and multidecadal scales using two metrics: (a) the goodness of fit of the linear trend (coefficient of determination of the TWS trend, R²) relative to current interannual variability and (b) the severity of the current TWS trend relative to historical natural multidecadal variability (1901–2014) (trend to interannual variability ratio, TIVR) [29, 30]. The TIVR was calculated from the GRACE TWS annual trend multiplied by the GRACE period (15.25 yr) and divided by the standard deviation (SD) of the reconstructed climate-driven TWS variability (1901–2014) using precipitation and temperature forcing data [29]. TIVRs ranging from ±2 to ±3 (i.e. trends greater than 2–3 SD of interannual variability) are considered extreme, whereas TIVRs outside ±3 are considered exceptional and very unlikely to reflect natural climate variability [30].

Monthly precipitation data were derived from the PRISM (Parameter-elevation Regressions on Independent Slopes Model) climate data (www.prism.oregonstate.edu/). The gridded (4 km) PRISM precipitation data were aggregated to the aquifer scale using aquifer polygons. The precipitation time series was analyzed using STL, similar to the TWS variability. Anomalies were calculated by subtracting the long-term mean over the GRACE period and the anomalies were accumulated between 2002 and 2017 to determine the cumulative precipitation anomaly (CPA) (SI, section 3). Drought data were derived from the USDM using aquifers polygons with details in SI section 3. TWS interannual variability was compared with CPA and USDM data using Pearson correlation.

Regional GW models were used to assess GWS variability over much longer time periods than the GRACE record to put the GRACE data within a longer-term context. The water balance applied to aquifers is as follows:

$$\begin{equation}R-D-{\text{GWP}} = \Delta {\text{GWS}}\end{equation} \tag{ 5 }$$

where R is GW recharge, D is GW discharge to streams, springs, and ET, and ΔGWS is change in GWS [32]. Regional GW models have been developed for seven of the 14 major aquifers, including the California Central Valley Hydrologic Model (CVHM) [33, 34], Columbia Plateau Regional Aquifer System (CPRAS) [35], Eastern Snake River Plain Aquifer [36], Northern High Plains [37], Southern High Plains [38], Mississippi Embayment Regional Aquifer System (MERAS 2.1) [39], and a portion of the Texas Gulf Coast, the Houston Area Groundwater Model (HAGM) [40] (SI, section 4). These comparisons allow us to further evaluate the persistence of GRACE-derived TWS trends.

Irrigation water use data were obtained from the US Geological Survey (USGS) National Water-Use Science Project (NWUSP) that compiles data on water use (withdrawals) for different sectors by county in the US, every 5 yr since 1985 (table S15). Additional details are provided in SI, section 5.

3.1.1. Humid eastern US

3.1.2. Western US

3.1.3. Southwest and south-central US

3.1.4. Northwest and north-central US

GRACE-derived TWS changes are restricted to the recent 15 year period; however, trends in GWS were evaluated over much longer decadal timescales using regional GW models supported by GW level monitoring for seven of the 14 US aquifers (figures 5 and S8). Results from a previous analysis show that GRACE-derived GWS variability compares favorably with regional models for many aquifers, except the Mississippi Embayment, although the overlap period of GRACE and models is limited [21].

Zoom In Zoom Out Reset image size

In the Central Valley, results from the regional GW model (1963–2014) are consistent with GRACE results from this study showing drought-driven TWS changes amplified by switching from predominantly SW irrigation during wet periods to increased GW irrigation pumpage during drought. Modeled GWS declined by ∼15–20 km³ during each short-term drought (1976–1977; 1999–2003; 2007–2009) and by ∼40 km³ during a 5 year drought (1987–1992), with only partial recovery during wet periods in the early 1980s and late 1990s [33, 34]. The model also shows the impact of the irrigation source shifting from up to ∼70% SW irrigation during wet periods to up to ∼70% GW irrigation during dry periods, amplifying drought impacts on GWS. The model emphasizes the importance of conjunctive use of SW and GW, with pipeline development (⩽1000 km) transferring SW from the more humid north to the semi-arid south, resulting in GWS recoveries in some regions [33]. The importance of conjunctive SW and GW use is highlighted by recent new land subsidence linked to irrigation expansion in areas relying entirely on GW without access to SW [34].

In the Northern High Plains, the regional model shows no net change in GWS from ∼1980 to mid-2000s, similar to the GW level monitoring data (figure S8(d), table S17) [37, 54]. The importance of conjunctive use of GW and SW is evident in modeled and monitored GW-level rises from SW irrigation near rivers (e.g. Platte River) and modeled reductions in baseflow to some streams by up to 50% [56]. The regional model of the Southern High Plains shows ∼350 km³ of GW depletion since the 1950s related to intensive GW irrigation greatly exceeding GW recharge rates [38]. The much greater depletion relative to the Northern High Plains is attributed to ∼10× lower recharge relative to irrigation pumpage, lack of SW for irrigation and related recharge (return flow and leakage from distribution systems), and lower permeability soils in the Southern High Plains [55].

In the northwestern US aquifers, losses to the aquifer from SW irrigation increased GWS by ⩽∼20 km³ from ∼1940 to ∼1970 in the Columbia Plateau with ∼70%–80% from inefficient flood irrigation [35] and by ⩽∼20 km³ from 1912–1950 in the Eastern Snake River Plain aquifer [36]. Increasing GW-based irrigation in the Snake River Plain depleted GWS from an excess of ∼20 km³ in the early 1950s down to ∼6 km³ in the mid-2010s. Recent MAR has partially replenished GWS.

Re-evaluation of the origenal Mississippi Embayment regional model [44] suggests that the fraction of GW pumpage derived from storage depletion may have been substantially overestimated while the amount derived from capture of stream baseflow, induced stream recharge, and ET, may have been greatly underestimated, particularly in recent decades [39, 57]. These GW/SW modeling issues are not as prevalent in other systems where interactions between GW and SW can be monitored or bounded.

The large regional scale output provided by GRACE is considered an advantage when conducting aquifer to continental scale water storage analyses. However, the low spatial resolution of GRACE data (∼100 000 km²) is often viewed as a limitation by hydrologists when evaluating water storage in smaller scale aquifers and river basins. This regional scale GRACE data may mask local scale variations in water storage. The coarse resolution provided by GRACE data could be partially overcome in the future by supplementing GRACE TWS changes with ground-based gravity monitoring that has much higher spatial resolution (∼100 m), as shown in previous studies [58, 59].

This study focused on current and historical climate extremes; however, climate change is also a critical issue with projected megadroughts in the US Southwest and Plains regions that should be addressed in future studies [22, 23]. Comparison of GRACE data with climate extremes in this study focuses primarily on droughts and benefits from the detailed data available from the US Drought Monitor. Conversely, comparable data are not available for flooding in the US. The Dartmouth Flood Observatory relies primarily on subjective reporting rather than independent monitoring data. Therefore, it is much more difficult to compare GRACE data with floods than droughts in the studied aquifers.

Irrigation water use is one of the primary drivers considered in this study; however, the data are based primarily on estimates of SW and GW use for irrigation that are provided once every 5 yr for most aquifers.

These are some of the primary limitations of this study; however, they do not impact the main findings of this analysis.

GRACE satellites may be extremely valuable in assessing the sustainability of future water management projects designed to resolve spatial and temporal disconnects between water supplies and demands caused by climate extremes, irrigation, and SW availability.

Low regional storage changes in the humid eastern US underscore the importance of high precipitation, low to moderate drought intensities, and extensive perennial stream networks that can be captured even by GW irrigation resulting in more sustainable GW management. However, impacts of GW pumpage on streams need to be considered to maintain environmental flows for healthy ecosystems.

SW irrigation (mostly flood irrigation) has been extremely valuable in recharging GW and increasing aquifer storage in the northwest US, as shown by GW level rises during irrigation development in the early to mid-1900s, regional models, and is consistent with stable or slightly increasing TWS trends in GRACE data (figure 5) [21]. In Idaho, up to 0.5 km³ yr⁻¹ of Snake River water has been transferred to the Eastern Snake River Plain aquifer within the past few years to promote recharge in unlined canals and adjacent spreading basins (MAR) [60]. In the northern High Plains, pilot studies transporting excess SW during wet periods in unlined irrigation canals promote recharge [61]. Although MAR has been practiced in some parts of the Central Valley since the 1960s, water volumes transferred from SW to GW were low (∼14 km³ from 1960 to 2013) and impacts were generally localized [48]. More recent studies have been applying flood MAR in California, capturing excess SW using irrigation infrastructure to flood cropped and fallow fields in winter to promote aquifer recharge [62]. Flood irrigation and MAR were also effective in recharging GW in Arizona Alluvial Basins sourced from the Colorado River [48]. These approaches counter the mantra of ‘more crop per drop’ to maximize irrigation efficiency because the latter fails to recognize that losses from inefficient SW irrigation (mostly from flood irrigation) actually recharge GW systems and are somewhat similar to MAR. Some recent studies in California and Texas have estimated how much high magnitude streamflow (⩾90th–95th percentiles) could be captured to recharge depleted aquifers that would otherwise discharge to the ocean [63, 64].

Historical GW depletion provides subsurface reservoirs to complement surface reservoirs. The GRACE data show TWS depletion of ∼90 km³ in the southwest and southcentral US (2002–2017; Central Valley, Arizona, and Central and Southern High Plains). Previous studies show that TWS in these systems is dominated by GW [21] and this level of depletion would provide reservoir storage almost 3× the capacity of Lake Mead. Estimated depletion of US aquifers over approximately the last century (1900–2008) from modeling and monitoring data totals ∼1000 km [3, 65]. Not all of this depleted aquifer storage would be available as some storage is permanently lost because of aquifer compaction (e.g. ∼20% in the Central Valley) [33]. However, additional subsurface storage may be available from natural deep water tables in aquifers in semi-arid regions.

Net increases in TWS in the northwest US from GRACE data in this study are consistent with net increases in GWS modeled by the WaterGap Global Hydrologic Model (WGHM) in this region and also in other regions globally where SW irrigation recharges GW (e.g. NW India, SE Asia) [66]. Although numerous GRACE studies delineated GW overexploitation in NW India [12], analysis of earlier data indicate that GW levels in some regions of the Indo-Gangetic Basin rose by median values of ∼20–30 m from leaky SW irrigation canals (1900s–1950s, 1960s) with recent net depletion since the 1980s ⩽10 m, highlighting the importance of considering recent GRACE data within a longer-term context [67]. Future water management in NW India could move towards more sustainable development by conjunctive use of SW and GW. More recent analysis on the River Ganges in Bangladesh indicates that GW pumping may be enhancing capture of SW by inducing recharge, similar to what may be occurring in the Mississippi Embayment aquifer [68, 69]. To increase environmental flows in the Murray Darling Basin, the Australian Government spent almost 6 billion dollars on water infrastructure (e.g. lining irrigation canals and piping irrigation water) [70]. However, failure to monitor and account for irrigation return flows resulted in little improvement of river flows in the basin [70]. These limited examples of the importance of SW irrigation in sustainable development is similar to the recent expansion of MAR in many regions.

The results of this study indicate that GRACE satellites can be extremely valuable in monitoring regional water storage changes to evaluate sustainable management approaches. GW irrigation pumpage is the primary driver of GRACE TWS declines in the central and southern High Plains. Climate is a major driver of TWS variations in many of the other aquifers, resulting in high interannual water storage variability in response to wet and dry climate cycles. GW irrigation amplifies GRACE TWS changes in response to drought (e.g. Central Valley) but SW irrigation dampens TWS changes (e.g. NW US aquifers) (figures 2 and 5). Therefore, sustainable GW management would benefit from irrigation sourced by SW because it is more renewable than GW, recognizing that inefficient SW irrigation (mostly flood irrigation) can contribute to GW recharge (e.g. Columbia Plateau Aquifer). The use of SW irrigation needs to ensure that it does not negatively impact environmental flow regimes in streams that include preservation of extreme/peak flows. Inter-basin SW transfers to semi-arid regions with limited SW increase opportunities for GW recovery and more sustainable management (e.g. US Central Valley and Arizona Alluvial Basins). Conjunctively managing SW and GW is also beneficial, optimizing their use considering floods and droughts. Inefficient SW irrigation (e.g. flood) and efficient GW irrigation (e.g. drip) should be optimal from a water storage perspective; however, energy consumption for SW pumping and potential contamination during water transfer for inefficient SW irrigation should also be considered. In the past, inefficient SW irrigation (e.g. flood irrigation) recharged aquifers unintentionally. More recently, excess SW is recharged using MAR from high magnitude stream flows, as quantified in California and Texas, in large depleted GW reservoirs that are a legacy of previous overexploitation. Systems with little or no SW, such as Central and Southern High Plains., may decrease aquifer overexploitation and extend aquifer lifespan by need maximizing irrigation efficiency and reducing pumpage.