Decadal climate variability and the spatial organization of deep hydrological drought

Ana P Barros; Jared L Hodes; Malarvizhi Arulraj

doi:10.1088/1748-9326/aa81de

1. Introduction

Baseflow can be the sole source of streamflow during extended periods without precipitation (Priest 2004, Peters et al 2005, Van Loon and Laaha 2015), and thus it is an essential element of the hydro-ecological resilience of natural streams and rivers to meteorological drought. Besides precipitation variability, the space-time variability of baseflow reflects watershed topography and hydro-geomorphic properties of river networks (e.g. Vogel and Kroll 1992, Wolock et al 2004), as well as hydrogeological controls of recharge and groundwater connectivity from regional (deep subsurface flow) to local (shallow subsurface flow) scales (e.g. Fan 2015, Van Loon and Laaha 2015, Schaller and Fan 2009, Bloomfield et al 2009). Whereas deep baseflow is normally associated with slow groundwater flow, Brun and Barros (2014, hereafter BB14) reported evidence of remote groundwater connectivity in the US South Atlantic (USSA) region in the quick rise and large differences in water-well levels west and east of the Blue Ridge Front in response to heavy hurricane rainfall (their figure 16), which they attributed to fast preferential flow at regional scale through networks of connected fractures and karst terrain (e.g. Trapp and Horn 1997, Swain et al 2004, their table 1).

Figure 1. Refer to the following caption and surrounding text. — **Figure 1.** (a) USGS streamflow gauge sites (n = 126) with at least 54 years of data (1962–2015 water years), divided by physiographic regions as defined by USGS. Stations in different physiographic regions are shown in different colors. (b) Spatial distribution of Blue Ridge stations with station 42 marking the location of the Mills River USGS gauge. Note that the stations are numbered from 1 to 126 by longitude (west to east, see table S1).
Download figure:
Standard image High-resolution image

The objectives of this study are (1) to characterize the temporal variability of baseflow at scales relevant for detecting regional impacts of climate variability on hydrological drought (Faye and Mayer 1990, Trapp and Horn 1997, Lloyd and Lyke 1995, Kienzle 2006, Dai et al 2010, Krakauer and Temimi 2011), and (2) to elucidate the spatial organization of baseflow variations with regard to known physical properties and boundaries of groundwater systems in the USSA. The motivation is to demonstrate that baseflow variability captures nonlinear interactions among climate and hydrological processes mediated by topography and geology at regional scale, and it is therefore a robust indicator of severe (long-range) hydrological drought.

Groundwater basins are understood here as the hydraulic and hydrologic continua linking recharge areas to aquifers and rivers across geologic and topographic boundaries (figure S1(a)). Baseflow at streamgauge locations includes the subsurface flow from the shallow saturated zone and from deep regional-scale groundwater systems (Schaller and Fan 2009, their figure 1). Common metrics to quantify baseflow contributions to streamflow include the BFI (Baseflow Index, that is the ratio of the long-term mean baseflow over the total streamflow), and two hydrograph recession parameters: the shape parameter b that describes the changes in streamflow Q with time after rainfall stops $( - \frac{{{\rm d}Q}}{{{\rm d}t}} \propto Q^b )$ , and the time-scale τ relating Q and basin storage S $( - \frac{{{\rm d}Q}}{{{\rm d}S}} \propto \frac{{Q^b }}{\tau })$ (e.g. Krakauer and Temimi 2011, Ahiablame et al 2013, Beck et al 2013 among many others). The BFI captures the long-term mean spatial behavior of baseflow, but not temporal variability. Applying differential methods to long records consisting of many hydrographs is challenging due to complex signatures of surface-subsurface interactions, as well as nonlinear baseflow response and groundwater storage described by time-varying b and τ (e.g. White and Sloto 1990, Vogel and Kroll 1992, Szilagyi and Parlange 1998). Alternatively, digital filters have been successfully applied to extract baseflow systematically from streamflow records at long time-scales (Nathan and McMahon 1990, Arnold and Allen 1999, Eckhardt 2008, Murphy et al 2009). A two-step digital filter based on White and Sloto (1990, hereafter WS90) is used in this study to extract baseflow time-series from streamflow observations. The space-time variability of baseflow is subsequently examined using empirical orthogonal functions (EOF), wavelets and wavelet coherence analysis methods. The study is conducted in the USSA because of the wealth of available data, and because comprehensive understanding of groundwater systems in the region is necessary to examine the physical basis of predictive relationships between climate modes of variability and baseflow.

The manuscript is organized as follows. The study area and the streamflow data are described below. A brief description of the baseflow separation algorithm followed by a detailed demonstration are provided as part 1 of supplementary data (SD) available at stacks.iop.org/ERL/12/104005/mmedia. Part 2 of SD is an overview of time-series analysis methodology used in this study (EOF, wavelet and wavelet coherence analyses) including an illustrative application. Results are analyzed and interpreted, followed by conclusions and discussion.

2. Data

To examine how baseflow varies among physiographic regions, all United States Geological Survey (USGS) streamgauges in the region of study defined by (34N–37N, 74W–88W) were selected. Gauges with more than 3% missing data over their length of record (most missing data occurred in the early history of these gauges) were truncated from the beginning of record until missing data accounted for only 3% of the total data. Finally, 126 stations with record length greater than 20 000 days (54.7 water years—October through September) beginning on 10/1/1961 or earlier, and ending on 9/30/2015 or later, were retained (figure 1(a), table S1 in SD).

The selected streamgauges are distributed by five physiographic regions including the Coastal Plain (CP), Piedmont (PD), Blue Ridge (BR), Valley and Ridge (VR), and Appalachian Plateau (APL) Provinces (figure 1(a)). The VR is characterized by elongated mountain ridges with layered sedimentary rock intensively faulted and folded at high elevations, and easily eroded layers of limestone, dolomite, and shale in the valleys. In contrast, sandstone and limestone rocks in the APL are nearly flat-lying or gently tilted and warped (Heath 1980). The BR mountain belt consists mostly of igneous and metamorphic rocks with small areas of sedimentary rocks and high density of alluvial fans in the inner mountain valleys due to widespread landslide activity. The PD spans from the Appalachian Mountains (APM) foothills to the CP lowlands at the Fall Line east of which sand and limestone aquifers overlay crystalline rocks (figure 1 and figure 2). The shallow subsurface flow contribution as a fraction of total streamflow is expected to be greater in the steep and moderate topographic relief of the BR due to the rapid response of interflow (Faye and Mayer 1990, Priest 2004, Tao and Barros 2013). Recharge in the BR and PD is from precipitation that enters aquifers through porous regolith, and some of the water reaches bedrock and enters fractures in the crystalline rocks (Trapp and Horn 1997, Lloyd and Lyke 1995). Whereas fractures provide a faster mode of water transport especially in regions with steep topographic relief and for heavy rainfall events, long-term baseflow yield is very small, because deep storage capacity is low and such bedrock is less permeable (depending on the nature of the fracture system; Wendland et al 2009).

Figure 2. Refer to the following caption and surrounding text. — **Figure 2.** Schematic representation of surface-groundwater interactions in the Blue Ridge and Inner Piedmont (top panel, modified after figure 92 in Trapp and Horn 1997, credit: US Geological Survey), and in the Coastal Plain (bottom panel, modified after figure 2 in Trapp and Meisler 1992, credit: US Geological Survey).
Download figure:
Standard image High-resolution image

The conceptual illustrations of hydrogeological controls of surface-groundwater interactions depicted in figure 2 contrast the dominant modes of water circulation in the BR and inner PD (top panel) against the CP (bottom panel). In the inner PD (terrain elevations > 200 m), average lateral transmissivity estimates of 32 m day⁻² reported by Daniel et al (1997) correspond to 10–20 years travel times in the saturated regolith assuming uniform aquifer thickness and hydraulic conductivity over the region (figure 1(a), 150–300 km width). In the outer PD and CP, high-frequency (fast) baseflow response is explained by shallow semi-consolidated and unconsolidated sediments consisting of silt, clay, and sand, with some gravel and lignite, underlain by consolidated beds including limestone and sandstone (Trapp and Horn 1997, refer to figure 8 from King and Beikman 1974).

3. Results

3.1. Baseflow variability

The scale averaged variance of annual baseflow for each streamgauge station (e.g. SD part 1, figures S1−S5) was normalized by the station's maximum variance to produce values in the [0, 1] range for each and every station, and thus facilitate detection of temporal variability co-occurrence patterns. Figure 3(a) shows time-series of normalized 1–20 years mean variance for all 126 gauges organized by physiographic region (color code) and longitude (west to east, table S1) to highlight periods of concurrent high variance across USSA watersheds. Time-scales associated with statistically significant (95% confidence level) peaks in the global wavelet power (GWP) spectra are identified in figure 3(b). Variability across physiographic regions along the west-east transect in BB14 (figure S6(a) can be examined in SD (figure S6(b)).

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

Marked differences before and after the 1990–1995 period at BR stations (figures 3(a) and S6(b)) are in step with the transition from negative (cold) to positive (warm) phases of the AMO (Atlantic Multidecadal Oscillation, figure S6(c)). High variance from west to east in the mid 2000s is explained by increase in drought-busting land-falling hurricanes and tropical cyclones (TCs), and in particular two major storms (Ivan and Frances, BB14) that crossed the Southern Appalachians in 2004 and propagated northeastward along the APM's topographic crest (figures S7(a)–(b)). Wet years in 2004 and 2005 (figure S7(b)) were followed by a period of severe drought in 2006–2007 Accordingly, every station other than the 74–93 block shows increased variability (i.e. higher variance) in the early 2000s centered around 2004–2005. Nonconforming stations (74–93) are in the rainfall gaps between TC tracks in 2004 (figure S7(a)). Baseflow and precipitation are tightly linked for BR stations 45 and under to the south as shallow upland soils transition to deep alluvial deposits (and high soil moisture storage) in the inner mountain valleys. To the north, however, the normalized variance of the easternmost BR stations (figure 1(b), 50–71 in table S1) is closer to PD stations (figure 3(b)) due to smoother topography and different geology. Note the similar evolution of the 1–20 years scale averaged variance for stations 43–74 in the VR, BR, and south PD regions. These stations are distributed along a summit-to-sea transect that runs from the APM to the Fall Line with deep surficial aquifers to the south and shallow surficial aquifers to the north and northeast (Heath 1980, Aucott 1996, Campbell et al 2010), and this variance behavior is consistent with long-range groundwater connectivity (i.e. regional inter-basin flow) through a complex network of fractures in the crystalline rock aquifer (e.g. see figures 36, 37 and 38 in Heath 1980).

Physiographic and geologic controls become more evident by contrasting the time-scales of baseflow variability at CP and PD stations (figure 3(b)). CP stations exhibit common behavior with widespread variability peaking in the 1960–1970s, or in the late 1990s–early 2000s; most PD stations exhibit strict bimodal variance with peaks in the 1970s and 2000s. This difference can be attributed to distinct subsurface responses in the same climate regime (e.g. AMO warm phase) linked to basin hydrogeology that impact recharge and surface-groundwater interactions. The baseflow fraction of streamflow at low elevations in the PD and CP is small on average, though high rainfall in wet periods (1960s and 1970s) results in significant enhancement of baseflow and total streamflow (figure 3(b)). The shift in baseflow variance in the mid 1970s (after the minimum of the AMO cold phase) is concurrent with macroscale shifts in precipitation variability in the Southeast US (Wang et al 2010, and Li et al 2011).

Baseflow variability in the PD and inner CP shows multiple peaks in the 10–25 years range suggesting multiscale complex interactions between climate forcing and the regional groundwater system (figures 3(a)–(b)). Only the outer CP stations, where a significant fraction of the baseflow originates from shallow surficial aquifers (figure 2), show a strong annual cycle as well as peaks at 3–5 year scales (independently of basin area). Regional flow is constrained by confining beds of reduced lateral hydraulic conductivity in the Coastal Plain (figure 2), thus effectively separating high-frequency (local) response (< 6 years) from regional-scale decadal (long-range) variability (figure 3(b)). The impact of dam construction is unequivocal (e.g. stations 14, 16, 19, 65, 67, 100 and 110) in reducing the amplitude of floods and smearing the variance (figure 3(a)).

PL, VR, and BR stations exhibit strong variability at multi-decadal scales (10–30 years) (figure 3(b)). The variance peak at the 5–7 year scales for BR stations only and in the inner mountain region alone is attributed to easterly low level convergence of moist air from land-atmosphere interactions in the PD and CP to the southeast that is key to the local diurnal cycle of rainfall (Wilson and Barros 2017 and Wilson and Barros 2015).

Figure 4. Refer to the following caption and surrounding text. — **Figure 4.** Right column: EOF analysis of normalized baseflow of all stations after removal of large basins and stations affected by dams. Spatial patterns are presented in (a) EOF1; (c) EOF2; and (e) EOF3, where the radius of the circle is proportional to the local EOF value (black positive, blue negative). Fraction of total variance explained by each principal component is indicated in the title of the corresponding spatial patterns. Left column: wavelet analysis of Principal Components. Temporal patterns are presented in (b) PC1; (d) PC2; and (f) PC3. Thin black lines isolate the cone of influence unaffected by edge effects due to the limited record length. Thick black lines indicate the 95% confidence level.
Download figure:
Standard image High-resolution image

Figure 5. Refer to the following caption and surrounding text. — **Figure 5.** Wavelet coherence (WTC) of the first three PCs of the normalized baseflow and AMO (top row) and Niño 3.4 (middle row) at annual time-scales: (a) PC1-AMO; (b) PC2-AMO; PC3-AMO; (d) PC1 N iño 3.4; (e) PC2- Niño 3.4; (f) PC3- Niño 3.4; and WTC of AMO and Niño 3.4 at annual (g) and seasonal scales, respectively (h) December-January-February (DJF) and (i) June-July-August and September (JJAS). Thin black lines mark the cone of influence unaffected by edge effects due to the limited record length. Thick black lines indicate the 95% confidence level. Right pointing arrows indicate that the time-series are in phase; left pointing arrows indicate out of phase; downward pointing arrows mean that the first time series leads the second by 90 degrees; and finally, upward pointing arrows indicate that the second time series leads the first.
Download figure:
Standard image High-resolution image

3.2. Principal components of baseflow

Stations at the outlet of watersheds that exceed the small and medium-size criterion (area < 5000 km²) and stations where streamflow variability is modified by dam operations upstream were removed for EOF analysis (marked in table S1). The first three spatial EOFs (see Methods in SD−part 2) of the annual normalized baseflow after removing stations affected by dams and large basins are shown in figures 4(a), (c), (e), while figures 4(b), (d), (f) show results from wavelet analysis of PC1, PC2 and PC3. PC1 follows closely the time-series of 2 year precipitation anomalies (figure S8, r² ~ 50%), and EOF1 (figure 4(a)) values are of the same sign (positive) with larger values on the BR and close to the topographic divide (APM crest) consistent with the regional climatology of orographic precipitation. In contrast, EOF2 and EOF3 show spatial organization according to regional topography and geology (Trapp and Horn 1997, Hack 1982): EOF2 (figure 4(c)) separates regions underlain by limestone and sandstone aquifers³ west of the BR (positive coefficients) from regions underlain by crystalline rock aquifers to the east (negative coefficients), and EOF3 (figure 4(e)) isolates complex topography along the Blue Ridge front and foothills and in the inner PD (positive coefficients; figure S9 and figure S12) from low elevation regions to the east and west underlain by surficial unconsolidated aquifers (negative coefficients). Statistically significant non-stationarity is present in PC1 with intense variability at 7–12 years from the mid-1990s onward (figure 4(b)). PC2 and PC3 display statistically significant high-frequency variability, respectively in the 1–2 (inter-annual) and 2–5 years intervals (figures 4(d), (f)). PC2 also displays intense variability at 8–12 years concurrent with PC1, whereas weak non-stationarity at long time-scales (centered around 16 years), albeit present, is not statistically significant in the case of PC3.

Attribution of baseflow variability to climate teleconnections requires an integrative hypothesis linking precipitation and hydrologic processes, specifically recharge and groundwater dynamics. The impact of both ENSO and AMO on North American precipitation, and thus streamflow, is well documented in the literature (Enfield et al 2001, Tootle et al 2005, Dong et al 2006, Trenberth and Shea 2006, Mo et al 2009, Hu and Feng 2012, Nogueira et al 2013). Specifically, warmer Atlantic temperatures during the positive phase of the AMO lead to an increased number of tropical storms and hurricanes (McCabe et al 2004, Curtis 2008, Maxwell et al 2013, Labotka et al 2015), although the co-occurrence of warm AMO and cold Pacific Decadal Oscillation (PDO) has been linked to drought (McCabe et al 2004, Curtis 2008, Mo et al 2009). Wavelet analysis of the AMO SST time-series since 1856 exhibits non-stationarity switching from positive to negative SST phases at multi-decadal time-scales interspersed with decadal scale variability trends (peaks of the global wavelet spectrum at approximately 62–70 years, 26–33 and 9.3 years, not shown). The beginning date of the period of analysis in this study (1961) coincides with a transition from a positive to a negative phase of the AMO with a declining trend (cooling) till the mid 1970s followed by a positive trend (warming) after 1980 through the next phase transition to a positive phase in the early 1990s (figure S6(c)).

The impact of ENSO in SE US climate is stronger in the winter: in the positive phase of ENSO (El Niño), winters tend to be cool and wet; during the negative phase (La Niña) winters tend to be warm and dry. SE summers during strong to very strong El Niño events are generally warm and dry with few Atlantic hurricanes. Hurricane activity peaked in the 2000s (2004 and 2005) during weak El Niño activity. Li et al (2011 and 2013) did not find strong dependence of summertime precipitation variability on ENSO in contrast with previous studies (e.g. McCabe et al 2004) likely due to relying on regional-scale averages that removed physiographic contrasts from their analysis.

From 1961 through 2015, PC2 trends positively at the same rate as the warming trend of the AMO (not shown), whereas PC3 mirrors the AMO over the same period. There is significant variability of streamflow and subsurface flow components at time-scales between 3 and 8 years linked to ENSO. PC2 and PC3 are in phase with ENSO before 1993 (using SST anomalies over the Niño 3.4 region as the ENSO index, hereafter referred to as Niño 3.4), and intermittently in- and out-of-phase thereafter. Wavelet decomposition of PC2 reveals non-stationarity at 2–4 and 8–10 years in the 1970s followed by a break in the 1980s and intensification at time-scales from 1–2 and from 8–12 years in the 1990s and early 2000s. The changes in inter-annual variability from the 1970s to the 1990s in PC2, and the amplification at 2–5 years scales in PC3, are in line with the increase in the variability of regional-scale summertime (June-July-August, JJA) precipitation in the Southeast US identified by Li et al (2013).

3.3. Wavelet coherence analysis

Wavelet coherence analysis (WTC) was conducted at annual and seasonal scales using water-year and DJF (December-January-February) and JJAS (June-July-August-September) AMO and Niño 3.4 averages for each year. WTC results of PC1, PC2 and PC3 with the AMO are shown in figures 5(a)−(c) and figure S9, and with Niño3.4 in figures 5(d)−(f) and figure S10. In particular, the arrows in figure 5 (and figures. S9 and S10) indicate that the AMO and Niño 3.4 lead PC2 and PC3 at short time-scales (< 6 years).

PC1 and PC2, the two principal components linked to regional-scale APM topography and broad hydrogeological controls (limestone and sandstone aquifers westward of the BR, fractured crystalline aquifer eastward), both exhibit statistically significant multi-year and decadal scale variability after 1990. Non-stationarity before and after 1990 is consistent with Atlantic climate dynamics in the transition from negative to positive phases of the North Atlantic Oscillation (NAO) ahead of the AMO (e.g. McCarthy et al 2015). The WTCs of PC1 and PC2 with the AMO (figures 5(a) and (b), S9(a) and (b)) and ENSO (Niño 3.4, figures 5(d) and (e), S10(a) and (b)) indicate that PC1 and PC2 variability post 1990 is slightly behind and, or in phase with the AMO at 6–8 years and with ENSO at 8–12 years This close relationship on multi-year time-scales across the region of study suggests deep (and remote) response to changes in precipitation and recharge, and thus strong surface-groundwater interactions and long-range groundwater connectivity (i.e. regional-scale inter-basin flow). PC3 does not show statistically significant co-variability at decadal scales with either AMO or ENSO (figures 5(c) and (f)); this changes at sub-decadal scales as discussed next.

Figure 6. Refer to the following caption and surrounding text. — **Figure 6.** Mean Baseflow Fraction (MBF) at each of the 126 gauge stations according to physiographic region organized by longitude (a) and streamgauge elevation (b). Note the contrast between MBF's spatial organization (shaded clusters) as a function of longitude east of the continental divide (left panel) and as a function of elevation (right panel). CP-Coastal Plain (blue); OPD−Outer Piedmont (light green); IPD−Inner Piedmont (dark green); APM—Appalachian Mountains (maroon).
Download figure:
Standard image High-resolution image

Indeed, at short time-scales (<6 years), WTC analysis reveals complex seasonal non-stationary phase relationships over the last four decades: (1) AMO leads PC2 and PC3 throughout the period of record; (2) ENSO leads PC3 in JJAS in the 1980s; and (3) ENSO leads PC2 in DJF in the 1990s (figure S10(e)). Such behavior emerges from nonlinear interactions among ENSO and AMO processes captured by WTC analysis of AMO with Niño 3.4 on an annual (figure 5(g)) and seasonal basis (figure 5(h), DJF; JJAS, figure 5(i)). The WTC results show that the magnitude and phase of the AMO-ENSO relationship for DJF vis-a-vis JJAS flips in 1990, which is depicted by PC2 for the AMO (figure S9) and by PC3 for ENSO (figure S10). In the context of severe drought, this behavior is indicative of regional susceptibility to inter-annual climate variability, and in particular frontal systems in the cold season (DJF) and TC frequency and intensity in the warm season (JJAS).

In summary, PC1 and PC2 capture decadal-scale variability in baseflow driven by precipitation changes connected to the AMO and ENSO at the same scales. High variability in PC2 and PC3 in the 2–6 years range is attributed to nonlinear interactions between the AMO and ENSO that in turn impact seasonal precipitation. Ultimately, baseflow variability proper depends on the overlap of the spatial patterns of precipitation and recharge areas in the context of regional hydrogeology.

3.4. Spatial EOF patterns and hydrogeology

The spatial organization of the annual mean baseflow fraction (MBF) of streamflow is examined in figures 6(a) and (b) with respect to longitude and elevation, respectively. In the context of EOF2, PL, VR, and BR stations are separated from the PD and the CP exposing two distinct modes of groundwater connectivity respectively west and east of the topographic divide as discussed earlier (sections 3.1 and 3.2, figure 4(c)). A secondary transition of EOF2 coefficients from positive to negative occurs in the outer PD in North Carolina (~82 W, figure 4(c)). This sub-region is characterized by tapering off topography in the transition zone with significant changes in the depth of the Coastal Plain surficial aquifer (deeper to the south and shallower north and eastward), and in the groundwater flow direction from easterly or perpendicular to the coast in the south-central CP to southeasterly or oblique to the coast northward (Campbell et al 2010, Galleci and Lautier 2010, Aucott 1996, Heath 1980). In the context of EOF3, the spatial organization of the MBF conforms to the hydrogeological controls elucidated by EOF3 at sub-decadal scale (figure 2, figure 4(e), and figure S11): high MBF in complex terrain watersheds with shallow soils (and small water storage capacity) at high elevations; low MBF at low elevations elsewhere in regions with surficial aquifers west (alluvium in the valleys) and east (unconsolidated sediments) of the APM.

The coefficient of variation of MBF within specific physiographic regions varies between 0.5 and 1.0 as a function of longitude (figure 6(a)) reflecting USSA topography and geology heterogeneity. Overlaps among stations in different physiographic regions in figure 6(a) are an artifact of latitudinal changes in the orientation and lay-out of key landscape features (e.g. APM swerving north and eastward). The spatial organization of MBF with elevation (figure 6(b), figure S12) is less ambiguous with rare overlaps from summit (BR, east of the APM topographic divide) to sea (CP, Atlantic Ocean) with well delimited boundaries among physiographic regions. Note the significant differences between the inner and outer Piedmont (IPD and OPD, respectively) reflecting the transition from complex terrain in the IPD to the largely smooth topography and emergent unconsolidated surficial aquifers in the OPD along the Fall Line (figure 2, figure S12). Nevertheless, because the MBF is the average of the annual baseflow fraction over the period of record, it does not capture the high values identified by BB14 in the AP and VR wells in response to extreme rainfall events.

Finally, the outcome of systematic EOF analysis of annual baseflow normalized by the local maxima replicates regional traits reported previously: EOF2 is in agreement with the hydrologic-landscape regions delineated by Wolock et al (2004), and EOF3 is consistent with the spatial BFI distribution derived by Wolock (2003), including the minimum in the outer Piedmont where sandy soils are dominant. The spatial patterns of EOF2 and EOF3 capture therefore long-range and deep (decadal scale) and local (seasonal to inter-annual scales) variability in baseflow dynamics tied to the variability of warm season (JJAS) precipitation, respectively and the regional hydrogeology.

4. Conclusions and discussion

High spatial variability of deep hydrological drought, here defined by decadal-scale decreases in regional groundwater flow contributions to streamflow, is apparent from the analysis of historic water levels in monitoring wells, but it is difficult to model because hydrogeological features and connectivity associated with specific aquifers are not known generally (e.g. Campbell et al 2010). Combining information based on physiographic region including climate, topography and underlying geology is necessary to interpret long-range baseflow response, and to examine nonlinear interactions among climate forcing and hydrologic drought (Bachman et al 1998, Bloomfield et al 2009, Lowman and Barros 2016). Here, baseflow was first extracted using a two-pass filter algorithm applied to streamflow records from gauges at the outlet of small and medium-size USSA watersheds (area < 5000 km²), followed by space-time modal analysis using EOF and wavelet methods.

EOF analysis of annual and seasonal baseflow normalized by the local maxima allowed us to identify large-scale climate (AMO and ENSO), topography (slope, elevation, landform), and hydrogeology controls of baseflow (e.g. aquifer type; regolith structure; basement boundaries). EOF1 and PC1 capture the decadal- to multi-decadal scale variability of the annual cycle of precipitation including orographic enhancement in the Blue Ridge; EOF2 and PC2 capture multi-year and decadal scale climate variability east and west of the topographic divide that can be linked to the AMO, and separates the regions underlain by basement crystalline rock with open fractures to the east from regions underlain by sandstone and limestone aquifers to west; EOF 3 and PC3 capture the impact of inter-annual and multi-year variability of ENSO and nonlinear interactions between ENSO and the AMO contrasting low MBF fraction from surficial unconsolidated aquifers in the CP and OPD regions against high MBF fraction in the complex terrain of the BR and IPD.

Inter-basin transfers and the effective areal extent of the groundwater source basin that is associated with a given watershed vary from one year to the next according to the space-time variability of precipitation forcing, and thus regional recharge. Baseflow exhibits nonstationary behavior at temporal scales ranging from 1–20 years, with most variability associated with extremely wet and, or dry periods at inter-annual and sub-decadal time-scales (< 6 years) attributed to interactions among multiscale climate phenomena (i.e. ENSO and AMO) governing warm season precipitation east of the APM. Specifically, inter-annual and sub-decadal variability becomes significant only after the 1990–1995 transition from the cool to the warm phase of the AMO, and decadal and multi-decadal variability intensifies with the warming trend of the AMO after 1980. Together with nonstationary ENSO behavior (e.g. Wolter and Timlin 2011, their figure 3), this association with the AMO provides insight into long-range (near century-scale) variations of baseflow (persistent dry and wet conditions) and the historical record of deep (severe) droughts⁴ (Ting et al 2009, Enfield and Cid-Serrano 2010).

Long-range streamflow memory has been explained by nonstationary climate forcing and by the cumulative spatial aggregation of streamflow through large-scale river networks (e.g. Potter 1976, Tootle et al 2005, Mudelsee 2007). A contribution of this work is to explain long-range streamflow memory as a multiscale process chain that propagates nonstationary climate forcing to regional groundwater systems and deep baseflow linked to decadal and multidecadal hydrological drought in small and medium-size watersheds. Further, interpretation of the findings of this study in the light of existing knowledge of USSA groundwater systems indicates that there is a physical basis for long-range predictability of hydrological drought from baseflow without region-specific constraints. Therefore, the data-driven methodology employed here can be systematically applied to streamflow observations anywhere. Because slow and intermediate components of baseflow include inter-basin groundwater flow, this framework also facilitates indirect monitoring of relative changes in natural groundwater resources at the watershed scale using streamflow observations, which are more widely available worldwide than groundwater wells.

Despite incomplete understanding of the mechanisms of ocean-atmosphere interactions associated with the AMO (e.g. Buckley and Marshall 2016), links between Atlantic basin hurricane activity, specifically the number of Tropical Cyclones and AMO multi-decadal variability, are well established (Caron et al 2015, Enfield et al 2001). Note that nonlinear impacts of increased precipitation from singular events such as intensification of tropical cyclones linked to anthropogenic forcing (Sobel et al 2016) on regional groundwater flow are uncertain, and should be largely tied to the specific storm tracks. Watersheds in the USSA region are particularly vulnerable to meteorological drought during the cold phase of the AMO because of the progressive reduction of the slow components of baseflow into streams and rivers, eventually resulting in the lack of natural capacity to meet minimum hydro-ecological instream flow requirements, let alone water supply needs (e.g. Petts 2009). Presently, the cooling trend of the AMO signals the next transition from warm to cold phase in the upcoming decade (McCarthy et al 2015), thus establishing conditions favorable to intensify deep hydrologic drought in the 2030s, albeit this impact should be modulated at local and intermediate scales by ENSO-related inter-annual variability. One immediate benefit of this predictive understanding is to increase the lead-times for developing robust and spatially adapted drought preparedness and mitigation plans.

Acknowledgments

This research was supported with funding from NSF grant CNH-1313799. We are thankful to two anonymous reviewers for their constructive comments and suggestions.

Author e-mails

Author affiliations

Dates

Peer review information