Long Term (1998–2019) Changes in Water Quality Parameters as a Function of Freshwater Inflow in a River–Bay Continuum

Highlights

• Water quality parameters were studied in the Delaware Inland Bays watersheds.
• Freshwater inflow (FWI) had the greatest effect on dissolved N in the Inland Bays.
• Dissolved P depended on the combined effects of FWI and metabolic processes.
• Dissolved N and P were higher than the respective standard during the growing seasons of submerged aquatic vegetation.
• Coastal developments and changing land use had effects on load transport.

Abstract

Freshwater inflow is important in transporting nutrients to a bay. We hypothesized that freshwater inflow was transporting dissolved nitrogen and phosphorus to the Inland Bays. We analyzed long term (1998–2019) water quality data collected from Indian River, Indian River Bay, Lewes-Rehoboth Canal, Little Assawoman Bay, and Rehoboth Bay watersheds. Freshwater inflow altered nitrite+nitrate (N-NO_{2_3}) concentrations in all but Lewes-Rehoboth Canal watershed, whereas phosphate (P-PO₄) concentrations in all watersheds were altered by freshwater inflow and metabolic processes in the water. The average N-NO_{2_3} and P-PO₄ were higher than the standard (0.14 and 0.01 mg/L for N-NO_{2_3}+N-NH₃ and P-PO₄, respectively) for growing seasons (March–October) i.e., 0.83 + 0.14 and 0.09 mg/L in Indian River; 0.79 + 0.10 and 0.06 mg/L in Indian River Bay; 0.21 + 0.15 and 0.09 mg/L in Lewes-Rehoboth Canal; 0.49 + 0.10 and 0.11 mg/L in Little Assawoman Bay; 1.0 + 0.08 and 0.06 mg/L in Rehoboth Bay. Average total suspended solids in the Indian River (33), Indian River Bay (22), and Lewes-Rehoboth Canal (31) were higher than the standard concentrations, i.e., 20 mg/L for the Inland Bays. With the evidence of higher dissolved nutrients and low dissolved oxygen concentrations, need for nutrient load reduction and water quality monitoring are paramount for the sustainable management of Inland Bays.

1. Introduction

In recent decades, humans have changed land use patterns resulting in more storm and wastewater discharges [1]. Increased development and higher intensity farming practices add nutrient sources in a watershed, these nutrients when drained to nearby water affect water quality condition. Furthermore, overfertilization in watersheds, to fulfill human’s demand, affects nutrient transports to the nearby water and promotes algal bloom [1,2]. Watershed transport of chemicals and particulates depend on changes in amount of flow and land usage. Streams with periodic high discharge showed increasing transport of dissolved and particulate forms of inorganic and organic nutrients [3] The high pulses of peak surface inflow that transport nutrients can affect primary productivity in adjoining coastal waters [4]. Additionally, turbidity pulses in coastal areas have negative impacts on aquatic ecosystems, such as coral reefs, seagrass meadows, and shellfish beds [5,6,7]. In one of the Northeastern US studies, urban sprawl was related to the low dissolved oxygen concentration in Hudson River and Raritan Bay due to the transport of higher organic load [8]. Higher nutrient loads from urban and agricultural runoff have negative impacts on abundance of submerged aquatic vegetation in Chesapeake Bay sub estuaries [9].

Delaware’s Inland Bays consist of Rehoboth, Indian River, and Little Assawoman Bays (Figure 1). Watersheds that drain to the Delaware Inland Bays are no exception to effects of land use changes, experiencing eutrophication due to excessive nutrient loadings [10,11,12]. Unsustainable agricultural practices that use excess amounts of nutrients were the main reason for eutrophication in Delaware Inland Bays [11,13]. Although the application of best agricultural management practices has started in recent years, there has been an increase in coastal developments and impervious surfaces over the last several decades, potentially leading to higher nutrient inputs to Delaware Inland Bays. Another study at a sub-watershed scale in Delaware identified the area as a significant source of nitrogen to the Inland Bays throughout the year, with seasonal variation from different point sources [14]. Wastewater was the dominant source of phosphorus to the Rehoboth Bay, DE, before 2002; however, since 2002, phosphorus loadings from wastewater were significantly reduced due to the technical improvements in wastewater plants [15].

Delaware was ranked second nationally in agricultural sales per farm, with broiler production in the state accounting for 75% of the value of agricultural production [16] (https://www.usda.gov/media/blog/2019/06/21/delaware-small-state-big-agriculture, accessed on 31 August 2021 at 11.42 a.m.). Sussex County, where the Inland Bays are located, was identified as the largest broiler producer in the country (USDA, 2017 [16]). With such domination in agricultural production, wash off fertilizers and organic materials affect nutrient concentrations. Inland Bays’ waters were categorized as polluted water in terms of nutrients and classified as dissolved oxygen impacted water [17]. One of the long-term visions for States, under Clean Water Act section 303(d), is to prioritize watersheds for restoration and protection. To do so, more scientific studies on pollutants transport and management plans are needed from each watershed. The Delaware Inland Bays Comprehensive Conservation and Management Plan (CCMP) was established to improve the environmental conditions within the Bays and their watersheds [18]. For the protection of Inland Bays, the CCMP recommended to oversee water quality changes and pollution load reduction over time. The total maximum daily load (TMDL) allocation put forward in 1998 was one of the efforts to continue comprehensive conservation management of Inland Bays. This study’s efforts to understand changes in water quality parameters is one of the CCMP recommendations to monitor water quality condition in the Bays. We used long-term water quality data to identify watershed transport and to identify the multivariate relationship between water quality parameters in the adjoining Inland Bays. We hypothesized that nitrogen and phosphorus transport to the Inland Bays are impacted by the inflow.

2. Site Descriptions

Delaware’s Inland Bays lie in Sussex County (Figure 1). Most of the headwater drainage areas of the bays are agricultural land, even though coastal developments have been higher in recent decades. Sussex County’s population has increased by nearly 19% from April 2010 to July 2019, compared to about a 9% increase in the Delaware State (https://www.census.gov/quickfacts/fact/table/sussexcountydelaware, DE/PST045219, accessed on 2 July 2021 at 11:28 a.m.). The Inland Bays are bar-built estuaries and separated from the Atlantic Ocean by barrier islands or sandbars. Typically, these well-mixed estuaries are built by river inflow or oceanic current and usually have low water volume throughout the year. The average water levels in these Inland Bays are less than 2 meters [19]. Chemicals and materials input to the Inland Bays are controlled by freshwater inflow and tidal input.

The Inland Bays are composed of Rehoboth Bay (RB) and Indian River Bay (IRB), and Little Assawoman Bay (LAB). Rehoboth Bay is the northern most among the three bays and is connected to Delaware Bay by Lewes-Rehoboth Canal (L-R Canal), to the south Rehoboth Bay is connected to Indian River Bay. Most of the freshwater inflows to the Rehoboth Bay are from small creeks and overland flows. The Indian River (IR) drains about 95% of the freshwater inflow to Indian River Bay [20,21] and is the major inlet for the Delaware Inland Bays. In addition, Indian River Bay also gets freshwater from small creeks and overland flows. Indian River Bay is connected to the Atlantic Ocean by Indian River Inlet. Little Assawoman Bay is the southern most of the three Inland Bays. It receives freshwater inflow by small creeks, i.e., Dirickson Creek and Miller Creek. Assawoman Canal connects Indian River Bay water with the Little Assawoman Bay. Similar to Indian River Bay and Rehoboth Bay, Little Assawoman Bay is separated from Atlantic Ocean by a barrier island and is connected to the Assawoman Bay in Maryland, USA.

3. Methods

3.1. Data

The present study analyzed long term (1998–2019) water quality data collected at sites located in the Inland Bays watersheds. These monthly or bi-monthly water quality data were collected by Delaware Department of Natural Resources and Environmental Control (DNREC) and stored in the US EPA’s storage and retrieval Data Warehouse (STORET) (https://www.waterqualitydata.us/portal/, accessed on 10 July 2020). Delaware has monitored surface waters since the 1950s. Sampling sites, procedures, frequency, and analyses of various water quality parameters are part of the surface water monitoring program. USGS and EPA manuals for the sample collection and analyses follow. Eleven different water quality parameters, i.e., nitrite+nitrate (N-NO_{2_3}; measured in mg/L), ammonia (N-NH₃; measured in mg/L), phosphate (P-PO₄; measured in mg/L), total suspended solids (TSS; measured in mg/L), salinity (Sal; measured in parts per thousand), pH, water temperature (Temp_water; measured in °C), dissolved oxygen (DO; measured in mg/L), chlorophyll_a (chl_a; measure in µg/L), chloride (Cl; measured in mg/L), and Secchi depth (Secchi; measured in inches), were used for this study. For the analysis, dissolved nutrients and salinity were used as indicators to changes in environmental flow to the Inland Bays, while DO and water temperature were used as indicators for changes in metabolic processes in the Bays. Use of chloride concentration was solely to add one extra component to understand salinity changes, as chloride is one of the most abundant ions in saline water and alteration of which also indicates changes in salinity.

DNREC has been sampling tidal and non-tidal stations in the Indian River, Indian River Bay, Rehoboth Bay, Little Assawoman Bay, and Lewes-Rehoboth Canal watersheds. Stations that have a tidal effect, which was identified by changes in salinity, were chosen for this study. Stations in Indian River watersheds are varied from fresh to brackish water, and to eliminate data bias we used salinity as an indicator. Stations with historic salinity data of greater than 1 ppt were included in the analysis. Stations in the four watersheds except Lewes-Rehoboth Canal are located along horizontal salinity gradients. For the analysis, we assumed each watershed as an independent system and changes in each watershed transport has effects on each sampling locations within that watershed or system. In total, 14 stations were chosen from Delaware’s stream monitoring network, in addition discharge data was downloaded from USGS 01484525 station at Millsboro Pond, and continuous DO data was downloaded from USGS stations at Millsboro Pond (USGS 01484525) and Massey Ditch (USGS 01484680) (Figure 1).

3.2. Statistical Analysis

Normality distribution of the multivariate data is one of the important assumptions underlying multivariate analysis. In the present study, except pH, as it is on log scale, the variables used for the analysis were log transformed to meet normal distribution of residuals. For the principal component (PCA), analysis of variance (ANOVA), and regression analyses all the log transformed data were standardized to mean of 0 and standard deviation of 1, using the PROC STANDARD module in SAS software [22,23,24]. The long-term monthly water quality dataset used for the analysis had some missing values, which were replaced by the average by the standardized procedure mentioned above. The multi-variable correlation coefficient (r) and p-value were evaluated to identify the significance of correlation between variables.

Principal component analysis (PCA) is a variable reduction technique and used to reduce a large number of variables to a small set of loads by preserving most of the variance in the original data set. The PCA technique assumes that total variance of the variable is the common variance. PCA was used on the standardized data set. The PROC FACTOR module on the correlation matrix, with Varimax rotation in SAS software, was used to perform PCA analysis.

Analysis of variance (ANOVA) was performed using the PROC GLM procedure with the TUKEY option to identify variability in water quality parameters among the five watersheds.

The PROC REG procedure was used to perform linear regression analysis. Linear regression analysis was performed to identify dissolved nutrients’ dependencies on other water quality parameters.

3.3. Land Use Changes

Land use data was retrieved from Delaware FirstMap data portal (https://firstmap.delaware.gov/, accessed on 1 August 2021) for the years 1997, 2007, and 2017 and analyzed in Esri ArcGIS Pro version 2.8. The intersect analysis tool was used to isolate land use data specifically within the Inland Bays Basin and acreage was calculated based on these data using the calculate geometry attributes data management tool. Land use categories were reclassified from the state’s modified Anderson Classification System (https://www.arcgis.com/home/item.html?id=4c21a2b79352453a9a8446195302dea7, accessed on 1 August 2021) to agriculture, developed, and natural in order to simplify the dataset and highlight the overall land use changes over time. Data was summarized using the summary statistics analysis tool based on the simplified land use categories and tabular data was exported and used to plot change over time.

4. Results

Flow measured at USGS Station, Millsboro Pond and salinity at tidal stations are negatively related. Salinity variability in the Inland Bays tidal stations is controlled by freshwater inflow sources and is statistically significant (p < 0.0001 between flow vs. yearly average salinity in Inland bays). Dissolved nitrite+nitrate (N-NO_{2_3}) in the five watersheds is significantly correlated with flow (

Table 1. The first part of the table depicted yearly average discharge measured at USGS Station Millsboro Pond—01484525 in. Abbreviation: Q = discharge (m³/s); short form of year shown. The second part of the table depicted correlation coefficient (r) between flow measured at USGS station, Millsboro Pond, dissolved nitrite+nitrate (N-NO_{2_3}), and dissolved phosphate as phosphorus (P-PO₄) in the Inland Bays watersheds. Correlation coefficient (r) and p-value (p) are presented.

Year	‘98	‘99	‘00	‘01	‘02	‘03	‘04	‘05	‘06	‘07	‘08	‘09	‘10	‘11	‘12	‘13	‘14	‘15	‘16	‘17	‘18	‘19
Q	3.5	1.9	2.6	2.3	1.9	5.2	2.9	3.0	2.8	2.0	1.6	3.8	2.9	1.3	1.6	3.7	2.5	2.2	2.9	2.3	3.9	2.6
Watershed Basins										Flow N-NO2_3 (r; p)							Flow P-PO4 (r; p)
Indian River										0.49; 0.005							0.12; 0.54
Indian River Bay										0.39; 0.01							0.06; 0.76
Lewes-Rehoboth Canal										−0.01; 0.002							−0.06; 0.95
Little Assawoman										0.50; 0.0002							−0.11; 0.91
Rehoboth Bay										0.46; 0.001							0.16; 0.41

), whereas dissolved phosphorus (P-PO₄) is not correlated with flow (Table 1). Salinity variability was clearly visible in the five watersheds with the lowest average in Rehoboth Bay (i.e., 17) and the highest average in Lewes-Rehoboth Canal (i.e., 24.7) (

Table 2. Summary statistics for all five watersheds. Units of measurement were inches for Secchi depth; ppt for salinity; degree C for temperature; mg/L for dissolved oxygen, nitrite+nitrate, ammonia-N, phosphate-P, total suspended solids, chloride; and µg/L for Chlorophyll-a.

	Indian River		Indian River Bay		Lewes-Rehoboth Canal		Little Assawoman		Rehoboth Bay
Variable	Mean	Std Dev	Mean	Std Dev	Mean	Std Dev	Mean	Std Dev	Mean	Std Dev
Secchi Depth	4.1	10.1	5.6	16.2	5.1	13.2	3.4	8.7	5.8	14.8
Chlorophyll-a	59.2	130.3	10.3	11.6	11.9	8.4	31.5	70.9	22.3	48.1
Chloride	9819	3522	12594	6949	14182	2346	9759	5496	9398	6657
Salinity	17.7	5.5	23.1	10.5	24.7	3.5	17.8	8.8	17.18	11.4
Water Temperature	19.0	8.5	15.1	7.6	15.4	8.0	15.9	8.5	16.6	7.9
Dissolved Oxygen	8.1	2.2	8.1	2.3	7.1	2.6	8.1	2.5	8.2	2.1
Nitrite+Nitrate	0.83	0.80	0.79	1.68	0.21	0.44	0.49	0.92	1.0	1.1
Ammonia-N	0.14	0.15	0.10	0.13	0.15	0.17	0.10	0.10	0.08	0.1
Phosphate-P	0.11	0.10	0.06	0.06	0.09	0.05	0.11	0.14	0.07	0.06
TSS	33.3	27.5	22.5	31.5	31.3	27.5	20.2	18.9	20.1	19.2
pH	7.6	0.4	7.6	0.5	7.5	0.35	7.9	0.3	7.08	0.7

). Salinity variability in each of the five watersheds is depicted by horizontal error bars (Figure 2). The average N-NO_{2_3} ranged from 0.21 to 1.0 mg/L, the lowest in Lewes-Rehoboth Canal and the highest in Rehoboth Bay (Figure 2A). Average chlorophyll-a (chl-a) was the highest in Indian River (59 µg/L) and the lowest in Indian River Bay (10.3 µg/L) (Figure S1 and Table 2). Chl-a and N-NO_{2_3} in the five watersheds were inversely related with salinity, whereas no such relationship was identified with dissolved ammonia (N-NH₃) and P-PO₄ (Figure S2 and Figure 2B). Average dissolved ammonia measured as N (N-NH₃) was the highest (0.15 mg/L) in Lewes-Rehoboth Canal (Figure S2 and Table 2). N-NH₃ distribution along the salinity gradient of five watersheds salinity was similar to total suspended solids (TSS) with the lowest average concentration in Rehoboth (Figures S2 and S3). Average dissolved oxygen (DO) in all five watersheds was higher than 7 mg/L (Figure S4), though low DO (<4 mg/L) was identified in USGS continuous monitoring sites in Indian River and Indian River Bay (Figures S5 and S6).

The first principal component (PC1) explained 38.5%, 48.6%, 48.6%, and 35.1% of variability, respectively, in the Indian River, Indian River Bay, Rehoboth Bay, and Little Assawoman Bay data set (Figure 3A–D). N-NO_{2_3} and salinity are inversely correlated that indicates that PC1 represents freshwater inflow. It is because salinity decreases with the increase in freshwater inflow which brings the dissolved form of nitrogen from a watershed. The second principal component (PC2) explained 30.9% and 25.81% in Indian River and Rehoboth Bay, respectively, while the third principal component (PC3) explained 25.68%, and 31.56% of variability, respectively, in the Indian River Bay and Little Assawoman Bay data set (Figure 3A–D). The PC2 (in Indian River and Rehoboth Bay) and PC3 (in Indian River Bay and Little Assawoman Bay) had DO and temperature inversely correlated that indicates PC2 and PC3 in respective watersheds as metabolic processes. With the increase in temperature biota will become active; thus, consuming oxygen resulted in low DO. In all watersheds, P-PO₄ had a negative value in PC2 that indicates possible release of P during higher temperatures.

The PCA analysis with the Lewes-Rehoboth Canal data set was different to the other four watersheds. The PC1 explained 42.57% of variation in the Lewes-Rehoboth Canal data set and it had DO and temperature inversely correlated (Figure 3E). The low PC1 scores indicates high freshwater inflow during respective years (Figure 4). In four of the five watersheds, freshwater inflow (represented by PC1) had the greatest effect on dissolved nitrogen. With the exclusion of Lewes-Rehoboth Canal dissolved nutrient variables score data (as we do not see effects of freshwater inflow in the Canal), in most of the years Indian River had the lowest PC1 scores, while mostly the data from Indian River Bay and Little Assawoman had the highest PC1 scores (Figure 4). Rehoboth Bay scores are mostly in the middle of the plot (Figure 4). Theis distribution on PC1 (which indicates freshwater index) indicates freshwater inflow has the most impact on Indian River Bay and Little Assawoman Bay; however, the distribution also indicates all four watersheds’ dissolved nitrogen are affected by freshwater inflow.

Water quality parameters among the five watersheds were tested to identify the differences in watersheds. DO, chlorophyll-a, pH, water temperature, TSS, salinity, dissolved phosphorus (P-PO₄), and dissolved nitrogen (nitrite+nitrate; N-NO_{2_3} plus dissolved ammonia—N-NH₃) were significantly different among five watersheds (

Table 3. Result of the ANOVA summarizing p values for null hypothesis and variance component analysis explained by watershed as the predictor variable. Abbreviation: Chl_a = chlorophyll a; Cl = chloride; Sal = salinity; TempWater = water temperature; DO = dissolved oxygen; N-NO_{2_3} = nitrite+nitrate; P-PO₄ = dissolved phosphate as P; Secchi = Secchi depth; N-NH₃ = dissolved ammonia as N; TSS = total suspended solids.

Variable	DF	p-Value	Variance Component (%)
Secchi	4	<0.001	88
Chl_a	4	<0.001	44
Cl	4	<0.001	34
Sal	4	<0.001	33
TempWater	4	<0.001	21
DO	4	<0.001	12
N-NO2_3	4	<0.001	95
N-NH3	4	<0.001	95
P-PO4	4	<0.001	87
TSS	4	<0.001	28
pH	4	<0.001	3

). Dissolved nitrogen had the highest percentage (95%) of variance explained by watershed, while dissolved phosphorus had 87% of variance explained by watersheds (Table 3). Mean differences between log transformed variables were identified using Tukey test. Nitrite+nitrate concentration was the highest in Rehoboth Bay while P-PO4 was the highest in the Little Assawoman Bay (Table S1). Total suspended solid concentration was the highest in Indian River (33.3 mg/L) and Lewes-Rehoboth Canal (31.3 mg/L). Dissolved oxygen was not significantly different between Indian River and Rehoboth Bay and was comparatively higher than the other three watersheds (i.e., Lewes-Rehoboth Canal, Little Assawoman, and Indian River Bay) (Table S1). Except in Indian River, chloride concentration differences between watersheds were similar to salinity concentration (Table S1).

Suspended materials, especially TSS and chlorophyll-a (with higher than 25% variability), were the significant predictor variables for P-PO₄ in all five watersheds (Supplementary Material—Equations (S₁–S₉)). TSS and chlorophyll-a were negatively correlated to N-NO_{2_3} in most of the watersheds (

Table 4. Correlation between dissolved phosphorus (P-PO₄) and nitrite+nitrate (N-NO_{2_3}) with chlorophyll-a (chl-a) and total suspended solids. All values have p < 0.05.

Watersheds	Chl_a and P-PO4	TSS and P-PO4	Chl_a and N-NO2_3	TSS and N-NO2_3
Indian River	0.61	0.49	−0.22	−0.28
Indian River Bay	0.30	0.27	−0.17	−0.14
Rehoboth Bay	0.56	0.38	−0.18	−0.18
Little Assawoman Bay	0.59	0.36	−0.07	−0.05
Lewes-Rehoboth Canal	0.08	0.50	0.05	0.04

). Lewes-Rehoboth Canal P-PO₄ concentration was significantly positively correlated with TSS, while in all other watersheds P-PO₄ was significantly positively correlated with TSS and chlorophyll-a (Table 4).

Notable changes occurred in land use during the 20 year period between 1997 and 2017 (Figure S7). Agriculture and natural areas decreased, although the change in natural land use was much less pronounced. The increase in developed lands was mostly seen in coordination with the loss of agricultural land. Since 1997, there was 17% loss in agricultural land in the Inland Bays watersheds, whereas developed land increased by 46%. Most of the increase in developed land occurred between 1997 to 2007, with a 26% increase and a corresponding 12% loss in agricultural land.

5. Discussion

This study encompassed periodic changes in dissolved nutrients, dissolved oxygen, salinity, total suspended solids, and algal growth (proxy measurement by chlorophyll-a concentration) in the Inland Bays. Despite the predominant land use being agriculture, the Inland Bays watersheds are not exempt from the pressures of development. Overland flow draining agricultural land transported nitrogen and sediment loads to the receiving bays and estuaries. These conditions are exacerbated by increased coastal development [8], which adds more nutrient loads to the system. Overland flow was the primary source of dissolved nitrogen in the Inland Bays, while dissolved phosphorus concentration was dependent on flow and metabolic processes in the Bays. The transported loads affect dissolved oxygen concentration with rapid algal bloom [25,26]. Average dissolved oxygen was above 7 mg/L in all bays; however, there were periods of low DO events that could have had detrimental effects on biota. Locals and citizen water quality monitoring groups have reported evidence of fish kill in the Inland Bays. Notably, the USGS Millsboro Pond outlet station had periods of low DO during the summer (lower than 4 mg/L) and harmful algal blooms were identified during those periods.

Freshwater inflow has an important role in altering salinity regimes [27,28]. The salinity of bays in four of the five watersheds responded with the hydrologic flow, except Lewes-Rehoboth Canal which received very low flow volume. The four peak flow events lowered salinity of the Inland Bays in all watersheds and had potential of fluctuating TSS. Water quality monitoring performed by the State of Delaware, university researchers, and citizen monitoring during 1998 identified that waters of Indian River, Indian River Bay, and Rehoboth Bay were highly enriched with nitrogen and phosphorus that resulted in eutrophication [16]. The TMDL developed during 1998 put thresholds on 85% reduction in nitrogen and 65% reduction in phosphorus loads from the major tributaries lie in the higher amount of reduction requirement areas. The required reduction levels for the lower reduction areas were 40% reduction of nitrogen load and 40% reduction of phosphorus loads draining into Inland Bays. Furthermore, in the tidal portion of Inland Bays, the applicable standards during the submerged aquatic vegetation (SAV) growing season (March to October) were established at 0.14 mg/L or below for dissolved nitrogen (i.e., dissolved nitrite+nitrate+ammonia), 0.01 mg/L or below for dissolved phosphorus (DIP), 5 mg/L or above daily average for DO, and 20 mg/L or below for TSS. In the present study, we identified that the average dissolved nitrogen in all watersheds was higher than the applicable standard (Figure 2A, Table 2, and Figure S2), similarly average dissolved phosphorus in all five watersheds was higher than the applicable standard (Figure 2B). The average dissolved nitrogen in the Rehoboth Bay, Indian River, Indian River Bay, Little Assawoman, and Lewes-Rehoboth Canal were approximately 8, 7, 6, 4, and 2 times higher than the standard for SAV growing seasons. Similarly, average DIP was approximately 11, 11, 9, 7, and 6 times higher than the applicable standard in the Little Assawoman, Indian River, Lewes-Rehoboth Canal, Rehoboth Bay, and Indian River Bay, respectively. The average TSS in the Lewes-Rehoboth Canal (31 mg/L), Indian River (33 mg/L), and Indian River Bay (22/L) were higher than the applicable standard (Table 2 and Figure S3). The higher TSS concentration in the Lewes-Rehoboth Canal was attributed to the wind-driven resuspended sediments, whereas combining effects of inflow and wind were influential in sediment resuspension in the Indian River and Indian River Bay.

Comparatively, the Indian River had higher DIN and DIP concentrations than other watersheds, and thus it was reasonable to observe a higher average chlorophyll-a concentration (59 µg/L) (Table 2 and Figure S1). Relatively, lower average chlorophyll-a (12 µg/L) in Lewes-Rehoboth Canal was associated with higher suspended sediment blocking sunlight and low DIN concentration, while lower average chlorophyll-a in Indian River Bay (10 µg/L) could be due to the combined effect of low DIP and high TSS. Interestingly, even though chlorophyll-a was variable, average dissolved oxygen in all five watersheds was above 7 mg/L (long-term DO obtained from monthly discrete grab samplings) and meeting applicable standard described in the 1998 DNREC TMDL report. However, there were cases with dissolved oxygen in USGS continuous (15 min interval) monitoring sites at Millsboro Pond, Indian River and Massey Ditch, Massey Landing (see Figure 1 for site location) below 5 mg/L, with few monthly cycles of less than 2 mg/L (Figures S5 and S6). Even though there were fewer (since 2007 in Millsboro Pond and since 2011 in Massey Ditch) continuous DO data when comparing to the long term (since 1998) data, the frequent low DO concentrations (i.e., <2 mg/L) indicate consistent nutrient pollution, especially in the Indian River. The low DO in the Indian River could be the reason for a higher number of fish killed in the region, especially during summer season (identified by photos and personal communication with citizen monitoring groups). This prompts the need for continuous DO monitoring in the Indian River and Indian River Bay.

Watershed variability for dissolved nitrogen, i.e., N-NO_{2_3} (95%) and N-NH₃ (95%), was the highest compared to other water quality parameters. This indicates that dissolved nitrogen varied the most among five watersheds. One of the main reasons for N-NO_{2_3} and N-NH₃ variability was due to the variable amount of inflow draining into the respective five systems. River inflow was significantly positively correlated with N-NO_{2_3} concentrations in four of the five watersheds, except Lewes-Rehoboth Canal—where no relationship was observed, and freshwater inflow had positive impacts on changing N-NO_{2_3} concentrations in Inland Bays (Figure 3). Freshwater inflow had the greatest effect on fluctuating N-NO_{2_3} concentration in Indian River, also comparatively, freshwater inflow supplied the least amount of dissolve nitrogen to the Indian River Bay and Little Assawoman Bay, as Indian River drains about 95% of freshwater (Figure 4). Comparatively, Rehoboth Bay had high and low effects of freshwater inflow, we argue that it could be because of the lack of sampling stations along the salinity gradients.

In contrast to nitrite and nitrate, transport of dissolved P-PO₄ concentration was comparably less dependent on freshwater inflow in all five watersheds (Table 1). A study conducted during 1998–2002 identified soluble reactive and dissolved organic phosphorus transported from the watershed of Rehoboth Bay [13,14]. Since 2002, with the technical improvement to wastewater plants, the direct discharge of phosphorus load was reduced to the Rehoboth Bay and other Inland Bays [14]. These improvements to wastewater plants in the Inland Bays reduced dissolved phosphorus. The productivity need of dissolved phosphorus for the faunal community was compensated by metabolic processes. Furthermore, in all five watersheds, P-PO₄ was inversely related to dissolved oxygen and positively correlated with water temperature, which indicates release of P-PO₄ during warm temperature due to bacterial activities (Figure 3). Transported agricultural soil was the major P source along the salinity gradient of Love Creek, Delaware (one of the tributaries to Rehoboth Bay) [29], which supports our findings that dissolved phosphorus in the Inland Bays was from the release during metabolic processes. Dissolved phosphorus was positively correlated with chlorophyll-a and TSS concentrations when suspended fractions were less than 50 mg/L [30] and Paudel, Barnegat Bay, New Jersey [unpublished data]. In the Inland Bays, with average TSS of less than 35 mg/L, P-PO₄ was positively correlated with chlorophyll-a and TSS. Labile P is the most bioavailable sediment-bound P [31,32,33], and its release from suspended fractions by microbial activities could be the reason for the positive correlation. A comprehensive biogeochemistry study is recommended to understand phosphorus release mechanism in the Inland Bays sediment.

The present study’s findings of freshwater inflow as the major source of dissolve nitrite+nitrate to Inland Bays is plausible as the watersheds have high poultry production (poultry production is about 570 million chickens in 2020 data, extracted from Delmarva Chicken Association, https://www.dcachicken.com/facts/facts-figures.cfm accessed on 31 August 2021) and with the evidence of homogenous atmospheric nitrate distribution [34] (i.e., precipitation sample analyzed between Lewes, DE and Indian River, DE sites). In the present study, we identified that the N-NH₃ concentration in Lewes-Rehoboth Canal was higher than the other four watersheds (Figure S2), which is associated with release due to the metabolic processes. The Canal is currently receiving discharge from the Lewes wastewater treatment plant and until 2018 it received discharge from the Rehoboth wastewater treatment plant. These discharges are seasonal and higher during summer, which may also have contributed to the higher ammonia concentration in the Lewes-Rehoboth Canal. Metabolic process was identified to account for 42.5% of variation in Lewes-Rehoboth Canal, and the positive correlation between N-NH₃ and water temperature indicates that the release of NH₃-N from the treatment plant was higher during warmer seasons (Figure 3E). The significantly higher airborne ammonium concentration, identified by Scudlark [34], in Indian River, DE compared to the Lewes, DE sites, also supports our conclusion of metabolic release of N-NH₃ in Lewes-Rehoboth Canal.

6. Conclusions

Watershed transport is the primary source of dissolved nitrogen and could be one of the major contributors of eutrophication in the Inland Bays. We recommend the study of sediment processes to understand phosphorus release in the Inland Bays. Evidence of very low dissolved oxygen and associated fish kill in the Indian River Bay demonstrate the need for continuous DO monitoring in the Inland Bays. A higher than applicable standard concentration for SAV growing seasons of dissolved nitrogen and phosphorus highlighted the need for load reduction for the sustainable management of healthy Inland Bays. Even though TMDL is in place for the nutrients’ load reduction in the Inland Bays, this study provides the need for a TMDL implementation strategy and better management of overland nutrient flow.