Simulated mussel mortality thresholds as a function of mussel biomass and nutrient loading

Mesocosm setup

Four 140 L, flow-through mesocosms (Fig. 2) continuously received untreated Iowa River water during the 107-d experiment, which culminated in an intensive 10-d water chemistry sampling period. Two mesocosms contained mussels collected from the Iowa River and two were without mussels (control). Twelve adult A. plicata and 13 adult L. cardium were placed in one mesocosm and 13 A. plicata and 12 L. cardium were placed in another mesocosm. This approximates a density of 70 mussels m⁻², which although high, is still a realistic density in some reaches of the UMR (Newton et al., 2011). Across both mesocosms, shell length (±1 standard deviation) was 95 ± 20 mm in A. plicata and 120 ± 25 mm in L. cardium. Initially, all mesocosms contained 8 cm of purchased sand substrate, but particulate deposition from the river water altered this composition over time. A gravity-fed, constant head system provided a controllable flow rate between 9 and 55 L h⁻¹. The flow rate during the 10-d intensive sampling period was 8.5 L h⁻¹ (16 h hydraulic residence time). Complete mixing in each mesocosm was provided by 1,500 L h⁻¹ submersible pumps, and two 1,000-watt solar simulators provided illumination on a 12:12 h light–dark cycle. Additional details regarding the mussel mesocosm system are available elsewhere (Bril et al., 2014).

Mesocosm sampling and analyses

Data from a 10-d intensive sampling period (days 97–107 of the 107-d experiment) were used for model evaluation. We intentionally delayed the start of the intensive sampling by 97 days so that the mussels could acclimate and bacteria responsible for nitrification and denitrification could establish. Electronic water chemistry sensors (model DS5; Hach Chemical Company, Loveland, CO, USA) were used to measure highly time-resolved (30-min) water chemistry data in the overlying water of each mesocosm and in the influent head tank. The sensors measured chlorophyll a (chl-a), NH${}_4^{+}$, NO${}_3^{-}$, pH, and temperature. Custom-made flow measurement devices with magnetic reed switches were used to quantify influent flow. Photosynthetically active radiation (PAR) sensors (model SQ-120; Apogee Instruments, Logan, Utah) were used to measure solar irradiance at the substrate and water surface of each mesocosm. All measurements obtained by the sensors were collected and stored using two data loggers. The model inputs for influent river temperature, food, NH${}_4^{+}$, NO${}_2^{-}$, NO${}_3^{-}$ and org N (Fig. 3) were measured values from within the river water head tank during the 10-d sampling period.

Discrete water chemistry samples were collected and analyzed at five time points during the 10-d sampling period from the overlying water and porewater of each mesocosm and from the influent head tank. The discrete samples were analyzed for chl-a, NH${}_4^{+}$, NO${}_2^{-}$, NO${}_3^{-}$, org N, and total N. Chl-a was measured by fluorescence. Measured chl-a concentrations (µg L⁻¹) were converted to “food” biomass (mg L⁻¹) based on literature values for phytoplankton chl-a content (Kasprzak et al., 2008). The fraction of nitrogen in food biomass (mg-N L⁻¹) was calculated using the empirical formula C₁₀₆H₂₆₃O₁₁₀N₁₆P (Chapra, 1997). NH${}_4^{+}$ was determined using the Salicylate Method, and NO${}_3^{-}$ was determined using the Dimethylphenol Method (APHA, 1996). NO${}_2^{-}$ was measured using the Diazotization Method, and total N was measured using the Persulfate Digestion Method (APHA, 1996). Sample measurements for org N were estimated by subtracting the sum of NH${}_4^{+}$, NO${}_3^{-}$, and NO${}_2^{-}$ from the total N measurements. A more detailed description of the mesocosm sampling and analysis setup is available (Bril et al., 2014).

Model calibration and sensitivity analyses

Seven days of the 107-d experiment were intensively sampled and previously reported (Bril et al., 2014) for food, NH${}_4^{+}$, NO${}_2^{-}$, NO${}_3^{-}$, org N, and temperature; these values were used as model calibration inputs. Linear interpolation between discrete samples was used where 30-min measurements were unavailable (org N, NO${}_2^{-}$, and total N), and ranges for unmeasured model variables (e.g., nitrification rate, denitrification rate) were obtained from the literature (Table 1). The model, created in Stella (version 8.0, ISEE Systems, Inc., Lebanon, New Hampshire), was initially calibrated using the no-mussel control data, then refined using data from mesocosms containing mussels to properly parameterize clearance and excretion rates (Bayne, Hawkins & Navarro, 1987; Englund & Heino, 1994; Haag, 2012). The optimized values used in the model calibration are given in Table 1. The optimized calibration values were determined by comparing model outputs to sensor and discrete sample measurements and then minimizing normalized mean error and maximizing R² values (Table 2).

Sensitivity analyses were conducted to identify the most important variables contributing to net system dynamic concentration response. A single variable sensitivity analysis was completed by adjusting the model variables based on a range of literature values (Table 1). When such information was unavailable, the value of the variable used in model calibration was adjusted by ±50%. Ten sensitivity model runs were completed for each variable using values obtained by sampling the range of literature values (or ±50% adjustments) at 10 equal intervals. The sensitivity analysis was considered for the normalized sensitivity coefficient (Fasham, Ducklow & McKelvie , 1990) (NSC): (1)${\displaystyle \mathrm{NSC}=\left(\frac{\frac{\phi }{{\phi }_0}}{\frac{P}{P_0}}\right)\frac{P_0}{{\phi }_0}}$ where, φ = mean value of a parameter (e.g., NH${}_4^{+}$, NO${}_3^{-}$) over the simulation period for the sensitivity run (mg-N L⁻¹), φ_o = mean value of a parameter over the simulation period for the calibrated model (mg-N L⁻¹), P = value of model variable in sensitivity run, and P_o = value of model variable in calibrated model. The NSC values for each sensitivity run were averaged to determine a net NSC for each model variable.

Mussel mortality threshold simulations

Based on 28-day laboratory toxicity tests with juvenile fat mucket mussels (Lampsilis siliquoidea), Wang et al. (2011) reported that 100% mortality occurred at 2.08 mg L⁻¹ total ammonia nitrogen (TAN). Given the pH (8.2) and temperature (20 °C) of that study, of the 2.08 mg L⁻¹TAN, about 1.9 mg-N L⁻¹ would be in the NH${}_4^{+}$ form and about 0.18 mg-N L⁻¹ would be in the NH₃ form. Given that our models were developed at a similar pH (8.2) and temperature (24 °C) to the Wang et al. (2011) study, we selected 2.0 mg-N L⁻¹ NH${}_4^{+}$ in porewater as a surrogate mortality threshold for Lampsilis mussels. Furthermore, the US Environmental Protection Agency (EPA) determined species mean chronic values of NH₃ for Lampsilis siliquoidea and L. fasciola to calculate a geometric mean chronic NH₃ value of 2.1 mg-N L⁻¹ for the genus Lampsilis (US Environmental Protection Agency, 2013).

The average measured porewater concentrations for NH${}_4^{+}$, NO${}_3^{-}$, NO${}_2^{-}$, org N, and food during the 10-d evaluation period (3.9, 0.2, 0.06, 5, and 0.1 mg-N L⁻¹, respectively) were used as initial conditions for porewater in the model. The average overlying water concentrations for the same variables were 0.05, 5, 0.05, 2.8, and 0.1 mg-N L⁻¹, respectively, and the “river water” inputs for 90-d model simulations were initially set to these values. The mussel density in our mesocosms was converted to estimated biomass (g m⁻²) using the allometric function, M = aL^b, where M is tissue dry mass (g) and L is length (mm) and with values for “a” and “b” for A. plicata taken from the literature (Newton et al., 2011). The resulting mass of 6.0 g mussel⁻¹ was multiplied by 35 mussels m⁻² (half the population) to determine an estimated biomass of 210 g m⁻² for A. plicata. In the absence of allometric data for L. cardium, the tissue dry mass was assumed to be 10 g mussel⁻¹(167% of A. plicata), and when multiplied by 35 mussels m⁻² resulted in a biomass of 350 g m⁻². Adding these values gave a maximum biomass of 560 g m⁻²,  which was used as the upper bound for the simulations. To simulate changes in porewater NH${}_4^{+}$ concentration as a function of mussel biomass and food availability, mussel biomass was varied at zero, 140, 280, 420 and 560 g m⁻² while the N content of food was varied at zero, 0.1, 1, 5 and 10 mg-N L⁻¹.

Model Evaluation

For the river water head tank (pH 8.2), a combination of sensor data (temperature, NO${}_3^{-}$, “food,” and NH${}_4^{+}$) and interpolated discrete data (org N and NO${}_2^{-}$) were collected and used as input to the numerical model on a 30 min time step (Fig. 3). For overlying water in mesocosms, the “food” sensor data were converted to a 25-d moving average (Fig. 4A) to condition the inherently noisy signal to enable visual comparison to the model output. The discrete sample results for NO${}_2^{-}$ concentrations in the overlying water were similar in magnitude, but did not agree closely with the model output (Fig. 4B). The model output for NH${}_4^{+}$ and NO${}_3^{-}$ concentrations (Figs. 4C and 4D) compared well with the sensor measurements. Overall, the model was capable of outputting results that accurately predicted the concentrations, and most of the dynamics, of the major N species at a 30 min time interval for the 10-d evaluation period.

The model was evaluated quantitatively using the standard deviation (SD) of the measured data variable compared to the root mean square error (RMSE) of the model output. If the RMSE was less than half the SD, the model output for that variable was deemed “accurate” (Singh et al., 2005; Moriasi et al., 2007). For comparative purposes, values for the mean bias, mean error, normalized mean bias, normalized mean error, and R² are reported along with the SD and RMSE for food, NH${}_4^{+}$, NO${}_2^{-}$, NO${}_3^{-},$ org N, and total N for the 7-d model calibration and 10-d evaluation periods (Table 2). The RMSE to SD ratio was ≤0.5 for the sensor-measured data for food, NH${}_4^{+}$ and NO${}_3^{-}$ for the 10-d evaluation period. The model evaluation based on discrete sample data yielded mixed results with RMSE to SD ratios of 0.55, 0.60, and 0.52 for NO${}_3^{-}$, NO${}_2^{-}$ and total N, respectively. The RMSE to SD ratios for food, NH${}_4^{+}$, and org N were 4.0, 1.4, and 0.86 for the discrete sample data, respectively. The lower accuracy determinations based on discrete sample data were likely a function of the small sample sizes, as compared to sensor measurements, and the low concentrations of food and NH${}_4^{+}$ which challenged the analytical limits of quantitation for these variables.

Sensitivity analysis

The modeled nitrogen species were collectively most sensitive to changes in temperature, hydraulic retention time, and mussel biomass (Table 3). Temperature was expected to be an influential variable since the majority of the first-order rate expressions are temperature dependent. Hydraulic retention time was also expected to be influential since the influent river water has a major impact on mesocosm water chemistry in a continuous-flow system. Mussel biomass was an unexpectedly sensitive model variable. However, given the influence of mussels on food, NH${}_4^{+}$, NO${}_2^{-}$, and NO${}_3^{-}$ concentrations shown in our previous work (Bril et al., 2014), this result, in hindsight, should have been anticipated.

Mussel mortality threshold simulations

At the highest simulated mussel biomass (555 g m⁻²) and the lowest simulated influent water food concentration (0.1 mg-N L⁻¹), the porewater NH${}_4^{+}$ concentration after a 2,160 h timespan in the absence of mussels, was 0.5 mg-N L⁻¹ compared to 2.3 mg-N L⁻¹ in the presence of mussels (Fig. 5). The food concentration in mesocosms without mussels was visibly higher than in mescocosms with mussels while NH${}_4^{+}$ and NO${}_2^{-}$ concentrations in overlying water were lower in the absence of mussels. Mortality threshold contours were estimated by varying mussel biomass and N concentration in the model (Fig. 6). Even when the simulated overlying water food availability was low, the mortality threshold was reached at a mussel biomass of about 480 g m⁻². At a food concentration of 10 mg-N L⁻¹ the mortality threshold was reached at a biomass of about 250 g m⁻².

In eastern Iowa, the median total N concentration in rivers and streams is commonly >10 mg-N L⁻¹ (Kalkhoff et al., 2000), which can place juvenile freshwater mussels at particular risk to ammonia toxicity. Minnesota has a draft criterion for aquatic life of 4.9 mg-N L⁻¹ total N, which was exceeded in 68% of samples collected in a study of Iowa waters between 2004 and 2008 (Garrett, 2012). The US EPA national recommended final acute ambient water quality criterion (AWQC) for protecting freshwater organisms from potential effects of ammonia is 17 mg-N L⁻¹ and the final chronic AWQC for ammonia is 1.9 mg-N L⁻¹ at pH 7.0 and 20 °C (US Environmental Protection Agency, 2013). At a total N concentration of 10 mg L⁻¹, our model predicts the mortality threshold to be reached when mussel biomass is about 400 g m⁻². However, the maximum total N concentration measured between 2004 and 2008 was 37.8 mg-N L⁻¹ (Garrett, 2012). Our model suggests the mortality threshold for juvenile Lampsilis could be exceeded at low mussel biomass if even a short exposure occurs at such a high total N concentration.

Reflecting on the relationships between nutrients and freshwater mussels conceptualized by Strayer (2014), we concur that high nutrient loads (particularly N in the agricultural Midwest) are a threat to the well-being of mussels. Conversely, our model predicts a somewhat linear mortality threshold relationship as mussel biomass and total N are varied, whereas Strayer stated this relationship would probably be non-linear and non-monotonic. In agreement with Strayer, our model suggests that mussels spatially focus nutrients from the overlying water to the sediments as evidenced by elevated porewater NH${}_4^{+}$ in mescosms with mussels. However, our previous work (Bril et al., 2014), and the model developed here, show elevated concentrations of NO${}_2^{-}$ and NO${}_3^{-}$ in overlying waters as an indirect consequence of mussel activity. This still represents a spatial focusing of nutrients by mussels, but the impact is not seen in the sediment alone.