LetterThe following article is Open access

Exploring the InSAR signature associated with river-sourced recharge in California's San Joaquin Valley

, , and

Published 11 July 2024 © 2024 The Author(s). Published by IOP Publishing Ltd
, , Citation Wesley R Neely et al 2024 Environ. Res. Lett. 19 074072DOI 10.1088/1748-9326/ad5855

Download Article PDF
DownloadArticle ePub

You need an eReader or compatible software to experience the benefits of the ePub3 file format.

1748-9326/19/7/074072

Abstract

Recharge is a critical component for understanding aquifer systems and the sustainable management of groundwater resources, yet this process is challenging to measure at policy-relevant spatiotemporal scales. Building upon previous research, we tested the hypothesis that InSAR can be used to observe river-sourced recharge if the underlying recharge pathways are associated with sufficient clay content. Our analysis leveraged the decomposition of InSAR time series with interpretations of 3D resistivity models derived from airborne electromagnetic (AEM) surveys. We focused our analysis on two study sites where high density AEM data were available and river-sourced recharge is determined to have occurred during wet years: (1) near Fresno, California and (2) near Visalia, California. Sediment type and hydrogeological structure from AEM supported our hypothesis with the InSAR signature attributed to river-sourced recharge occurring only in the study site with semi-confined to confined conditions and relatively high fraction of interbedded clay within recharge pathways. The timing to peak amplitude, the key feature we wanted to isolate in the InSAR data, near Visalia was interpreted as a pressure pulse front associated with river-sourced recharge propagating into the San Joaquin Valley. This study further validated the potential of InSAR, coupled with AEM data, to map and monitor river-sourced recharge in aquifer systems. As InSAR data become more accessible, this approach holds promise for broader applications in groundwater science and management worldwide.

Export citation and abstractBibTeXRIS

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 license. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

1. Introduction

Groundwater currently supplies roughly a third of global freshwater demand [1] and is projected to be increasingly relied upon as changes in climate threaten the accessibility of surface-water resources [2]. In California's Central Valley, groundwater utilisation is even higher, satisfying up to two-thirds of water needs during very dry years [3]. In response to increasingly frequent and persistent droughts where groundwater replenishment is limited, many aquifer systems in this region have been overexploited resulting in significant land subsidence [37], water quality degradation [8, 9], permanent storage loss [10, 11], and wells running dry [12]. This has motivated state-wide legislature to ensure sustainable management of groundwater supplies by 2040 through the Sustainable Groundwater Management Act (SGMA) [13]. To achieve sustainability goals, there is a critical need to understand and monitor how recharge occurs in these aquifer systems.

2. Background and motivation

The aquifer systems of California's Central Valley are composed of coarse-grained sand and gravel interbedded with fine-grained clay and silt (constituting ∼60% of sediment volume) [14, 15]. Although complex and spatially heterogenous, these aquifer systems are often conceptually described as a single unconfined to semi-confined aquifer in the northern third of the Central Valley—the Sacramento Valley—and an upper unconfined to semi-confined aquifer overlying a lower confined aquifer, separated by a laterally continuous clay unit (the Corcoran Clay), in much of the southern two-thirds of the Central Valley—the San Joaquin Valley [3, 16]. Natural recharge in the Central Valley is currently presumed to occur through a combination of (1) direct infiltration from mountain-sourced rivers and streams, (2) lateral subsurface flow paths from mountain bedrock aquifers, and (3) infiltration of 'local' water from precipitation [17, 18]. However, the timing, volume, and mechanisms for natural recharge at local to regional scales are either poorly understood or available data/models for these processes are inadequate for monitoring efforts at sufficient spatial and temporal resolutions. Here, we explored the potential of utilising surface displacement observations from interferometric synthetic aperture radar (InSAR) to identify locations and timing of recharge at high spatial resolutions.

Recharge results in an increase of groundwater volume within the aquifer-system. The volume of groundwater stored in an aquifer is commonly referred to as storage. Any storage change in the sediments of an aquifer system can be expressed as a change in saturation (the filling/draining of pore space), or as poroelastic deformation due to a change in pore pressure, or as a combination of both [11]. In unconfined systems, storage changes are primarily expressed as changes in saturation with little to no associated poroelastic deformation [11]. Poroelastic deformation, which is governed by the law of effective stress [19], in unconsolidated sedimentary basins such as the San Joaquin Valley, is associated with the expansion or compaction of confined aquifer systems [4]. When pore pressure falls below its historical low, this type of deformation can result in permanent (inelastic) subsidence. Otherwise, deformation can be reversible (elastic), expressed as land surface subsidence or uplift. For semi-confined systems, both a change in saturation and poroelastic deformation can occur with the relative components depending on the degree of confinement. Adopting this framework, Lees and Knight [11] demonstrated the ability to monitor storage changes expressed through poroelastic deformation in the San Joaquin Valley by using surface displacements measured by InSAR. Lees and Knight [11] substantiate the potential of geodetic methods for improving our understanding of recharge.

Neely et al [20] found patterns of seasonal surface displacements, using InSAR time series, in the San Joaquin Valley—with particular attention given to the region between the San Joaquin River and Deer Creek—that they interpreted as revealing recharge. Specifically, they found that the estimated timing to peak amplitude, conceptually shown in figure 1, during the wet 2017 water year (WY2017; 1 October 2016–30 September 2017) displayed patterns generally characterised by earlier uplift near the eastern valley edge with progressively later uplift towards the valley center; we refer to this as the 'InSAR signature'. Though they acknowledge that storage-related deformation will result from a combination of both anthropogenic (e.g. groundwater pumping and managed recharge) and natural (e.g. precipitation) sources, the timing and spatial extents of these patterns were most consistent with river-sourced recharge as the driving component during a wet year. They suggested that this InSAR signature could be used to map recharge from Sierra Nevada sourced rivers and streams, and inferred that these recharge pathways would likely be associated with sufficiently confined permeable sediments. Without further corroborating hydrogeological information, the claims in Neely et al [20] remained largely speculative.

Figure 1. Refer to the following caption and surrounding text.

Figure 1. Example of seasonal decomposition applied in Neely et al [20] and here. These diagrams show an idealised time series of deformation at a single location (pixel) which can be modeled by a linear annual term () and a sinusoidal seasonal term () where is time, is the annual rate, is the peak amplitude, and is timing to peak amplitude in a given year.

Standard image High-resolution image

In this study, we built upon the idea that the characteristic InSAR signature could observe river-sourced recharge, with the timing to peak amplitude mapping the path of the pressure pulse front associated with the change in storage expressed by poroelastic deformation. We first replicated the analysis performed by Neely et al [20], using updated InSAR data, for the study region between the San Joaquin River and Deer Creek shown in figure 2(a). Using resistivity data from airborne electromagnetic (AEM) surveys, we derived information about the subsurface hydrogeologic structure and sediment type. Multiple AEM studies in the Central Valley have suggested that higher resistivity regions interpreted as coarse-grained dominated sediments provide potential recharge pathways [2123]. Using both InSAR and AEM, we tested the following hypothesis: There will be an InSAR signature due to river-sourced recharge only if there is high clay content associated with recharge pathways. While recharge pathways need to be permeable, clay is essential for generating detectable poroelastic deformation in the InSAR data as (1) the presence of clay can create a level of confinement in the aquifer system sufficient to produce large pore pressure changes in response to recharge and (2) the material properties of clay itself—higher specific storage values than coarse-grained sediments—allow for a greater poroelastic deformation response [24].

Figure 2. Refer to the following caption and surrounding text.

Figure 2. Maps of our study region and data coverage over the San Joaquin Valley. Panel (a) shows the study region between the San Joaquin River and Deer Creek with notable rivers and streams represented by light blue lines. Gray circles mark the locations of Fresno, California (to the north) and Visalia, California (to the south). Our two study sites are indicated by the black ellipses and referred to as the Fresno site and Visalia site. Panel (b) shows the average rate of displacement (; mm yr−1) over 1 October 2015–30 September 2022 observed by InSAR. Yellow, orange, red, and purple values indicate subsidence while blue values indicate uplift. Black dots correspond to AEM sounding locations. The red outline represents the extent of our 3D resistivity model.

Standard image High-resolution image

We focused on two study sites for this analysis; the locations are shown in figure 2(a). While local precipitation has been recognised as a dominant component to groundwater recharge in the San Joaquin Valley, geochemical analysis of groundwater near the Kings and Kaweah rivers, the major rivers associated with our study sites, showed distinct isotopic signatures, interpreted as river-sourced in origin [18]. These spatial patterns of the isotopic signatures are described as plume-like rather than broadly mirroring the smoothly varying precipitation patterns across the region [18]. Further, increases in groundwater levels [25] at these sites allowed us to assume that river-sourced recharge was occurring during wet years. However, there was an InSAR signature at one site but not at the other. Both sites were the focus of dense AEM surveys [2123] enabling us to examine the underlying hydrogeological structure and sediment type critical for testing our hypothesis. One was near Fresno, California where the Kings River is considered one of the primary local water sources for groundwater replenishment [26], but no InSAR signature was observed previously [20]. The other was the site near Visalia, California, located between Cross Creek and Mills Creek which connect to the Kaweah River during wet years, where Neely et al [20] observed the most prominent InSAR signature attributed to river-sourced recharge.

3. Data

3.1. InSAR displacement observations

InSAR time series of vertical displacement over our study region (figure 2(b)), were generated by TRE Altamira and made publicly available by the California Department of Water Resources [27]. We used the vertical displacement time series provided as input into our InSAR time series decomposition (section 4.1). While a full description of processing and calibration/validation workflows are outlined in supporting data documentation [27], we provide a summary here. Data spanning 1 January 2015–1 January 2023 from the European Space Agency's Sentinel-1 missions was processed using the SqueeSAR method [28] over groundwater basins in California. The resulting line-of-sight (LOS) displacement estimates were then calibrated using a network of 536 continuous Global Navigation Satellite Systems (GNSS) stations. Projecting the 3D GNSS measurements into the satellite 1-D LOS, the average displacement velocity (linear trend) of each station was compared against the displacement velocity of the average time series determined from InSAR pixels within a 100 m radius. At regional scales, pairwise differences between GNSS and InSAR velocities were used to estimate and remove a first order plane from all InSAR pixels. InSAR time series were then shifted by a common time series of residuals, with respect to the GNSS LOS time series, to force the relative InSAR measurements into the absolute reference frame of the GNSS network. This has the additional benefit of minimising non-geophysical signals in the time series, such as tropospheric delays or changes in scattering properties [29, 30]. To create a consistent dataset from multiple orbital geometries and tracks, full resolution InSAR LOS results were subsampled to a common spatial grid (100 × 100 m). Similarly, the time series for each pixel on the common spatial grid are interpolated and resampled to a common temporal grid with an approximate 7 d interval, resulting in 481 measurement times. Where multiple orbital geometries (ascending and descending) were available, these re-gridded LOS results were combined to estimate the component of vertical displacement. Otherwise, the LOS results were projected into vertical direction using the satellite look angle.

Using an independent set of 181 continuous GNSS stations, the Root Mean Square Error (RMSE) of time series for each station and calibrated InSAR dataset resulted in a state-wide RMSE value of 8.86 mm with individual stations ranging from 29.37 to 1.33 mm. The quantitative validation analysis determined a vertical accuracy assessment of 18 mm (95% confidence interval) for the InSAR time series [28].

3.2. AEM resistivity data

AEM data were previously acquired in the study region using SkyTEM 508, 312 and 304 systems [2123, 31], and inverted to obtain electrical resistivity models of the subsurface. Shown in figure 2(b) are the flightlines along which the AEM data, corresponding to the resistivity models used in this study, were acquired. While data were collected continuously, measurements were stacked at ∼35 m-intervals (referred to as AEM soundings) along the flightlines.

For all three previous surveys, a spatially constrained inversion was used to obtain models of the subsurface resistivity from the AEM data. This approach prioritised a resistivity model that fit the data while favoring smooth transitions of resistivity values in vertical and lateral directions [22, 32]. A resistivity model is composed of 1D resistivity profiles at each AEM sounding. While measurement errors and inversion decisions introduce uncertainty in the resistivity profiles, the penalisation of resistivity variations not supported by the data during the inversion process provide strong confidence in the large-scale structure imaged [22]. These 1D resistivity profiles consisted of 30–40 vertical cells. The thickness of the cells ranged from 1–3 m at the surface, then increased with depth with a constant geometric factor (e.g. 1.07) to compensate for the degrading resolution of the AEM data [33]. The depth of investigation (DOI), which describes the maximum depth at which there is measurable change in the AEM data due to changes in subsurface resistivity, was on average 300 m below ground surface (mbgs), ranging from 150–400 mbgs, in the three surveys.

4. Methods

4.1. InSAR time series decomposition

To analyse the InSAR data, we adopted the curve fitting approach of Neely et al [20] using linear and sinusoidal terms to approximate the InSAR time series:

where is the modeled displacement time series (mm), is time (year), is the linear term estimating the rate of displacement (mm yr−1), is the peak amplitude (mm) for the seasonal term, is the timing to peak amplitude for the seasonal term (year, where = 0 represents a timing in the peak amplitude term on 1 October, the start of the California water year and suggesting groundwater storage peaks around this time), and is a constant displacement shift (mm). Terms in equation (1), conceptually shown in figure 1, were estimated using least-squares minimisation techniques applied to the vertical displacement time series and computed for each water year. The average to the year-by-year standard deviation to peak amplitude and timing to peak amplitude are shown in figure S1. Examples of model fits to the time series are shown in figures S2–6.

Neely et al [20] showed that these terms differed between a dry year (WY2016) and a wet year (WY2017), with the determination of dry vs. wet based on drought conditions and water availability relative to the average [34, 35]. We applied the same methodology to observe the behavior on a year-by-year basis. While the previous study [20] was limited to two water years, we assessed the consistency of terms for 7 years. These include 5 dry years (WY2016, WY2018, WY2020–2022) and 2 wet years (WY2017 and WY2019).

To relate displacement values to recharge via storage changes e.g. [11], we assumed InSAR is observing primarily poroelastic deformation in our study region [4, 3638]. As pixels with low peak amplitude values are expected to give higher uncertainties for the timing to peak amplitudes, we excluded values for pixels that had values less than the accuracy of the InSAR data (18 mm peak-to-trough) [27]. For visualisation purposes, we filled in regions of missing data using a median filter (5 km diameter) on our derived maps of rate of displacement, peak amplitude, and timing to peak amplitude.

4.2. 3D resistivity model

We used the available 1D resistivity profiles to create a laterally continuous 3D model of subsurface resistivity for the study region. We first aligned each sounding to a common reference frame by vertically shifting each by their ground elevation. As the 1D resistivity profiles have variable thickness of cells increasing with depth, we then resampled resistivity values at 5 m intervals in the vertical direction with depths below the DOI excluded from the analysis. A nearest-neighbor interpolation scheme with a 12 km search radius was applied to each vertical interval using Generic Mapping Tools [39] to create 400 m × 400 m resolution horizontal grids.

4.3. Relating resistivity to sediment type

To interpret sediment type from resistivity in our study region, we applied the method presented in Knight et al [21]. Using sediment type reported in driller's logs located within 100 m of a 1D resistivity profile (349 well logs), we determined the distribution of resistivity values that corresponded to coarse-grained dominated sediments and fine-grained dominated sediments (hereafter referred to as coarse- and fine-dominated sediments respectively). We limited this analysis to the saturated zone, using water table levels reported at the time of AEM data acquisition (10–20 mbgs at the Fresno site and 20–40 mbgs at the Visalia site). We presumed that the resistivity distributions in the relatively small region above the water table would broaden (due to the saturation state varying from dry to saturated) and shift to higher resistivity values as shown in a detailed study of the impact of saturation on AEM-derived resistivity values [40].

5. Results and discussion

5.1. InSAR year-to-year analyses

The InSAR time series decomposition using a curve fitting approach resulted in three parameter terms, (1) the average rate of displacement, , (2) the peak amplitude, , and (3) the timing to peak amplitude, , on a year-by-year basis. Although the focus of this study was to understand how InSAR data might capture river-sourced recharge during the wet years, we also analysed the signals observed in the dry years for completeness.

5.1.1. Dry water years

Our curve-fitting approach for the dry water years (WY2016, WY2018, WY2020–2022) displayed similar spatial patterns in the linear rate of displacement, peak amplitude, and timing to peak amplitude with respect to each other dry year. The rate of displacement term featured a large subsidence bowl with surface elevations falling more than 250 mm yr−1 (purple in figure 3 column 1). The magnitude and spatial pattern of subsidence, often attributed to groundwater pumping, was consistent with previous interferometric studies conducted over the San Joaquin Valley during drought [7, 20, 41, 42]. We found that in regions of higher subsidence rates, peak amplitude was generally higher (up to 80 mm peak-to-trough motion and indicated by yellow in figure 3 column 2). The timing to peak amplitude occurred between late winter (pink in figure 3 column 3) and early summer (yellow-orange in figure 3 column 3).

Figure 3. Refer to the following caption and surrounding text.
Standard image High-resolution image
Figure 3. Refer to the following caption and surrounding text.

Figure 3. Maps of seasonal curve fitting analysis across our study region for each water year and the 7 year span. Column 1 shows the linear rate, , of displacement (mm yr−1). Column 2 shows the peak amplitude, , of the sinusoidal fit (mm). Column 3 shows the timing to peak amplitude, , with low peak amplitude (<9 mm) pixels excluded. Water year is 1 October–30 September (e.g. WY2016: 1 October–30 September 2016). Rows labeled in red text are considered dry while rows labeled in blue text are determined to be wet. The last row labeled with black text shows results using the 7 full water years available. Our two study sites are indicated by the black ellipses with the Fresno site to the north and the Visalia site to the south.

Standard image High-resolution image

For these dry years, we did not observe the general spatial pattern in the timing—progressing from east to west—that had been interpreted as river-sourced recharge in the study by Neely et al [20]. Rather, the timing to peak amplitude for these years were in-phase with expected seasonal oscillations of groundwater storage which, on average, peaks around April for this region [16, 37, 43]. The expected seasonal oscillations are aligned with local precipitation primarily falling between winter through early spring [16, 43] and with lowest groundwater elevations due to groundwater pumping occurring in autumn. The overall drier conditions and increased reliance on groundwater, as evidenced by decreases in groundwater elevations [25], over the study site during these dry years limited our ability to assess InSAR signatures related to river-sourced recharge.

5.1.2. Wet water years

Given our hypothesis, we focused our analysis on the wet water years, when we expected to see river-sourced recharge. The wet water years (WY2017 and WY2019), when compared to each other and to the results in Neely et al [20] for WY2017, had consistent spatial patterns in the linear rate of displacement, peak amplitude, and timing to peak amplitude. The rate of displacement for both years showed broad subsidence with values up to 150 mm yr−1 (dark orange in figure 3 column 1). While these wet years exhibited increases in groundwater levels (upwards of 20 m for both study sites) [25], subsidence due to the residual compaction of clay, can occur on decadal to centurial time scales [44]. We found peak amplitude values up to 70 mm peak-to-trough (indicated by yellow-orange in figure 3 column 2) at various locations, including the Visalia site. In locations where peak amplitude values were greater than the reported accuracy of the InSAR displacement time series (18 mm), we had confidence in our ability to quantify timing to peak amplitude—the critical feature we wanted to isolate in the InSAR data. Further providing confidence in our analysis, the average standard deviation for the timing to peak amplitude term, shown in Figure S1b, was less than 30 d for our Visalia site. As the peak amplitude values at our Fresno site were less than 18 mm, there was effectively no InSAR signature at this site that could be considered as related to recharge.

The timing to peak amplitude at the Visalia site exhibited a distinct spatial pattern in both WY2017 and WY2019. Near the eastern edge of the valley, peak amplitudes occurred in late winter (pinks and reds in figure 3 column 3). Moving towards the valley center, peak amplitude occurred at progressively later times in the year (orange—blue in figure 3 column 3). Taking a transect across this feature, the progression of the timing to peak amplitude is exemplified in figure 4 with the time series and models shown in figures S2–6. This InSAR signature was focused around Cross Creek and Mills Creek; on the order of 10 s km in scale. This differs from the broad, regional behavior in the timing to peak amplitude observed in the dry years. Although the most prominent example of this spatial pattern in timing to peak amplitude was at the Visalia site, this pattern was present elsewhere in the study region; for example along Tule River and Deer Creek. Due to the absence of high density AEM resistivity data, we did not assess these other locations of interest.

Figure 4. Refer to the following caption and surrounding text.

Figure 4. Comparison of the timing to peak amplitude (left) at points across the WY2017 InSAR signal at the Visalia site (right). The left panel shows the normalised modeled time series of the seasonal component for 5 points. Time series are detrended and normalised by their amplitude. Point 1 is to the northeast relative to Point 5. Time series show a clear progression of the timing to peak amplitude with peak seasonal uplift occurring earlier for Point 1 (March) and Point 5 occuring later (June). These points coincided with the transect (black line) used for the cross section shown in figure 6(b). The base map shows the timing to peak amplitude for WY2017 over at the Visalia site (black ellipse).

Standard image High-resolution image

Given the wet conditions [34, 35], locally focused seasonal deformation patterns, proximity to rivers and streams, and evidence of increased groundwater levels in the area [25], we interpreted the InSAR signature near Visalia to be the expression of river-sourced recharge.

5.2. Resistivity ranges for sediment/rock types at the sites

The AEM resistivity values found at the study sites ranged from 3–1000 Ohm·m. The highest resistivity values (300–1000 Ohm·m) in the region were interpreted as the bedrock of the Sierra Nevada, as suggested by earlier studies [21, 23].

Using the relationship derived with the methodology presented in section 4.3, we found the distribution of resistivity values corresponding to fine-dominated sediments ranged from 11–18 Ohm·m while for coarse-dominated sediments, the distribution ranged from 25–50 Ohm·m, as shown in figure 5. Resistivity values below 18 Ohm·m were interpreted as clay-rich sediments, as clay (rather than silt) was the primary sediment type within the fine-dominated classification of the driller's logs; resistivity values above 25 Ohm·m were interpreted as sand and gravel. Resistivity values between 18 and 25 Ohm·m were interpreted as mixtures of clay, sand, and gravel, with the amount of sand and gravel increasing as resistivity increases.

Figure 5. Refer to the following caption and surrounding text.

Figure 5. Esitmated relationships between resistivity (from airborne electromagnetic data) and sediment type (in-situ driller's logs). Shown are the resistivity distributions for fine-dominated (blue) and coarse-dominated (red) sediments. For visualisation, the plot is limited in the x-axis to a resistivity range of 6–100 Ohm·m.

Standard image High-resolution image

5.3. Detailed analysis of study sites

At the Fresno site, the presence of the perennial King's River to the east is a known source for recharge [26]. Increases in local groundwater elevations during wet years [25] and previous isotope studies [18] supported our assumption that river-source recharge occurred. When looking at the resistivity structure of this study site, shown in figure 6(a), we interpreted three key features. On the eastern edge of the valley, the high resistivity values (red in colour) corresponded to the granitic bedrock of the Sierra Nevada. Moving into the valley, the base of the resistivity model showed a low resistivity package (blue—green in colour), interpreted as clay-rich sediments, with the top of the package reaching near the ground surface at the eastern edge and dipping to the greater depths to the west. This clay-rich sediment package would suggest low permeability and be unlikely to provide a recharge pathway into the deeper aquifer. In the remainder of the displayed section there is a package, extending from the ground surface to the top of the clay-rich sediments, dominated by the resistivity values corresponding to sand and gravel (orange—red in colour) which likely provide pathways for recharge. The absence of evidence for significant clay within this package is consistent with reported unconfined conditions in the aquifer system across the entirety of the Kings subbasin [26].

Figure 6. Refer to the following caption and surrounding text.

Figure 6. Resistivity model profiles across our two study sites. Panel (a) shows the resistivity structure (in Ohm·m) underlying the Fresno site with no observable InSAR signature. Red colours indicate high resistivity interpreted to be bedrock—when along the eastern edge of the valley—and sand and gravel. Blue and green colours indicate the low resistivity range interpreted to be clay-rich sediments. The colours of orange and yellow indicate a resistivity range the corresponds to coarse-dominated sediments but with the clay content increasing in moving from orange to yellow. Depth profiles of resistivity have a 50x vertical exaggeration. The surface map is the timing to peak amplitude for WY2017 shown in figure 3 overlayed on Google Earth imagery. Panel (b) shows the resistivity structure underlying the Visalia site with an observable InSAR signature. Interpretations of groundwater recharge pathways are displayed as dashed black arrows. Inset maps show the surface traces (black dashed lines) of the cross-sectional view for each study site. These profiles were selected to intersect with the locations of highest AEM flightline density within each study site.

Standard image High-resolution image

While there is evidence of river-sourced recharge in the Fresno site [18, 25], the low clay content associated with the permeable pathways limited the potential for poroelastic deformation to occur or be observed. This, combined with the absence of an InSAR signature, supported our hypothesis that clay is essential for observing an InSAR signature.

At the Visalia site, where an InSAR signature interpreted as being related to river-sourced recharge is present, there are the ephemeral Cross Creek and Mills Creek that connect to the Kaweah River system during wet years. As with the Fresno site, there is a distinct isotopic signature attributed to river-sourced recharge and local groundwater elevations rose during WY2017 and WY2019, supporting assumptions that river-sourced recharge occurred [18, 25]. The underlying resistivity model shown in figure 6(b) has a more complex structure than the Fresno site. Along the eastern edge, we observe high resistivity values corresponding to the Sierra Nevada bedrock (red in colour). This is juxtaposed by a low resistivity package (blue—green in colour) interpreted as clay-rich sediments. The top section of the resistivity model, overlying the clay-rich sediments, has resistivity values primarily in the range of ∼18–35 Ohm·m (yellow—orange in colour) with thickness of <10 m to the east and thicken to the west. In contrast with the interpreted sand and gravel package near the surface of the Fresno site, the top sections of the model at the Visalia site displayed lower resistivity values. We interpreted this as sand and gravel with relativity high amounts of interbedded clay. Intersecting this package in the west, and underlying the InSAR signature, is a ∼30 m thick lower resistivity unit (blue—green in colour) dipping the southwest with depth ranging from 50 to 150 m. We interpreted this unit as a clay-rich confining layer, likely corresponding to the Corcoran Clay.

We interpreted, shown on figure 6(b), the sand and gravel as recharge pathways for river-sourced water to migrate from the surface to greater depths in the aquifer system. The high clay content from interbedded clays and confining layer creates the semi-confined and confined conditions needed for observable poroelastic deformation responses. As such, an InSAR signature is observed during times of recharge. Our hypothesis of the essential role high clay content has in determining whether or not we can observe an InSAR signature due to river-sourced recharge is supported by the underlying hydrogeological structure and sediment type in the Visalia site.

5.4. Considerations of limitations and assumptions

This study explored the InSAR signature associated with river-sourced recharge in California's San Joaquin Valley. Interpretations of hydrogeologic structure and sediment type from AEM resistivity data supported claims in Neely et al [20] and our hypothesis necessitating there to be relatively high clay content associated with recharge pathways. Below are considerations for the extension this work to other aquifers systems and InSAR signatures.

We made the reasonable assumption, as demonstrated by a number of studies with focus in the San Joaquin Valley [4, 3638], that the vertical displacement measurements are predominately attributed to poroelastic deformation. This may not be the case for other aquifer systems. For example, elastic loading due to changes in the mass of ice and snow in the Sierra Nevada have been shown to be the largest process driving seasonal changes in surface displacement in the Sacramento Valley to the north [36]. Thus, InSAR signatures found in other aquifer systems, where poroelastic deformation is not the dominant driver of seasonal surface displacements, may not reflect recharge or may need to be interpreted differently.

This study followed Neely et al [20] by fitting the time series with a sinusoidal function with an annual period. While this is simple and easily implemented, it cannot provide a sense of the dominant oscillation for deformation. For example, a time series where only seasonal groundwater extraction occurred would still result in a model that suggests both seasonal subsidence and uplift of equal magnitudes. Thus, this modeling approach describes the relative seasonal displacements. However, increases in groundwater elevations for both study sites (up to 20 m) during the wet water years [25] as well as true uplift observed in the InSAR time series (shown in figures S2–6) provided confidence that recharge occurred rather than a halting in groundwater storage declines. The use of more sophisticated methods such as polynomial fitting and principal component decomposition may improve the estimation for the timing to peak amplitude and be useful for investigating non-seasonal recharge events. Further, we focused on addressing the question of the necessary hydrogeological conditions in order to observe an InSAR signature in the San Joaquin Valley. Future detailed analysis of these InSAR signatures with the other components of the InSAR time series, the annual rate and peak amplitude, may provide deeper insight into the driving forces of deformation and the underlying heterogeneity of the aquifer system.

This investigation was limited by the availability of dense AEM surveys. Extension to other aquifer systems would likely face similar constraints. Serendipitously, our two sites demonstrated ideal endmembers of clay content and confining conditions with the Visalia site coinciding with an InSAR signature. The work accomplished here provides the foundational framework for the interpretation of InSAR signatures attributed to recharge. Collection of additional AEM data over other regions and aquifer systems, with and without InSAR signatures, may prove invaluable for further understanding the hydrogeological controls on the surface expression of river-source recharge and recharge in general.

6. Conclusions

Surface displacement observations from InSAR provide a rich dataset, sensitive to storage changes expressed by poroelastic deformation. We have demonstrated that under the correct hydrogeological conditions, InSAR may be used to map and monitor river-sourced recharge.

In California, interferometric displacement data are becoming increasingly accessible. Data from the European Space Agency's Sentinel-1 missions and NASA's upcoming NISAR mission promise to provide a continuation of freely available data with sampling at relatively short time scales (1–2 weeks) for the foreseeable future. Recognising the importance of land subsidence monitoring and to aid groundwater sustainability plan development and implementation, InSAR displacement time series are provided by the California Department of Water Resources and updated quarterly as a part of the SGMA technical assistance program with users able to easily view the data through the SGMA data viewer [25, 27]. As more data become available and accessible, these displacement data may become an even more important in applications to groundwater science and management.

The focus of this study was evaluating the InSAR signature—and lack of signature, along the eastern edge of the San Joaquin Valley where recharge occurs with water delivered by the mountain-sourced rivers and streams. An obvious extension of this study is to use the two data types which we have employed to design field campaigns at locations throughout the world, with suitable hydrogeologic conditions, to map and monitor other forms of both natural and managed recharge. Determining the locations and timing of recharge is essential to improve the ability to predict the impact of climate change on groundwater aquifers.

Acknowledgments

This research was supported by funding to R Knight from the Gordon and Betty Moore Foundation; Grant Number GBMF6189.

Data availability statement

The data that support the findings of this study are openly available at the following URL/DOI: https://purl.stanford.edu/vj301kx0276 [45]; https://purl.stanford.edu/yc041mz9691 [46]; https://purl.stanford.edu/yf065fy8317 [47].

undefined