Description of Study Area

The study area comprises the upper White River hydrologic unit (05120201) (hereafter referred to as the “upper White River Basin”), located predominantly in central and east-central Indiana (fig. 1). The upper White River Basin drains an area of approximately 2,718 square miles (mi²) and contains all or part of 16 counties that, in 2020, included 7 of the 20 most populated cities in Indiana (Indianapolis, Carmel, Fishers, Muncie, Noblesville, Greenwood, and Anderson) (STATS Indiana, 2022). Indianapolis was the most populated city in Indiana, with a 2020 population estimated at more than 887,600 and a population of more than 2,111,000 in the Indianapolis metro area that includes Indianapolis, Carmel, and Anderson, Indiana (STATS Indiana, 2022).

There are more than 2,180 miles (mi) of streams (Tedesco and others, 2005) and four water-supply reservoirs (Eagle Creek, Geist, Morse, and Prairie Creek Reservoirs) in the upper White River Basin. Although Eagle Creek Reservoir has been used for water supply for the City of Indianapolis since 1976, the reservoir was developed primarily for flood control (Tedesco and others, 2005).

The drainage area upstream from the study gages is predominantly within the Tipton Till Plain physiographic region, with a small northerly part of the drainage area lying in the Bluffton Till Plain physiographic region. The Till Plains were formed by continental glaciation during the last Ice Age and consequently have low topographic relief. The Till Plains have a covering of glacial till and outwash material composed of 100 to 200 feet of silty clay till interspersed with thin (5 to 10 feet) layers of sand and gravel (Tedesco and others, 2011). Till Plain soils have attracted widespread agricultural land use but tend to have low infiltration rates so that tile drains are commonly installed to improve drainage and reduce surface runoff.

Based on an analysis of land-cover data from the 2019 National Land Cover Database (Dewitz and U.S. Geological Survey, 2021), the dominant land cover in the upper White River Basin was agricultural (constituting approximately 55.6 percent of the basin) followed by developed (29.5 percent) and forested (13.3 percent) land covers. In general, agricultural land cover in the upper White River Basin decreased as a percentage of the drainage area from the headwaters in Randolph County to downstream locations on the White River. Those decreases in agricultural land cover were predominantly offset by increases in developed land covers.

Central Indiana has a humid-continental climate characterized by distinct summer and winter seasons, large annual temperature changes, and highly variable weather patterns (Tedesco and others, 2005). For the period 1991–2020, the average annual temperature and precipitation for Indianapolis was 53.7 degrees Fahrenheit and 42.62 inches, respectively (National Oceanic and Atmospheric Administration, 2022). According to Widhalm and others (2018), Indiana has warmed, and average precipitation has increased since 1895, with more precipitation falling in heavy downpours.

Purpose and Scope

The purpose of this study is to estimate annual mean concentrations and fluxes of nutrients (nitrate plus nitrite, total Kjeldahl nitrogen, and total phosphorus) and total suspended solids (TSS) at three streamgage locations on the upper White River, and to assess temporal trends in concentration, flux, streamflow, and precipitation. The estimates and trend results in this report update and extend results from a previous study (Koltun, 2019) using three additional years of data and newer estimation methods. Nutrient and TSS concentration data from several sampling locations in the upper White River Basin were used to look for spatial patterns in their sampling distributions. In addition, a variety of anthropogenic factors that could influence the concentrations and flux of nutrients and TSS were examined to better understand their spatial and temporal variations. The analyses are based on water-quality data collected in water years 1992–2020 by Federal, State, and local agencies; and streamflow data collected by the USGS.

Streamflow and Water-Quality Data Sources

The following agencies collected water-quality data used in the analyses:

The data (Koltun, 2023) were from sampling locations where concentrations of one or more of the following constituents were sampled:

Water-quality data collected at the sampling sites were paired with streamflow data measured at the study gages to estimate annual mean concentrations and fluxes, and to evaluate temporal trends in concentration and flux. Nutrient and sediment data are instantaneous concentrations but were treated as daily mean concentrations in the analyses. Table 1 lists water-quality sampling site locations, the agencies whose data were used for those sites, and the study gages that were paired with the sampling sites. Water-quality results were paired with daily mean streamflow data measured at one of the study gages if (1) there were no known nutrient or sediment inputs between the study gage and sampling site that would likely cause the concentrations measured at the sampling sites to be unrepresentative of concentrations at the study gage, and (2) there was less than 10 percent difference between the drainage areas at the sampling site and the study gage. The only parameter and site combination that might not have met those criteria was TSS for the study gage near Centerton, Indiana. IDEM’s site WWU160-0004 was the only sampling site near the Centerton study gage where TSS data were collected. There were three facilities, all non-publicly owned treatment works that discharged to the White River between the gage and the IDEM sampling site, which reported TSS concentrations as part of their permit requirements. Two of the facilities were classified as minor dischargers. The third facility was classified as a major discharger for non-contact cooling water; however, the facility was required to sample for TSS in storm-water runoff. It is not known whether discharges from those facilities resulted in appreciable changes in concentrations measured at the sampling site relative to concentrations at the gage; however, concentrations of TSS measured at the sampling site were assumed to be representative of concentrations at the gage.

Ideally, the same analytical period would be used for all sites for a given water-quality constituent to help facilitate comparisons between sites. To do so would have severely limited the analytical period for some constituents. For most constituents, there were suitable data to facilitate analysis of the 29-year period extending from water years 1992 to 2020 (calendar years 1991 to 2020); however, because of unavailability of sample data for part of that period, shorter analytical periods were necessary for total Kjeldahl nitrogen at the Muncie and Centerton gages and TSS at the Centerton gage (table 2).

Data were retrieved on July 15, 2021, from the National Water Quality Monitoring Council’s Water Quality Portal (https://www.waterqualitydata.us) for evaluating how distributions of sampled nutrient and TSS concentrations varied across sites within the upper White River Basin. The Water Quality Portal includes data from the USGS National Water Information System and from the U.S. Environmental Protection Agency’s Water Quality Exchange (formerly the STORET database) which in turn contains data from more than 400 State, Federal, Tribal, and local agencies. Data from the Water Quality Portal were from the Indiana IDEM, Indiana STORET, or the USGS Indiana Water Science Center. Data retrieved from the Water Quality Portal were filtered for water years 1992–2020, and there was some overlap with the dataset used to estimate concentrations and fluxes, and to evaluate trends. Data from multiple collecting agencies were combined if more than one agency collected water-quality data for a particular constituent at the same location.

Most agencies collect dip or grab samples (samples collected by filling an open container held beneath the surface of the water), whereas other agencies (such as the USGS) frequently collect samples using depth and (or) width integrating techniques and isokinetic samplers intended to produce concentrations that are more representative of a flow-weighted mean (Edwards and Glysson, 1999; U.S. Geological Survey, 2006). To use all available water-quality data, analytical results from samples collected by all agencies were treated as flow-weighted mean concentrations.

Data Quality

Many agencies perform internal quality-control (QC) checks of the water-quality data they produce. Some agencies (for example, the USGS and CEG) document the results of their QC checks in their databases. The water-quality data from the various agencies initially were assumed to be valid unless there was indication that a result was rejected because of QC (or other) issues.

Many constituents (such as sediment and total phosphorus) transported in runoff tend to increase in concentration with increasing streamflow; however, concentrations of constituents not strongly associated with runoff might decrease because of dilution. In either case, a pattern of either increasing or decreasing concentration with increasing streamflow is common. The relation between streamflow and pollutant concentrations typically is not strictly monotonic; concentrations for a given streamflow vary over time because of a variety of factors. Even with that variability, concentrations associated with a given streamflow typically tend to scatter within a limited range. Values that plot far away from the majority of the data in scatter plots are referred to as “outliers.”

Transport plots are scatter plots in which streamflow is on the x-axis and concentration (or flux) is on the y-axis (Glysson, 1987). Rudimentary quality-control screening of the data was done by creating transport plots for each constituent and examining them for outliers (values that plot far away from the rest of the data). Outliers were subsequently evaluated to determine whether they were likely erroneous and should be omitted from the dataset. That evaluation resulted in removal of water-quality results for constituents associated with station 03351000 (the gage near Nora) for a single day (January 27, 1994) and one TSS result associated with station 03354000 (the gage near Centerton) on May 13, 1996. All other water-quality results were retained. Constituent-specific arithmetic mean concentrations were substituted for individually measured concentrations if more than one sample was collected on the same day.

Concentration and Flux Estimation and Trends

Weighted Regressions on Time, Discharge and Season (WRTDS) as implemented in the Exploration and Graphics for RivEr Trends (EGRET) package (version 3.0.6; Hirsch and DeCicco, 2015) was used to estimate and assess temporal trends in annual mean concentrations and fluxes of nutrients and TSS. The WRTDS Kalman filter method (WRTDS-K) was used to estimate concentrations and flux. WRTDS-K makes the estimates of concentration equal to the measured values on sampled days and uses a first-order autoregressive model to account for the autocorrelation structure of model residuals for computing estimates for unsampled days (Zhang and Hirsch, 2019). WRTDS-K provides generally better daily estimates of concentration and flux (Zhang and Hirsch, 2019) than the base WRTDS method used in the earlier Upper White River analyses by Koltun (2019). The WRTDS-K analyses done for this study used the default values for the number of iterations and the lag-one autocorrelation value (rho). The WRTDS-K method was used only to obtain estimates of annual mean concentration and flux. The base WRTDS estimation method was used for evaluating trends.

The EGRET package was used to evaluate whether streamflow characteristics changed over time. EGRET analyses of long-term changes in streamflow characteristics are based on time-series smoothing methods pioneered by Cleveland (1979) and Cleveland and Devlin (1988). EGRET performs locally weighted scatterplot smoothing (LOWESS) on annual streamflow statistics (relevant to low, high, and mean streamflows) to produce plots that show patterns of change over time spans of about a decade or more (see Hirsch and others [2010] for details on the LOWESS method). LOWESS plots were also used to evaluate trends in annual precipitation totals.

Quantile-Kendall plots were prepared for each study gage. A Quantile-Kendall plot is a scatter plot of streamflow non-exceedance probability versus the trend slope for streamflows associated with that probability calculated with a Mann-Kendall trend test (Mann, 1945), expressed in percentage change per year. Larger non-exceedance probabilities are associated with larger streamflows and smaller nonexceedance probabilities are associated with smaller streamflows. Each point in the plot is color coded according to the p-value for a test of the null hypothesis that the trend slope is zero. All hypothesis tests in this study were based on an alpha level of 0.05.

In addition to LOWESS-based assessments of trend in streamflow statistics, two tests for nonstationarity in streamflow were performed. Nonstationarity is a condition in which the probability distribution of some process changes with time. Pettitt tests (Pettitt, 1979) were used to identify step trends in the natural logarithms of water-year annual streamflow statistics. Step trends differ from gradual trends in that they are abrupt, therefore step trends are not easily addressed with tools such as flow normalization in WRTDS (Hirsch and others, 2010). The null hypothesis for the Pettitt test is that there are no step trends in the time series (in other words, the observations are independent and identically distributed); whereas the alternative hypothesis is that there is at least one step trend in the time series. The Mann-Kendall test (Mann, 1945) was used to test for the presence of more gradual monotonic (unidirectional) trends in the natural logarithms of water-year annual streamflow statistics and calendar year precipitation totals. The Mann-Kendall trend test does not assume the shape of the trend (for example, whether the shape is linear), but instead tests whether values tend to increase or decrease overall. The null hypothesis for the Mann-Kendall test is no trend in the series. Sen’s slope (Sen, 1968) (the median of the slopes of all lines through pairs of points) computed the median linear rate of change in the natural logarithms of water-year annual streamflow statistics. The Pettitt and Mann-Kendall tests and Sen’s slope computations were done using the trend package (Pohlert, 2018) in R (R Core Team, 2017).

The WRTDS-K model calculates the expected value of concentration as a function of streamflow and time. The form of each WRTDS-K model (Hirsch and others, 2010) is:

where The equation for the WRTDS base method is identical to the equation for the WRTDS-K method but omits the σ_i term that accounts for the autocorrelation structure of the model residuals. The sine and cosine terms in equation 1 help account for seasonality. WRTDS uses a weighted regression technique where the weights are determined as a function of the “distance” between the estimation point (defined by Q_i and t_i) and the sample points. The measure of distance is defined in three dimensions: ln(Q), t, and season (proximity to the same time of year).

WRTDS accommodates water-quality results that are left censored or interval censored. Left-censored results are less than some value and interval-censored results are less than one value, but greater than another value. There were some left-censored results in the water-quality data reported by the various agencies, but none that were interval censored. Reporting of censored data differed by agency. For example, the Indiana Department of Environmental Management reported a laboratory result and a laboratory reporting limit. If the result they obtained was less than the laboratory reporting limit, then they reported the result as -1. By comparison, the USGS reported a left-censored value with a “less than” remark code (<) and a result value equal to the laboratory reporting limit. Each agency’s convention for reporting censored data was considered when creating input files in the format required for left-censored data in WRTDS.

The concentration and flux of a constituent can be strongly influenced by the time history of associated streamflow conditions. It can be unclear whether a year-to-year decrease in flux of a constituent is primarily because of changes in streamflow conditions or decreases in concentrations of the constituent associated with a given range of streamflows. For accurate information about temporal changes to a watershed’s constituent transport characteristics, it is crucial to remove the effects of year-to-year variation in streamflow without removing the influences associated with seasonal and long-term trends in streamflow. The method used to accomplish this in WRTDS is referred to by Hirsch and DeCicco (2015) as “flow normalization.” The flow-normalized concentration and flow-normalized flux on a given day are calculated as:

and where The flow-normalized concentration on a day t of the record is the mean of the concentrations estimated from the daily mean streamflows that occurred on that calendar day-of-year over the entire record. The corresponding flow-normalized flux is the mean of the product of the estimated flow-normalized concentrations, the daily mean streamflows (from which the flow-normalized concentrations were estimated), and a units conversion factor. For a given calendar day, flow-normalized values in sequential years are calculated using the same set of daily streamflows (the streamflows that occurred on that calendar day in each year of record). The concentrations estimated for those streamflows change from year to year because they are computed on the basis of observations whose weights are changing as a function of their proximity to the time of estimation. Daily values are averaged over the year to compute annual flow-normalized mean concentration and flux.

For complex statistical methods (such as the smoothing procedure applied in WRTDS), it generally is not feasible to calculate the uncertainty of results using simple mathematical expressions (such as those that apply to ordinary least-squares regression). Instead, bootstrapping is a common approach to describing the uncertainty of more complex analyses (Hirsch and others, 2015). Consequently, trends in flow-normalized concentrations and flux were evaluated by use of the WRTDS Bootstrap Test (WBT) contained in the EGRETci R package (version 2.0.4; Hirsch and others, 2015). The WBT uses a random sampling procedure with replacement to create multiple subsets of measured concentrations from the original dataset. Flow-normalized annual concentrations and fluxes are then calculated from each subset. Results from model iterations are used to estimate (1) the uncertainty of the flow-normalized annual values of concentration and flux and (2) the level of significance of changes in flow-normalized annual values between selected water years. Hirsch and others (2015) adopted the definitions shown in table 3 to describe the degree of statistical support that the dataset and WBT results provide regarding the likelihood associated with a direction of change over time. Because trends can be upward or downward, table 3 lists alternate descriptions depending on whether describing the likelihood of an upward or downward trend. For example, if an upward trend is “highly likely” that also means that a downward trend is “highly unlikely.”

In the WRTDS analysis, a one-way level of significance of 0.05 was chosen to identify a statistically significant trend. Consequently, a trend categorized as “highly likely” is statistically significant. Analyses discussed in this report were completed using RStudio version 1.1.463, R version 3.6.1, EGRET version 3.0.6 and EGRETci version 2.0.4. The R code used to perform the WRTDS analyses can be downloaded from ScienceBase (Koltun, 2023).

Assessment of Spatial Patterns in Sampling Distributions of Nutrients and Total Suspended Solids

Spatial patterns in the sampling distributions of nutrient and TSS concentrations were assessed by comparing boxplots of concentrations measured at selected sites in the upper White River Basin. Data retrieved from the Water Quality Portal were filtered to include only those data for sites that had 50 or more observations of total Kjeldahl nitrogen, total phosphorus, or TSS, and 10 or more observations of nitrate plus nitrite. The minimum threshold of 50 observations was a tradeoff between having (1) sufficient observations to produce representative sample distributions and (2) a large and spatially distributed set of sites to compare. The lower minimum threshold for sites with nitrate plus nitrite data was because no site had more than 11 observations and only observations from July, 2019 to July, 2020 were in the Water Quality Portal. Because of the smaller number of observations of nitrate plus nitrite, there is increased potential for sample distributions to be uncharacteristic of their parent populations (because of sampling error) than the other constituents; therefore there is more uncertainty when comparing distribution characteristics between sites.

Evaluation of Potentially Influential Anthropogenic Factors

Data were gathered from a variety of sources on anthropogenic factors that could influence the concentrations and fluxes of nutrients and TSS in the upper White River. The data were processed to compute time series of summary results for the upper White River hydrologic unit or sub-drainages corresponding to the intervening drainage areas at the study gages. Plots illustrating how the factors varied as a function of time were constructed to facilitate visual comparisons with water-quality trend results.

Temporal Trends in Streamflow and Precipitation

Temporal trends in streamflows between water years 1978 and 2020 at the study gages were assessed using EGRET and using Mann-Kendall and Pettitt tests. Only streamflow records after the 1977 water year were used for trend analyses to avoid changes in streamflows associated with the adoption of Eagle Creek Reservoir as a water supply for the City of Indianapolis in 1976. Without exception, the annual maximum, mean, and median daily streamflows and the annual minimum 7-day average streamflows (hereafter referred to as the “annual 7-day minimum streamflows”) showed indications of temporal trend in the EGRET analyses (figs. 2–4).

The positive tau (Kendall's nonparametric correlation coefficient) values from the Mann-Kendall analyses indicated that the streamflow statistics at all study gages generally trended upward over time between water years 1978 and 2020 (table 4). However, only trends for annual mean and maximum daily streamflows were statistically significant for all study gages (table 4).

Like the Mann-Kendall tests, the quantile-Kendall plots (figs. 5–7) indicated upward temporal trends in streamflow and which parts of the streamflow regime in which non-zero trend slopes are most likely. Very low streamflows (non-exceedance probabilities less than about 0.05) and mid-range and larger streamflows (non-exceedance probabilities greater than about 0.5) showed the greatest likelihood of non-zero trend slopes at the study gages at Muncie and near Nora. Similar results were found for the study gage near Centerton; however, there were weaker indications of significant trend in a small part of the streamflow quantiles greater than about 0.5. At the three study gages, streamflows in the third quartile of non-exceedance probabilities (between 0.5 and 0.75) had positive trend slopes with most quantiles showing high likelihood of temporal trend (p-values less than or equal to 0.05).

The Pettitt tests indicated that a step trend (abrupt change) occurred in annual mean daily streamflows at each of the study gages between water years 1978 and 2020 (table 5), and that the step trends most likely occurred at each gage around water year 2001. Other statistically significant step trends were indicated for the annual maximum daily streamflows at the gages near Nora and Centerton, and the annual median daily streamflow at the gage at Muncie. No statistically significant step trend was indicated for the annual 7-day minimum streamflows at the study gages (table 5); however, the Pettitt test is less effective at detecting changes in extremes than changes in means or medians (Mallakpour and Villarini, 2016). For statistically significant step trends, the water year change point was consistently 2001. The step trend in streamflows adds uncertainty to the interpretation of trends in the WRTDS analysis near the year of the step trend (in this case, water year 2001). That uncertainty decreases as you move further away in time from the year of the step trend.

To assess whether trends in streamflow might be related to changes in precipitation, a Mann-Kendall analysis was done using annual precipitation measured at the Indianapolis International Airport (network ID GHCND: USW00093819) between calendar years 1932 and 2020. That analysis indicated a weak but statistically significant upward trend (tau=0.20, p=0.005) in annual precipitation with a Sen’s slope of 0.089 inches per year. Consequently, though not definitive, this weak trend supports a hypothesis that upward trends in streamflows resulted (at least in part) from the upward trend in precipitation. LOWESS smooth lines fit through the annual precipitation totals measured at the Indianapolis International Airport (fig. 8) and annual mean streamflows measured at the White River near Nora gage (03351000) (fig. 9) indicate that precipitation and streamflow began trending upward more rapidly in the mid to late 1970s.

Evidence of the step trend in annual mean daily streamflow identified in the Pettitt tests for the gage near Nora (table 5) can be seen in figure 9. Annual mean daily streamflows less than about 1,000 cubic feet per second (ft³/s) occurred in approximately 39 percent of the years during water years 1931–2000. However, during water years 2001–20, the annual mean daily streamflow exceeded 1,000 ft³/s in all but one year.

Estimated Mean Concentrations and Flux of Sediment and Nutrients

Water-quality data were collected for a wide range of streamflows and in all months. The medians and interquartile ranges of streamflows at the study gages on the days when samples were collected were similar to those for the entire analytical period at the gage; however, streamflows at the study gages did not range as high as streamflows for the analytical period. Streamflows on sample dates included nearly the entire range of streamflows for the gage near Nora; however, there were fewer samples from high flow days at the Muncie and near Centerton gages. The absence of sample data at very high flows adds uncertainty because computation of concentrations and fluxes on those days requires extrapolation.

Water year annual mean daily concentrations and fluxes estimated with WRTDS (with and without flow normalization) are in tables 6–9 for TSS, total phosphorus, nitrate plus nitrite, and total Kjeldahl nitrogen. Confidence limits for the estimates have been omitted from tables 6–9; however, the annual mean concentration and flux estimates (from the WRTDS-K method) and their 90 percent confidence limits are available for download (Koltun, 2023). Estimates from the WRTDS-K method tended to be of equal or smaller magnitude than estimates from the base WRTDS method. For example, cumulative loads of total phosphorus for water years 1992–2020 computed with the WRTDS-K method at the three study gages ranged from about 90.1 to 96.9 percent of loads for the same period computed with the WRTDS base method.

To facilitate comparisons between gage locations, analytical-period and mean-annual loads and yields (load divided by drainage area) are summarized for the longest periods of concurrent record of each constituent at the three study gages (tables 10–11). Analytical-period and mean-annual loads of each of the constituents increased sequentially from the most upstream gage at Muncie to the most downstream gage near Centerton. Yields did not consistently follow the same pattern as loads. In this case, the largest analytical-period and mean-annual yields of TSS and total phosphorus were at the most upstream gage at Muncie, but the largest yields of nitrate plus nitrite were at the most downstream gage near Centerton. It is common for yields of nonpoint constituents in streams to decrease in a downstream direction because stream channel and land surface gradients typically decrease, and there is more opportunity for constituent losses (Koltun, 2019). The two most upstream gages had about equal total Kjeldahl nitrogen yields that were each about 13 percent larger than the yield of total Kjeldahl nitrogen at the most downstream gage.

Changes and Trends in Flow-Normalized Concentration and Flux

Table 12 lists WBT results that show the magnitude, direction, and likelihood of change in flow-normalized concentrations and flux of TSS, total phosphorus, nitrate plus nitrite, and total Kjeldahl nitrogen between water years 1997 and 2019. A similar analysis previously examined changes between 1997 and 2017 (Koltun, 2019). Differences in results between the two analyses reflect the addition of 2 years of data as well as using a more recent end year (2019) to compare against 1997. With two exceptions, concentration and flux changes that were statistically significant in Koltun (2019) showed the same statistically significant change direction. The first exception was nitrate plus nitrite concentration at the White River at Muncie, Indiana, which showed a downward change; however, the likelihood of change decreased from “highly likely” to “likely.” The second exception was for the White River near Centerton, Indiana, which previously showed a “highly likely” downward change in TSS concentration but now shows an upward change classified as “as likely as not.” Because TSS concentration data were not collected at or near the gage near Centerton between water years 2008 and 2017, the previous study’s change results for that constituent were based on data collected only between water years 1992 and 2007. Results from the current study include sample results from water years 2018–20 because sampling for TSS was restarted at the gage near Centerton in water year 2018. This shift in change direction demonstrates the degree to which perceived trends in concentration and flux are a function of the endpoints chosen for comparison. To provide a more temporally continuous picture of trends, plots were constructed showing flow-normalized annual mean concentrations and fluxes of nutrients and TSS as a function of time (figs. 10–13). The flow-normalization process attempts to remove the effects of interannual variation in streamflow on annual mean concentrations and fluxes so that trends driven by changes in the relation between streamflow and concentration are more apparent, whereas the WRTDS-K estimate includes the effects of changes in the streamflow-concentration relation and year-to-year changes in streamflow. Consequently, in a given year, the annual flow-normalized values can differ appreciably from annual mean concentrations and fluxes estimated with the WRTDS-K method (figs. 10–13). Uncertainties associated with the flow-normalized annual mean concentration and flux estimates vary over time; however, confidence intervals for the estimates have been omitted from the plots to facilitate comparisons of temporal trends between study gages.

Flow-normalized annual fluxes increased at the study gages from upstream to downstream (downstream order) during all years in the study period (figs. 10–13). In the period between 2010 and 2020, flow-normalized nutrient fluxes at the gage near Centerton showed some of the most rapid decreases of the three study gages. Flow-normalized nutrient fluxes at the gage near Nora also decreased (albeit less rapidly) over the same period; however, they increased some at the Muncie gage. Flow-normalized TSS fluxes at the gages at Muncie and near Nora (fig. 10) also decreased between 2010 and 2020. At the gage near Centerton, flow-normalized TSS flux was not reported for water years 2008–17 because of lack of sample data; however, flow-normalized TSS fluxes computed for water years 2018–20 are larger than the values computed for water years 1992–2017 and more than 4 times larger than corresponding water-year flux at the gage at Nora (fig. 10). This increase might be an artifact of the large sampled data gap (2008–2017), therefore the 2018–20 flow-normalized flux results for TSS at the gage near Centerton and the related change-analysis results might be misleading. Because of the weighting schemes used in WRTDS, TSS trend and flux results for 2018–20 for the gage near Centerton are more strongly influenced by data collected during 2018–20 than the trend and flux results for the period preceding the data gap (1992–2007). Hirsch and DeCicco (2015) discuss issues that can arise when using WRTDS with datasets that have large data gaps.

Annual flow-normalized concentrations did not increase or decrease consistently in downstream order at the study gages over the study period (figs. 14–17). Flow-normalized concentrations of nutrients were always larger for the gages near Nora and Centerton than for the gage at Muncie. Flow-normalized concentrations of TSS were not consistently larger at one site than at the others. The gage near Nora was the only site where flow-normalized concentrations of a constituent (TSS) consistently decreased over the analytical period.

The WRTDS-K and flow-normalized annual mean concentrations are time-weighted means and consequently are influenced by the most frequent hydrologic conditions, which tend to be closer to median streamflow conditions than the larger magnitude mean-streamflow or high-streamflow conditions. Fluxes, on the other hand, tend to be more heavily influenced by conditions that occur during high-streamflow periods. These computational differences can cause counterintuitive outcomes in which concentrations are decreasing over time while fluxes are increasing, or in which concentrations are increasing over time while fluxes are decreasing.

Spatial Patterns in Sampling Distributions of Nutrients and Total Suspended Solids

Figures 18–21 show boxplots depicting the distributions of sampled concentrations of nutrients and TSS at 20 mainstem, tributary, and distributary water-quality sampling locations in the upper White River Basin (table 13; fig. 22). Sampling locations designated A–I are on the mainstem of the upper White River, and sampling locations designated J–T are on other waterways in the upper White River Basin. If a boxplot for a sampling location is not in a figure for a given constituent, the sampling location did not meet the minimum observation criteria for that constituent.

One pattern in figures 18–19 and 21 is that the largest median concentrations of TSS, total phosphorus, and total Kjeldahl nitrogen were associated with mainstem upper White River sampling locations downstream from Indianapolis (sites G–I). Another pattern is that median total phosphorus, and total Kjeldahl nitrogen concentrations at the upper White River site immediately downstream from Muncie (site C) were elevated relative to bracketing upstream/downstream mainstem sampling locations. Location C also had one of the highest median nitrate plus nitrite concentrations and the highest concentrations of total phosphorus, total Kjeldahl nitrogen, and nitrate plus nitrite of all sampling locations.

The median concentrations of all constituents at the Indianapolis Waterway Canal sampling location (K) were larger than at the other non-mainstem sites. In contrast, the Little Buck Creek sampling location (S) had the lowest median concentrations of total phosphorus, and total Kjeldahl nitrogen compared to other locations. Except for the Indianapolis Waterway Canal sampling location (K), median concentrations of nutrients and TSS at non-mainstem upper White River sampling locations generally had magnitudes in the same range or smaller than the mainstem upper White River locations. The canal was originally intended to be part of a transportation corridor, which was never completed because of insufficient funding. The canal was later repurposed to convey water from the White River to a downtown water treatment plant.

Most of the Water Quality Portal sampling locations had no associated streamflow data; consequently, interpretation of concentration data alone is impeded by a lack of hydrologic context. For example, boxplots of concentrations measured during predominantly low to medium streamflows might look very different if concentrations were measured for the full range of streamflows. This variable relation between streamflow concentration was illustrated by Koltun (2019) with the “U” or check-mark-shaped transport plots for total phosphorus, nitrate plus nitrite, and total Kjeldahl nitrogen for the study gages near Nora and near Centerton. Without associated flow data, there is no way to determine whether differences in the distributions of concentrations measured at different sites might be because of sampling strategy versus other factors.

Population

The relationship between human populations and water resources is complex. Human populations can affect a variety of hydrologic and environmental characteristics through the way they modify and (or) manage the land and water resources. Census data were obtained for 1990–2020 for the Indianapolis and Muncie metro areas and the City of Anderson, Indiana (STATS Indiana, 2022), the three most populated city/metro areas in the upper White River Basin. Figure 23 shows that the population of the Indianapolis metro area increased rapidly at a uniform rate between 1990 and 2020, while populations in the City of Anderson (which is part of the Indianapolis metro area) and the Muncie metro area decreased slightly. Between 1990 and 2020, the population increase in the Indianapolis metro area far exceeded the population decrease in the Muncie metro area.

Land Cover

Changing land cover could change the hydrology and water chemistry of streams draining the affected areas. Changes in hydrology and water chemistry can result from (1) how quickly and how much water runs off land surfaces in each land-cover category, as well as (2) the quantity and wash-off/transport characteristics of contaminants associated with the various land covers.

Land-cover characteristics for the period 2001–19 were computed from data in the National Land Cover Database (Dewitz and U.S. Geological Survey, 2021) for the intervening land areas draining to the study gages. The percentages of the upper White River Basin in the intervening areas draining to the study gages at Muncie, near Nora, and near Centerton are about 8.9, 36.0, and 45.1 percent, respectively. About 10 percent of the upper White River Basin drainage is downstream from the gage near Centerton. Because the gage at Muncie is the most upstream study gage, its intervening land area includes the entire headwater drainage upstream from the gage. Land-cover characteristics in the drainage to the gage at Muncie changed little during the 2001–19 period (fig. 24). More than 73 percent of the area draining to the Muncie gage was classified as cultivated crops, with deciduous forest as the second largest land-cover category at about 8.5 percent of the area. Some land-cover characteristics in the drainage between the gage at Muncie and the gage near Nora changed slightly during the 2001–19 period (figs. 25–26). Cultivated crops decreased from about 68.2 percent to 66.8 percent of the area, whereas developed land covers cumulatively increased from about 18.6 to 20.6 percent of the area over the same period.

Land-cover characteristics in the drainage between the gage near Nora and the gage near Centerton (fig. 26) were different in character than those for the other two gages and showed noticeable changes in some land-cover categories during the 2001–19 period. In general, the intervening drainage between the gages near Nora and Centerton had smaller percentages of cultivated crops and larger percentages of developed land cover than the intervening areas of the other two study gages. Cultivated crops in the intervening drainage between the gages near Nora and Centerton decreased from about 40.1 to 36.6 percent of the area during the 2001–19 period and developed land covers cumulatively increased from about 37.5 and 42.3 percent of the area. Pasture and hay also decreased from about 8.4 to 7.2 percent of the area.

Crops and Operational Tillage

In a landscape with predominantly agricultural land uses, the delivery of nutrients and other contaminants to streams and lakes might be influenced by the amount and type of crops that are grown, use of winter cover crops, crop rotation, and tillage practices. For example, a recent analysis of data from midwestern streams indicated that yields and flow-weighted mean concentrations of nitrate from areas of soybean cultivation were larger than from areas of corn cultivation on an equal percentage-of-area basis (Piske, 2019). Rotational planting of corn and soybeans is common in the Midwest, which lets farmers use less nitrogen fertilizer for growing corn, as compared to continuous corn production. The use of less fertilizer, in turn, reduces the risk of having large amounts of land-applied fertilizer available for transport by runoff during storms. Winter cover crops and no-till/reduced-till practices reduce erosion by wind and water (Snapp and others, 2005) by binding and (or) covering the soil, which also helps prevent nutrients and other contaminants from running off agricultural lands.

Crop cover percentages for the period 2007–20 were computed from the U.S. Department of Agriculture (2021) Cropland Data Layers for the intervening land areas draining to the study gages. In most years during 2007–20, soybeans, followed by corn, constituted the largest crop percentages in the intervening drainage areas to the study gages (figs. 27–29). Other agricultural land covers constituted much smaller percentages except for in the intervening drainage area to the gage near Centerton, where grass/pasture percentages nearly equaled corn and soybeans until 2013, and steadily decreased afterward (fig. 29). Grass/pasture percentages steadily decreased over the entire 2007–20 period in the intervening drainage areas to the gages at Muncie and near Nora. The lines representing corn and soybean percentages intertwined over time in almost a mirror-image fashion (figs. 27–29), reflecting the influence of corn-soybean crop rotation. Central Indiana received nearly 19 inches of precipitation in the spring 2019, exceeding normal precipitation by almost 6 inches (Hill, 2019), which resulted in larger areas than typical left fallow/idle in 2019 (figs. 27–29).

Operational Tillage Information System data (Conservation Technology Information Center, 2021) on tillage practices and winter crop cover were obtained for the upper White River Basin. Figures 30–32 show percentages of crops within the upper White River Basin that were no-tilled (plant residue cover of more than 30 percent for all crops other than corn, for which the residue cover is more than 50 percent), reduced tilled (plant residue cover of 16–30 percent for all crops other than corn, for which there is a reduced tillage category with a larger residue of 31–50 percent) and conventionally tilled (plant residue cover 0–15 percent), respectively, prior to planting the next crop. There were considerable year-to-year variations in the percentages of previously planted crops in the three tillage categories and no discernable consistent upward or downward trends. In most instances, the annual percentages of each previously planted crops that were no-tilled were larger than the corresponding percentages conventionally tilled. Crops previously planted in corn or soybeans generally had larger percentages of reduced till than no-till treatment whereas the opposite was true for small grain crops.

The percentages of corn, soybeans, and small grain crops within the upper White River Basin for which a cover crop was planted during the following winter were plotted as a function of time (figs. 33–35). The percentage of winter cover crops that were planted following planting of corn varied irregularly over time between 2005 and 2018, and the percentage was typically small (totaling less than 1.5 percent) (fig. 33).

The percentage of winter cover crops planted following planting soybeans was larger than for corn, with winter commodity crops (for example, winter wheat) typically being the most common (fig. 34). Percentages of winter commodity cover crops planted following soybeans have been trending downward since about 2008, other cover crops (cereal rye for example) have been trending upward since about 2013, and percentages of perennial cover crops have remained small (typically less than 0.1 percent) and mostly unchanged over time (fig. 34).

The percentage of winter cover crops established following planting small grain crops was larger than for corn or soybeans, with cover crops other than winter commodity or perennial being the most common (fig. 35). The percentages of winter commodity and perennial cover crops remained mostly steady over the period 2005–18; however, there was an upward trend in the percentage of other winter cover crops planted. For example, the average percentage of other winter cover crops planted following small grain crops during 2005–11 was 6.8 but that average increased to about 18.7 percent during 2012–18. Small grain crops averaged only about 1 percent of the area of the upper White River Basin planted in row crops during 2012–18, so changes in winter cover crop practices are unlikely to have appreciably influenced nutrient and TSS transport at a hydrologic unit scale.

Fertilizer Application

Land-applied commercial fertilizers and manure might become a source of water pollution if they are used in excess, they are washed off the land surface, or they enter an aquifer before they are taken up by plants. Fertilizer application estimates were derived from county-level data reported by Falcone (2021a). Falcone (2021b) described the process by which the county-level estimates were made at 5-year increments between 1987 and 2017.

Nutrient mass estimates from commercial fertilizers and manure were summed for 13 of the 16 counties constituting the upper White River Basin. Brown, Monroe, and Owen counties were omitted from the sums because the upper White River Basin only includes small parts of those counties. Figure 36 shows the estimated mass of nutrients applied as farm and non-farm commercial fertilizers (includes commercial fertilizers applied to lawns and other non-farm urban and suburban green spaces) or manure as a function of time during the period 1987–2017. Farm-applied commercial fertilizer was the largest source of nitrogen and phosphorus. Manure was the next largest source followed by non-farm commercial fertilizer. The estimated masses of farm-applied nitrogen and phosphorus from commercial fertilizers were larger than from manure (the second largest source of applied nitrogen and phosphorus) on an annual basis by factors of at least 5.4 and 3.2, respectively.

The estimated masses of nitrogen and phosphorus applied as manure seem to have gradually trended downward between 1987 and 2017, whereas nitrogen and phosphorus masses from non-farm commercial fertilizers gradually trended upward (fig. 36). The remaining fertilizer sources did not show indication of temporal trend. The estimated mass of nitrogen from farm-applied commercial fertilizers did not change monotonically between 1987 and 2017; however, the most recent three estimates (for 2007, 2012, and 2017) were larger than estimates for earlier years except 1992.

Upgrades to Wastewater Treatment Systems and Delivery Processes

Nutrients and TSS are discharged in treated municipal wastewater and combined sewer overflows (CSOs). These discharges might constitute a substantial part of the nutrient and TSS flux in the receiving stream under some flow conditions. The mass of nutrients and TSS contributed depends on the level of treatment of the wastewater and the volume of wastewater discharged. Consequently, even though CSO discharge volumes are small compared to the volume of treated wastewater, the mass of nutrients and TSS discharged from CSOs are large because the CSO discharge contains a mixture of untreated sewage and stormwater. Improving wastewater treatment and eliminating or reducing CSOs are two ways of reducing the mass of nutrients and TSS delivered to the river.

The two largest population centers in the upper White River Basin (Indianapolis and Muncie, Indiana) provided lists of changes to their wastewater treatment systems or delivery processes between 1992 and 2020 that likely had an appreciable effect on nutrient and (or) TSS loadings to the White River. The earliest reported change to wastewater treatment in Indianapolis was in November 2002 and in Muncie in September 2016 (fig. 37). Changes identified as most impactful for reducing nutrient and (or) TSS loadings in Muncie occurred in 2018 when CSO with a discharge between 100 and 200 million gallons per year was separated so that untreated sewage was no longer discharged to the river (Rick Conrad, Director, Muncie Sanitary District Bureau of Water Quality, written commun., 2021). Changes identified as most impactful for reducing nutrient and (or) TSS loadings in Indianapolis occurred after 2011, including (1) the treatment improvements in December 2012 and April 2017 and (2) the completion of the Deep Rock and Eagle Creek tunnels in March 2018 (Olivia Hawbaker, Project Manager, Citizens Energy Group, written commun., 2021).