SWAT+ model

The Soil and Water Assessment Tool Plus (SWAT+) ecohydrological model was used to simulate the hydrologic balance and resulting streamflows for the MARB system, which is based on predecessor SWAT codes [39–43]. The SWAT/SWAT+ models together represent extensive algorithm evolution over a period of three decades [43–46], and have been applied worldwide for a broad suite of water resource problems and environmental conditions for study regions that range from field-scale drainage areas to transnational river systems [41, 46–53]. SWAT applications have also been performed for all of the major MARB water resource regions as noted by Schilling et al. [23].

The SWAT+ simulations incorporated in this study have been performed within the National Agroecosystem Model (NAM), which provides the capability of executing the model for different areal resolutions including specific crop fields [54–57]. The utilization of SWAT+ within NAM provides enhanced flexibility in depicting various watershed elements (e.g., rivers, ponds, wetlands, rice paddies, crop fields, urban areas), and the hydrologic and pollutant flows between watershed elements. The NAM system is currently being calibrated and refined for the entire contiguous U.S., using a structure that consists of 2,121 individual SWAT+ watershed simulations [57]. These watersheds coincide with the boundaries of 8-digit watershed hydrologic unit codes (HUCs), which are part of the hierarchical system of basins and watersheds that have been delineated by collaborating federal agencies [23, 58]. The NAM structure also accommodates other HUC levels, such as providing output for the 86,016 12-digit watersheds across the contiguous U.S. [57]. The NAM MARB structure is comprised of 845 8-digit (HUC8s) and 34,944 12-digit watersheds (HUC12s) [36].

The smallest computation units incorporated in NAM are represented by hydrologic response units (HRUs) in SWAT+, which currently consist of roughly 4.6 million agricultural fields and 2.3 million non-agricultural areas [57]. SWAT+ is executed for each of these agricultural and non-agricultural HRUs, which respectively average 24 ha and 323 ha in size [57]. A complete set of management practices, such as planting, tillage, fertilizing, grazing, harvesting and irrigating, along with conservation practices and subsurface tile drainage, are accounted for the agricultural HRUs as applicable. However, output is generated at the 12-digit watershed level at present rather than for specific HRUs.

Intensive testing of SWAT+ within NAM is being conducted in two primary phases for the MARB region and the broader continental United States (CONUS). The first phase is focused on evaluating model components on a “soft data” basis, which consists of establishing the ability of SWAT+ to realistically replicate fundamental hydrologic characteristics for different stream systems (e.g., surface flow versus baseflow, subsurface tile drainage flow, groundwater effects), crop growth and yields, and other known attributes of specific watersheds. Calibration and validation of SWAT+ predicted corn and soybean yields has been reported for the entire CONUS [56], and similar testing is ongoing for other crops such as wheat and rice. Execution of different groundwater functions, including groundwater pumping, tile drainage, and groundwater-surface water interactions, has been performed with the SWAT+ gwflow component for six 8-digit watersheds representative of different U.S. regions [55]. The model’s ability to represent overall water balance, including discharge to precipitation ratios, has been confirmed for the MARB (Fig 2) and other CONUS regions. This is the stage of SWAT+ testing within NAM that represents the streamflow data used in this study.

Open in a separate window

Fig 2

Comparison of discharge-precipitation ratios measured at USGS stream gauging sites with ratios determined with SWAT+ simulations.

The second phase of NAM testing has been initiated which consists of calibrating flow duration curves (FDCs) and estimating sediment and nutrient loads across CONUS. Reliable results have been obtained with comparing FDCs versus the water balance that was obtained on the basis of the soft calibration process. A more rigorous testing step currently involves comparing FDCs with USGS gauge sites representative of major stream systems. Average annual suspended sediment, total nitrogen and total phosphorus loads estimated with SWAT+ are also being compared with corresponding loads determined with Load Estimator (LOADEST; https://water.usgs.gov/software/loadest/) at USGS gauges throughout CONUS. It is stressed that the testing of SWAT+ within NAM is an evolving process that incorporates refinements as improved monitoring and other data becomes available.

Definition of streamflow metrics

The calibrated SWAT models simulated streamflow at each HUC12 subwatershed from January 1, 2000, to December 31, 2018—a period of 19 years (6,940 days). The model outputs included daily estimates of the flow’s runoff, groundwater, and subsurface drainage components, which were summed to obtain total streamflow for each HUC12. Flow values were normalized by their tributary area to account for the variability in subwatershed size. Thus, modeled time-series of streamflow (mm/day) were generated in all ~35,000 HUC12s within the MARB. We calculated several statistics for each time-series, most of which are regularly used in hydrologic applications. These metrics were identified as potentially being informative of a basin’s flood risk, thus warranting exploration throughout the MARB. They are described herein and summarized in S1 Table.

The first group contained standard parametric statistics. The arithmetic mean is the average daily streamflow in a HUC12 and is defined as

where Q_i is the daily streamflow on day i and 6940 is the number of days in our 19-year simulated timeframe. The arithmetic mean describes the watershed’s overall water volume and can be interpreted as a measure of relative wetness. The standard deviation describes the dispersion in the daily flow values. It is defined as

Skewness quantifies the asymmetry of flow values. Daily streamflow in free-flowing rivers throughout the MARB consistently exhibits positive skewness [59], where high-flow events are often orders of magnitude larger than low-flow conditions. Skewness is defined as

Several standard nonparametric statistics were also found, including the 50^th (median), 90^th, and 99^th percentiles and the maximum daily flow.

Hydrologists and geomorphologists have identified parameters that describe the typical flow conditions under which streams exceed their bankfull capacity [60–62]. The ratio between the 99^th percentile and the mean daily flow has been used to approximate bankfull flow thresholds at sites with daily hydrologic observations but limited geomorphologic data [63, 64]. This statistic (hereafter referred to as ratio) provided a way to quantify bankfull discharge throughout the MARB. In this study, it was defined as

Stream flashiness is a hydrologic characteristic that describes the day-to-day volatility of streamflow values. The most common metric for flashiness is the Richards-Baker flashiness index (RB; [65]). For a given timeframe, the RB statistic calculates the daily changes in streamflow and normalizes this by the river’s total discharge. It is formally defined as

The standardized streamflow index (SSI) is often used to describe anomalous streamflow conditions [66, 67]. It quantifies the number of standard deviations from the mean for a given flow value. While its typical application is for drought assessment [68], the SSI can also assess abnormally wet periods [69]. The SSI was calculated on a daily basis in all HUC12s. SSI calculations require that that time-series are transformed into a normal distribution [70]. Lognormal distributions can adequately model daily streamflow in the MARB [59], so each SWAT flow output was log-transformed to normalize the datasets. Days with 0.0 mm of streamflow were set to 0.01 mm prior to transformation—a value corresponding to the minimum nonzero modeled flow. The daily SSI values were calculated using the following equation

where q_i is the log-transformed value on day i, $\overline{\mathrm{q}}$ is the mean of the log-transformed flows, and σ is the standard deviation of the log-transformed flows. For each HUC12, we found the maximum SSI value (maxSSI) and number of days where SSI was greater than 3 (SSI_3; i.e., the number of days the flow was more than 3 standard deviations above the mean).

Finally, we created an additional statistic that quantified the percentage of discharge occurring in the wettest days. A high proportion of flow discharged during only a few days of record indicates a watershed prone to intermittent high flow events and potential flooding. We termed this metric Top Days (TD) and defined it as

where k is the number of maximum flow days to include in the calculation, e.g., for k = 4, the TD metric sums the flows from 4 wettest days and divides this by the basin’s total flow. Therefore, TD is always a percentage bounded between 0.0 and 1.0. TD metrics were calculated using k values of 1, 4, 21, 69, and 139, corresponding to the 100, 99.9, 99.5, 99, and 98 flow percentiles over the record of 6940 days.

Each of these statistics can be calculated by utilizing daily streamflow records. They characterize various hydrologic regimes in all of the HUC12s that may help prioritize flood risk. It was unknown whether these parameters exhibit similar spatial tendencies—e.g. if watersheds with large maximum flows also contained the highest skewness values. In such an instance, using both the max and skew metrics to prioritize watersheds would prove redundant, as much of the same hydrologic behavior is being captured. To investigate the potential of informational redundancy among the metrics, we calculated their correlations and created spatial maps that color-coded the values by watershed. The statistics for each HUC12 are included in the S1 File.

Spatial patterns in the MARB

Spatial patterns of streamflow metrics evaluated in this study varied among the HUC12 watersheds within the MARB (Fig 4). Watersheds with higher mean flows were primarily located in the eastern half of the MARB (Fig 4A). Similarly, higher median and maximum daily flows were situated in the southeast and eastern portions of the MARB, respectively (Fig 4B and 4C, respectively). Watersheds with the highest 90^th and 99^th percentile flows were likewise located in the eastern portion of the MARB. In contrast, streamflow values associated with metrics such as skewness, SSI, and various TD classifications were much higher in western MARB watersheds (Fig 4C, 4D, and S1 File). The RB index and ratio metric did not show obvious regional patterns, with higher values more evenly distributed in HUC12s across the MARB (Fig 4F).

Open in a separate window

Fig 4

Maps of the HUC12 subwatersheds within the MARB color-coded by various flood metrics: A) mean, B) median, C) maximum daily flow, D) skewness, E) maximum SSI, and F) RB flashiness index.

The boundary between the eastern and western regions follows a line of longitude of approximately 100°W (Fig 1). East of this line, streamflow values were higher for metrics describing greater overall flows (mean, median, standard deviation, and percentiles) whereas west of the 100°W line, streamflow metrics favoring extreme values were highest (SSI, skewness, TD). A correlation matrix among the various streamflow metrics supports these spatial patterns (Fig 5). Spearman correlation coefficients of HUC12 streamflow values were highly positively correlated among the mean, standard deviation, median, and percentile rankings (r = 0.63–0.99). Similar positive correlations were observed among the SSI, skewness, and top day metrics (r = 0.72–0.99). However, the two general classes of streamflow metrics were negatively correlated with each other (Fig 5). The RB index was not well correlated with higher average flows (i.e., means) but had a higher degree of correlation with extreme flow metrics (i.e., SSI, TD).

Open in a separate window

Fig 5

Matrix containing the Spearman correlation coefficients among various flood metrics.

Flood damages

Expected annual losses (EAL) from flooding including damages to buildings and agricultural land were aggregated at the HUC12 level for the MARB. Monetary losses for buildings greatly exceeded expected annual losses for agricultural lands (Table 2). In watersheds where damages were reported, mean annual losses for agricultural lands ($3,815) were substantially lower than mean annual damages reported for buildings ($25,232), but the standard deviations associated with the losses were very high for both. Indeed, the minimum reported damage was $0, and there were no damages reported up to the 25^th percentile of HUC12 watersheds, and only $31 and $419 of losses were reported for agricultural lands and buildings, respectively, for the 50^th percentile (Table 2). EAL were highly skewed for both damage measures, with maximum building losses exceeding $16 million for one HUC12 watershed. Because building losses were so much greater than agricultural losses, their behavior dominated total annual losses.

EAL for agricultural lands and buildings were weakly correlated within the population of HUC12 watersheds of the MARB (Spearman correlation coefficient of 0.45; Fig 6). A cubic root transformation for both parameters was needed to account for the high degree of skewness and high number of zero values in the datasets.

Open in a separate window

Fig 6

Scatterplot showing the relation between agricultural losses and building losses.

Both datasets have been transformed using cubic roots.

Comparison of flood metrics

Hydrologic results showed a clear spatial difference within the MARB among two different classes of metrics, those describing higher overall flows such as the mean, median, and percentiles, and those characterizing extreme flows such as SSI, skewness, and TD. The 100°W line of longitude provides an approximate line of demarcation between the two flood-metric types (see Fig 4). In the eastern half of the MARB, greater annual precipitation leads to greater mean and median annual streamflows and greater flows associated with various statistical percentiles (i.e., 90^th, 99^th, maximum). Flooding is common in the Mississippi and Ohio River basins (e.g., [80–82] and is often driven by climate perturbations [8, 83]. Among these “normal-expected” flood metrics, there was a high degree of correlation, ranging from 0.63 to 0.94 (Fig 4), indicating these metrics were largely reporting on the same statistical pattern. The mean flow metric had the highest degree of correlation with other “normal-expected” variables.

In contrast, west of the 100°W line, conditions in the High Plains and Mountain West are considerably drier, and average daily streamflows are considerably lower. Despite dry conditions, sudden severe storms, rapid snowmelt, and other hydrologic events in mountainous areas can trigger extreme runoff in the western region and generate streamflows that greatly exceed normal conditions. Extreme flooding has been observed in many western locations (e.g., Big Thompson and Rapid City flash floods; [84]; Platte River flooding; [85]).

Metrics describing these rarer floods include streamflow skewness and SSI, a metric that describes the number of standard deviations from the mean for a given flow value [67]. We also introduced a new metric to capture rare and extreme floods that we called “top days.” With this metric, we quantified the percentage of total discharge occurring in the wettest days associated with the 100, 99.9, 99.5, 99, and 98 flow percentiles over the period of interest (in this case, 6940 days). Although spatial patterns of 1- and 4-day metrics are similar to those shown by skewness and SSI, the addition of more days within the various TD classifications gradually extends the distribution of higher values to the east as more water is included in wetter areas (see S1 File). However, in general, the streamflow metrics describing rare but extreme floods appear to produce similar spatial distribution within the MARB.

The ratio metric, defined as the ratio between the 99^th percentile and the mean daily flow, was highly correlated with the “rare-extreme” flood metrics (r = 0.66–0.97). This metric is typically used to approximate bankfull flow thresholds [63, 64], and the correlation of this metric with other extreme value metrics indicates that it is capturing the bankfull capacity of larger events [60–61]. The RB flashiness index was less correlated with the extreme event metrics (0.29–0.53) and was more uniformly distributed across the MARB. Further, the area of the basin highlighted by a higher RB index was unlike other spatial patterns (favoring the lower high plains and lower Mississippi areas) and thus was not particularly useful in our analysis.

The spatial overlap of flood metrics reveals that many of them broadly describe the same type of hydrologic behavior. Most statistics evaluated in this study seem well-categorized using two spatial groups: “normal-expected” and “rare-extreme.” Using multiple metrics within a group to identify flood risk would likely yield similar results due to the co-linearity of many values in subwatersheds. Studies exploring hydrologic indices have also found intercorrelation [86, 87] and noted that the resulting statistical redundancy might lead to ineffective decision-making when applying related metrics simultaneously [88].

Flood metrics and flood risk

Among the flood metrics evaluated in this study, those focused on higher overall flows were more correlated to EAL than the extreme value metrics (Fig 7). Although the correlation with agricultural land flood losses was not particularly strong (0.21–0.26), the correlation to building losses was higher (0.41–0.46), and the overall correlation was 0.43–0.48. In contrast, the extreme value metrics were more weakly, negatively correlated to EAL.

The better correlation of EAL to higher overall “mean” flows is consistent with geographic patterns of human settlement and agricultural production occurring in the eastern portion of the MARB (east of 100°W). Human societies evolved to be alongside rivers [89] and towns, cities, and larger metropolitan areas adjacent to rivers are prone to flooding and building losses. Recent major riverine floods include those in eastern Iowa [90], Memphis [91], and New Orleans [92], but most river floods do not receive scientific scrutiny. Instead, flood damages are often used to track impacts [28, 30]. While building losses within the MARB are concentrated in HUC12s that include major rivers, many watersheds report few, if any, losses. Indeed, half of the HUC12 watersheds in the MARB had building losses less than $419 (Table 2). However, when building losses were reported in a few outlier watersheds, damage amounts reached upwards of $16 million. Hence, a cubic root transformation was needed to account for the abundance of zero values and the high degree of skewness in the data.

Agricultural losses were also concentrated in the eastern part of the MARB, primarily in Iowa, Illinois, Indiana, the eastern Dakotas (glaciated regions), and along the floodplains of major rivers. The Upper Midwest regions dominated by corn and soybean production are continually threatened by excess precipitation [93]. For example, nearly $2 billion in federal relief funds were distributed for the Midwest Flood of 1993 [94]. Regions with higher mean flows correlate with areas of rain-fed agriculture in the Upper Midwest and production along the floodplains of major rivers. Hence, metrics that reflect higher overall flows (mean, standard deviations, and various flow percentiles) are better correlated with agricultural damages.

Although extreme flow metrics (SSI, skewness, TD) were negatively correlated with EAL, the negative relation mainly indicates that HUC12 watersheds with greater risk to extreme flows have fewer agricultural lands and buildings damaged by flooding. Extreme floods are primarily associated with dryer and mountainous western regions (Fig 3). When extreme flooding occurs in the area, much of the floodwater is dissipated in floodplains or larger valleys where there is less infrastructure to be damaged, and the land cover is dominated by perennial vegetation. Because of this, flood metrics associated with these western areas do not represent the remainder of the MARB particularly well.

A goal of this study was to identify a metric that could be used to prioritize flood risks across ~35,000 HUC12 watersheds in the MARB where future NI practices could be installed for flood risk reductions [23]. It is important to recognize that EAL values associated with flood damages are inappropriate flood risk indicators. EAL values were based on a fixed 23-year period (1996–2019), which does not necessarily indicate future risk. A small minority of watersheds witnessed catastrophic river flooding floods during that period [90–92], resulting in a subset of HUC12s with extremely high damages. Conversely, numerous HUC12s (~25%) saw no damages from 1996–2019. The EAL values appear too prone to the random fluctuations associated with flooding to be used as a viable tool for identifying long-term risk [95].

At this point, we suggest that the mean flow metric should be used as the best indicator for flood risk. The mean flow metric correlates with building and agricultural losses occurring in HUC12 watersheds, and it represents a statistical parameter that is easily determined, well understood by users and widely comparable across the basin. Other measures of flows describing higher overall water yields, such as the median, maximum, or various flow percentiles, were similarly correlated to EAL, but the statistical differentiation and spatial distribution of these metrics were more limiting across the MARB. Unlike the median, all mean values of streamflow were nonzero. Risk identification based on the median would assign zero flood risk to a handful of watersheds, but this is prevented by using the mean. In a previous assessment of combined flood and nitrate-nitrogen risks across the MARB, flood risk was assumed to be proportional to mean flows [23]. With this study, we provided additional consideration of other flood metrics and, in the end, returned to the original assumption.

USDA is an equal opportunity provider and employer.