Continental-scale prediction of hydrologic signatures and processes

Understanding how dominant hydrologic processes and their drivers vary across diverse continental-scale landscapes is critical for hydrologic modeling and water management applications. Our research addresses this question by synthesizing large-sample watershed datasets, Caravan and GAGES-II, and developing random forest models to identify patterns in hydrologic function. We assessed dominant processes by examining hydrologic signatures – summary indicators of watershed function derived from hydroclimatic time series and random forest models across 14 146 gauged United States watersheds. The results reveal clear continental-scale gradients in hydrologic processes, including baseflow, overland flow, storage, and water balance losses. Our map of dominant processes highlights, for example, the transition from baseflow to fast responses and back to baseflow along the elevation gradient from the Appalachian spine, through the Piedmont, to the Eastern Coastal Plain; a distinct outer ring around the Great Lakes region; and sharp contrasts between coastal and inland processes in the West. Variable importance analysis from random forest models show that processes in the western U.S. are primarily controlled by climate, whereas in the eastern U.S., soil, geology, and topography play larger roles, with distinct human influences apparent in urban areas. Our approach of estimating dominant processes and their drivers facilitates extending process knowledge from research watersheds to the continental scale, assessing current hydrological understanding, and evaluating hydrological model structures.

1 Introduction

Estimating the contributions of different hydrologic processes to streamflow generation at a continental scale is essential for flood forecasting and water resources management. Optimal management strategies, including the design of grey and green infrastructure, differ depending on which processes dominate hydrological response (Oswald et al., 2023; Thompson et al., 2020), which vary substantially by regional environmental conditions (Blöschl, 2006; Paola et al., 2006; Penna, 2024). Understanding how water is partitioned, stored, and transported through different parts of the terrestrial systems is a fundamental question in the hydrologic sciences (Brooks et al., 2015). To simulate a diverse set of processes at large-scale, a new generation of hydrologic models with flexible and heterogeneous structures has emerged (Clark et al., 2015; Frame et al., 2025; Johnson et al., 2023). However, despite these technological advances, we still lack an estimate of dominant hydrologic processes controlling streamflow generation at continental scales (McMillan et al., 2025; Reinecke et al., 2025). Developing this understanding is a critical step toward unified hydrologic theory (Sivapalan, 2005) and can provide a blueprint for robust model development and informed decision making.

Previous efforts to map multiple hydrologic processes at continental scales are scarce. Most large-scale studies have focused on one process. For example, Buchanan et al. (2018) assessed the likelihood of infiltration excess flow occurrence by comparing whether rainfall intensity exceeds saturated hydraulic conductivity, finding that saturation excess dominates across the contiguous United States (U.S.), while infiltration excess is regionally likely in the central U.S. Similarly, studies on baseflow indices have shown their strong dependence on climatic and soil properties (Beck et al., 2013; Xie et al., 2024), and Fang and Shen (2017) quantified the runoff-storage connectivity through correlations between anomalies in streamflow gauge and satellite water storage observations, highlighting large-scale interactions among groundwater table, soil thickness, topography, and snow. In contrast, studies that examined multiple processes have been typically focused on a single or small groups of watersheds. A study in Alaska shows that the use of multiple streamflow statistics can help distinguish and assign hydrologic regions (Barnhart et al., 2022). One of the few studies to holistically investigate storm runoff generation at continental scales (Wu et al., 2021) demonstrates that large-sample analysis aligns with previous conceptual understanding of runoff mechanisms across landscapes (Dunne, 1978) while highlighting the roles of precipitation volume and geology. Model-aided studies have simulated global patterns of multiple indices: water partitioning into green and blue water, streamflow response elasticity to rainfall, and streamflow flashiness (Ji et al., 2025), U.S.-wide indices for water balance seasonality (Berghuijs et al., 2014). Another model-based approach has involved inferring hydrologic processes through parameter sensitivity analysis (Hay et al., 2023). These synthesis studies present promising descriptions of spatial patterns and directions for future progress toward a holistic understanding of runoff generation mechanisms, which still remains elusive.

Much of the research for generalizing watershed function has focused on summarizing flow regimes (Dettinger and Diaz, 2000; Lane et al., 2017; Lee et al., 2015; Lins, 1997) and predicting shifts in flow regime under future climate (Brunner et al., 2020; Hodgkins et al., 2024). Many studies cluster streamflow gauges using flow indices that target general (Almagro et al., 2024; Ariano and Ali, 2025; Mosley, 1981; Wu et al., 2021), intermittent (Sauquet et al., 2021), or seasonal streamflow patterns (Dhungel et al., 2016; Haines et al., 1988; Kennard et al., 2010). However, most of these studies aim to define the similarity of flow regimes rather than the underlying runoff generation processes. Furthermore, the results from clustering approaches are constrained to gauged locations and lack spatial coherence, making it challenging to extrapolate to ungauged watersheds.

To estimate watershed processes in ungauged locations, hydrologists have conventionally used maps derived from physiographic datasets. For example, in the United States context, the Environmental Protection Agency's Ecoregions (Omernik, 1987, 2004), an ecosystem classification based on the physical and biotic characteristics, is a common reference when discussing hydrologic processes (Falcone et al., 2010). Other classifications include the United States Geological Survey's Water Resources Regions (Seaber et al., 1987) based on streamflow networks, Hydrologic Landscape Regions (Santhi et al., 2008; Winter, 2001; Wolock, 2003a) based on physiographic and climatic datasets, and the United States Department of Agriculture's Hydrologic Soil Groups (Web Soil Survey, 2025) based on soil surveys. Nevertheless, regionalization based on physiographic data often fails to capture the full variability of watershed function (Ali et al., 2012; Oudin et al., 2008) because hydrologic processes can differ even among physiographically similar watersheds (McMillan et al., 2014). Capturing watershed processes at a continental scale calls for a scalable method to draw information from hydroclimatic datasets. To date, no studies have attempted to develop comprehensive maps of runoff generation processes based on streamflow observations that can effectively capture watersheds' functions.

Hydrologic signatures are metrics that quantify hydrologically-relevant dynamics, and offer a promising way to infer watershed processes with minimal data requirements (McMillan, 2021). Hydrologic signature calculations require only widely-available datasets, such as streamflow and precipitation, and can be related to various watershed processes, such as runoff generation and water storage dynamics (McMillan, 2020; Wlostowski et al., 2021). Using hydrologic signatures, expert knowledge, and landscape characteristics, Fenicia and McDonnell (2022) inferred dominant runoff processes and developed perceptual models at the regional scale; and Pechlivanidis and Arheimer (2015) mapped process differences at the national scale in India. Hydrologic signatures can capture the functional streamflow responses to climatic forcings and can discriminate different processes across landscapes (Araki et al., 2022; Gnann et al., 2020, 2021b; Janssen and Ameli, 2021). This enables a signature-based exploration of the relationship between landscape form and function (Bracken et al., 2013; Sivapalan, 2005).

Watershed attributes describe the physical characteristics of watersheds, which can be used to identify the drivers of hydrologic processes and to transfer hydrological knowledge to ungauged locations (Tarasova et al., 2023). The link between watershed attributes and signatures of streamflow response can be explored via machine learning approaches on large watershed samples. Regional and global applications include studies in the U.S. (Addor et al., 2018; Janssen and Ameli, 2021; Wu et al., 2021), Australia (Trancoso et al., 2017), Zimbabwe (Mazvimavi et al., 2005), Brazil (Almagro et al., 2024), Europe (Rudlang et al., 2025; Kuentz et al., 2017), and globally (Beck et al., 2015). Across all studies, climate emerged as the primary control on signatures. Non-climatic factors (i.e., landscape attributes), such as soil, geology, vegetation cover, and topography, had weak or limited predictive power. However, substantial evidence from field-based studies shows that landscape forms are a primary control of watershed function (Angermann et al., 2017; Fan et al., 2019; Jackisch et al., 2017; Jefferson et al., 2010; Lohse and Dietrich, 2005; Pfister et al., 2017; Zimmer and Gannon, 2018).

Weak predictive power of non-climatic drivers can be attributed to lack of high-resolution, accurate landscape attributes that describe regionally important processes (Gnann et al., 2021b; Tarasova et al., 2023). For example, wetlands are key regulators of low flows in the U.S. (Worland et al., 2018) and have been left out of previous studies (Addor et al., 2018). Similarly, weathering and glaciation have primary impacts on baseflow storage and generation (Neff et al., 2005; Tague and Grant, 2004), but rock permeability and porosity predictors did not clearly capture the relationship (Wu et al., 2021). Coarse spatial resolution, or limited quality and consistency of global datasets may reduce their predictive power (do Nascimento et al., 2025; Beck et al., 2015). Additionally, large-sample studies across broad climatic gradients may be obscuring the influences of landscape attributes. Regional analysis can mitigate this effect and elucidate the non-climatic drivers; for example, regional random forest models have revealed physiographic and anthropogenic controls on flow regimes (Almagro et al., 2024; Hammond et al., 2021). However, smaller regional sample sizes may limit prediction accuracy if datasets only provide tens of watersheds per region (Willard et al., 2024).

Lastly, the quality of signatures can compromise data-driven model performance and interpretation for process understanding. Examples include the sensitivity of flow duration curve slope to measurement errors (McMillan et al., 2017), the sensitivity of signatures to rating curve uncertainties (Westerberg et al., 2016), lack of process representativeness (McMillan et al., 2022), and inaccurate parameterization of storm separation algorithms (McMillan et al., 2023). Minimizing the impact of signature uncertainty is important for differentiating different regional watershed functionalities (Westerberg et al., 2016).

This study presents the first hydrologic processes map for the contiguous United States (CONUS). We synthesized hydrologic signatures as process indicators, going beyond pattern identification from single signatures. We hypothesize that signature combinations can represent six key hydrologic processes (McMillan, 2020; McMillan et al., 2022): baseflow and storage, water balance and seasonal flow variability, and saturation and infiltration excess overland flow. Using random forest models, we demonstrate the explanatory power of landscape metrics to predict hydrologic signatures and their regional variations, and thus the underlying processes, across CONUS.

We address the limitations of previous studies in predicting hydrologic signatures. First, we improved the quality of non-climatic attributes by: (i) incorporating new geological and wetland landscape attributes that have demonstrated strong connections to baseflow processes (Holt and McMillan, 2025); and (ii) utilizing watershed attributes from GAGES-II datasets (Falcone, 2011), derived from survey-based and higher-resolution products. Second, we interpret random forests using Shapley values (Shapley, 1953) following Husic et al. (2025), as well as permutation importance values within a regional model-building approach, following Hammond et al. (2021), which extends prior work to elucidate the regional contributions of non-climatic, landscape attributes to hydrologic processes. Furthermore, our work assessed 14 146 U.S. watersheds and was trained on 10 261 watersheds, nearly ten times more sample watersheds than previous studies; we leverage the Caravan and GAGES-II – the most extensive open-source large-sample datasets currently available (Falcone, 2011; Kratzert et al., 2023). Third, we utilize a set of hydrologic signatures proven robust across large-sample watershed studies and have a clear connection to critical-zone processes (McMillan et al., 2022), with their parameters further tuned to local storm characteristics. With these improvements, we expand watershed coverage and uncover more detailed spatial patterns of watershed processes than previously possible, using widely-available hydroclimatic datasets and physiographic attributes.

2 Data

We used two primary sources of streamgages and watershed attribute data to expand the number of samples: Caravan v1.5 (Kratzert et al., 2023, 2024) and U.S. Geological Survey GAGES-II (Falcone, 2011; Falcone et al., 2010). See Fig. 1 for the spatial distribution of the study watersheds. Caravan is an open-source dataset of global watersheds; its CONUS subset consists of 9234 watersheds sourced from CAMELS-US (Addor et al., 2017) and HYSETS (Arsenault et al., 2020). GAGES-II is a geospatial dataset of 9067 watersheds in the United States, selected for their quality to characterize natural and altered flow regimes.

Table 1 Hydrologic signatures used for building process hypotheses. The signature descriptions are adapted from McMillan et al. (2022).

2.1 Hydroclimatic dataset

We calculated hydrologic signatures listed in Table 1 using daily hydroclimatic timeseries data from watersheds within the contiguous United States (CONUS). For Caravan watersheds, we used U.S. Geological Survey (USGS) streamflow measurements paired with daily ERA5-Land forcings provided. For the GAGES-II watersheds, we obtained the USGS streamflow records (U.S. Geological Survey, 2025) using the dataRetrieval package (DeCicco et al., 2018) and gridMET forcings from Wieczorek et al. (2023). For calculating infiltration excess overland flow signatures of Wu et al. (2021; “RC_Pint”), we used the hourly precipitation from the North American Land Data Assimilation System 2 (NLDAS-2; Xia et al., 2012) provided through CAMELSH: a Large-Sample Hourly Hydrometeorological Dataset and Attributes at Watershed-Scale for CONUS (Tran, 2025; Tran et al., 2025).

2.2 Watershed attributes

We combined watershed attributes from three sources: (1) Caravan, (2) GAGES-II, and (3) geologic age and wetland attributes (Holt and McMillan, 2025). We added average geologic age and isolated wetland fraction metrics because of their strong link to baseflow processes, which were missing from previous large-sample analyses (Holt and McMillan, 2025). From the Caravan and Holt and McMillan (2025) attribute sets, we excluded binary or categorical attributes, monthly climate variables, uninformative attributes for the CONUS context (e.g., permafrost extent, gross domestic product), and highly correlated attributes (Spearman's rho >0.8 or $<-\mathrm{0.8}$; see Sect. S1 in the Supplement). Where available, Caravan attributes were substituted with GAGES-II attributes, as described in Sect. 3.2 and Table S1 in the Supplement. Table 2 lists the 23 attributes used in the random forest analysis. The purpose of merging Caravan and GAGES-II dataset is to maximize the sample size of watersheds and better capture regional hydrologic variability (see Table S1 and Fig. S1 in the Supplement for the watershed coverage across datasets).

Table 2 Landscape attributes used in training the random forest model. Descriptions are adapted from Falcone (2011), Falcone et al. (2010), Holt and McMillan (2025), Kratzert et al. (2023), Linke et al. (2019). For predictions, when certain attributes are unavailable, equivalent attributes are substituted (e.g., Caravan equivalents are used when predicting signatures for watershed samples available only in Caravan). The combinations are detailed in Table S1. An asterisk (*) in the unit column indicates that the landscape attribute unit from GAGES-II was converted to the Caravan equivalent (Fig. S14 shows the comparison). “n/a” stands for not applicable.

3 Method

We analyzed 14 146 gauged U.S. watersheds; our map of processes was based on observational data from 10 261 gauged sites and extended using random forest predictions to an additional 3885 watersheds. See Table S1, Figs. S1 and S2 for an overview of the datasets used and the workflow.

3.1 Calculating hydrologic signatures

A total of 12 signatures (four baseflow and groundwater signatures, four water balance and seasonality-related signatures, and four overland flow signatures) were used to characterize hydrologic dynamics (see Table 1). The signatures were selected based on their reliability in representing processes (McMillan et al., 2022). We calculated signatures using the TOSSH toolbox (Gnann et al., 2021a) and tuned the parameters for event separation for each hydroclimatic region (see Tables S2, S4).

We filtered out watersheds from our signature calculations based on quality criteria for watershed area and snow used by previous studies, and on the timeseries length needed for signatures to stabilize. First, we removed watersheds from our analysis with uncertain topographic boundaries, showing high discrepancies (>25 %) in the estimated drainage area between GAGES-II and Caravan datasets. Errors of <20 % are possible due to differences in watershed delineation tools or missing small tributaries (Ray, 2018). Second, for overland flow signature analysis, we excluded snow-dominated watersheds (>20 % snow fraction of total precipitation; a >30 % criterion were used in McMillan et al., 2022, and Wu et al., 2021); this is because our overland flow signatures can be heavily influenced by periods with no flow response due to snow or frozen conditions. Third, we excluded watersheds with less than 5 years of streamflow observation record, and those with over 30 % missing daily data over the period where streamflow was recorded (yielding at least three years of available data). Studies suggest that temporal hydrologic variability is adequately captured with 3 to 5 years of data (Refsgaard and Storm, 1990; Klemeš, 1986; Merz and Blöschl, 2009).

3.2 Training random forest models and predicting hydrologic signatures

We developed random forest models to examine potential drivers of hydrologic processes. Random forest models have been widely used for this task (Eng and Wolock, 2022; Lapides et al., 2023; Zipper et al., 2021) for their interpretability, relatively low computational demands, and robustness to multi-collinearity (Addor et al., 2018). For each signature, we constructed a random forest model to predict its values based on watershed attributes using the caret R package (Kuhn, 2008; R Core Team, 2024). Each model used 500 trees with the optimal number of features randomly resampled at each split, selected by minimizing root mean squared error (RMSE) through 10-fold cross-validation.

Only the signatures calculated from quality-controlled streamflow records (as described in Sect. 3.1) were used for training. Training on all Caravan watersheds yielded R²<0.4 for many of the signatures, so we limited the training samples to the 4748 Caravan watersheds with streamflow gauge IDs overlapping with GAGES-II to attain model performance comparable to previous studies (see Sect. S2). When a gauged watershed was present in both datasets, we prioritized CAMELS over HYSETS, and Caravan over GAGES-II. This is to ensure the broader applicability of our method across different countries, as Caravan is available at global scale. This yielded a total 14 146 watersheds for signature analysis (overview in Table S1).

We then used the trained model to predict hydrologic signatures for 3885 watersheds where observations did not meet data quality standards as described in Sect. 3.1. Previous study (do Nascimento et al., 2025) and our preliminary experiments (Sect. S2) showed improved model performance when the watershed attributes were derived from higher-resolution datasets based on detailed field surveys, such as in GAGES-II and Holt and McMillan (2025). Therefore, we used GAGES-II attributes and when unavailable, used the coarser resolution Caravan attributes (see Table S1).

Figure 1 Method used to obtain hydrologic signatures. Signatures are derived either from observed data (“Observed”: Caravan samples, n=7465; GAGES-II samples, n=2807; total n=10 261) or predicted using random forest models (“RF”; n=3885). Predicted samples are categorized as: “Caravan + GAGES-II overlap” (present in both the Caravan and GAGES-II datasets; n=618), “Caravan only” (exclusive to Caravan; n=2424), and “GAGES-II only” (exclusive to GAGES-II; n=843). State boundaries are indicated by grey lines.

3.3 Interpretation of hydrologic signatures as process descriptors

We combined signatures calculated from observed streamflow data and predicted with random forest models to develop a comprehensive map of processes for watersheds across the U.S. (Fig. 1). A bivariate space of hydrologic signatures was used to infer process dominance. For each selected process, we used the two signatures most strongly related to the process inferred from previous work (Bolotin and McMillan, 2024; McMillan, 2020; McMillan et al., 2022; Wu et al., 2021). Each signature was categorized based on the quantiles of signatures, from low (0 %–25 %), mid-low (25 %–50 %), mid-high (50 %–75 %), to high (75 %–100 %). When both of the two target signatures had mid-high (50 %–75 %) or high (75 %–100 %) values, we interpreted this as indicative of process dominance. This bivariate matrix can highlight the complexity of hydrologic processes where two signatures do not necessarily show the same trends.

The process hypotheses tested are described in Table 1 and cover six major hydrologic processes: baseflow, watershed storage, water balance, seasonal variability, overland flow dominance, and its type. Baseflow sustains discharge across seasons and reflects groundwater connectivity, while overland flow drives stormflow and influences flood peaks. Storage governs buffering and recession dynamics, and water balance losses through evapotranspiration and deep percolation determine how much precipitation is converted to streamflow. Together, these processes span the continuum from slow to fast hydrologic response and integrate both vertical and lateral fluxes, making them essential for hydrologic theory and modelling (Bergström, 1992; Kirchner, 2009; Berghuijs et al., 2014).

3.4 Interpretation of process drivers using Shapley values

We quantified feature importance using Shapley values (Shapley, 1953), which provide a robust and consistent measure to interpret random forest models (Lundberg et al., 2018). Shapley values represent the average marginal contribution of a feature (i.e., a landscape attribute) to a prediction, given the effects from all combinations of the considered features. Shapley values allow for local and global interpretation of machine learning model predictions, helping to uncover site-specific and generalizable linkages between hydrology and landscape features (Husic et al., 2025). We used the interpretable machine learning (iml) R package (Molnar et al., 2018) to calculate Shapley values over the training data.

To evaluate the regional effects of watershed attributes, we computed summary statistics on Shapley values. Shapley values are site-specific: ${\mathit{\varphi }}_x^{(y,i)}$ is the Shapley value calculated for an attribute x for a signature y at location i. Summing the Shapley values across watershed attributes x at a single location gives the deviation of the predicted signature value y_i at location i from the mean signature value across all sites. To compare effects of a landscape attribute x across sites, we normalize Shapley values by the total absolute contribution from all attributes at a site i; this gives a metric for the relative contribution of an attribute x to signature y at site i as:

$$\begin{array}{}\text{(1)}& R_x^{(y,i)}=\mathrm{|}{\mathit{\varphi }}_x^{(y,i)}\mathrm{|}/{\sum }_{x\in A}\mathrm{|}{\mathit{\varphi }}_x^{(y,i)}\mathrm{|}\end{array}$$

where A is the set of all watershed attributes, and |⋅| denotes the absolute value. To investigate which types of landscape characteristics are influential, we classified the watershed attributes into five categories (see Table 2), namely, topography, land-cover, soil & geology, human alteration, and climate.

Then, the average relative contribution of category k for signature y at location i, ${\stackrel{\mathrm{‾}}{R}}_k^{(y,i)}$, is calculated as:

$$\begin{array}{}\text{(2)}& {\stackrel{\mathrm{‾}}{R}}_k^{(y,i)}={\displaystyle \frac{\mathrm{1}}{K}}{\sum }_{x\in C_k}{R}_x^{(y,i)},\end{array}$$

where C_k is the set of watershed attributes belonging to category k, and K is the number of categories (in our case, five).

3.5 Interpretation of process drivers using permutation importance

To further evaluate locally important watershed attributes, we computed permutation importance, which measures the change in model performance when a feature (i.e., a landscape attribute) is removed. Prior work has shown that permutation importance derived from random forest models trained on regional samples is more effective than a continental approach for identifying physiographic, landscape controls on hydrologic responses, as it allows assessment under consistent climate conditions (Almagro et al., 2024; Hammond et al., 2021; Holt and McMillan, 2025). Therefore, we calculated permutation importance as the average changes in mean squared error (MSE), normalized by its standard deviation using the caret R package (Kuhn, 2008), from random forest models trained on regional watershed samples. Six climate regions were defined using a Gaussian mixture model in Scikit-learn (Pedregosa et al., 2011) based on relevant Caravan, GAGES-II, and Hammond et al. (2023) climate attributes (Table S4), and separate random forest models were trained for each region. Figure S4 shows the identified climate regions

Figure 2 Map of dominant processes estimated based on our hypothesis (defined in Table 1 and Sect. 3.3). Note that when baseflow and overland flow both occur, their colors are overlaid to give purple hues. “Unclassified” means a watershed is deemed neither baseflow- nor overland-low-dominant. “Region boundary” indicates the areas described in Sect. 4.2.1–4.2.3 (East and South (bottom right), Midwest and Central, West and Southwest (left in the figure)).

Figure 3 Hydrologic signatures of each process hypothesis, shown in bivariate maps (a–e). See the legend at the bottom for explanation. The high-process quantile from (a) is used to infer “Baseflow” in Fig. 2; from (b) to infer “High storage capacity”; from (c) to infer “Water balance losses”; and from (e) to infer “Overland flow”. Panel (f) shows the differences between the two signatures related to infiltration-excess (IE) flow and saturation-excess (SE) flow (i.e., values of IE Correlation (RC_Pint) minus SE Correlation (PC_Pvol). In the overland flow panels (e) and (f), watersheds dominated by snow (i.e., where more than 20 % of annual total precipitation falls as snow) are not shown. For the overland flow type pane (f), watersheds are not shown when the correlations between the event runoff coefficient and both rainfall characteristics (i.e., storm rainfall volume and maximum intensity) are negative. For maps of each signature value, see Figs. S5 and S6. See Table 1 for signature names' abbreviations.

4 Results

4.1 Mapping dominant processes across the contiguous U.S.

Figures 2 and 3 show the maps of dominant processes derived from the hypotheses outlined in Table 1. Figure 2 presents the signature of each process hypothesis in a bivariate map. Figure 3 provides a summary, displaying the four primary hydrologic processes when it is deemed dominant (i.e. both signatures are in the mid-high (50 %–75 %) or high (75 %–100 %) quantiles). Together, these maps highlight distinct regional patterns in hydrologic processes across the study area. The following sections examine these patterns in greater detail by region: the East and South. (Sect. 4.2.1), the Midwest and Central (Sect. 4.2.2), and the West and Southwest (Sect. 4.2.3). See Fig. 4 for the geographic features used to describe regional hydrologic processes.

Figure 4 Map of the contiguous United States showing (i) areas described in Sect. 4.2.1–4.2.3 (East and South, Midwest and Central, West and Southwest; bolded brown line) (ii) geographical boundaries used for the USGS National Water Availability Assessment (Qi and Mason, 2023; Stets et al., 2025; Van Metre et al., 2020) (beige line) (iii) topographic and geological features named in the text (pink annotations).

4.2 Spatial patterns of hydrologic processes inferred from signatures

4.2.1 Region 1: East and South

This humid region has moderate to high precipitation (1000–1500 mm yr⁻¹; calculated based on the 10th and 90th percentiles of sample watershed attributes), with low precipitation seasonality except in Florida. Temperatures vary widely from snow-dominated areas in the NorthEast to subtropical areas in Florida, with mean annual temperature ranging from 7–19 °C (Fig. S7). The landscape is old with deeply weathered soils and characterized by predominantly low-lying elevation (mean watershed elevation ranges between 40–600 m), though there is a primary elevation gradient from the Appalachian Mountains and Piedmont to the Eastern coastal plains, with peaks exceeding 1000 m (Fig. S11). In Fig. 3, signature values show that these climate and landscape conditions produce slowly-varying, baseflow-dominated flow regimes and mid-quantile signature values showing a lack of hydrologic extremes. Runoff ratios (Total Runoff Ratio and Event Runoff Ratio; Fig. 3c) are moderate or high and seasonal variability in flow and recessions is moderate to low. Storage capacity (Average Storage) is overall moderate, but recession shapes (Recession Parameter b) are variable (Fig. 3b). Evidence for overland flow is weak with saturation excess prevailing when it occurs (Fig. 3e,f).

The gradient along the geographical transect from the Appalachian spine to the Eastern coastal plain is apparent in several processes. The Appalachians have strong baseflow influence, shown by high baseflow index and slow recessions (Fig. 3a). Nonlinear recessions (high Recession Parameter b; Fig. 3b) indicate multiple groundwater reservoirs supplying baseflow. In contrast, the Piedmont has lower baseflows and fast recessions, relating to lower storage. The Eastern coastal plain, especially towards the South, has high baseflow and moderate to slow recessions (Fig. 3a). Linear recessions suggest a single dominant groundwater reservoir supplying baseflow in this sandy, coastal plain aquifer (Fig. 3b). Lower runoff ratios in the coastal plains indicate losses to deep groundwater including offshore discharge, especially in Florida's karst area (Figs. 3c, S9). The karst area stands out for its high dynamic storage and seasonality in recessions. Saturation excess dominates overland flow in the Coastal plain (Fig. 3f), although evidence for overland flow is weak (Fig. 3e) in contrast to a previous study (Wieczorek and LaMotte, 2010) that suggests the Florida panhandle has the highest fraction of saturation excess overland flow in the US.

In inland areas such as the valleys of the Tennessee-Missouri region, baseflow is moderate and recessions are relatively fast (Fig. 3a). The Gulf Coast region has lower baseflow and faster, linear recessions. Infiltration excess flow largely occurs in the narrow ocean margin of the Gulf Coast region but does not extend far inland (Fig. 3f). Exceptions to the area's fast runoff occur in the Ozark Mountains and the west of the Mississippi embayment where limited areas of high baseflow and slow recessions occur.

4.2.2 Region 2: MidWest and Central

The landscape of the Midwest and Central region is dominated by the gradient from recently-glaciated, sandy, forested watersheds of the Great Lakes region, to the poorly-drained, clay-rich but highly developed for agriculture and populated region of the Souris-Red-Rainy and Midwest regions. Across the Midwest and Central area, mean watershed elevation ranges from 200 to 700 m, and mean annual precipitation varies from 500 to 1000 mm. Moving west into the Central and Northern High Plain regions, elevation gradually increases, precipitation decreases, and population density decreases (Figs. S11, S7, S10). The region experiences mean annual temperatures between 6 to 13 °C. The absence of major topographic barriers results in a continental climate characterized by intense thunderstorms in summer and heavy snowfall in winter.

Signature values show that storage capacity is moderate throughout the Midwest (Fig. 3b). Storage in this region is provided by a moderate snowpack and high depth to bedrock (Fig. S8). Most of the region was previously glaciated, leaving a thick layer of glacial drift. The soil texture is graded from coarse and sandy around the Great Lakes to clay-rich further South and West, forming a distinctive outer ring around the Great Lakes region (Miller and White, 1998; Fig. S9). Following this gradient, there is no significant evidence for overland flow around the Great Lakes, changing to stronger evidence further South-West (Fig. 3e, f). Some occurrence of infiltration excess is consistent with evidence of this process from Midwest agricultural watersheds (Abban et al., 2014; Davis et al., 2014; Wilson et al., 2012). Streamflow seasonality follows the same gradient (Fig. 3d), with low seasonality around the Great Lakes where sandy aquifers sustain discharge year-round, and higher seasonality further SouthWest (Miller and White, 1998; Fig. 3d). A second gradient occurs in the MidWest from West to East, following precipitation and aridity gradients (Fig. S7). In the west, high aridity leads to high water balance losses to ET and low runoff coefficients at the annual and event scale (Fig. 3c).

4.2.3 Region 3: West and Southwest

The landscape of the West and Southwest region is dominated by the mountain ranges of the Coastal Ranges, Cascades, Sierra Nevada and Rocky Mountains, with mean watershed elevation ranging from 400 to over 2700 m. Dense populations in the coastal cities give way to sparsely populated inland areas. The climate exhibits strong gradients. The Pacific Northwest and Sierra Nevada mountain ranges receive substantial amount of precipitation than interior, with mean annual precipitation ranging from 460 to over 2100 mm yr⁻¹ across the region. The region shows a north-south temperature gradient with coastal moderation. Mean annual temperature ranges from 2 °C in northern and high mountain areas to over 20 °C in inland southern desert regions (Fig. S7). Precipitation patterns follow Mediterranean or semi-arid climates characterized by winter precipitation peaks and dry summers.

High baseflows with slow recession are prevalent across most of the Western region, where deep snowpacks drive sustained baseflow processes (Figs. 3a, S8; Barnhart et al., 2016; Tague and Grant, 2009). Inland areas tend to have faster recessions while retaining high baseflows, while coastal areas – where snow is rare – have lower baseflow while retaining slow recessions. The Southwest desert contrasts with the rest of the region, having low baseflows and fast recessions typical of the arid or semi-arid climate with water tables far below the land surface (Goodrich et al., 1997). Storage capacity and retention follow the same gradient from high in the Pacific Northwest to low in the South-East, but the high storage region is more constrained to the Rocky, Cascade and Sierra Nevada mountains (Fig. 3b). Water balance patterns contrast the pattern still further, with only the high mountains having high runoff ratios in contrast to low ratios throughout the remainder of the Western U.S. (Fig. 3c) Seasonal variability in processes is higher in the South (primarily California) where the seasonal Mediterranean climate pattern occurs with hot, dry summers and cool, wet winters (Figs. 3d, S6).

Processes in the coastal margin are markedly different from those inland. The moderating influence of the coast is strongly apparent in storage capacity (Fig. 3b): the northern Coast Ranges have lower average storage compared to high storage inland areas, while the southern coastal band has higher storage compared to low storage inland areas. Overland flows are strongly indicated all along the coast, but more weakly inland (Fig. 3e). Most overland flow favors saturation excess, although inland watersheds of the Southwest desert show areas of infiltration excess (Fig. 3f).

4.3 Inferred climate and landscape drivers of hydrologic processes

In this section, we interpret the random forest models to understand which aspects of climate and landscape are most important in controlling hydrologic processes in different regions of the U.S. We hypothesize that variable importance statistics from Shapley and permutation analysis reflect the relative importance of hydrologic process drivers. Random forest models performed reasonably well (R²>0.4) for most signatures (Fig. 5), consistent with previous studies using similar model setups (Addor et al., 2018; Beck et al., 2015; Bolotin and McMillan, 2024; Kuentz et al., 2017). Performance was higher for baseflow, water balance loss, and seasonality signatures, but lower for overland flow signatures. Figure S12 presents the regional model performances for each signature.

Figure 5 Ten-fold cross-validation performance of the random forest model trained on 4748 CONUS-samples, where gauge IDs overlapped with Caravan and GAGES-II. Bars show the average R² between observed and predicted signatures, with error bars representing the standard deviation. See Table 1 for signature definitions.

Figure 6 (a) The landscape attribute category that contributes most to hydrologic responses was identified based on the average relative contribution of each category, ${\stackrel{\mathrm{‾}}{R}}_k^{(y,i)}$ (derived from Shapley values; see Sect. 3.4). For each watershed, the most important category k was determined using the median of ${\stackrel{\mathrm{‾}}{R}}_k^{(y,i)}$ across all hydrologic signatures. Results are displayed for the watershed samples included in the random forest training. Numbers in the legend indicate the frequency that each category was identified as the most important. (b) Frequency of watershed attributes ranked among the top three most important variables in permutation importance (IncMSE %) across all signatures in six U.S. climate regions. The x axis indicates how many times each attribute appeared in the top three. See Sect. 3.4 and Table 2 for attribute names.

Figure 6 provides an overview of variable importance results: Figure 6a focuses on spatial patterns, showing the landscape attribute category that has the strongest contribution to predictions of signatures and processes for each watershed, calculated using aggregated Shapley values; Fig. 6b provides deeper insights into the ranking of landscape attributes, ordered by permutation importance, for predicting signatures in each region. Figure S13 complements Fig. 6a by showing the importance of landscape attribute categories in each region, based on permutation importance.

4.3.1 Region 1: East and South

In the East and South, a wide variety of landscape attribute categories dominate process predictions, including topography, soils and geology, climate and human alteration (Fig. 6). Climate attributes dominate in cooler areas in the Northeast and along the Appalachian spine, while topography attributes dominate on the Eastern coastal plain. Along the Gulf Coast, either climate or soils and geology may dominate. Human alteration attributes dominate clusters of watersheds around cities including New York, Philadelphia, Washington D.C., Raleigh and Atlanta.

Overall, and particularly for signatures relating to storage and water balance in the East and South Region (Total Runoff Ratio, Event Runoff Ratio, Recession seasonality, Average Storage, Recession Parameter b), the random forest models show that climate drivers are less important than in the rest of the U.S., and soils and geology, topography, and land cover drivers are more important (Figs. 6, S13). Human influence (population density) is a more important driver here than in other regions across most signatures, consistent with large areas of high population (Fig. S10). In addition to the major cities, highly developed areas of Western Florida show anomalous areas of low baseflow, as do developed Piedmont areas (Zimmer and Gannon, 2018).

In the NorthEast, across all signatures, the drivers that most often appeared in the top three controls of random forest performance were Silt fraction, Precipitation, Geologic age and Population density – representing the effects of geology, soils, climate and human development (Fig. 6b). Climate characteristics appear more often for signatures related to water balance and overland flow. In the South, Silt fraction, Aridity, Precipitation and Slope occur most often, representing gradients in elevation and soils from the Appalachians to the coastal plain and into Florida (Figs. S9, S13).

4.3.2 Region 2: MidWest and Central

In the Midwest and Central area, a wide variety of landscape attribute categories dominate process predictions, including topography, soils and geology, climate and human alteration, showing strong spatial patterns (Fig. 6). Soils and geology attributes dominate in the Great Lakes region, and in the arc of clay-rich soils in the High Plains and Midwest regions (Fig. S9). A mixture of climate and topography attributes dominate in the Souris-Red-Rainy region. Human alteration attributes dominate in clusters of watersheds around Chicago, Detroit and Cleveland.

Overall in the Midwest and Central area, the random forest models show that land cover and topography drivers are more important than in the rest of the U.S., while climate drivers are less important. Across all signatures, the drivers that most often appeared in the top three controls of random forest performance were Clay fraction, PET, Precipitation and Slope – representing the effects of soils, climate and topography (Fig. 6); this is consistent with the gradual variation in signatures' spatial patterns accompanying gradients in glacial drift and climate in this region, as discussed in Sect. 4.2.2. Despite the flat topography of the region, several topographic attributes appear in the top ten, perhaps reflecting the effect of unusual topographic features such as the driftless area. Land cover metrics (wetland, cropland, pasture) were secondary drivers, appearing for signatures related to storage and overland flow.

The impact of climate is spread between multiple drivers: PET, Precipitation, Low precipitation frequency and Aridity. Climate drivers in the Midwest and Central area show multiple distinct spatial patterns, with aridity and low precipitation metrics showing an east-west gradient, temperature and PET having a north-south gradient, and precipitation and seasonality having a Northwest-Southeast gradient (Figs. S7, S8). Thus, each part of the Midwest and Central area has a unique holistic climate combination. Climate patterns differ distinctly from the NorthEast-Southwest pattern of the soils and land cover.

4.3.3 Region 3: West and Southwest

In the West, climate attributes dominate process predictions across most watersheds in the Pacific Northwest and Mountain West (Figs. 6a, S7). Some mountain areas have dominant topographic attributes, and topography drivers are more important in the Southwest region compared to the wider U.S.. Climate properties that appear most often include Snow fraction, Precipitation, Aridity and Seasonality (Fig. 6b: regions Pacific Northwest, Southwest, Mountain West). These attributes describe the primary climatic features of the West and Southwest U.S., which are governed by precipitation and aridity gradients from North to South, and from coasts to inland (Fig. S7). Inland mountain chains influence flow regimes by providing spring snowmelt and mountain block recharge, among the many influences of topography on hydrologic processes (Gnann et al., 2025). These controls are demonstrated by the importance of snow fraction alongside topographic attributes, elevation and slope. Soil control on runoff process is seen by the importance of clay fraction in the Pacific Northwest, reflecting Oregon's common clay soils (Miller and White, 1998).

5 Discussion

This study creates comprehensive maps of hydrologic processes across the contiguous United States by using machine learning to analyze streamflow signatures and connecting these signatures to dominant watershed processes. The analysis from over 10 000 watersheds shows distinct regional patterns in estimated hydrologic processes and its potential drivers. These process maps provide novel information for selecting appropriate hydrologic models across large domains and help hydrologists anticipate how watersheds will respond to environmental changes such as altered climate or land use. In the following sections, we discuss how these maps provide new benchmarks (Sect. 5.1), inform hydrologic modelling (Sect. 5.2), and outline directions for future work (Sect. 5.3).

5.1 New benchmark maps of process understanding over large domains

Our results build on previous work to map hydrologic processes and drivers. Our map of baseflow process importance shows similar patterns to previous studies into baseflow and groundwater contribution to streamflow (Beck et al., 2013; Santhi et al., 2008; Xie et al., 2024). As with those studies, our approach of using observations and machine learning methods provides finer detail than can be estimated using statistical interpolation or by hydrologic or climate models. By combining multiple recent datasets, we increase the number of observations used in our analysis. In our study, we used >10 000 observed watershed data within CONUS, representing a substantial advancement compared to the >600 to >3000 observation samples used in previous studies (Addor et al., 2018; Beck et al., 2013, 2015; Janssen and Ameli, 2021; Wu et al., 2021). Our analysis therefore provides a new benchmark, offering the most comprehensive coverage and highest spatial characterization of hydrologic processes across the contiguous United States to date. While larger datasets have been analyzed elsewhere, for example, >8000 watersheds (Santhi et al., 2008), >23 000 watersheds (Xie et al., 2024), those efforts focused exclusively on baseflow index. Beck et al. (2013) found sometimes differing drivers of baseflow index and recession slope despite their close connection: by using bivariate plots, we could more clearly highlight regions where patterns of these two signatures diverge. Those areas include the Pacific Northwest coast with lower baseflow index but slow recessions, and the central high plains with high baseflow index but fast recessions.

Previous studies investigated patterns of overland flow generation across the U.S. using soil maps and rainfall intensity (Buchanan et al., 2018) streamflow signatures (Wu et al., 2021) and modeling approaches (Wolock, 2003b). Like us, Buchanan et al. (2018) and Wu et al. (2021) found infiltration excess runoff important throughout the high plains, and saturation excess in the valleys of the Tennessee-Missouri region, and a mixture of saturation and infiltration excess in the Southwestern U.S. Substantial overland flow occurs in Southwest chaparral systems (Valeron and Meixner, 2010), and although deep groundwater tables suggest infiltration excess, we found a mixture of mechanisms. This may reflect vegetation shifting the inferred overland flow mechanism toward saturation excess. Infiltration excess is inferred when overland flow is related to storm intensity rather than storm size. In arid and semi-arid catchments, vegetation can locally increase infiltration capacity and soil water retention, reducing the extent of infiltration excess overland flow (Stein et al., 2021). Additionally, where smaller storms are intercepted by canopies, signatures may incorrectly attribute the runoff to saturation excess rather than infiltration excess. However, our results are supported by global studies that show saturation excess is always more common than infiltration excess even in arid regions, as saturation excess is generated in riparian zones and topographic convergence areas where water tables are higher (McMillan et al., 2025).

By mapping and categorizing the primary drivers of runoff processes, we can untangle which physical characteristics drive the hydrologic response in each region. In the East and South, soil, geology, and topography emerged as primary drivers, which is consistent with regional hydrologic process knowledge. Topography is important in the Appalachian Piedmont, where wide and wet valley bottoms generate fast responses (Zimmer and Gannon, 2018). Soils are important along the Gulf Coast where clay-rich soils promote infiltration-excess overland flow (Miller, 1999; Fig. S9), producing mixed storage and water balance signatures despite deep bedrock (Fig. S8) and semi-consolidated sand aquifers; and on the Eastern Coastal Plain where sandy soils, seasonal flooding, and wetlands likely support a single dominant groundwater reservoir supplying baseflow (Fig. 3b; Holt and McMillan, 2025; Hupp, 2000). The machine learning approach is especially powerful for this purpose, as multiple landscape attributes often contribute simultaneously to the hydrologic response.

Our maps of primary drivers based on Shapley values extend previous work to analyze the drivers of hydrologic signatures. For example, Addor et al. (2018; their Fig. 3) show that climate (aridity, seasonality, snow fraction) is the primary driver across most signatures, with topography (elevation, slope) and land cover (forest, leaf area index) being secondary drivers. Figure 6a similarly shows climate and topography as dominant, but adds spatial information to show that, for example, climate is dominant in the mountainous western U.S., but soils and geology dominate the Midwest and much of the Northeastern U.S. Geological age, a recently-proposed attribute to summarize watershed geology, was often in the top random forest attributes (Holt and McMillan, 2025). This highlights the need and opportunity for development of new landscape attributes that characterize the subsurface, echoing the call by Tarasova et al. (2023) and do Nascimento et al. (2025).

In four of the six regions, soil texture, particularly silt or clay fraction, was identified as a recurring primary driver (Fig. 6b), though their roles differ by context. In the Northeast, silt dominates variable importance; silt is found in glacial till layer and supports high water storage and baseflow (Shanley et al., 2015) while facilitating subsurface stormflow under wet conditions (Detty and McGuire, 2010). In the South, despite silt being identified as a primary driver, clay is the dominant soil texture in many areas (Miller and White, 1998); in the Mississippi embayment, extensive confining units of clay and silt separate aquifers and control the groundwater flow (Renken, 1998; Clark et al., 2011). These two cases suggest that Shapley or permutation-based methods may not fully separate correlated variables due to their treatment of joint variable distributions, and high clay content may be implicitly captured through the absence of silt in regional analyses.

5.2 Informing model selection and evaluation

Our results support hydrological modeling by enabling hydrologists to check whether key processes in a watershed are well-represented by a candidate model prior to application. A wide range of hydrologic models with differing process representations, structures and complexities are available (Knoben et al., 2020). Hydrologists must make choices on whether to include simulations of additional processes such as snowpack or deep groundwater, and the complexity required such as including energy balance at the land surface. Our maps of hydrologic processes provide a pre-screening tool to match hydrological models with appropriate process representations to regions. This approach aims to reduce model structural errors by discouraging use of models ill-suited to the dominant processes (e.g., using a bucket model in overland flow-dominated regions).

Many previous studies have assessed preferred model structure in individual research watersheds, often using in-depth data analysis to ensure that modeled processes are consistent with observed processes (e.g. Hrachowitz et al., 2014; Kavetski and Fenicia, 2011). This study provides a method to support transparent model justification in applied studies without the resources to conduct model structure investigations, and to upscale model structure decisions to large domains. For example, if selecting models from the MARRMoT toolbox (Knoben et al., 2020), models for regions of dominant overland flow should include saturation excess and/or infiltration excess pathways, and models for regions of complex storage and retention should include multiple parallel groundwater reservoirs. The ability to choose appropriate models for thousands of watersheds is needed for new, flexible model frameworks such as the U.S. Next-Generation National Water Model Framework (Cosgrove et al., 2024; Johnson et al., 2023; Ogden et al., 2026). Our observation-based method complements previous large-domain model-based methods that use analysis of model sensitivities (Markstrom et al., 2016) and performance (Prieto et al., 2021; Spieler et al., 2020) Therefore, where hydrologists seek to evaluate models against process representation, this study offers an opportunity to enhance model benchmarking frameworks by adding process realism as a metric.

5.3 Limitations and future work

The hydrologic process maps produced by this study are limited to the contiguous U.S. Recent streamflow observation datasets offer the opportunity to extend this method to other regions or globally. Such datasets include the community Caravan dataset (Kratzert et al., 2023), and the international dataset of watersheds with limited human influences, Reference Observatory of Basins for International hydrological climate change detection (ROBIN; Turner et al., 2025). If extending the method globally, caution is advised with scaling, in order to represent different ranges of signature values in different regions. In this study, we plotted signature values as quantiles based on the U.S. distribution, but other countries may have very different signature distributions (McMillan et al., 2022). Therefore, watershed processes that are considered important in a U.S. context, may be considered less important in a global context. Further, some regions of the U.S. are excluded or poorly represented in the dominant process maps presented in this paper, due to a low spatial coverage of USGS stream gages. For example, there are significant gaps in the arid southwest where perennial streamflow is rare (Kiang et al., 2013; Krabbenhoft et al., 2022). In such regions there is a need for alternative process-mapping methods that do not rely on streamflow records.

Hydrological signatures in this study are long-term averages of the multi-year streamflow dynamics, which may not fully capture temporal variability in watershed processes. Future studies should account for long-term hydroclimatic changes (Hobeichi et al., 2022; Gudmundsson et al., 2025), as well as inter-annual variability (Vogel et al., 1994) and seasonal variations in watershed function (Payn et al., 2012; Gomi et al., 2008). Another complication is that hydrologic signatures are often confounded by multiple processes (McMillan, 2020, 2023), whether driven by natural flow dynamics or impaired by human activities. For example, water abstraction by reservoirs reduces downstream flow variability and increases water balance deficits (Salwey et al., 2023; Veldkamp et al., 2017), but changes in vegetation or climate could induce similar effects. Disentangling these impacts remains challenging without testing narrower hypotheses about watershed function, incorporating expert knowledge, or having detailed information about human interventions. In this study, we partially mitigated this issue by using multiple signatures to characterize processes, and by representing human alteration through population density, which showed strong explanatory power for the signatures. Nevertheless, considerable effort is still needed to isolate the combined impacts of multiple processes, as well as the effects of urban development and agricultural practices on flow dynamics (Grantham et al., 2022) for improving the large-scale application of signatures.

A limitation of this study that would become more apparent at a global scale is the quality of precipitation, streamflow, and attribute data. A previous study noted issues with limited quality and consistency of the global attribute data for soils and geology that reduced their predictive power (Beck et al., 2015). Continental scales necessitate the use of gridded precipitation products, but in areas with low density of observations these products may be insufficient to analyze localized, flashy processes such as infiltration excess flow (McMillan et al., 2023). In small, headwater watersheds, precipitation grid size may be large compared to watershed area, and headwaters are also underrepresented in streamflow observations (Golden et al., 2025). Additionally, errors in watershed boundary delineation would affect signatures that use drainage area to normalize flow, such as runoff ratio (Total Runoff Ratio, Event Runoff Ratio) and water balance (Average Storage). In snowy areas, signature values can be compromised because liquid water inputs to the watershed come from snowmelt rather than directly from precipitation. In our study, we excluded snow-dominated watersheds for signatures related to overland flow, as these require event-scale surface water input that are particularly affected by frozen or snowmelt conditions. Products such as NLDAS3 (Case et al., 2025) or surface water inputs considering rain-on-snow and snowmelt (Hammond, 2024; Hammond and Kampf, 2020) may provide future abilities to estimate overland flow processes in snow areas using estimates of hourly snow accumulation and melt. While our study used potential evapotranspiration (PET) information in only one signature (Average Storage), uncertainty in PET is a major issue of global datasets and needs to be addressed (Clerc-Schwarzenbach et al., 2024; Destouni and Zarei, 2024) before this approach can be expanded to a variety of (eco)hydrologic processes.

A further limitation is the extent to which continental scale maps of dominant processes can be validated. Large-domain signature datasets can be evaluated for data quality, for interpolation quality using cross-validation, and compared with previous datasets. However, it is more difficult to determine how accurately signatures relate to processes over large domains. Research watersheds offer “ground truth” points at which processes are already well understood (Penna, 2024). Previous studies used a handful of U.S. critical zone observatory watersheds for evaluation (McMillan et al., 2022). However, the large number of past and present research watersheds across the globe offer an interesting future opportunity for wider-scale validation of process mapping techniques (McMillan et al., 2025; Sebestyen et al., 2025). Similarly, validation of process drivers remains challenging. While Shapley values and permutation importance provide explanatory power for random forest models, they have some limitations. Both metrics characterize model interactions within a given dataset; therefore, the variety of processes covered in the dataset matters, and data or model uncertainties may propagate into the interpretations (Husic, 2025). Shapley values do not capture joint distributional effects among multiple interacting variables (Lundberg and Lee, 2017). Developing an explanatory framework that maximizes both model performance and interpretability remains an ongoing research area in hydrology (Robert Maier et al., 2024; Willard et al., 2024).

6 Conclusion

A fundamental question in hydrology is how hydrologic processes are organized over large scales, and how they are controlled by climate and landscape (Blöschl et al., 2019). In this study, we contribute towards answering this question by mapping hydrologic processes and their drivers across the contiguous U.S. Our approach used hydrologic signatures to describe streamflow dynamics, and connected these dynamics to dominant processes in the associated watersheds using established relationships between signatures and watershed processes. We analyzed 14 146 gauged U.S. watersheds; our map of processes was based on observational data from 10 261 gauged sites and extended using random forest predictions to an additional 3885 watersheds with insufficient record length or completeness. Our method enables knowledge transfer from gauged basins with well-established conceptual models to ungauged or poorly instrumented watersheds.

Our results comprise maps of hydrologic process importance across the contiguous U.S., including baseflow, overland flow, water storage, seasonal variation and water balance processes. Using interpretable machine learning methods, we create maps of process drivers that explain which climate and landscape attributes are dominant in controlling hydrologic processes in each watershed and each region. We find clear patterns at the continental scale in hydrologic processes, with infiltration excess overland flow dominating the high plains., saturation excess flow prevalent in the valleys of the Tennessee-Missouri region, and varying baseflow contributions across regions. The novelty of this study is in demonstrating that incorporating more detailed landscape attribute elucidates non-climate variables as dominant controls on hydrologic processes, even with a large sample across multiple climate regions. Specifically, the results showed that climate primarily controls hydrologic processes in the western U.S., while soils and geology dominate in the Great Lakes region, topography controls processes in the Southeast, and human influences are most important around large cities across the East.

Our findings extend and generalize process understanding from research watersheds to large domains, revealing regional heterogeneity within broader physiographic provinces that are often treated as hydrologically uniform. Hydrologic process maps provide essential support for new, large-domain model frameworks that must select model structure across thousands of watersheds. These maps enable hydrologists to select models that adequately represent the dominant processes of a watershed. Identification of dominant processes in each region further enables hydrologists to anticipate streamflow response to environmental change, by identifying which processes are most sensitive to shifts in driving variables. Such analysis has the potential to support scenario testing for future land use or climate, to guide selection of green and grey infrastructure compatible with dominant processes, and to inform risk assessments for regions prone to flash flooding, streamflow depletion, or altered seasonal flow regimes.