Predicting runoff in data-sparse drylands: a case study evaluating regionalization success for two parsimonious loss models

Reliable flash flood forecasting is essential for sustainable storm infrastructure planning and hazard mitigation in sensitive dryland landscapes. Loss or infiltration models are central to flash flood prediction in arid watersheds, where runoff is dominated by infiltration excess. Parsimonious models are critical for use in practice. Accurate prediction remains a challenge because calibration data are typically scarce or unavailable. Regionalization can help fill this data gap: model parameters from a physically similar donor basin with available calibration data are transferred to a target watershed with no or limited data. In this study, we tested whether regionalization success depends on model selection. Two simple loss models available in the widely-used HEC-HMS software – the SCS curve number (CN) and the linear-constant (LC) methods – were calibrated and validated based on long-term observations from the Walnut Gulch Experimental Watershed in Arizona. Optimized parameters derived from this donor basin were then transferred to the target basin, a partially urbanized watershed in New Mexico with a dense network of rain and flow gaging stations. The target basin has a limited observational record and therefore serves as a suitable test case for regionalization under data-limited conditions. Validation results from the donor basin demonstrate that the LC model outperformed the CN model for large storm events across four performance metrics. Performance for small events was not significantly different. When parameters were transferred to the target basin and tested for a rare, high-intensity storm event recorded in 2025, the LC model again performed better for most metrics. While regionalization can be a viable strategy to support stormwater infrastructure design in data-scarce, semi-arid urban watersheds, this study indicates that its success depends on model selection.

Data from the World Meteorological Organization (2021) reveals a rise in hydrologic disasters over the last two decades. Changes in streamflow, timing, and associated extremes are predicted to have a negative impact on many watersheds worldwide (Intergovernmental Panel on Climate Change, 2023; Soltani et al., 2023). From 2000 to 2019, flooding accounted for nearly half of all natural disasters, affecting 1.6 billion people (Soares et al., 2025). Flash floods represent the most frequent type of natural disaster and are associated with the highest global fatality rates among hydrologic hazards (Kumar, 2023). Arid and semi-arid regions are susceptible to flash flooding, leading to the destruction of infrastructure, harm to people, and severe economic damage (Alkaabi et al., 2025). This natural vulnerability is increasingly compounded by rapid urban development; approximately 50% of the nation’s most rapidly expanding metropolitan centers are situated within these arid and semi-arid regions (Schoener, 2022).

Managing the risk of flash flooding in urbanizing areas requires effective forecasting (Wang et al., 2023). Hydrologic models are crucial tools for flash flood prediction. As part of any hydrologic model, the loss or infiltration module determines partitioning of precipitation into infiltration (or loss) and runoff. Over the past century, various infiltration methods of varying complexity have been developed (Schoener et al., 2021). Simple or parsimonious models are widely used in practice. They require fewer parameters and are therefore easier to implement than more complex options. One of the most popular models for estimating rainfall loss and excess among practicing engineers and hydrologists is the one-parameter SCS curve number (CN) method (Chin, 2023). Other simple loss models, such as the linear-constant method, have been proposed and incorporated into widely used hydrologic modeling software (Rassa et al., 2024).

While parsimonious models require fewer input parameters and are thus easier to apply, model calibration and validation in data-sparse arid and semiarid regions still pose a major challenge (Chen et al., 2024). Here, regionalization approaches can offer a solution (Guo et al., 2021). The concept behind regionalization is to calibrate model parameters for a catchment with available rainfall and discharge measurements and then apply those parameters to an ungauged catchment in close spatial proximity and/or with similar physical characteristics (Mihret et al., 2025). Regionalization techniques have been successfully applied across diverse data-scarce environments worldwide (Kim and Kaluarachchi, 2008; Mihret et al., 2025; Oudin et al., 2008; Samuel et al., 2011; Goodarzi et al., 2021; Wufu et al., 2026). Despite advances in regionalization practice over the past two decades, limited comparisons are available for dryland regions. Moreover, the success of parameter regionalization depends on the model structure (Guo et al., 2021). In choosing a particular infiltration model, the modeler makes an important structural decision. To our knowledge, a comparison of regionalization performance across different parsimonious loss models has not been completed to date.

To address this important knowledge gap, our study poses the following research questions: (1) Can regionalization be successfully applied for loss model parameterization in an arid catchment? And (2), are there discernible differences in regionalization success between two parsimonious loss models? To answer these questions, we first calibrate and validate hydrologic models for three basins in the Walnut Gulch Experimental Watershed (WGEW) in Arizona, USA, where high-resolution rainfall-runoff data have been collected for decades. We then apply parameters to the Star Heights test watershed in New Mexico, USA, where the authors maintain a dense network of rain and flow gaging stations. We compare model performance for a rare, large-magnitude storm event recorded in 2025, along with several smaller events.

Results provide actionable guidance for practicing engineers and hydrologists facing model selection and parameterization decisions in data-sparse dryland watersheds.

Study areas

The Walnut Gulch Experimental Watershed (WGEW) in southeastern Arizona served as the donor basin for regionalization in this study. The target basin was the Star Heights test watershed in central New Mexico.

The WGEW (Figure 1) was established by the United States Department of Agriculture in 1953 and equipped with rain gauges and runoff-measuring flumes organized across a network of nested subwatersheds (De Graaff et al., 2019). The watershed is a tributary of the San Pedro River. Three non-nested subbasins within WGEW – Basins 4, 103, and 102 with drainage areas of approximately 2.3 km 2, 0.036 km 2, and 0.017 km 2, respectively, were selected for this study.

The Star Heights test watershed, located in Rio Rancho, New Mexico (Figure 2), is part of the Black Arroyo watershed and drains to the Rio Grande. Six non-nested basins ranging in size from 0.56 to 1.57 km 2 were analyzed as part of this study. The Star Heights area comprises a mix of open space and urban development, primarily single-family residential areas, with a combination of impervious surfaces and pervious yards featuring natural desert landscaping and open channels designed to convey stormwater runoff. A dense network of rain and flow gaging stations (Figure 2) was installed in spring 2025 by the Southern Sandoval County Arroyo Flood Control Authority in collaboration with the University of New Mexico to monitor precipitation and runoff. The gauge network is operated and maintained by the authors.

Table 1 compares donor and target basin characteristics. WGEW basins were selected strategically to provide the closest possible match to the Star Heights test watershed. The range of drainage areas of Star Heights basins (0.56–1.57 km 2) falls within the range of WGEW basins (0.02–2.30 km 2). Average slopes and longest flow paths calculated from digital elevation models are also comparable. Soils at WGEW consist primarily of sandy loams of hydrologic soil group B (Osterkamp, 2008). Field soil surveys at Star Heights, conducted by the local flood control authority (personal communication), confirm that soils are predominantly sandy loams of the same hydrologic soil group. One notable difference between donor and target basins is urbanization: while WGEW basins remain largely in their natural state, Star Heights basins feature varying degrees of imperviousness due to urbanization. We address this difference by using regionalization only to parameterize loss models for pervious areas and treating impervious surfaces independently (see section 2.5).

Most runoff events at both the New Mexico and Arizona sites are caused by the North American Monsoon, which delivers roughly two-thirds of the annual precipitation during the summer months (Goodrich et al., 2008). Average annual precipitation at WGEW is 290 mm; this value reflects the 20-year average (2006–2025) based on rain gauges used in this study. The Star Heights test watershed receives an average of 273 mm of rainfall per year (NOAA, 2026).

Precipitation and runoff data

High-resolution precipitation and runoff records for the WGEW were retrieved from the Southwest Watershed Research Center Data Access Project (US Department of Agriculture, 2026). The dataset spans from 2000 to 2025 and consists of discharge and precipitation depth measurements at a one-minute temporal resolution. The original records were provided in a breakpoint format, representing discrete changes in rainfall intensity and discharge over time; these were subsequently processed in Python to a consistent one-minute interval by filling temporal gaps and interpolating values to ensure uniform synchronization for model input. Discharge measurements were collected at three flumes located at the outlets of basins 4, 102, and 103. Precipitation measurements were collected by seven rain gauges (Figure 1).

Since the concern motivating regionalization in ungauged catchments is often flood prediction or infrastructure design, the top 100 runoff events (based on cumulative runoff depth) for all WGEW basins were selected for this study. This is equivalent to 53% of runoff-producing storms in basins 102 and 103 during the 26-year study period and 83% of runoff-producing storms in basin 4. Cumulative runoff depths ranged from 0.5 to 32.5 mm in basins 102 and 103, and from 0.01 to 11.4 mm in basin 4. Box plots in Figure 3 display distributions of cumulative rainfall for seven rain gauges located in or adjacent to WGEW basins 4, 102, and 103. Precipitation depths for basins 102 and 103 (RG 038 and 384, Figure 3) ranged from 6.9 mm to 69.2 mm. For basin 4, precipitation depths ranged from 2.0 mm to 79.4 mm.

Rain and flow data from the Star Heights study area were collected as part of this study for the period April 2025–August 2026. Eight tipping bucket rain gauges (Figure 2) (HOBO Data Loggers, Bourne, MA) measured precipitation in five-minute increments at a resolution of 0.25 mm. Flow depth at the outlet of six subbasins was monitored using non-vented pressure transducers (HOBO Water Level Logger, Bourne, MA). Raw pressure time series were corrected for atmospheric pressure fluctuations based on co-located atmospheric pressure sensors and converted to depth measurements. Discharge time series were then derived from depth measurements using rating curves developed for each stream gaging station. All Star Heights gaging stations are located in concrete-lined open channels with uniform cross-section geometry. For each station, a hydraulic model was developed using the US Army Corps of Engineers River Analysis System software (HEC-RAS). Each model was run with 100 flow rates spanning the expected flow range; paired depth-discharge data from each monitoring location were then used to develop rating curves.

On July 3, 2025, a rare storm event impacted the Star Heights test watershed. As much as 77 mm of rain fell in the basin over a period of approximately 2 hours, with maximum intensities exceeding 100 mm h −1. Figure 2 shows cumulative precipitation measured by rain gauges during the July 3, 2025, event. Shaded background colors represent cumulative precipitation estimates derived from Next Generation Weather Radar (NEXRAD Level3/DTA). Figure 4 compares intensity-duration frequency curves (solid lines) for the Star Heights study area with observed precipitation for the July 3 event for seven gaging stations (dashed lines). The comparison shows that, for all but two gaging stations, the storm data approached or exceeded the 100-year recurrence interval (or 1% annual exceedance probability) for a 2-h duration.

In addition to the large July 3 event, three smaller storms that occurred during 2025 (May 7, July 19, and August 25) as well as one event in 2026 (July 23), generated runoff in the Star Heights test watershed and were included in the analysis. Average total precipitation for each of the smaller storms was 5 mm, 22 mm, 17 mm, and 8 mm, respectively. Maximum precipitation intensities were 30 mm h −1, 64 mm h −1, 21 mm h −1, and 76 mm h −1, respectively.

Hydrologic model setup and forcing

The US Army Corps of Engineers HEC-HMS model is a widely adopted, open-source rainfall-runoff simulation tool (Rassa et al., 2024) capable of estimating surface runoff responses to precipitation events. It was selected for this study due to its status as a standard industry framework, which provides a robust and computationally efficient environment for simulating complex rainfall-runoff processes in both natural and urbanized watersheds (US Army Corps of Engineers, 2016). In this study, version 4.12 and 4.13 Beta 2 were used for modeling study basins. Hydrologic models were developed using a semi-distributed approach. Watershed boundaries and stream networks were delineated using a 1-m resolution digital elevation model (DEM). Elevation data for both study areas were obtained from the US Geological Survey1.

While WGEW basins 102 and 103 were modeled as single catchments, basin 4 and the Star Heights catchment were discretized into various subbasins (see Figures 1, 2). This higher resolution was necessary to account for spatial heterogeneity in precipitation. For the Star Heights catchment, subbasin delineation was also informed by locations of flow gaging stations and flood control ponds. Meteorological modeling was performed using inverse distance interpolation between gauge locations. Precipitation from two rain gauges was interpolated for WGEW basins 102 and 103 (RG 83, 384). Five rain gauges were interpolated for basin 4 (RG 27, 31, 71, 74, and 80), and eight gauges for the Star Heights catchment (RG 10–17).

Infiltration and excess

Two infiltration methods were evaluated and compared in this study: the SCS curve number (CN) method (Equation 1) and the linear-constant (LC) model. The CN method was developed by the Soil Conservation Service (US Department of Agriculture, 1954) and remains one of the most widely used models in practice for partitioning rainfall into infiltration and excess. Originally derived from data collected in agricultural watersheds in the eastern United States, the CN method has since been extended for application in arid and urban areas (Fennessey et al., 2001; Schoener, 2018). Cumulative runoff Q (mm) is estimated as:

Where P is the cumulative rainfall (mm), S is the potential maximum retention after runoff begins (mm), and I a is the initial abstraction (mm). The potential maximum retention, S, is derived from the dimensionless CN parameter:

The CN model assumes that no runoff occurs until the initial abstraction is fully satisfied. Initial abstraction is often estimated at 20 percent of S. This assumption may fail to capture the rapid infiltration-excess runoff response characteristic of high-intensity convective storms in semi-arid regions.

The linear-constant (LC) infiltration model (Equation 3) (Schoener et al., 2021) does not use an initial abstraction term and expresses infiltration potential as a function of cumulative infiltration:

Where f t (mm h −1) is the potential infiltration at time t; m is the infiltration decay coefficient (h −1); F t (mm) is the cumulative infiltration at time t; f 0 (mm h −1) is the initial infiltration capacity; F c (mm) is the initial deficit; and K eff (mm h −1) is the effective hydraulic conductivity.

Transform and routing

Rainfall excess was transformed into runoff hydrographs using the SCS synthetic unit hydrograph method. This method is widely used in small watershed hydrology to predict basin hydrologic response (Shatnawi and Ibrahim, 2022). Basin lag time, defined as the duration from the center of mass of effective rainfall to the hydrograph peak (Dandupati et al., 2024), was estimated as 60% of the time of concentration (T c). T c for each subbasin was estimated using the velocity method (US Department of Agriculture, 1986), which is based on the physical principle that travel time depends on the runoff flow length and the flow velocity. The velocity method accounts for the sum of travel times across three distinct flow regimes: sheet flow, shallow concentrated flow, and open channel flow (Fang et al., 2007; McCuen et al., 1984).

The Muskingum–Cunge method was applied for routing runoff through open channels. This routing model uses measurable channel characteristics (Salvati et al., 2024): channel length, slope, and cross-section geometry were derived from 1-m digital elevation models. Manning’s roughness coefficients were assigned based on surface material, with values ranging from 0.028 for sandy, alluvial channels to 0.013 for concrete-lined conveyances. In addition to open channels, the Star Heights test basin also includes four flood control ponds. Pond routing was accomplished using elevation-storage-discharge curves obtained from record drawings and pond design reports (personal communication).

In the Star Heights watershed, impervious surfaces, including paved roads and parking lots, driveways, and rooftops, were mapped based on orthoimagery and land use data obtained from the City of Rio Rancho (personal communication). Zero infiltration was assumed for impervious surfaces; all rain falling on the impervious portion of each Star Heights subbasin was converted to runoff using the “percent impervious” option in HEC-HMS.

Walnut Gulch model calibration

To minimize the discrepancy between observed and simulated hydrographs, automated optimization trials were conducted within HEC-HMS. Table 2 shows parameters adjusted during model calibration at WGEW. To ensure physical realism and prevent the models from converging on hydrologically implausible solutions, all calibration parameters were constrained within defined physical ranges (Odey and Cho, 2025).

Prior studies have shown that curve numbers in the WGEW determined using various approaches rarely drop below 75 (Duan et al., 2025; Simanton et al., 1996); the CN range was therefore set to 75–99. The initial deficit (F c) range was set equal to 1–70 mm, whereby the upper bound represents the maximum measured precipitation in the WGEW. The range for effective hydraulic conductivity (K eff) was informed by prior work (Coulter et al., 2026; Schoener et al., 2024) and reflects plausible values for conditions ranging from bare soil to 75% plant cover in desert landscapes comparable to the study area. The differential evolution global optimization algorithm was utilized during model calibration to find optimal values of CN, F c, K eff, and Lag for each storm event in WGEW basins 4, 103, and 102. Since its introduction in 1995, differential evolution has been widely adopted in optimization research, owing to its proven efficacy as a global search algorithm (Eltaeib and Mahmood, 2018). The objective function used in this study minimized the root mean square error (RMSE) between observed and simulated discharge. At present, the HEC-HMS model only allows single-objective optimization. RMSE was chosen in this study because it is sensitive to high flows, which are crucial for flood forecasting and infrastructure design.

In addition to curve number optimization described above, CN values were also calculated from measured cumulative rainfall-runoff pairs for two initial abstraction scenarios: using the standard assumption that initial abstraction is equal to 20% of the potential maximum retention (I a = 0.20 S), and applying the lower 5% initial abstraction ratio (I a = 0.05 S) recommended by various studies (Shi et al., 2009; Brandao et al., 2025; Garcia et al., 2026). For both scenarios, S can be calculated directly from rainfall-runoff data using Equations 4, 5, respectively:

Curve numbers can then be determined from S by applying Equation 2. Four loss model optimization scenarios were thus assessed in this study: optimized curve number and linear-constant scenarios (CN opt, LC opt) based on RMSE minimization, and two calculated curve number scenarios for initial abstraction ratios of 0.20 and 0.05, respectively (CN 0.20 and CN 0.05) using Equations 2, 4, 5.

Walnut Gulch model validation

To ensure the models were not over-fitted, a leave-one-out cross-validation protocol was employed for model validation. Leave-one-out cross-validation is particularly advantageous for studies with limited datasets, as it remains robust even when conventional partitioning methods are unreliable or insufficient (Geroldinger et al., 2023). The leave-one-out protocol is executed as follows: (1) a single event is excluded from the cohort; (2) model parameters are calibrated using the remaining events; (3) mean calibrated parameter values are applied to simulate the runoff for the excluded target event (Zhang et al., 2018). This cycle is repeated until every event within the cohort has served as an independent validation case. Since this method avoids the randomness inherent in traditional training/testing splits, the resulting point estimates of accuracy are constant for the dataset, providing a deterministic and reproducible measure of model performance (Wong, 2015).

Model performance metrics

Using the concept of similarity-based regionalization, mean optimized parameters from the WGEW basins were applied to parameterize loss models for pervious areas of the Star Heights study area. The performance of HEC-HMS simulations, both for WGEW validation trials and for the Star Heights watersheds, was quantified using a suite of performance metrics. These indicators were selected to evaluate the models’ ability to replicate the overall hydrograph shape and to assess errors in peak discharge and runoff depth.

The Nash-Sutcliffe Efficiency (NSE) (Equation 6) is commonly used as a goodness-of-fit statistic in hydrology. NSE values range from ∞ to 1, where 1 indicates a perfect match between observed and simulated data (Nash and Sutcliffe, 1970):

Where O i = observed value, S i = simulated value, Ō = mean of observed values, and n = number of observations. The root mean square error (RMSE) (Equation 7) was employed to measure the deviation between the simulated and observed values in the same units as runoff (Moriasi et al., 2007):

In addition to NSE and RMSE, we calculated absolute errors in runoff peak and depth. The accuracy of peak discharge prediction is essential for infrastructure sizing and flood hazard mapping, while runoff depth (or volume) is important for the design of flood control dams.

Statistical significance of differences in performance metrics between model scenarios at the Walnut Gulch study site was assessed using the non-parametric Wilcoxon signed-rank test (α = 0.05).

Model optimization results, Walnut Gulch basins

Optimized model parameters for Walnut Gulch basin 102 are plotted against cumulative precipitation for 100 observed storm events in Figure 5. Optimized linear-constant model parameters K eff and F c (Figures 5a,b) show no strong relationship with cumulative rainfall as indicated by low coefficients of determination (R 2 = 6 and 2%, respectively). Calculated and optimized curve numbers, on the other hand, display a stronger dependence on precipitation depth, whereby curve numbers decrease with increasing rainfall (Figures 5c–e). This decreasing trend is strongest for the CN 0.20 scenario, where the linear regression explains 50% of the variability in the data (R 2 = 50%).

Figure 6 shows box plots comparing model results from leave-one-out cross-validation in Walnut Gulch basin 102. Four performance metrics are included for the linear-constant method (LC opt) and three curve number scenarios (CN 0.05, CN 0.20, and CN opt). Mean values and standard deviations are included under Supplementary Tables S1, S2. Given the demonstrated dependence of curve number parameters on cumulative rainfall depth, results are presented separately for large rainfall events with cumulative precipitation >24 mm (top) and small events with a total rainfall depth between 0 and 24 mm (bottom). The 24 mm threshold corresponds to the mean precipitation value across all Walnut Gulch rain gauges for the study period. Significant differences between performance metrics are indicated with square brackets.

The comparison indicates that for large rainfall events, the LC models show significantly improved performance over all CN scenarios, as indicated by NSE, RMSE, absolute peak error, and absolute volume error. Among curve number scenarios, CN 0.05 generally performed best, and CN opt performed worst. The only exception was the absolute peak error metric, where CN 0.05 resulted in the largest error. For small rainfall events, the CN 0.05 scenario performed better than most other scenarios with respect to most performance metrics. No other systematic differences are evident from the comparison. Results for basin 103 (not included) are very similar to basin 102.

Results for basin 4 are shown in Figure 7 (see also Supplementary Table S3 for mean values and standard deviations). For large events, the LC and CN 0.05 scenarios showed similar performance. Both generally performed better than the other curve number scenarios. For small events, differences in model performance were generally not statistically significant.

Regionalization – knowledge transfer from WGEW to the Star Heights basins

Knowledge gained from Walnut Gulch model optimization was applied to parameterize infiltration models for six gauged subbasins in the Star Heights test watershed in central New Mexico. Table 3 presents the mean values of optimized parameters for the linear constant model (K eff, F c) and three curve number scenarios across all three WGEW test basins. Mean parameter values are calculated for all events (Table 3, top), and only for large precipitation events with cumulative rainfall exceeding 24 mm (Table 3, bottom). Both parameterization scenarios were applied and tested in Star Heights test basins for the large storm of July 3, 2025. Only loss parameters were transferred from the donor to the target basins. Transform and routing parameters for Star Heights basins were estimated based on watershed characteristics.

Figure 8 contains a map of the Star Heights watershed. Subbasins draining to each flow gaging station (circle) are shaded in different colors. Charts show precipitation intensity (inverted axis) and observed runoff (shaded area), along with flow simulated by four infiltration models for the parameterization scenario using the mean of all events (Table 3). Comparison of observed and simulated hydrographs reveals that all three curve number scenarios overestimated runoff at four of six gauges.

Box plots in Figure 9 (see also Supplementary Table S4) compare four performance metrics for two parameterization scenarios: (1) loss model parameters in target basins were set based on the average of all events from the donor watershed (Figure 9, top), and (2) target basins were parameterized based on the average of the largest events (>24 mm cumulative precipitation) only (Figure 9, bottom). The second scenario was included because optimized parameters at Walnut Gulch showed a clear dependence on storm magnitude for the curve number method. Figure 9 shows that the LC model consistently performed better than all three CN scenarios. Curve number model performance improved slightly when parameterization was based on large events. The CN 0.05 scenario yielded better results than the other tested CN options with respect to NSE and RMSE, but led to higher peak and volume errors. Figure 10 (see also Supplementary Table S5) shows box plots of performance metrics for the Star Heights basins across four smaller storm events observed during the 2025–2026 study period. Similar to findings from Walnut Gulch, model performance for small events was generally poor, as indicated by negative NSE values. Across all metrics, the LC model showed slightly better performance than the curve number scenarios.

Parameterizing infiltration models in data-sparse areas remains a major challenge for practitioners in engineering and hydrology (Xu et al., 2024). Application of handbook values or upscaling of laboratory or point-scale field measurements can introduce substantial uncertainty (Tügel et al., 2022). Regionalization is a commonly used method in ungauged watersheds, whereby optimized parameters from a gauged catchment are transferred to a comparable, data-sparse basin (Guo et al., 2021; Tügel et al., 2022). This paper evaluated the regionalization performance for two parsimonious loss models: the curve number (CN) and linear-constant (LC) method. Model parameters were first optimized using long-term precipitation and runoff measurements from the Walnut Gulch Experimental Watershed in Arizona. Optimized parameters were then applied to the Star Heights watershed in New Mexico using the concept of physical similarity. Model performance was assessed based on observations from a rare, large magnitude storm event that occurred in 2025, as well as four small magnitude storms.

Influence of loss model selection on regionalization success

Results from this case study indicate differences in regionalization success between the two tested loss models, with the LC model outperforming CN model scenarios across most evaluated performance metrics. We attribute this difference in performance at least in part to the curve number dependence on cumulative precipitation. Decreasing curve number trends with increasing cumulative rainfall depths have been documented in previous studies (Bertotto et al., 2021; Chin, 2021a; Hawkins et al., 1985; Panigrahi and Ramadas, 2025; Schoener et al., 2021). Since calculated or optimized curve numbers increase with decreasing rainfall depth (Figure 5), parameters for regionalization depend on the characteristics and distribution of storm events used during optimization. For the WGEW example, curve numbers are thus much higher when they are calculated as the mean of all events (see Table 3, CN = 82–91), compared to considering large rainfall events only (CN = 75–86). Comparison of performance metrics (Figures 9, 10) reveals that curve number scenarios parameterized based on large events performed better than those parameterized based on all events. This suggests that for successful regionalization of curve numbers, storm magnitude between optimization and application should be similar, as suggested by several authors (Chin, 2021a, 2021b; Panigrahi and Ramadas, 2025). This limits the transferability of curve numbers between locations. It also limits the application of regionalization for studies focused on large magnitude, rare events, such as the design of flood control infrastructure. The flood control design storm in the Star Heights watershed is the 1% chance event with a 6-h duration (personal communication). According to the National Weather Service precipitation frequency data server2, cumulative precipitation for that recurrence interval and duration is 59 mm. Out of the 100 events analyzed for WGEW basins 102 and 103, only three exceeded that threshold in a 26-year period. Gaging stations with long periods of record in dryland regions are exceedingly rare.

Findings from this study suggest that the linear-constant method is less susceptible to the impact of precipitation depth. Figure 5 shows no strong dependence between optimized parameters and cumulative storm rainfall. Consistent with prior findings (Rassa et al., 2024), the LC model demonstrated clear and statistically significant performance benefits for small study basins (102, 103) and large storm events (>24 mm) in the WGEW (see Figure 6). Model performance was not significantly different between the LC and CN scenarios for most small events. Moreover, both loss models performed better for larger events. Decreased dependence on event characteristics and overall improved performance suggest that the LC model is better suited for regionalization and thus provides a benefit over the CN method in ungauged catchments.

Uncertainty, limitations, and future work

Model forcing data and discharge measurement uncertainty

Uncertainty in precipitation input has been found to contribute considerably to model error (Bárdossy et al., 2022), particularly in arid and semiarid areas with high spatial and temporal rainfall heterogeneity (Bao et al., 2025). In WGEW basin 4, models did not predict any runoff for 43 out of 100 events on the low end of the magnitude spectrum. Heterogeneity in precipitation coverage and corresponding uncertainty in model forcing data likely contribute to this discrepancy (Meles et al., 2024). To assess the potential impact of precipitation heterogeneity on model results and conclusions for the target basin, a sensitivity analysis was conducted for the extreme storm event of July 3, 2025. To that end, the LC model was run by alternately omitting one precipitation gauge, and performance metrics were compared with the baseline model that used all eight gauges. Findings from this gauge leave-one-out interpolation assessment (Otieno et al., 2014) are included under Supplementary Figure S7. Results show generally low variability in metrics when individual gauges were omitted, with mean changes in error metrics of less than 10%. The exception is rain gauge 10RG, which, when excluded, caused marked deterioration across all performance indicators (Supplementary Figure S7).

In addition to uncertainty related to model forcing data, errors in observed discharge can contribute to overall model uncertainty (Westerberg et al., 2022). While the rating curves developed for this study to transform stage measurements into discharge hydrographs could not be validated against field measurements due to safety concerns under flash-flood conditions and logistical constraints, the stage at all gaging stations remained within the hydraulically modeled range. Field surveys conducted immediately following the July 3, 2025 event confirmed that the channel geometry remained unchanged relative to modeled conditions. At two of the six gaging stations (52F and 56F), high water marks left by floodwater were surveyed and found to be within 2% of peak depths measured by pressure transducers.

Structural uncertainty

Structural uncertainty can arise from the lack of knowledge regarding the representation of a hydrologic system (Beven, 2016). Poor model performance for small events across all test basins may stem from the simple structure of both loss models and the lumped simulation approach. Treating each subbasin as a lumped unit for infiltration and runoff calculation neglects the complexity of arid landscapes. In reality, dryland watersheds feature intricate mosaic patterns of vegetated and bare surfaces with vastly different infiltration and runoff behavior (Schoener et al., 2026). Runoff from less permeable areas within a subbasin can re-infiltrate in run-on areas (Howes and Abrahams, 2003) or permeable ephemeral channels (Schoener, 2022) before reaching the watershed outlet. Improving the representation of hydrologic processes in mosaic dryland watersheds is an important area of future research.

Parameter and regionalization uncertainty

Equifinality or non-uniqueness (Beven, 2006) – the phenomenon that more than one set of parameter combinations can yield similar outcomes – contributes to parameter uncertainty during regionalization. Although we chose parsimonious models with few parameters in this study, equifinality effects are nevertheless a concern, particularly between loss and transform parameters. Supplementary Figures S1–S6 plot optimized loss model parameters (K eff, F c, CN) against optimized lag times for three Walnut Gulch basins. Comparisons generally show no correlation, except for F c and lag in basins 102 and 103 (R 2 = 25 and 18%, respectively). Only loss parameters were transferred from the donor to the target basins. Lag times for the Star Heights test basins were estimated solely from watershed characteristics.

Similarity-based regionalization is based on the premise that physical basin characteristics between the donor basin (i.e., WGEW) and the target catchment (Star Heights) are similar (Guo et al., 2021). While this is true for soil, vegetation, and rainfall characteristics, the donor and target catchments in this study differ in that Star Heights is partially urbanized, while WGEW remains largely in its natural state. We accounted for this difference by modeling impervious surfaces associated with urbanization separately and using regionalization only to parameterize loss models for pervious areas (open space, natural desert landscaping) in the Star Heights watershed. Good model performance, particularly for the LC model, provides strong validation for this approach.

In this study, we tested the regionalization performance of two parsimonious loss or infiltration models available in the widely-used HEC-HMS software: the curve number (CN) and the linear-constant (LC) methods. We first calibrated and validated the parameters for both models using long-term rainfall and runoff observations from three basins in the Walnut Gulch Experimental Watershed (WGEW), the donor basin. We then applied optimized loss parameter values to the target basin, a partially urbanized, densely instrumented watershed in central New Mexico. Regionalization performance was compared for an extreme event observed in 2025, as well as four smaller-magnitude storms. Three main conclusions emerge from the results of this study: (1) Validation based on long-term data from the WGEW shows that the LC model significantly outperforms the CN model for large storm events with cumulative precipitation exceeding 24 mm. For small events, performance was in general not significantly different based on four performance metrics (NSE, RMSE, absolute peak, and absolute volume error).

(2) Optimized curve numbers showed a dependence on cumulative precipitation depth, while optimized LC model parameters did not.

(3) When optimized parameters from the donor basin (WGEW) were applied to the target watershed (Star Heights), the LC model outperformed the CN model for most performance metrics.

This case study suggests that regionalization can be a viable strategy to support stormwater infrastructure design in data-scarce, semi-arid urban watersheds; however, its success appears to depend at least in part on the selected loss model.