Critical Aquifer Overdraft Accelerates Degradation of Groundwater Quality in California's Central Valley During Drought

Drought‐induced pumpage has precipitated dramatic groundwater‐level declines in California's Central Valley over the past 30 yr, but the impacts of aquifer overdraft on water quality are poorly understood. This study coupled over 160,000 measurements of nitrate from ∼6,000 public‐supply wells with a 30 yr reconstruction of groundwater levels throughout the Central Valley to evaluate dynamic relations between aquifer exploitation and resource quality. We find that long‐term rates of groundwater‐level decline and water‐quality degradation in critically overdrafted basins accelerate by respective factors of 2–3 and 3–5 during drought, followed by brief reversals during wetter periods. Episodic water‐quality degradation can occur during drought where increased pumpage draws shallow, contaminated groundwater down to depth zones tapped by long‐screened production wells. These data show, for the first time, a direct linkage between climate‐mediated aquifer pumpage and groundwater quality on a regional scale.

groundwater-level decline to observed trends of degrading groundwater quality related to legacy agricultural land use (Burow et al., 2013;Jurgens et al., 2020).
Approximately 3 million residents of the Central Valley depend on public-supply wells (PSWs) for drinking water (Belitz et al., 2015), which are clustered around major urban centers within the largely agricultural landscape (Figure 1b). Typical PSWs are long-screened and open to deeper portions of aquifer systems than domestic wells (Voss et al., 2019), making them comparatively resilient to supply failure during drought. However, water-quality degradation is a major and potentially costly issue for drinking-water purveyors in the Central Valley where legacy agricultural recharge with elevated concentrations of nitrate, fumigants, salinity, and uranium has penetrated to depths commonly exploited for public drinking-water supply (Burow et al., 2008Hansen et al., 2018;Jurgens et al., 2010). Migration of agricultural recharge to depth throughout the Central Valley has been driven by intensive irrigation and pumpage from production wells, which, in aggregate, have increased downward vertical seepage rates sixfold from pre-development conditions (Faunt, 2009;Williamson et al., 1989). Although seasonal variations in pumpage have been shown to cause transient incursions of nitrate-laden groundwater to PSWs from shallower depth zones (Bexfield & Jurgens, 2014;Yager & Heywood, 2014), long-term effects of overdraft on groundwater quality have not been systematically evaluated.
The purpose of this study was to statistically evaluate dynamic relations between groundwater levels and water quality on an unprecedented spatial and temporal scale. To this end, we leveraged over 320,000 measurements of groundwater level and nitrate collected from separate monitoring networks throughout the Central Valley over the past 30 yr. The period of interest (1990-2018) encompasses the last significant 3 of 10 pluvial (1993)(1994)(1995)(1996)(1997)(1998)(1999) and the past three California droughts (2001-2002, 2007-2009, and 2012-2016), which represent periods of anomalously low soil moisture and water availability predicted to recur with greater frequency in the 21st century (DWR, 2020;Ullrich et al., 2018). We hypothesize degrading water quality from production wells becomes more prevalent in the Central Valley during extended droughts where intensive pumpage causes the rapid decline of groundwater levels and increased vertical seepage rates.

Study Area and Geostatistical Approach
The Central Valley is an asymmetric structural trough containing up to 15 vertical km of layered continental and marine sediments, the upper 30-300 m (100-1,000 ft) of which are commonly utilized for public drinking-water supply (Voss et al., 2019). The Sacramento Valley (SAC) subregion occupies the northern third of the Central Valley and has a wetter climate and higher surface-water flows than the SJV subregion in the southern two thirds, where the effects of pumping on groundwater levels have been more pronounced ( Figure 1a). Average annual precipitation ranges from 130 mm in the southern SJV to 660 mm in the northern SAC (Faunt, 2009). The SJV can be further subdivided into the San Joaquin Basin and the internally drained Tulare Basin ( Figure 1a). Pumpage-induced downward groundwater flow is superimposed over regional flow systems that recharge at the valley's margins and move laterally to discharge zones along its axial trough (Faunt, 2009).
This study focuses on 6,031 PSWs throughout the Central Valley that had available water-quality data during the period of study (Figure 1b). We calculated both overall (by well) and spatially weighted (by area) statistics to evaluate groundwater quality. Spatial aggregation of groundwater-quality data can be used to calculate areal proportions of aquifers where groundwater quality exceeds a specified concentration threshold and prevents zones of high well density (e.g., cities) from overrepresenting aquifer conditions within a broader study region (Belitz et al., , 2015. For this study, we projected a regularized hexagonal grid over the Central Valley so cells were evenly distributed and not oriented to or against regional hydraulic gradients. We optimized grid resolution to the scale of major urban centers (∼300 km 2 ; Text S1).
We used nitrate as a proxy for overall groundwater quality in the Central Valley because it is required by law to be sampled from active PSWs on an annual basis and often co-occurs with a host of other contaminants that are monitored with less frequency (e.g., fumigants, salinity, and uranium; Burow et al., 2019;Rosen et al., 2019). To assess ambient conditions in aquifers used for public drinking-water supply, we compiled nitrate measurements for over 160,000 raw (untreated and unblended) groundwater samples collected from PSWs for regulatory-compliance monitoring (SWRCB-DDW, 2019a) and the California State Water Board's Groundwater Ambient Monitoring and Assessment (GAMA) program . Nitrate measurements for individual wells were seasonally averaged by "water-level year" to ensure wells with higher sampling frequency were not overrepresented in our analysis (Text S2). Similar to the water year, we defined the water-level year as beginning November 1 of the previous calendar year and extending through October, which is the month when groundwater levels are typically lowest in the Central Valley following the summer pumping season (Faunt, 2009). We obtained groundwater-level measurements from the California Department of Water Resources (DWR) periodic groundwater level measurement data set (DWR, 2021b) and interpolated groundwater-level data to the location of study PSWs to assess the general direction and magnitude of hydrologic change at water-quality monitoring sites (discussed further in Section 2.2).
We assessed water-quality degradation using a vulnerability threshold of 5 mg of nitrate as nitrogen per liter (mg-N L −1 ), which is one half the Federal maximum contaminant level (MCL; 10 mg-N L −1 ) and the level above which PSWs are required by law to increase nitrate monitoring frequency from annual to quarterly schedules (SWRCB-DDW, 2019b). Nitrate concentrations in excess of the MCL may be underrepresented in regulatory monitoring data because violating wells are taken off-line if downstream mitigation via treatment or blending is not possible. Thus, a lower vulnerability threshold more accurately assesses the leading edge of groundwater-quality degradation. Analysis of the regulatory data showed the MCL was exceeded within 1, 5, and 10 yr of the first exceedance of the vulnerability threshold, respectively, for 15%, 47%, and 74% of all wells where the vulnerability threshold was exceeded prior to the MCL.
We calculated overall and spatially weighted detection frequencies of wells with "high-nitrate" (>5 mg-N L −1 ) groundwater by water-level year and subregion. Detection frequencies are expressed as decimal fractions (count of wells with nitrate >5 mg-N L −1 divided by count of sampled wells) and 95% confidence intervals for overall detection frequencies were quantified using the Jeffreys interval for the binomial distribution . Spatially weighted detection frequencies represent the area-weighted average of all cell-scale detection frequencies (Text S1). We evaluated interannual trends for detection frequencies and corresponding groundwater-level anomalies using Kendall's Tau and Theil-Sen slopes (Helsel et al., 2020).
In addition to detection frequencies, we assessed relative "year-to-year" changes in nitrate concentrations at individual PSWs by differencing the nitrate concentration for a given year with that of the year prior for wells sampled in consecutive years. We used interpolated groundwater-level data to calculate a corresponding year-to-year groundwater level change for each differenced nitrate value. Aggregate year-to-year changes within a given cell were summarized by the mean of all differenced nitrate values (to be more sensitive to change magnitude) and the median of all differenced water level values (to be more sensitive to change direction) for each water-level year. Relations among nominal variables assigned to cells based on median groundwater-level change (rise or decline) and mean nitrate change (increase or decrease) were evaluated using contingency analysis (Helsel et al., 2020).

Groundwater-Level Modeling
We analyzed over 160,000 depth to groundwater measurements to create modeled groundwater-depth surfaces on an annual timestep for the Central Valley from 1990 to 2018. Our modeling objective was to simulate the relative direction and magnitude of interannual groundwater-level changes in unconfined and semi-confined aquifers used for public drinking water supply. Numerous discontinuous silt and clay layers interbedded with coarse-grained sediments in Central Valley aquifers can create locally confining conditions where head responses to seasonal pumpage can exceed 30 m (100 ft; Faunt, 2009;Williamson et al., 1989). To minimize seasonal effects on our analysis, we only modeled "spring" water levels (median depth to groundwater from January through April) when aquifers are mostly recovered from intensive summer pumpage (DWR, 2015, Appendix E). Additionally, we filtered out deep wells (depth to perforations or median water level >213 m [700 ft]) and those with high seasonal fluctuations (median seasonal difference >15 m [50 ft]) that are more likely to be influenced by locally confining conditions and are less representative of ambient pressure conditions in the surrounding area (Text S3, Figure S1). Although PSWs in the Central Valley are screened through heterogeneous sediments under a range of confining conditions, modeled groundwater levels were considered diagnostic of ambient pressure responses to hydrologic stressors in a given area and were not intended to reproduce head values for specific depth zones or geologic strata.
We used a combination of multiple linear regression (MLR) and ordinary kriging to impute gaps in time-series data and interpolate groundwater-depth surfaces in space, respectively. The imputation approach of Evans et al. (2020) was implemented with a novel, iterative algorithm to fill data gaps at target wells with data for >50% of study years using MLR on the top five most correlated groundwater-level records in the subregion. This resulted in complete annual time series for spring groundwater levels at 3,239 locations throughout the Central Valley from 1990 to 2018 (Levy, 2021;Text S3). We interpolated depth to groundwater for individual water-level years at a spatial resolution of 25 km 2 (township scale) over the Central Valley using ordinary kriging (Pauloo et al., 2020;Pebesma, 2004). The imputation routine allowed for the same control points to be used for interpolation at each timestep, which optimized quantification of relative groundwater-level change over areas encompassing heterogeneous aquifer systems with strong vertical head gradients. Modeled values accurately approximated depth to groundwater for independent validation data with a mean absolute error of 5.0 m (16 ft) and were highly correlated with observed interannual fluctuations in validation time series (Pearson's r > 0.95; Text S3, Figures S2 and S3). Groundwater-level anomalies were calculated at 6,031 study PSW locations by subtracting the mean from modeled time series.

Regional Linkages Between Groundwater Levels and Nitrate
Interannual changes in groundwater levels and nitrate at PSWs were responsive to regional drought cycles in the SJV from 1990 to 2018 (Figure 2). During the last major pluvial (1993)(1994)(1995)(1996)(1997)(1998)(1999), which followed a major regional drought from 1987 to 1992, groundwater levels at PSWs in the SJV rose to a median elevation of 3.7 m (12 ft) above average for the period of record (Figure 2b). This was followed by a long-term trend of decline in the median groundwater-level anomaly from 2000 to 2018 (Kendall's Tau = 0.95, p-value < 0.001, Theil-Sen slope = 0.6 m yr −1 [1.9 ft yr −1 ]) punctuated by three droughts of increasing duration and severity, during the latter two of which the rate of decline increased by a factor of 2-3 (Figures 2b and S4a). Groundwater levels reached their lowest elevations by the end of the 2012-2016 drought when the median and interquartile range for depth to groundwater were 7.3 m (24 ft) and 3.4-11.9 m (11-39 ft) greater than average for the period of record, respectively. ncdc.noaa.gov/cag/) averaged by water-level year. Periods of moderate-to-extreme California drought (PDSI < −2) and are shaded pink. (b) Groundwater-level anomaly (deviation of depth to groundwater from the mean) in feet (ft) and meters (m). Solid line represents the median anomaly modeled for study publicsupply wells in the SJV, and dotted lines represent the interquartile range. (c) Annual detection frequencies of high-nitrate (>5 mg-N L −1 ) groundwater for study wells in the SJV (solid line and points) with the 95% confidence interval (95% conf. intvl.) for the binomial distribution (dotted lines). Dashed line is the median depth to perforations for wells with high-nitrate samples normalized to that value in 2002. (d and e) Analyses for Sacramento Valley as described for panels (b and c). (f) Relative changes for cell detection frequencies of high-nitrate groundwater during contrasting drought and recovery periods.
Detection frequency of high-nitrate (>5 mg-N L −1 ) groundwater in the SJV showed an episodic increase during drought superimposed on a long-term upwards trend, similar to patterns of groundwater-level decline (Figure 2c). During the pluvial of 1993-1999, however, overall detection frequencies dropped from above 0.25 to below 0.20. Water-quality statistics from 1990 to 2001 should be interpreted with caution because the number of wells reporting nitrate data in the Central Valley increased by a factor of 5 during this period before leveling off at a rate of ∼3,500 wells per year after 2001 due to changes in regulatory reporting requirements (Tables S1-S4; Jurgens et al., 2020). However, relative nitrate concentrations for specific wells with long-term records in the SJV also decreased during the wet period from 1993 to 1999, supporting statistical inferences from the broader data set.
Spatially weighted detection frequencies in the SJV had near-identical long-term trend (0.0021 yr −1 ) and drought (factor increase of 3-5) slopes but were on average 0.03 lower than overall detection frequencies ( Figure S4b). This indicates that although the spatial distribution of high-nitrate PSWs throughout the SJV is uneven (e.g., Nolan et al., 2014) drivers of dominant temporal trends were distributed and not localized to a handful of high-density well clusters (Figure 2f). Detection frequency of high-nitrate groundwater was highest during peak drought in 2015 with overall and spatially weighted values of 0.29 (about one third of all wells) and 0.24 (about one quarter of the aquifer area used for public supply), respectively.
In contrast, SAC showed minimal groundwater-level decline at PSWs from 2000 to 2018 (Kendall's Tau = 0.56, p-value = 0.012, Theil-Sen slope = 0.06 m yr −1 [0.2 ft yr −1 ]) and much lower detection frequencies of high-nitrate groundwater (∼0.12) with no long-term trend (Figures 2d and 2e). The small step increase in detection frequencies from pre-2000 levels in the SAC is likely due to an expansion of the population of wells reporting water-quality data during the former period as described above (Table S1). Groundwater storage loss in the SAC has been minimal compared to the SJV (Figure 1a) due in large part to wetter climate and abundant surface-water diversions, which are restricted in more southerly basins during drought (Faunt, 2009). Lower nitrate concentrations in the SAC have been attributed to dilution from greater rain and surface-water inputs, less agricultural land use proximal to PSWs, and higher attenuation by microbial action on finer-grained aquifer sediments than in the SJV where denitrification rates tend to be lower (Burow et al., 2013;Landon et al., 2011). Comparatively muted and variable responses of water quality to regional climate cycles in the SAC compared to the SJV suggest different processes control nitrate transport to PSWs among the two subregions.
We assessed percent area of the Central Valley with year-to-year cell changes for groundwater levels and nitrate to better understand the extent to which improving or degrading groundwater quality co-occurred with rising or declining water levels (Figure 3). Groundwater-level decline and rise were associated with nitrate increase and decrease, respectively, for aggregate cell changes from 2000 to 2018 (Pearson's χ 2 = 4.79, p-value = 0.029; Text S5, Figure S7). During inter-drought periods the total area of declining groundwater levels tended to be <50% and comprised roughly equivalent areas with increasing and decreasing nitrate concentrations. During droughts, however, the area of drawdown expanded to >75% and was dominated by areas of nitrate increase compared to decrease by almost twofold (Figure 3b).
The total area of increasing nitrate was >50% for most of the period of study, reflecting widespread trends of degrading groundwater quality in the Central Valley (Burow et al., 2013;Jurgens et al., 2020), and expanded up to 65% during drought (Figure 3a). However, the area of increasing nitrate contracted to <50% during brief recovery periods following drought peaks in 2004-2017. Although groundwater nitrate tended to decrease at the onset of recovery periods following drought, rising groundwater levels were also sometimes associated with increasing nitrate later in recovery periods-particularly in 2006. Despite heterogeneities expected on a regional scale, periods of more widespread nitrate increase were strikingly synchronous with regional drought cycles.

Implications for Sustainable Groundwater Management
Our results show that long-term trends of increasing groundwater-quality degradation in the Central Valley accelerate in response to aquifer overdraft during regional drought cycles. Central Valley aquifers are commonly age-stratified, with modern (recharged post-1950s), high-nitrate groundwater overlaying pre-modern, low-nitrate groundwater (Figure 4). In this setting, long-screened PSWs often produce binary mixtures of young and old water depending on the vertical position of the screened interval within the aquifer system (Jurgens et al., 2016). The modern fraction of produced waters tend to be several decades old and reflects histories of anthropogenic nitrate loading at the land surface (Burow et al., 2008(Burow et al., , 2013. Because long-screened PSWs integrate mixtures of groundwater recharged over broad timescales, episodic water-quality responses to drought are more likely to be driven by the effects of increased pumpage on aquifer pressure dynamics than contemporaneous recharge processes. Here, we show progressively deeper wells were affected by water-quality degradation during periods of long-term and episodic groundwater-level decline in the SJV (Figures 2c and S5), indicating effects of pumping stress on the vertical penetration of modern groundwater to depths tapped by PSWs.
We propose that rapid head loss during periods of intensive groundwater pumpage can shift the vertical groundwater-age structure of aquifers, increasing downward seepage rates and the proportion of younger, nitrate-laden discharge to production wells (Figure 4). This effect can be exacerbated by exhaustion of lateral flows through deeper,  lower-permeability strata during overdraft, as evidenced by land subsidence in extreme cases (Faunt et al., 2016). Intraborehole flow may also enhance downward vertical fluxes across confining layers under pumping stress (Faunt, 2009;Gailey, 2017). Unexpectedly, we found water quality improved at the onset of wet periods following droughts, which is most likely due to relaxation of vertical pumping gradients and recovery of lateral flows through deeper strata. Increased recharge during subsequent wet periods could either dilute or transport additional nitrate to shallow groundwater, but would likely only affect the quality of produced waters on multiannual to decadal timescales. Although the exact hydraulic mechanisms will vary depending on well construction and hydrogeologic setting, our results suggest sustainable management of pumping operations could have long-term co-benefits for both supply availability and quality in critically overdrafted basins of the Central Valley.
Groundwater-level decline is a symptom of unsustainable pumpage where aquifer outputs exceed inputs . However, related water-quality impacts are not one-size-fits-all and depend upon contamination sources and flow patterns within complex aquifer systems. Long-term aquifer development can alter geochemical conditions that affect contaminant fate and transport (Haugen et al., 2021). In some geologic settings, excessive pumpage can increase the discharge of deeper, older groundwater to production wells depending on the structure of layered aquifer systems (Izbicki et al., 2005). Relations between overdraft and water quality in the Central Valley are driven by high vertical-flux rates from intensive irrigation and pumpage, and are most applicable to aquifer systems where flow under pumping stress is predominantly vertical as opposed to lateral (e.g., Nunes et al., 2021). Although considerable heterogeneity exists on the scale of individual wells, our results evidence critical linkages between groundwater use and quality on a regional scale.

Data Availability Statement
Groundwater level and nitrate data used in this study are publicly available through DWR (2021b), SWRCB-DDW (2019a) and Jurgens et al. (2018). Construction data for PSWs were compiled from various sources (DWR, 2021a;SWRCB, 2021;Voss et al., 2019). Modeled groundwater levels are available in the data release of Levy (2021). study was funded by the California State Water Resources Control Board as part of its GAMA Program (Agreement 18-005-250). Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government.