Beyond total impervious area: a new lumped descriptor of basin-wide hydrologic connectivity for characterizing urban watersheds

Urbanization impacts on hydrologic response are typically indexed as a function of the fraction of total impervious area (TIA), i.e., the proportion of impervious areas in a basin. This implicitly assumes that changes in flood characteristics are somehow proportional to the extents of land development without considering that such impacts may vary widely depending on the location of the developed areas with respect to each other, the less developed land patches, the stream network, and the basin outlet. In other words, TIA is blind to the spatial arrangement of the different types of land patches within a basin and to the nuanced ways in which runoff volumes are variously generated over them and then subsequently retained or detained as they are routed towards the stream network and then the outlet. To overcome such limitations, we propose a new lumped index that measures the impacts of urbanization on basin response in terms of the emerging hydrologic connectivity, defined here as the distributed property that explains the ability of any hillslope location to quickly receive and transfer runoff to the stream network as driven by topographically induced runoff pathways and locally affected by the different land-use/land-cover types present in a watershed. This alternative, hydrologic-connectivity-based index of urbanization (HCIU) displays sensitivity to the spatial arrangement of both fully developed and less developed or undeveloped patches, each with different degrees of imperviousness, roughness, and other characteristics affecting their abilities to either generate or hydrologically retain or detain runoff, reflecting their distinct localized effects on hydrologic connectivity. The proposed HCIU can be readily obtained in a GIS environment from easily available raster geospatial data. We found that HCIU improves the predictive power of regional equations for peak flow in three large case-study homogeneous regions, when used in place of the traditional TIA.

1 Introduction

The ongoing expansion in land development across many regions of the world is a major driver of alterations in the hydrologic response of watersheds (Nirupama and Simonovic, 2007; Sillanpää and Koivusalo, 2015), with subsequent impacts on urban stream ecosystems (Walsh et al., 2005; Vietz et al., 2016a, b). Along with climate change, this trend in urbanization is expected to pose formidable challenges for water resources management in the years to come (Praskievicz and Chang, 2009; Bell et al., 2017; Zölch et al., 2017). Because developed land patches have quite different infiltration and interception capacities, as well as surface roughness characteristics, than undeveloped sectors, they have strong effects on stormwater runoff dynamics; this is why quantifications of the level of urbanization based on a basin's impervious area are widely used in the domains of engineering hydrology and urban river ecology (Bauer et al., 2007; Roy and Shuster, 2009; Lee et al., 2018; Gong et al., 2020; Liang et al., 2022). Urbanization descriptors specified as a fraction of the total watershed area (Bell et al., 2016) are extensively adopted in a range of stormwater management (Kong et al., 2017; Sultana et al., 2020) as well as flood risk assessment and mitigation practices (Suharyanto et al., 1997; Loperfido et al., 2014; Sohn et al., 2020). Popular applications of these metrics include the study of stormwater runoff dynamics (Meierdiercks et al., 2010; Fletcher et al., 2013; Yao et al., 2016), water-quality assessments (Fletcher et al., 2014; Lee et al., 2012; Li et al., 2021), and peak-flow prediction in ungauged (urban) basins (PUBs; Kennedy and Paretti, 2014; Southard, 2010), also involving urban planning and regulation (Smucker et al., 2016).

The total impervious area (TIA) of a basin has been historically the most widely adopted urbanization metric in hydrology (Shuster et al., 2005). Other representations of the level of urbanization are the directly connected impervious area (DCIA, i.e., the subset of TIA connected to the stream network through constructed drainage or other fully impervious pathways; Han and Burian, 2009; Sytsma et al., 2020) and the effective impervious area (EIA, i.e., an indirect estimate of DCIA; Boyd et al., 1993, 1994; Ebrahimian et al., 2016a, b, 2018). Despite their widespread use, these methods have intrinsic limitations, as highlighted by a growing body of research (Shuster et al., 2005; Law et al., 2009; Beck et al., 2016; Bell et al., 2016). One drawback relates to their inability to explicitly account for differences in the spatial arrangements of impervious patches (Shuster et al., 2005; Beck et al., 2016; Bell et al., 2016). However, basins with similar levels of land development can exhibit distinct hydrologic behaviors because the actual locations of the developed sectors within a watershed significantly impact surface runoff dynamics (Corbett et al., 1997; Pappas et al., 2008). Another issue with these kinds of descriptors is that they are not able to capture the complex distributed interactions between urbanized patches (with varying land-development intensities, depending on location) and undeveloped sectors with heterogenous hydrologic characteristics (Bell et al., 2016; Law et al., 2009) even though different spatial configurations will significantly affect a watershed's hydrologic response (Loperfido et al., 2014).

Despite the conceptual limitations of TIA and similar indices used to depict a watershed's impervious areas, these kinds of simple descriptors are needed to characterize the degree of urbanization in lumped hydrologic models, with their inherent trade-offs between basin representation detail and the spatial heterogeneity of the captured hydrologic processes (Hrachowitz and Clark, 2017). More recent urbanization metrics proposed by Yang et al. (2011) and Beck et al. (2016) try to incorporate information on the spatial distribution and geometric fragmentation of developed and undeveloped patches. These indices implicitly attempt to capture the degree of contiguity and interconnection of different land-use/land-cover (LULC) sectors based on their spatial density and granularity but overlook the effects of basin relief on runoff routing. However, spatial contiguity does not fully explain the hydrologic connectivity of patches with distinct LULC types as this property is ultimately determined by water pathways induced by topographic gradients. Following a conceptually different approach, Zhang and Shuster (2014) proposed the following two indices: (1) the average distance of impervious patches to the outlet, measured along topographic pathways, and (2) the mean number of pervious cells along those routes, regarded as a proxy for hydrologic disconnection. These metrics do consider the effects of topography but fall short in accounting for the heterogeneity across different LULC types and their varying effects on surface runoff dynamics as they adopt a binary, simplified pervious-versus-impervious LULC classification. Undeveloped areas located downstream of urbanized land patches can mitigate the adverse hydrologic impacts of the latter to varying degrees and through different mechanisms, depending on specific factors, such as soil infiltrability and vegetation type (and its effects on interception and roughness; Law et al., 2009). The fact that Zhang and Shuster's (2014) indices cannot account for this continuum of behaviors likely explains their poor correlation with simulated stormwater runoff volumes in their two case-study basins.

In summary, ongoing efforts in urban watershed characterization have brought to light the limitations of traditional impervious-area-based indicators like TIA as proxies for the effects of land development on stormwater hydrologic processes. Researchers have proposed alternative approaches that shift the focus towards hydrologic or geometric properties, considering either the degree of dispersion and granularity of patches with different LULC types or the interconnections between pervious and impervious patches as driven by topographic gradients. However, despite these advances, existing methods have yet to simultaneously address the synergistic impacts of topography and of the existence of a wide spectrum of LULC types with heterogeneous effects on stormwater runoff dynamics.

We propose a conceptual approach for deriving lumped urbanization metrics starting from the DEM and LULC map of a watershed by measuring the impacts of land development on hydrologic response in terms of its distributed effects on hydrologic connectivity, the spatially distributed property that explains the ability of any hillslope location to quickly receive and transfer runoff to the stream network, as driven by topographically induced runoff pathways and influenced by the spatial arrangement of LULC patches with different hydrologic properties, such as roughness and infiltrability (Hooke et al., 2021). Distinct LULC patches intensify or mitigate the hydrologic response of their contributing area to varying extents, depending on the hydrologic processes that take place under specific LULC conditions (e.g., high canopy and litter interception, as well as detention, due to high surface roughness in densely vegetated areas; increased runoff volumes and peaks because of negligible infiltration losses and shorter travel times due to low roughness in urbanized sectors; and enhanced retention due to higher infiltrability in undisturbed areas). Hence, basin cells with different land-use/land-cover types contribute to hydrologic connectivity in different ways. Conceptually, highly urbanized, fully impervious patches represent one extreme in the continuous spectrum of LULC potentials for generating runoff and increasing connectivity, with their smooth surfaces and absence of infiltration losses. In the other extreme instead are forested areas and other natural LULC types with high vegetation densities because of their ability to reduce runoff volumes and travel speeds through canopy and litter interception and temporary storage, high surface roughness, and soil infiltration.

Because it displays sensitivity to both the topographic structure and the heterogenous mosaic of LULC types of a basin, hydrologic connectivity provides a methodological framework for a conceptual yet quantitative and comprehensive assessment of the impacts of land development that considers not only the spatial arrangement of urban sectors with distinct land-development intensities, but also their interactions with undeveloped patches, with their range of flood-mitigating capabilities. Connectivity analyses have gained popularity in recent years within hydrology and geomorphology (Bordoni et al., 2018; Heckmann et al., 2018; Husic and Michalek, 2022; Martini et al., 2022). These methods define a connectivity index (Eq. 1) for each basin cell that can be interpreted as a measure of its potential for affecting runoff or sediment fluxes (Hooke et al., 2021), depending on its location within the watershed. However, a map of connectivity values, with possible local peaks induced by the presence of fully developed patches, does not represent a lumped measure of the impacts of land development per se, nor does it allow for a straightforward comparison of the effects of urbanization across a range of basins or of different urbanization scenarios within a given watershed; thus, further conceptualizations are required.

As the traditional definition of the connectivity index only accounts for topographically induced runoff pathways (Borselli et al., 2008), additional adjustments may be needed, depending on the level of urbanization and the scope of the analysis, to also include the effects of underground stormwater drainage infrastructure typically present in urban environments. Underground pipe flows may be regarded as an additional source of connectivity, which can alter and sometimes even reverse the connectivity induced by topography (e.g., when stormwaters are pumped against topographic gradients).

In this work, we derive a lumped metric of urbanization effects on hydrologic response, incorporating only topographically induced connectivity (i.e., neglecting any effects of underground storm sewer infrastructure), and test its performance as a predictor in regional peak-flow equations. Peak flows are among the hydrologic-response variables of greatest interest in urban flooding risk (Feng et al., 2021) and are the most important for design purposes (Vogel and Castellarin, 2017). While considering the additional source of connectivity introduced by the underground drainage network would be straightforward, as explained in the Discussion section, we could not account for it here because it was impossible to obtain stormwater sewer data for the hundreds of watersheds involved in our regional-scale analyses. However, for the scope of our investigation, which focuses on hydrologic response during severe flooding (peak flows with return periods from 2 to 500 years), considering only topographically induced connectivity should be acceptable. This approach allows us to capture the impacts of land development on the surface and near-surface phases of a basin's response and the effects of streams and watercourses, including the artificial ditches and canals that make up the so-called major drainage system of stormwater infrastructure (i.e., excluding the underground network, also known as the minor system; Martins et al., 2017). During severe flooding, it is surface dynamics that predominantly govern hydrologic response, as the underground stormwater infrastructure's capacity is typically exceeded.

We benchmark our hydrologic-connectivity-based index of urbanization (HCIU) against the traditional fraction of TIA by alternatively using one of these two metrics as a predictor in regional peak-flow equations for urbanized basins. Imperviousness descriptors expressed as a fraction of the total basin area (e.g., TIA, EIA, and DCIA) are still among the most popular approaches to quantify the effects of land development in lumped hydrologic and regional models (Bell et al., 2016; Yang et al., 2023). Among these, we chose TIA as a benchmark because HCIU and TIA both condense distributed surface basin information (i.e., LULC and the topographic structure and LULC only, respectively) into a lumped urbanization metric, making their comparison conceptually straightforward. On the other hand, EIA is an indirect estimate of the impacts of urbanization based on retrospective analyses of concurrent historic flow and precipitation data for the case-study watersheds (Ebrahimian et al., 2016b). In preliminary tests, we found much uncertainty in EIA values, possibly due to the challenges involved in reliably estimating precipitation depths across basins with varying sizes and, for the same watershed, across distinct storm events (depending, e.g., on the areal footprint and location of the storm relative to basin extent). We also discarded DCIA as its estimation would require knowing the configuration of the stormwater sewer network for each case-study basin, which was unfeasible, as mentioned above.

In the next section, we explain in detail our conceptual methodology for deriving HCIU. We then test this novel lumped urbanization metric against TIA as a predictor of peak flows in regional equations for three case-study regions with urbanized watersheds described in Sect. 3. We conclude by showing and discussing the performance of HCIU compared to TIA (in Sects. 4 and 5, respectively), highlighting strengths and weaknesses of the proposed conceptual framework and outlining possible future research directions.

2 Methodology

Of the many types of connectivity indices proposed in the literature (see, e.g., Bracken et al., 2013, for a comprehensive review), we adopt a formulation (Eq. 1) first introduced by Borselli et al. (2008), which measures the potential hydrologic connectivity based on a weighted topographic analysis (Bracken et al., 2013). We consider two alternatives for the weights as a function of either Manning's surface roughness coefficient, n (Eq. 2), or the curve number, CN (Eq. 3), of each basin pixel. Conceptually, both n and CN are distributed basin properties that consistently vary across LULC types with different surface runoff dynamics, e.g., due to their distinct water retention/detention capacity or infiltrability.

Below, we recall established formulations based on Eq. (1), including recommendations for the weighting factors. These methods provide a measure of the hydrologic connectivity at each basin cell, resulting in a connectivity map for the watershed. We then propose hydrologically driven criteria to obtain a lumped, hydrologic-connectivity-based index of urbanization (HCIU) able to summarize the effects of the spatial arrangement of the varied LULC patches in a watershed in terms of their distributed impacts on connectivity. Two alternative indices, HCIU(n) and HCIU(CN) (depending on whether n or CN is chosen as weighting coefficient) are derived in this work and tested against the traditional TIA as explanatory variables in predictive peak-flow equations.

2.1 Connectivity-index formulations

Borselli et al. (2008) proposed a widely used GIS-based index of connectivity to assess sediment erosion and transport, which was then modified by Cavalli et al. (2013), Persichillo et al. (2018), Zanandrea et al. (2019), Hooke et al. (2021), and Husic and Michalek (2022), among others, to focus on other basin dynamics, such as runoff generation or landslide occurrence. In general, irrespective of the formulation, computing the index of connectivity requires assigning a flow direction (by either the D8 or the D-Infinity algorithm; Hooke et al., 2021) and slope value, S, to each basin cell, from the DEM, as well as a weighting coefficient, W, that varies across formulations depending on some additional hydrologic properties of interest (e.g., potential for erosion or runoff generation).

The index of connectivity (IC_k) is estimated for each raster cell k of the basin hillslope component (i.e., excluding cells corresponding to the stream network) as the logarithm of the ratio of the upslope (D_up,k) and downslope (D_dn,k; Fig. 1a) components, as shown in Eq. (1) (Hooke et al., 2021).

$$\begin{array}{}\text{(1)}& {\mathrm{IC}}_k={\log }_{\mathrm{10}}\left({\displaystyle \frac{D_{\mathrm{up},k}}{D_{\mathrm{dn},k}}}\right)={\log }_{\mathrm{10}}\left({\displaystyle \frac{{\stackrel{\mathrm{‾}}{W}}_k{\stackrel{\mathrm{‾}}{S}}_k{\sqrt{A}}_k}{{\sum }_{i=k}^{n_k}\frac{d_i}{W_i{S}_i}}}\right)\end{array}$$

The upslope component, D_up,k, relates to cell k's upstream area (determined from the flow direction raster) and is proportional to its length scale, $\sqrt{A_k}$ (where A_k is the area draining to cell k); its average slope, ${\stackrel{\mathrm{‾}}{S}}_k$; and its average weighting coefficient, ${\stackrel{\mathrm{‾}}{W}}_k$. On the other hand, the downslope component, D_dn,k, accounts for the effects along the topographically determined flow path between cell k and the stream network, obtained as the summation of runoff travel distances d_i (weighted by their respective coefficient, W_i, and slope, S_i) across cells k, …, i, …, n_k, moving from cell k down to the pour point where the runoff pathway eventually meets the stream network.

Figure 1 Methodological steps for obtaining the hydrologic-connectivity-based index of urbanization (HCIU): (a) scheme for calculating Borselli et al.'s (2008) connectivity index at generic cell k; (b) create a virtual, totally impervious copy of the basin, with the same shape, topography, and stream network, but different LULC, i.e., fully developed at all cells; (c) separately calculate the raster maps of connectivity for both the actual basin and its totally impervious copy; (d) calculate the raster map of normalized connectivity for the basin by dividing, on a cell-by-cell basis, the connectivity of the actual basin by the connectivity of the totally impervious copy; (e) assign a weight w_k to each basin cell k depending on its distance to the outlet, as measured along the stream network, starting from the cell's pour point; (f) calculate HCIU as the weighted average of the normalized connectivities at each basin cell.

In the literature, the cell weighting factor, W, depends on the type of analysis. Some examples include (1) the RUSLE C factor (Renard et al., 1997), i.e., a measure of the potential for erosion, adopted in sediment transport studies (Borselli et al., 2008); (2) measures of topographic roughness, often used for morphologic characterization (Cavalli and Marchi, 2008) or landslide risk assessment (Husic and Michalek, 2022); and (3) quantifications of the hydrologic characteristics of different LULC types, typically expressed as a function of Manning's surface roughness coefficient, n, to study, e.g., anthropogenic effects on landscape and sediment transport changes (Persichillo et al., 2018), landslide occurrence (Zanandrea et al., 2019), or runoff generation (Hooke et al., 2021).

In general, independently of the choice of W, the index of connectivity, IC_k (Eq. 1), is higher for cells with a larger and/or steeper contributing area, reflecting that such conditions are associated with the generation of potentially larger runoff volumes that can concentrate faster at cell k. Indeed, D_up,k shows some similarity with the well-known topographic wetness index (TWI; Beven and Kirkby, 1979; Riihimäki et al., 2021), often regarded as a proxy for soil moisture (Riihimäki et al., 2021); both are proportional to the size of the upstream contributing area, but TWI only considers the local (at-the-cell) slope instead of the average slope of the area draining to the cell. In turn, longer travel distances to the stream network and runoff paths with milder slopes both increase the value of the downslope component, D_dn,k (in the denominator of IC_k), resulting in reduced connectivity values, to reflect the lower potential of distant (upstream) hillslope locations to readily contribute the upland runoff (or sediment) to the stream network.

2.2 Recommended weighting coefficients for deriving HCIU

Among the options discussed above, when deriving HCIU, we recommend choosing W values that primarily depend on the LULC type of each basin cell, considering both developed (urbanized) and more natural (e.g., barren, cropland, and forested) categories, possibly differentiating across distinct intensities of land development and dominant vegetation types for the developed and vegetated categories, respectively. In this way, the effects of pixels with different surface characteristics can be differentially weighted depending on their potentials for either generating and quickly transmitting surface runoff (e.g., in the case of developed cells) or retaining, detaining, or infiltrating water (e.g., in the case of cells with vegetated land cover), depending on the distinct hydrologic dynamics associated with different LULC types. In turn, this incorporates sensitivity to the presence and spatial arrangement of the wide spectrum of LULC patches, from natural to fully developed. Specifically, land development will locally increase connectivity at those locations that receive runoff water from urbanized pixels upstream as well as those sectors that contribute runoff volumes to the stream network through impervious water pathways, proportionally to the intensity of land development; in a similar way, cells with more natural LULC types, characterized by a lower range for the weights (i.e., reflecting higher capacity for interception, detention, retention, infiltration, or some other runoff-mitigation mechanism), will locally decrease connectivity.

A candidate LULC-sensitive expression for the weighting coefficient W is the n-dependent function given by Eq. (2) (Persichillo et al., 2018; Zanandrea et al., 2019; Hooke et al., 2021), which assigns larger weights to urbanized cells, with small roughness coefficients, compared to vegetated (undisturbed) cells, which typically display greater roughness (Liu and De Smedt, 2004; Hooke et al., 2021).

$$\begin{array}{}\text{(2)}& W\left(n\right)=\mathrm{1}-n\end{array}$$

Another candidate could be a weighting specified as a function of the curve number (CN; Rallison, 1980), which accounts for LULC and soil characteristics with synergistic effects on surface runoff dynamics. In this work, we propose the expression W(CN) given by Eq. (3).

$$\begin{array}{}\text{(3)}& W\left(\mathrm{CN}\right)=\mathrm{CN}/\mathrm{100}\end{array}$$

In what follows, we assess both weighting approaches (i.e., based on either n or CN) for deriving a hydrologic-connectivity-based index of urbanization, indicating the resulting metrics as either HCIU(n) or HCIU(CN), respectively. A table of Manning's surface roughness coefficients associated with different LULC types, adapted from Hooke et al. (2021) and Liu and De Smedt (2004), is reported in Appendix B. Similarly, CN values adapted from Wu et al. (2024) are tabulated in Appendix C for different combinations of LULC types and hydrologic soil groups (HSGs; Rallison, 1980; Ross et al., 2018).

2.3 From distributed connectivity to a lumped hydrologic-connectivity-based index of urbanization HCIU

We seek to obtain a meaningful, lumped basin urbanization metric that conceptually encapsulates the impacts of the spatial arrangement of land development and other land-use/land-cover types into a single number, starting from a basin's connectivity map. To achieve this, we propose to first determine a relative measure of the effects of urbanization on the connectivity at each basin cell by normalizing with respect to a benchmark. This is done by computing the ratio between the connectivity of each cell in the actual basin and the connectivity at the same cell but for a virtual, totally developed copy of the watershed (see Fig. 1b and c), with the same shape, relief, and stream network but fully urbanized conditions. It is clear that connectivity in this virtual, “fully paved” basin takes the highest possible value at each cell for a given watershed's shape, topography, and stream network. To normalize consistently positive quantities, we suggest a change to the traditional connectivity given by Eq. (1) by computing the hydrologic connectivity index HCI_k of cell k simply as the ratio of the upslope over the downslope component (Eq. 4), without considering the logarithmic transform so that connectivity is always positive. HCI_k differs from the traditional IC_k (Eq. 1) only in the scale and sign but maintains all other properties of the original formulation, including the sensitivity to topographic characteristics and to the spatial arrangement of patches with different LULC types.

$$\begin{array}{}\text{(4)}& {\mathrm{HCI}}_k={\displaystyle \frac{D_{\mathrm{up},k}}{D_{\mathrm{dn},k}}}={\displaystyle \frac{{\stackrel{\mathrm{‾}}{W}}_k{\stackrel{\mathrm{‾}}{S}}_k\sqrt{A_k}}{{\sum }_{i=k}^{n_k}\frac{d_i}{W_i{S}_i}}}\end{array}$$

The imposed change in land cover in the totally developed benchmark basin involves forcing a Manning roughness value n of 0.02 (or a CN value of 99) and a consequent W(n) of 0.98 (or 0.99 in the case of W(CN)), at each hillslope cell. To distinguish the weighting coefficients of actual basin cells from the corresponding ones in its virtual, fully developed copy, we indicate the latter as W_imp in what follows (for impervious). The resulting connectivity map for the totally impervious basin (Fig. 1c) is a benchmark for the localized effects of the varied land-use/land-cover characteristics in the actual basin. The connectivity values for the totally developed watershed represent a theoretical upper limit for the level of connectivity that could be achieved at each pixel in the study watershed for its given fixed shape, topography, and stream network. Using the proposed, non-negative formulation of the connectivity index given by Eq. (4) ensures that the normalization is operated on non-negative values, so that the resulting normalized variable is in the interval (0, 1]. Equation (5) provides the expression of the normalized connectivity index, ${\widehat{\mathrm{HCI}}}_k$, for a reference cell k (see also Fig. 1d).

$$\begin{array}{}\text{(5)}& \begin{array}{rl}{\widehat{\mathrm{HCI}}}_k& ={\displaystyle \frac{{\mathrm{HCI}}_k}{{\mathrm{HCI}}_{\mathrm{imp},k}}}={\displaystyle \frac{{\stackrel{\mathrm{‾}}{W}}_k}{{\stackrel{\mathrm{‾}}{W}}_{k,\mathrm{imp}}}}\cdot {\displaystyle \frac{{\sum }_{i=k}^{n_k}\frac{d_i}{W_{i,\mathrm{imp}}{S}_i}}{{\sum }_{i=k}^{n_k}\frac{d_i}{W_i{S}_i}}}\\ & ={\displaystyle \frac{{\stackrel{\mathrm{‾}}{W}}_k}{W_{\mathrm{imp}}}}\cdot {\displaystyle \frac{\frac{\mathrm{1}}{W_{\mathrm{imp}}}{\sum }_{i=k}^{n_k}\frac{d_i}{S_i}}{{\sum }_{i=k}^{n_k}\frac{d_i}{W_i{S}_i}}}\end{array}\end{array}$$

In Eq. (5), ${\stackrel{\mathrm{‾}}{W}}_{k,\mathrm{imp}}$ and W_i,imp both refer to the totally impervious benchmark basin, indicating the mean weighting coefficient of the drainage area upstream of cell k, and the coefficient of the generic cell i downstream of k, respectively. Since all upstream cells will have a constant coefficient value W_imp in the impervious basin, ${\stackrel{\mathrm{‾}}{W}}_{k,\mathrm{imp}}$ is equal to W_imp for any k. Similarly, $W_{i,\mathrm{imp}}=W_{\mathrm{imp}}$ for any i.

The normalized connectivity index, ${\widehat{\mathrm{HCI}}}_k$, in Eq. (5) is equal to the product of the ratio of the average weighting coefficients from the upslope components, ${\stackrel{\mathrm{‾}}{W}}_k/W_{\mathrm{imp}}$, and the ratio of the weighted distances of reference cell k to its pour point along the stream network, $\frac{\mathrm{1}}{W_{\mathrm{imp}}}{\sum }_{i=k}^{n_k}\frac{d_i}{S_i}/{\sum }_{i=k}^{n_k}\frac{d_i}{W_i{S}_i}$, measured in the totally impervious and the actual basin, respectively. For both factors, the numerator is always smaller than or equal to the denominator. In what follows, to emphasize the distinction between HCI_k and ${\widehat{\mathrm{HCI}}}_k$ whenever these two are compared, we refer to the former as absolute connectivity instead of simply as connectivity, while the latter is termed normalized connectivity.

With the normalized connectivity map of a watershed, we can derive HCIU for the basin as a weighted average of the normalized connectivity indices at each pixel. We consider a weighted average instead of a straightforward arithmetic mean because there is still one aspect that the connectivity index does not account for by construction, i.e., the distance of each cell's pour point to the basin's outlet. As shown in Fig. 1a, the at-a-cell connectivity index considers the flow path from each cell k to its pour point along the stream network. However, different pour points along the network are located at varying distances from the basin outlet; hence, hillslope cells with similar levels of connectivity (e.g., because they have comparable drainage areas and characteristics and are located at similar distances from some segment of the stream network) but at different locations with respect to the basin's main channel and its outlet will display strong variability in their potential for quickly contributing runoff to the outlet because of the different timing. To account for these differences, we consider, for each hillslope cell k, its corresponding along-the-stream-network distance, d_SN,k, to the outlet, measured starting from the pour point, SN_k, of cell k (see Fig. 1e). We then assign a weight, w_k, to each hillslope cell k in a way such that pixels whose pour point is closer to the outlet receive a larger weight than those which are further away from it. In this work, we propose the weighting function $f_{\mathrm{w},\mathrm{HCIU}}\left(d_{\mathrm{SN},k}\right)$ given by Eq. (6), which assigns weights in the interval [0.5, 1]. Those basin cells pouring at stream locations with the minimum distance d_min to the outlet receive a weight of 1, while cells that drain to stream locations at the maximum distance d_max to the outlet get a weight of 0.5. All other basin cells with distances d_SN,k between d_min and d_max receive intermediate weights that vary linearly with d_SN,k in the interval (0.5, 1). It is worth noting that d_max varies across basins, while d_min, the distance between that channel pixel adjacent to the outlet and the outlet, is clearly equal to 1 pixel for any basin.

$$\begin{array}{}\text{(6)}& w_k=f_{\mathrm{w},\mathrm{HCIU}}\left(d_{\mathrm{SN},k}\right)=\mathrm{1}-\mathrm{0.5}{\displaystyle \frac{d_{\mathrm{SN},k}-d_{\min }}{d_{\max }-d_{\min }}}\end{array}$$

In Eq. (6), the weight assigned to the generic cell k is indicated as w_k. Once a weight is assigned to each cell, depending on its along-the-stream-network distance, d_SN,k, to the outlet, the lumped hydrologic-connectivity-based index of urbanization (HCIU) is obtained as the weighted average of the normalized connectivity indices of each basin cell using Eq. (7), where the normalized connectivity of the generic cell k is indicated as ${\widehat{\mathrm{HCI}}}_k$ (see also Fig. 1f).

$$\begin{array}{}\text{(7)}& \mathrm{HCIU}={\displaystyle \frac{{\sum }_k{w}_k{\widehat{\mathrm{HCI}}}_k}{{\sum }_k{w}_k}}\end{array}$$

In summary, the proposed methodology provides a lumped metric (HCIU) that is able to conceptually capture the varied hydrologic effects arising from the spatial arrangement of different LULC patches, both natural and developed, depending on their relative location with respect to each other, the stream network, and the basin outlet. First, hillslope-to-stream connectivities, weighted depending on the hydrologic effects of distinct LULC types, are normalized with respect to a fully impervious benchmark (Fig. 1a, b, c, and d), which allows for the comparison between the effects of heterogeneous levels of urbanization both across and within basins. Then, HCIU is obtained as a weighted average of normalized connectivities across the entire watershed, assigning different weights to each pixel depending on the along-the-stream-network distance of that cell's pour point to the basin outlet (Fig. 1e and f).

The proposed two-step formulation – where the flow paths of hillslope cells to the pour points along the stream network and then the distances of those pour points to the basin outlet are considered separately – is different from other established, outlet-focused applications of the connectivity index, such as the IC_outlet distributed metric proposed by Cavalli et al. (2013). The latter is calculated following Borselli et al. (2008; with some adaptations to the weighting coefficient and the flow direction algorithm) but considering flow paths all the way to the outlet (hence considering both overland flows and subsequent channelized flows within the same path) instead of flow paths to the closest stream link, following only hillslope surfaces. The two main components of a basin's hydrologic response, i.e., overland and channel flow, generally involve quite different temporal scales because of the different orders of magnitude in roughness and water depths. The IC_outlet metric is able to capture these differences as IC_outlet raster maps typically exhibit the highest connectivity values along the watershed stream network (comparable only to connectivities in the hillslope sectors closest to the outlet), followed by connectivities in zero-order valleys or hollows adjacent to channels (Cavalli et al., 2013). On the other hand, the focus of our methodology is on the hydrologic effects of land development, which mostly influences the surface and near-surface components of the basin response by locally decreasing infiltration and increasing runoff speeds. Considering only the hillslope-to-stream connectivity in our first step allows us to enhance the method's sensitivity to the effects of land development on hydrologic response by focusing on how runoff interacts with the distinct LULC patches encountered along the hillslope path, which control (i.e., enhance or mitigate) the connectivity. Once runoff reaches the stream network the effects of travel distance along the stream network must still be accounted for, but this is performed in the separate, second step, considering a narrower range for the weights. This ensures that HCIU displays adequate sensitivity to urbanized sectors that are adjacent to the stream network, but at reaches located far upstream of the outlet.

Breaking down the calculations for HCIU in two parts (the hillslope-to-stream and then stream-to-outlet flow paths) also presents a practical advantage, particularly for large-scale implementation of the index. To ensure broad applicability of the proposed methodology, we need to be able to quickly compute HCIU for any basin (in a region, country, province, state, etc.) as selected by the final user. If we were to use a cell-to-outlet scheme, such as the IC_outlet metric, we would need to recompute everything from scratch every time a user chooses a different basin (i.e., a different outlet location along the stream network). Splitting the computations from cell to pour point and then pour point to outlet offers the opportunity to precompute “static” raster maps (i.e., independent of outlet location) of connectivity and normalized connectivity for all the pixels over large areas. In this way, later, when a user selects a specific outlet location, the final computation of HCIU only involves the much quicker weighted averaging of the precomputed, at-a-cell, normalized connectivities, only considering those cells within the selected basin and their along-the-stream-network distances to its outlet.

3 Data and case studies for testing HCIU

3.1 Hydrologically homogenous regions

To test the proposed connectivity-based urbanization index HCIU against the traditional TIA fraction as a predictive variable in regional peak-flow equations, we considered three distinct hydrologically homogeneous regions, as determined by the U.S. Geological Survey (Southard, 2010; Austin, 2014; Feaster et al., 2014; see Fig. 2). One homogenous region encompasses all major metropolitan areas in Missouri (MO) and surroundings, including 34 urbanized watersheds (Southard, 2010). The second involves urban centers in Virginia (VA), with a total of 112 developed basins (Austin, 2014). The third homogenous region (Feaster et al., 2014) is the largest, spanning parts of three states – Georgia and North and South Carolina; it includes 79 urbanized watersheds, encompassing the Piedmont and part of the Ridge and Valley ecoregions defined by the US Environmental Protection Agency (USEPA), resulting in a long band of land, moving from Georgia to North Carolina, with consistent flood frequency characteristics (Feaster et al., 2014). For brevity, in what follows, we refer to this third homogeneous region as EPA ecoregion (or EPAE). The VA and EPAE case studies only include basins with at least 10 % of TIA (Austin, 2014; Feaster et al., 2014), while the MO study considers a lower threshold, with all basins above 5 % TIA except for one with only 2.33 % (Southard, 2010). TIA values for all case-study basins are reported in Appendix A.

Figure 2 Case-study regions (basin locations are reported in Appendix E; map created with ArcGIS Pro; ESRI, 2024).

For each basin in their case-study regions, Southard (2010), Austin (2014), and Feaster et al. (2014) extracted annual maxima series from instantaneous discharge records (typically with 15 min or hourly temporal resolution) and performed flood frequency analyses to estimate peak-flow values for a range of return periods (see Appendix A) following the US national guidelines provided in Bulletin 17B (Interagency Advisory Committee on Water Data, 1982). Southard (2010) and Austin (2014) obtained the fraction of TIA from the National Land Cover Database (NLCD) 2001 (Homer et al., 2020) for the MO and VA case studies, respectively. On the other hand, for the EPA ecoregion, Feaster et al. (2014) considered the 2006 version.

Regional peak-flow equations proposed for these regions all include basin area, A, and the percentage of TIA (simply referred to as TIA in what follows) as explanatory variables for predicting the magnitude Q_T of the flood with return period T. Southard (2010) and Feaster et al. (2014) adopted the functional form given by Eq. (8), with TIA as the generic urbanization metric U, while Austin (2014) considered a different form, where the peak flow per unit area ($Q_T/A$) is modeled as a function of A and TIA. For convenience and consistency, in this work, we systematically consider the simple linear model given by Eq. (8) to test the predictive power of the hydrologic-connectivity-based index of urbanization (HCIU) against the traditional TIA. We alternatively use TIA or HCIU as the generic urbanization metric U in Eq. (8), each time fitting the regression model on A–TIA–Q_T and A–HCIU–Q_T data, respectively. The explanatory power of these competing variables would change if we considered other functional dependencies, but from a qualitative standpoint, the superiority of one variable over the other to explain peak flows should not be affected by this change.

$$\begin{array}{}\text{(8)}& \log Q_T={\mathit{\beta }}_{\mathrm{0}}+{\mathit{\beta }}_{\mathrm{1}}\log A+{\mathit{\beta }}_{\mathrm{2}}U.\end{array}$$

Figure 3 plots all the A–TIA pairs for the three case studies, considering basin areas in logarithmic scale. The datasets span 4 orders of magnitude of basin sizes and a wide range of land-development conditions. A table with all the basins considered in the different case-study regions, and related information, is provided in Appendix A.

Figure 3 Case-study basins from the EPAE, MO, and VA homogeneous regions, characterized by their area (A) and percentage of total impervious area (TIA). Areas are plotted in logarithmic scale.

3.2 DEM, LULC, and HSG data

For consistency in the comparisons across the different case studies, we used the same DEM and LULC data across all homogenous regions. Because of its fine resolution and thorough coverage across the US territory, we selected the $\mathrm{1}/\mathrm{3}$rd arcsecond DEMs by the USGS (U.S. Geological Survey, 2023), while for LULC, we adopted the same NLCD maps (Homer et al., 2020) as in the original studies to obtain TIA, thus using the 2001 version for the MO (Southard, 2010) and VA (Austin, 2014) case studies and the 2006 version for the EPAE case study (Feaster et al., 2014). The original NLCD maps do not distinguish between needleleaf and broadleaf dominant species within the forest categories even though these two types of tree cover have very different runoff retention capabilities and should thus be modeled using different Manning roughness coefficients (Liu and De Smedt, 2004). This information was obtained from the global, 300 m resolution LULC maps produced by the European Space Agency (ESA, 2017) by overlapping information from the two sources. Making this distinction results in the expansion of the two original NLCD forested LULC categories evergreen and deciduous into four classes: deciduous needleleaf, deciduous broadleaf, evergreen needleleaf, and evergreen broadleaf forest. To generate the curve-number-based urbanization metric, HCIU(CN), for our case studies, we also needed information about the hydrologic soil group for each basin cell, for which we considered the HYSOGs250m global dataset (Ross et al., 2018), providing worldwide HSGs on a 250 m grid.

Information from the DEM, (expanded) LULC, and HSG maps was extracted for each case-study basin using the watershed boundaries obtained from the USGS data repository by Krstolic (2006) for the VA case study and from the StreamStats web application developed and maintained by the USGS (U.S. Geological Survey, 2019), for the MO and EPAE case studies. Then, LULC types and LULC–HSG pairs at each basin pixel were mapped into values of Manning's roughness coefficient, n, and curve number, CN, respectively, following the tables given in Appendices B (adapted from Liu and De Smedt, 2004; Hooke et al., 2021) and C (adapted from Wu et al., 2024), respectively.

For each basin, the stream network was obtained from the flow accumulation raster (using the D8 algorithm) by setting a minimum threshold for the number of upstream cells. Instead of following the traditional procedure of considering a constant threshold for the whole basin (Tarboton and Ames, 2001), we locally selected the threshold for each headwater to closely match its location as per the National Land Cover Database Plus High Resolution (NHDPlus HR; Moore et al., 2019), outlining a digital stream network as similar as possible to the official blue lines provided in that dataset. We found this approach a preferable alternative to DEM-burning procedures that would enforce the NHDPlus stream network onto the DEM by locally lowering the elevations and introducing preferential flow directions (Getirana et al., 2009) because our method is better at preserving the actual connectivity patterns between the hillslope component and the stream network. DEM-burning algorithms, on the other hand, could have led to issues such as the occurrence of parallel streams arising from the non-alignment of the original DEM channels and the corresponding network links coming from an external source (Lindsay, 2016), with the consequent risk of diffused runoff pathway disruptions or alterations before reaching the stream.

Figure 4 shows the frequency distribution of the Manning surface roughness coefficient, n (Fig. 4a, b, and c); the CN (Fig. 4d, e, and f); and the slope (Fig. 4g, h, and i), observed at both pixel and basin scales (considering basin averages in the latter case), for the three case-study regions. Based on the functional dependency of HCI_k outlined by Eq. (4), the connectivity-based urbanization index HCIU, in its two versions HCIU(n) and HCIU(CN), is expected to display sensitivity to those distributed basin characteristics. We note that the VA case study presents a wider spectrum of n values (Fig. 4c) compared to the other two (Fig. 4a and b), which is imputable to differences in the vegetated areas, with VA being more dominated by the presence of broadleaf tree species, both deciduous and evergreen. While VA also presents a higher frequency of low-CN-value pixels (Fig. 4f), compared to the other two case studies (Fig. 4d and e), the difference in the spread of CN values is not as noticeable as is the case with n. This is because, in contrast to our n values, the CN for different combinations of LULC types and hydrologic soil groups do not weigh broadleaf and needleleaf tree vegetated areas differently (Wu et al., 2024). The VA case study also differs from the other two in its relief characteristics, with a wider range of point and basin-averaged slopes (Fig. 4i compared with Fig. 4g and h), associated with basins from the Appalachian region.

Figure 4 Frequency distribution of the Manning surface roughness coefficient, n (a, b, c); curve number, CN (d, e, f); and slope, S (g, h, i), for the three case studies as observed at each basin pixel across all watersheds (filled bars) and as basin averages (empty bars). (j, k, l) Mix of the four NLCD developed LULC types (associated with different extents of impervious area, IA) for all basins in the three case studies, respectively, expressed as percentages of the total extent of developed areas.

Figure 4 also illustrates the mix of developed LULC types in the basins by showing the distributions (box plots) of the extents of the four developed NLCD categories in each watershed for the three homogenous regions (Fig. 4j, k, and l, respectively). Those categories include developed, open space; developed, low intensity; developed, medium intensity; and developed, high intensity, associated with ranges of impervious area of less than 20 %, 20 %–49 %, 50 %–79 %, and 80 % or more, respectively. For each watershed, Fig. 4 provides developed LULC extents as percentages of the number of basin pixels with a given developed LULC category (say, low intensity) with respect to the total number of developed LULC cells in that watershed (i.e., the sum of the numbers of developed, open-space; low-intensity; medium-intensity; and high-intensity cells).

The four urbanized LULC categories are associated with distinct levels of imperviousness. For instance, open space is the least developed LULC category, with each pixel having less than 20 % impervious area (IA), while highly urbanized cells, with impervious areas from 80 % to 100 %, fall in the high-intensity category. Because of these differences, moderately developed areas contribute less to the overall TIA of a basin compared to highly impervious ones. In other words, larger areas of developed, open-space and developed, low-intensity patches are needed compared to developed, medium-intensity and developed, high-intensity ones, to contribute the same proportion of TIA in a watershed.

It is evident from Fig. 4j, k, and l that each region has its own characteristic mix of basin urbanized areas, even though there is a clear dominance of moderately developed LULC categories (i.e., open-space and low-intensity) over the more urbanized types, across all the case studies. VA (Fig. 4l) is the region with the smallest relative extents of highly urbanized areas (with respect to the total extent of land development), with most basins having less than 10 % high-intensity pixels (i.e., with IA between 80 % and 100 %) and wider (relative) extents of open-space areas (i.e., with IA less than 20 %) compared to basins from the other two case studies (Fig. 4j and k). MO, on the other hand (Fig. 4k), contains many watersheds with more concentrated urbanized areas, as indicated by the higher proportion of both medium- and high-intensity areas and low-urbanization environments more dominated by low-intensity (i.e., with IA between 20 % and 50 %) than open-space areas (i.e., with IA below 20 %). The EPAE region (Fig. 4j) displays intermediate conditions between those observed for VA (Fig. 4l) and MO (Fig. 4k), in terms of both the more urbanized (i.e., medium- and high-intensity) and the less urbanized LULC (i.e., open-space and low-intensity) categories.

4 Results

4.1 Interpretation of the intermediate raster data products

To get a qualitative understanding of how the new connectivity-based urbanization index works and examine differences between the n-based and CN-based formulations, Fig. 5 shows the intermediate raster products generated by applying the proposed methodology at one of the MO watersheds. Basin 06894000 (Fig. 5a) is characterized by quite a heterogeneous mosaic of land-cover patches, with many parks and forested areas as well as neighborhoods and parking lots. A zoom on a portion of the basin, in the circular window in Fig. 5a (obtained from OpenStreetMap; OpenStreetMap contributors, 2015), offers a more detailed glimpse into the variety of urbanized and vegetated areas.

Figure 5 (a) Basin 06894000 (MO); (b) raster map of slope, S; (c) raster map of the absolute hydrologic connectivity index, HCI(n), obtained as a function of the Manning roughness coefficient, (n); (d) raster map of the normalized connectivity index, $\widehat{\mathrm{HCI}}\left(n\right)$, obtained as a function of the Manning roughness coefficient (n); (e) raster map of the absolute hydrologic connectivity index, HCI(CN), obtained as a function of the curve number (CN); (f) raster map of the normalized connectivity index, $\widehat{\mathrm{HCI}}$(CN), obtained as a function of the curve number (CN). Background map in panel (a) was retrieved from OpenStreetMap (OpenStreetMap contributors, 2015). © OpenStreetMap contributors 2015. Distributed under the Open Data Commons Open Database License (ODbL) v1.0.

Figure 5c and e show the HCI raster maps for basin 06894000, obtained from the n- and CN-based formulations, respectively. The pairs of circular windows in Fig. 5c and e both depict the same enlarged portion of absolute connectivity maps obtained for the actual basin (HCI) and the totally developed virtual copy (HCI_imp), respectively. In all cases, absolute connectivity increases for cells located closer to the stream network. Slope also controls both formulations of HCI, with higher gradients of absolute connectivity observed in steeper valleys compared to flatter riverine zones (see Fig. 5b). This is clear when looking at the floodplain for the basin's main channel (the darker buffer area around the main channel in Fig. 5b), which is characterized by lower absolute connectivities than many of the riverine zones of its tributaries, with steeper slopes.

The raster maps of normalized connectivity $\widehat{\mathrm{HCI}}$ (Fig. 5d and f), on the other hand, display low sensitivity to topographic gradients, as there is no correlation with local spatial patterns in slopes (Fig. 5b). For instance, the floodplain of the main channel displays higher normalized connectivity ($\widehat{\mathrm{HCI}}\left(n\right)$) than many of the tributary valleys even though the latter show higher absolute connectivity (HCI(n); compare Fig. 5c and d). This is because the n-based normalized connectivity is most sensitive to differences in land-cover, as captured by W(n) (see also Eq. 5). This is easy to notice from the $\widehat{\mathrm{HCI}}\left(n\right)$ raster (Fig. 5d), where cells with lower normalized connectivities coarsely correspond to green patches in the map of Fig. 5a, while cells with high values of $\widehat{\mathrm{HCI}}\left(n\right)$ are typically associated with developed areas of the basin.

While spatial patterns in $\widehat{\mathrm{HCI}}\left(n\right)$ are easy to interpret retrospectively because n only reflects differences in land-cover types, spatial patterns in $\widehat{\mathrm{HCI}}$(CN) (Fig. 5f) correlate less with the spatial arrangement of land-cover types. This is expected as W(CN) reflects differences in not only LULC types, but also hydrologic soil groups. Depending on the specific soil group, some pervious cells with lower soil infiltrability can get CN values as high as those of urbanized cells (Wu et al., 2024), which may result in similar normalized connectivity values for some of the pervious patches compared to impervious sectors.

In summary, the two types of normalized connectivity raster maps, $\widehat{\mathrm{HCI}}\left(n\right)$ and $\widehat{\mathrm{HCI}}\left(\mathrm{CN}\right)$ (Fig. 5d and f, respectively), may look very different from each other even though the absolute connectivity maps, HCI(n) and HCI(CN) (Fig. 5c and e, respectively), present similar spatial patterns. This is because local topography and the proximity to the stream network are both strong controls for the absolute connectivity but not necessarily so for the normalized connectivity.

4.2 Conceptual differences between HCIU and TIA

The new proposed urbanization metric expresses the relative degree of hydrologic connectivity – with respect to the fully developed benchmark – that arises from the spatial arrangement of the spectrum of developed and undisturbed LULC patches present in a basin, with different hydrologic characteristics, contingent upon the topographic structure of the hillslopes and stream network. As such, it is expected to index the effects of land development differently compared to TIA, which is simply the proportion of impervious surfaces. This means that different basins with similar percentages of TIA will display variability in their associated HCIU values because of the heterogeneous effects on connectivity (as captured by W) of distinct spatial arrangements of the urban sectors (with different land-development intensities) and the undeveloped patches in the watershed. Figure 6 shows, in two separate subplots (a and b), the variability in HCIU obtained for our case studies for various levels of TIA, considering HCIU(n) and HCIU(CN), respectively. We include example basins with distinct spatial arrangements of natural and developed areas (Fig. 6c to i) to clearly illustrate the sensitivity of HCIU to LULC configurations. For the sake of simplicity, for this visualization, we adopt a simplified, four-category LULC classification: highly impervious, moderately impervious, moderately pervious, and highly pervious. We label developed pixels with impervious areas equal to or above 50 % (i.e., the NLCD medium and high-intensity classes) as highly impervious, while those with imperviousness below 50 % (i.e., the NLCD open-space and low-intensity classes) are classified as moderately impervious. On the other hand, the moderately pervious category encompasses those undeveloped cells with moderate values of Manning's coefficient, n, equal to or below 0.3 (see Appendix B), such as herbaceous, hay/pasture, and barren land LULC, while the highly pervious category covers the remaining undeveloped LULC types, including forested and other densely vegetated areas.

Figure 6 Case-study basins from the EPAE, MO, and VA homogenous regions in the TIA–HCIU plane, considering (a) HCIU(n) and (b) HCIU(CN), respectively. Seven specific basins (c–i) are used as examples to visualize the sensitivity of HCIU to different spatial arrangements of LULC patches even for basins with similar TIA; a simplified, four-category LULC classification is adopted as indicated in the legend.

In general, a wide range of HCIU values is associated with any given TIA across all case studies; this is expected since similar relative extents of impervious areas can have different effects on hydrologic connectivity depending on their spatial arrangement with respect to each other and also with respect to the less developed and undeveloped patches (Fig. 6a and b). However, we also note general positive trends, both within and across all the case studies considered, suggesting that basins with an increasingly larger TIA are characterized by overall higher HCIU values, as should also be expected.

Figure 6a and b also show that HCIU(n) and HCIU(CN) have different ranges of variability. Both have a lower bound greater than 0; this is intrinsic to the formulation of HCIU itself as a (weighted) average of the relative degree of connectivity of the watershed. The absolute connectivity of a basin with totally pervious land-cover type, which can be regarded as a lower-bound yet realistic scenario, is small (or very small) compared to the connectivity of the corresponding totally impervious virtual basin, but it is not zero. Therefore, the fact that HCIU's lower bound is greater than zero comes from physical (or at least conceptual) considerations, depending on the range of variability in the weighting coefficient, W, i.e., specifically on its minimum value (see Eq. 5). Since HCIU(n) and HCIU(CN) use different weighting coefficients, W(n) and W(CN), with different lower limits, their lower bounds will differ.

Figure 6a and b display another interesting characteristic of HCIU, i.e., its sensitivity to the heterogeneity in the undeveloped LULC types in a basin. Compared to the other two case studies, VA watersheds have overall smaller HCIU values for similar levels of TIA. This is because VA basins are characterized by pervious land-cover types with generally greater flood-mitigating capabilities (larger n values and smaller CN values), as indicated by the wider spread of Manning's coefficients in the range of small values (Fig. 4c), as well as the stronger positive skewness in the frequency distribution of CN values (Fig. 4f) compared to the other two case studies (Fig. 4a–b, and d–e, respectively). As a result, VA basins have lower HCIU values compared to those from the MO and EPAE case studies with similar relative extents of developed areas because the pervious fraction in VA basins is more effective in reducing the hydrologic connectivity of the urbanized patches. This effect is more evident with the n-based formulation for HCIU, where a clear separation also exists for the smallest TIA values. This is modulated by the diverse ways in which the effects on hydrologic connectivity of distinct undeveloped LULC types are differentially weighted by W(n) and W(CN).

Distinct spatial arrangements of pervious and impervious patches (see examples in Fig. 6c to i) get HCIU(n) and HCIU(CN) values (Fig. 6a and b, respectively) in agreement with our expectations, irrespective of the weighting approach (based on either n or CN). For instance, focusing on watersheds with similar TIA levels (around 30 %), basin 01645975 (Fig. 6h), with predominant highly pervious land cover and most urbanized areas concentrated in a limited area far upstream of its outlet, has significantly lower HCIU values (on both the HCIU(n) and HCIU(CN) scales) than basins 01657000 (Fig. 6f) and 02203800 (Fig. 6c), both of which display spatial arrangements of impervious areas that would be expected to cause stronger hydrologic effects.

Basin 01657000 (Fig. 6f) has a proportion of pervious and impervious LULC patches similar to basin 01645975 (Fig. 6h), but its urbanized areas (including a considerable extent of highly impervious patches) are all concentrated downstream, where they can contribute impervious runoff more effectively to the outlet, hence having a stronger impact on the overall hydrologic connectivity of the watershed. Basin 02203800 (Fig. 6c) shows even worse conditions, with a significantly larger proportion of moderately impervious areas than pervious patches, which impact a wider extent of the watershed, including areas near its outlet. These LULC conditions, even though they are associated with locally lower levels of imperviousness compared with highly impervious but spatially concentrated areas of basin 01657000 (Fig. 6f), systematically increase distributed connectivity, hence decreasing the overall response time of the watershed. Consistently, basin 02203800 (Fig. 6c) gets higher HCIU(n) and HCIU(CN) values than basin 01657000 (Fig. 6f).

These comparisons indicate that large developed areas, even with locally low levels of imperviousness, can have stronger effects on HCIU compared to highly urbanized but spatially concentrated patches. Because of this, basins that get the highest HCIU values are those where developed patches are dominant and uniformly spread. Basins 07005000 (Fig. 6d) and 01644291 (Fig. 6e), both in the higher end of the TIA–HCIU plots (Fig. 6a and b), exemplify this case, with their surfaces entirely covered by a mix of moderately and highly impervious land cover and negligible presence of moderately or highly pervious patches.

The stronger effects of large extents of moderately impervious patches on HCIU, compared with highly urbanized but concentrated areas, hold true also when comparing watersheds with different TIA levels. For instance, basin 02204070 (Fig. 6g) has a TIA of 20 %, lower than the 30 % TIA of 01657000 (Fig. 6f) and 01645975 (Fig. 6h); yet the HCIU value of the former is higher than the HCIU values of the latter two for both n- and CN-based formulations. This is, again, because the spatial arrangement of its impervious patches covers a wider (although less dense) spatial extent, resulting in a more widespread increase in connectivity.

The examples above indicate that HCIU displays sensitivity to the spatial arrangement of the heterogeneous LULC patches in a watershed irrespective of basin size and shape. However, HCIU can also capture similarities in these spatial arrangements across basins with different shapes and sizes. For instance, basins 02203800 (Fig. 6c) and 02336360 (Fig. 6i), with areas of 108 and 70 km², respectively, both display a preponderance of moderately impervious LULC types, with few, more intensely developed patches located far upstream with respect to the outlet and a finely dispersed pervious fraction that tends to condense into larger pervious patches moving towards the outlet. Because of these similarities in the LULC spatial configurations, these two basins get similar HCIU values, when using both HCIU(n) and HCIU(CN), despite their differences in size and shape.

Two highly urbanized VA basins, 01644290 and 01644291 (Fig. 6e), have unusually higher HCIU(CN) values compared to all others (Fig. 6b), because of their peculiar LULC conditions. Both are very small in size (< 0.20 km²) and are almost totally characterized by highly impervious LULC (developed, high-intensity NLCD type). Furthermore, the few, small undeveloped patches are unable to mitigate the connectivity of the developed areas since they are located far upstream, close to the water divide (see Fig. 6e for basin 01644291). Both aspects lead to LULC conditions that are almost identical to those for the virtual, totally impervious scenario considered for normalizing absolute connectivity, hence justifying the large HCIU(CN) values for these two watersheds, which are the highest among all basins and case studies. These two watersheds are outliers on the set of HCIU(CN) values for VA (Fig. 6b), while they are the highest on the TIA–HCIU(n) plane (Fig. 6a) but still clustered with other VA basins with high TIA.

4.3 Performance of HCIU in regional peak-flow equations

We tested the predictive power of HCIU as an alternative to TIA, to be used as urbanization metric in the development of regional peak-flow equations for the three case studies. Figure 7 shows the performance of the regional model (Eq. 8) fitted to basin area A, HCIU, and a range of flood quantiles Q_T (with a return period T of 2, 5, 10, 25, 50, 100, and 500 years; see Appendix A), comparing it to the benchmark model with the same functional dependency, but with TIA as the urbanization metric. The adjusted R² is considered the error metric.

Figure 7 Comparison of the performance of regional peak-flow equations calibrated on (1) A, TIA, and Q_T data (i.e., $Q_T\sim A+\mathrm{TIA}$ benchmark model); (2) A, HCIU(n), and Q_T data ($Q_T\sim A+\mathrm{HCIU}\left(n\right)$ model); and (3) A, HCIU(CN), and Q_T data ($Q_T\sim A+\mathrm{HCIU}\left(\mathrm{CN}\right)$ model) for the (a) EPAE, (b) MO, and (c) VA case studies. Quantile values Q_T associated with a return period T of 2, 5, 10, 25, 50, 100, and 500 years are considered (data in Appendix A).

Our results indicate that not only is HCIU a strong peak-flow predictor in combination with A, but it also systematically outperforms TIA when the n-based formulation is considered. Improvement is strongest for the EPAE and MO cases studies, while for VA, where the performance of the benchmark model (with A and TIA) was already the highest among case studies, only a marginal gain is obtained. On the other hand, the CN-based formulation for HCIU seems to be overall less robust, with varying behaviors depending on the specific case study. For the EPAE homogeneous region, HCIU(CN) outperforms the benchmark (but not its n-based counterpart) when fitting extreme flood values with return periods between 10 and 100 years, while it slightly underperforms for other flood quantiles. On the other hand, for the MO case study, HCIU(CN) displays a noticeably worse performance compared to the benchmark, except for the two largest flood quantiles (Q₁₀₀ and Q₅₀₀). VA is the only region where the CN-based HCIU performs similarly to the n-based HCIU, systematically outperforming the benchmark, even though only marginally, as all models perform well. In the Discussion section we hypothesize about what might explain the lower performance of HCIU(CN) compared to HCIU(n).

Figure 7 also shows that, regardless of the urbanization metric considered, model performance decreases with increasing return period (T) of the flood quantile. This may be due to a combination of two aspects. First, for such extreme events, differences across distinct LULC patches in the basin become increasingly negligible from a hydrologic perspective (Ogden et al., 2011) as wetter pervious patches infiltrate a smaller fraction of precipitation, causing any land-cover descriptor to lose predictive power. The second reason lies in the inevitable uncertainties associated with the estimation of such extreme quantiles from short flow records (Klemeš, 2000), which means that the models are fitted to highly uncertain data points. Because of these reasons, we suggest that improvements in prediction accuracy of urbanization metrics should be only pursued in the range of more frequent floods (say, below the 100-year return period). A generalized decrease in model performance is also observed when moving from intermediate to smaller quantiles, except for TIA in VA. This trend may be due to the increasing influence of the minor drainage system on the hydrologic response during smaller events, overshadowing surface runoff dynamics. However, both TIA and HCIU primarily focus on aspects related to surface runoff.

The main hydrologic application of regional peak-flow equations fitted to data from gauged basins is to extrapolate the relationships to other, ungauged basins. Therefore, a more informative way of testing the two competing urbanization metrics is to consider a k-fold validation framework, where the full dataset is split into a training set and a test set for fitting and evaluating the model, respectively. Following an approach similar to Dell'Aira et al. (2022), we produced a distribution of test errors (shown in Fig. 8) by repeatedly fitting the peak-flow equation (Eq. 8) and evaluating model performance using separate subsets of the full dataset. For each case study, we alternately considered a 66 % : 33 %, 75 % : 25 %, and 80 % : 20 % proportion for the training- and test-set sizes, respectively, associated with three, four, and five alternative training–test subset splits, respectively, in a way such that all the obtained test sets together span the full dataset without having duplicate basins. For each of these proportions, we repeated the k-fold procedure 10 times, considering 10 different random samplings to populate the folds, to avoid any potential bias associated with a single sampling. This procedure resulted in 120 blind assessments of the performance of the peak-flow equations, each time fitted and tested on different subsets for each case study.

Figure 8 Box plots of test errors for (1) the benchmark $Q_T\sim A+\mathrm{TIA}$, (2) $Q_T\sim A+\mathrm{HCIU}\left(n\right)$, and (3) $Q_T\sim A+\mathrm{HCIU}\left(\mathrm{CN}\right)$ models, when they are fitted and blind-tested on distinct basin subsets, for the (a) EPAE and (b) MO homogenous regions, respectively (see this same plots for VA in Appendix D). For each box plot, filled bars represent error values between the first and third quartiles, the upper (lower) whisker extends from the third (first) quartile by adding (subtracting) 1.5 times the interquartile range, and any outliers beyond the whiskers are marked by circles.

The box plots of errors in Fig. 8 confirm that HCIU(n) is the more robust urbanization metric for peak-flow regression equations. For all flood quantiles, HCIU(n) improves model performance (compared to the $Q_T\sim A+\mathrm{TIA}$ benchmark) both in average terms (as indicated by the systematically higher medians and the higher location of the box plots) and with respect to the error spread, which is narrower for all EPAE (Fig. 8a) and MO (Fig. 8b) quantiles, with the only exception being Q₂ for EPAE (Fig. 8a).

On the other hand, HCIU(CN) exhibits again a more heterogeneous performance. In the EPAE case study (Fig. 8a), the error box plots for $Q_T\sim A+\mathrm{HCIU}\left(\mathrm{CN}\right)$ are marginally better than those for the benchmark $Q_T\sim A+\mathrm{TIA}$ across all quantiles except the two smallest (i.e., Q₂, and Q₅), with slightly lower distributions of adjusted R² (in agreement with what was observed in Fig. 7). A different picture is obtained for the MO case (Fig. 8b), where HCIU(CN) significantly underperforms the $Q_T\sim A+\mathrm{TIA}$ benchmark for most quantiles, indicating that HCIU(CN) is not a robust urbanization metric for the MO region.

The same box plots for VA, moved to Appendix D since they do not add any relevant additional information, indicate comparable results as those already observed when testing the new urbanization metric (with both n- and CN-based weighting approaches) on the full dataset (Fig. 7c). Specifically, both HCIU(n) and HCIU(CN) lead to improvements in model performance as their related error box plots are shifted upward compared to the benchmark.

5 Discussion

Our results indicate that the proposed conceptual framework for deriving hydrologic-connectivity-based urbanization metrics does produce lumped basin descriptors that successfully encapsulate information about the flood-enhancing impacts of urban sectors with different land-development intensities (and thus with different imperviousness levels), considering their spatial arrangement, but is also sensitive to the spectrum of mitigation effects afforded by undeveloped (pervious) patches, depending on their relative location within the watershed. This is obtained by differentially weighting the effects of different LULC types through widely used and accepted conceptual criteria such as Manning's surface roughness coefficients and curve number values, hence capturing heterogeneity in the capacity of LULC patches to both generate and retain/detain runoff.

The resulting numerical value of HCIU for a single watershed represents a measure of the proportion of hydrologic connectivity arising from the specific mosaic of land-cover patches in that watershed relative to the maximum theoretical connectivity of that basin if its surface were completely paved. While this definition may not be as straightforward as that for TIA, which is simply a proportion of impervious areas, it results in a conceptually more comprehensive and hydrologically driven representation of the distributed impacts of urbanization on surface runoff dynamics in a lumped basin descriptor. Like TIA, HCIU can be used to characterize and compare different basins to either simply determine which basin is more impacted by land development or develop regional models. HCIU can also be utilized for planning to compare the expected hydrologic effects of different scenarios for land development in a given watershed.

We suggest that HCIU should also increase our explanatory power when predicting other event-related variables such as lag times and times of concentration. HCIU is indeed sensitive to not only the presence and spatial arrangement of LULC patches with different hydrologic characteristics, but also those locations where flows tend to concentrate, locally decreasing surface runoff travel times, as conceptually reflected in the upslope component, D_up,k (Eq. 4). HCIU also considers the distance of these surface runoff hotspots, where stormwater tends to concentrate and travel faster, to the stream network, as reflected by the downslope component D_dn,k. This in turn determines how easily those locations with accumulating flows will contribute to the overall basin response. Ultimately, HCIU conceptually summarizes in a single number the effects of all potential runoff travel paths occurring on the basin surface moving towards the stream, including interactions among converging surface flow paths, following a hydrologically driven approach. Because other response variables, such as lag time and time of concentration, are emergent basin properties arising from the interactions of all individual travel paths, their correlations with HCIU or other connectivity-based descriptors should be investigated in future research. The HCIU approach could also be further tested against TIA for predicting other hydrologic variables that are affected by land development, such as water-quality indicators. Depending on the specific application, it may be necessary to make some adaptations to account for additional sources of connectivity, induced by, e.g., the underground stormwater sewer infrastructure.

In the next subsection, we discuss in more detail the possible reasons for the poor performance of the CN-based HCIU. We then highlight advantages and disadvantages of the proposed methodology and point out challenges and possible future research directions to deploy the connectivity-based analysis framework to a range of different hydrologic applications.

5.1 Considerations about the low performance of HCIU(CN)

As noted in Sect. 4.1, the CN-based formulation of the normalized connectivity index and the resulting HCIU(CN) (obtained as the weighted average of the former) are more complex than the n-based formulation as they attempt to account for the combined effects of different coexisting land-use/land-cover and soil types. Because of this, the choice of the CN-based weighting coefficients and the interpretation of the resulting raster maps of normalized connectivity require particular attention.

For our analysis, we adapted CN values associated with different land-use/land-cover types and hydrologic soil groups from the work by Wu et al. (2024) based on hydrologic modeling. According to their classification, urban cells with varying degrees of land development may be assigned CN values over a wide range depending on the associated soil characteristics. In some cases, moderately urbanized cells may have a weight W that is not too dissimilar from those of some cells with natural land-cover types but low soil infiltrability, meaning that both types of cells are expected to generate similar amounts of runoff. However, it is clear that undeveloped patches would still mitigate hydrologic impacts better than urban areas, mostly through enhanced detention due to their higher roughness, delaying runoff contributions.

Interpretations for the overall low performance of HCIU(CN) may hence include the variability in stormwater hydrology dynamics due to storm intensity and soil moisture conditions. Depending on these aspects, dominant control mechanisms for hydrologic connectivity and disconnection, such as soil infiltrability, rainfall interception by vegetation and litter, and surface roughness, as well as their interplay, may significantly change (Saffarpour et al., 2016; Zölch et al., 2017). For instance, the potential for generating runoff in moderately urbanized cells may be mainly governed by the type of soil when regular events are considered – like in the case of the validation of the CN classification proposed by Wu et al. (2024). On the other hand, when the scope of the analysis focuses on more extreme events (as in our case, with peak-flow equations), the effects of land development may become preponderant (e.g., by decreasing response lag-times) since soils are more likely to reach saturated conditions. All of this in turn suggests the need to suitably adjust the CN values of urbanized cells to better differentiate them from pervious pixels even when the latter display low infiltrability.

Another likely explanation is that CN also depends on antecedent soil moisture conditions (ASMCs), with different sets of values associated with distinct ASMC categories (i.e., dry, average, and wet; Wu et al., 2024). In our work, we considered the average scenario. However, because of the interactions between event intensity and soil moisture conditions, different sets of CN values, associated with different ASMCs, should probably be considered instead, depending on the predicted event's return period. In the case of small, more frequent events (where we found the worst performances for EPAE and MO; Fig. 7a and b), HCIU(CN) might benefit from using lower CN values to reflect drier antecedent conditions. On the other hand, for greater return-period events, soils are saturated faster, and higher CN values, related to wet ASMCs, might improve the predictive power of HCIU(CN).

In summary, based on the overall low performance observed with HCIU(CN), the adopted weighting approach based on average ASMC curve numbers may not be the best for our proposed methodological framework, at least when the analysis aims at predicting extreme peak flows. A poor tuning of the CN-based weighting coefficients, which do not clearly distinguish undeveloped patches with low infiltrability from developed areas with highly permeable soils, may also explain why the HCIU(CN) seems less robust compared to the simpler n-based version. Further research is needed to study the sensitivity of HCIU to different CN-based weighting approaches and fine-tune the coefficients to maximize its predictive power.

5.2 Advantages and limitations of the proposed approach and future research directions

The proposed conceptual framework builds on and takes advantage of our qualitative understanding of some of the varied and complex dynamics affecting a basin's hydrologic response. The main strength of the proposed HCIU is that, in contrast with the traditional percentage of total impervious area or other metrics that index the impervious fraction of a basin, it explicitly takes into account the location and spatial arrangement of all types of LULC patches in relation with basin topography, the stream network, and the outlet in a continuum from the highly developed (which tend to generate more runoff faster) to the undisturbed forested patches (with the strongest mitigating effects on hydrologic response). This represents an advantage over other advanced urbanization metrics in the literature. For instance, while the methodologies adopted by Yang et al. (2011) and Beck et al. (2016) both frame the hydrologic interconnections across some main LULC categories simply based on geometric aspects such as the density, spatial adjacency, and granularity of distinct LULC patches, HCIU captures the connectivity across a spectrum of surface patches as driven by basin topography. Or, in contrast with the two metrics proposed by Zhang and Shuster (2014), which also consider the effects of basin relief but adopt a binary, pervious-or-impervious LULC representation, HCIU considers the continuous heterogeneity of hydrologic characteristics in the mosaic of LULC patches over a watershed.

This is achieved by alternatively using Manning's surface roughness coefficients or curve numbers, which are both well-established conceptualizations that quantify different LULC types' abilities to either facilitate the generation and transmission of runoff or retain, detain, or dampen water volumes, dependent on their surface roughness, extent (and type) of vegetated areas, and infiltrability. Highly impervious patches with smooth surfaces and limited losses represent one extreme in the continuous spectrum of LULC potentials for generating runoff, resulting in large peak flows; in the other extreme are forested areas and other LULC types with natural soils, presence of litter, and high vegetation densities because of their ability to diminish runoff volumes and travel speeds through interception (in the canopy and litter), temporary canopy storage, high surface roughness, and variable infiltration capacity.

Besides its comprehensive formulation that conceptually considers the synergistic effects of both topography and LULC, another advantage of HCIU is that it can be virtually derived for any basin in the US or the world from widely available data, i.e., DEM and LULC (or CN) maps. As a result, the proposed urbanization metric may easily find systematic application in peak-flow prediction models and, depending upon preliminary testing, also in a range of other water resources management fields that require quantifying the level of urbanization in a basin, such as water-quality assessments and stormwater infrastructure design. We considered two alternative approaches to weigh the effects of different LULC types (and HSGs in the case of CN-based weights) on hydrologic connectivity. However, modelers may also try other weighting criteria depending on their specific needs. To foster the dissemination of our new metric, a link to a Python program to calculate the HCIU of any basin is provided in the “Code availability” section (https://doi.org/10.5281/zenodo.14457110, Dell'Aira, 2024).

A further advantage of the proposed approach is that it yields simple and easy-to-interpret yet conceptually comprehensive assessments of the hydrologic impacts of urbanization across different basins. It can also be used to compare different land-development scenarios for the same watershed, aiding stakeholders in making urban planning decisions, or evaluating possible future LULC changes due to, e.g., the distributed implementation of candidate stormwater control measures.

A current limitation of HCIU is that it only frames connectivity patterns driven by basin relief (topography) even though the potential presence of stormwater drainage infrastructure would be another important source of connectivity in any urbanizing basin. Specifically, in its current version, the methodology does not capture the effects of underground stormwater sewer networks (also referred to as the minor system of stormwater infrastructure; Martins et al., 2017), although these are typically present in urban environments, especially in highly developed areas. However, stormwater drainage infrastructure usually includes not only underground pipe networks, but also surface flow pathways and canals, which make up the so-called major system. The major system is critical for handling larger, less frequent storm events. When calculating HCIU, major drainage system sections connected to natural channels are treated as part of the stream network (assuming that excess flow from the major system is poured directly into the stream network). This means that the connectivity of hillslope cells draining to the major system is calculated referring to the pour points along the major system. The contributions of these hillslope cells are then weighted based on the along-the-stream-network distance to the outlet, measured starting from the major-system pour point and following both the major-system and any subsequent natural-stream-network links downstream, when averaging the (normalized) connectivities to compute HCIU. This approach captures the effects of the stream network and major drainage systems, which have a stronger influence on the hydrologic response to extreme rainfall events compared to the minor system, whose capacity is typically overwhelmed by large runoff volumes. Consequently, the proposed HCIU should be a more reliable predictor of hydrologic-response variables under severe flooding conditions, as also suggested by the increase in model performance moving from small to intermediate peak quantiles, observed in Fig. 7.

On the other hand, when the analysis focuses on basin response to regular storms (e.g., in water-quality studies), the effects of the minor system should not be neglected as the underground network may be able to handle most of the (smaller) runoff volumes. Another scenario where it is highly recommended to explicitly consider underground connectivity is when dealing with heavily urbanized watersheds, typically characterized by the presence of extensive drainage infrastructure. In these basins, detention tanks and sections of the minor system pumping stormwater against topographic gradients may completely change the connectivity determined by topographically driven surface runoff pathways. Our results are for basins with heterogeneous LULC characteristics, where urbanized sectors with varying development rates are mixed with natural LULC patches, typically displaying a distribution of land-development intensities more skewed towards lower values (Fig. 4j, k, and l), as is common for residential areas. Among the three studied regions, only VA included watersheds with TIA above 50 %, but all of those were small in size. Because our dataset may not be representative of large, highly urbanized basins, for these cases (e.g., in countries where cities present generally higher land-development intensities compared to the US) we recommend considering the effects of the minor underground stormwater drainage infrastructure as well when deriving HCIU. If stormwater sewer data are not available for the study region, and HCIU is estimated only considering topographically induced connectivity, some preliminary testing of its predictive power on gauged basins should be required (e.g., using the validation approach depicted in Fig. 8) before using the index for systematically generating peak flows in ungauged, highly urbanized watersheds.

In principle, the proposed methodology allows for considering the effects of underground stormwater infrastructure as well; this would be achieved by first identifying pixels that are connected to each other and the stream network through stormwater sewer links and then computing the connectivity index in a way such that connectivity patterns due to the stormwater sewer infrastructure override those from topography, whenever applicable. As we have found, however, information about underground stormwater drainage networks is not easily available, which may impede considering this additional source of connectivity, particularly in the case of multi-basin studies. Depending on data availability, an expanded version of HCIU that also accounts for the effects of stormwater sewer infrastructure as an additional source of connectivity should be another topic for future research.

6 Conclusion

We proposed and tested a new, hydrologic-connectivity-based index of urbanization (HCIU) that can be obtained in a GIS framework from the digital elevation model and a land-use/land-cover (or curve number) map of a basin. We showed that, compared to the traditional fraction of total impervious area (TIA), HCIU helps capture more information about the impacts of land development on hydrologic response. TIA only indexes the proportion of impervious patches in a basin, while our new metric explicitly accounts for the spatial arrangement of the different land patches found in a watershed, both natural and developed. This is obtained by considering the spectrum of localized effects of distinct land-use/land-cover types on the hydrologic connectivity of surface runoff pathways.

The methodology builds on the well-established connectivity index, which has already found wide application in several hydrologic and geomorphic problems. Our specific interpretation of the connectivity index, in the framework of our approach, considers it as that distributed property explaining the ability of any hillslope location to quickly receive and transfer runoff to the stream network, depending on how topographically induced runoff pathways interact with urbanized sectors (with possible local differences in land-development intensity) and the undeveloped, typically more pervious patches in a watershed. We considered two alternative widely used conceptual descriptors for quantifying the potential for runoff generation of different land-use/land-cover types, i.e., Manning's surface roughness coefficients and curve numbers. Depending on these metrics, the contributions to hydrologic connectivity of basin pixels with distinct surface characteristics are weighted differently. We found that weighting factors specified as a function of Manning's surface roughness coefficients result in more robust HCIU metrics, compared to curve-number-based weighting methods, when explaining urban peak flows. However, we do not exclude that fine-tuning the latter might improve model performance.

Irrespective of the weighting criterion, the procedure for obtaining a lumped metric for the effects of urbanization on hydrologic response, starting from a connectivity map of the basin of interest, involves the following steps: first, we define a normalized connectivity map, with respect to the maximum connectivity scenario associated with a virtual, fully developed copy of the original watershed; then, we calculate HCIU as a weighted average of the normalized connectivities of all basin pixels that are not part of the stream network, depending on the distance from their pour point to the watershed outlet, following the stream network.

We have shown that our new urbanization metric improves the predictive power of existing peak-flow regional equations for three comprehensive case studies. Further research is required to test HCIU's explanatory power for other hydrologic-response (i.e., flood-related) variables, such as lag time and time of concentration, as well as other hydrologic variables of interest that have traditionally displayed correlation with TIA, such as water-quality indicators. Depending on the scope of the analysis, an expanded version of the current formulation to account for the additional source of connectivity introduced by underground storm sewer infrastructure may be necessary; in highly urbanized watersheds, the latter may be a stronger control over basin response than topographically induced connectivity, especially in the case of less intense, more frequent events.

Besides its direct application as a metric of urbanization effects on basin response, HCIU's sensitivity to the spatial arrangement of more developed and less developed (or undeveloped) sectors may provide a novel framework to facilitate comparisons of the hydrologic impacts caused by basin changes (e.g., due to urbanization or the introduction of stormwater runoff control measures), offering a valuable tool to stakeholders for informed urban planning decisions. More research is needed to study the benefits and the range of applicability of the hydrologic-connectivity-based index of urbanization.

Appendix A

Table A1 Case-study basins. Area, A; fraction of total impervious area, TIA; and flood quantiles Q₂ to Q₅₀₀ were retrieved from Southard (2010), Austin (2014), and Feaster et al. (2014) for the MO, VA, and EPAE case studies, respectively. We calculated HCIU(n) and HCIU(CN) following the methodology proposed in this work.

Appendix B

Table B1 Manning's roughness coefficients associated with the different NLCD land-use/land-cover types (adapted from Liu and De Smedt, 2004; Hooke et al., 2021).

Appendix C

Table C1 Curve numbers associated with NLCD land-use/land-cover types and hydrologic soil groups (adapted from Wu et al., 2024).

Appendix D

Figure D1 Box plots of test errors for (1) benchmark $Q_T\sim A+\mathrm{TIA}$, (2) $Q_T\sim A+\mathrm{HCIU}\left(n\right)$, and (3) $Q_T\sim A+\mathrm{HCIU}\left(\mathrm{CN}\right)$ models, when they are fitted and blind-tested on distinct basin subsets, for the VA homogenous region. For each box plot, filled bars represent error values between the first and the third quartiles, the upper (lower) whisker extends from the third (first) quartile by adding (subtracting) 1.5 times the interquartile range, and any outliers beyond the whiskers are marked by circles.

Appendix E

Table E1 Outlet coordinates, referring to the World Geodetic System (WGS84).