Introduction
Nitrogen (N) is essential for life but can be present in the environment in excess of growth requirements due to human activities. N is a common point-source pollutant to aquatic systems from wastewater treatment plants (WWTPs). Nitrate (
) and total ammonia nitrogen (TAN; where TAN includes both ammonia (NH
3) and ammonium (
)) are the two inorganic N forms that determine the critical loads beyond which aquatic ecosystems experience eutrophication or acidification (
Posch et al. 2001;
Schindler et al. 2006). The fate of these inorganic N species is a key determinant in the health of ecosystems and the services they provide to humans. TAN can be both a fertilizer of and detriment to aquatic life. At elevated concentrations, NH
3 is toxic to aquatic life (
Canadian Council of Ministers of the Environment 2010). Similarly, elevated concentrations of
degrade water quality by harming aquatic life (
Canadian Council of Ministers of the Environment 2012), and those above drinking water limits can lead to adverse health effects in people (
Iwanyshyn et al. 2008). Consequently, understanding the environmental fate of TAN and
discharged to surface waters is important for managing of human-disturbed aquatic ecosystems.
Many processes remove N from aquatic ecosystems. By understanding the relative contributions of each process and the factors that affect their rates, the environmental fate of N loading to aquatic ecosystems can be predicted (
Iwanyshyn et al. 2008). Successful nutrient mitigation strategies in larger aquatic ecosystems rely on using smaller, tractable ecosystems as realistic and replicatable systems (
Schindler 1998;
Dodds and Welch 2000;
Withers and Lord 2002;
Webster et al. 2003;
Sharpley et al. 2009). The concept of nutrient spiralling in streams was developed to describe the cycling and transport of nutrients in small lotic ecosystem (
Newbold et al. 1981,
1982,
1983) and is based on downstream changes in nutrient concentrations. Nutrient spiralling techniques have been improved by adding
15N-enriched compounds as a tracer and following them through different pools (e.g.,
Mulholland et al. 2000,
2004,
2008;
Tank et al. 2000;
Earl et al. 2006;
Hall et al. 2009). In a similar fashion, low nutrient streams can be spiked with nutrients and changes in the nutrient pulse can be used to understand ecosystem metabolism of nutrients (e.g.,
Davis and Minshall 1999;
Hall and Tank 2003). These studies are often restricted to short lengths of streams where the hydrology can be well characterized and to smaller systems in general. The understanding of nutrient spiralling in large impacted rivers is often confounded by a heterogeneous river morphology, frequent run-of-the-river dams, groundwater, and multiple nutrient inputs, and it consequently relies on the intensive work conducted in these smaller systems supplemented by sampling campaigns of both concentration and stable isotopes of N species. Further, observed values are a cumulative result of a plethora of contemporaneous N cycling processes with rates that change in relative importance with distance from inputs and time of day. Disentangling the relative rates of these processes in large rivers is greatly aided by the additional information supplied by stable isotopes and the development of numerical models (
Denk et al. 2017).
Novel technical developments in the analysis of stable isotopes have allowed for improved assessment of nitrogen cycling in rivers including the use of the differences in δ
15N-N
2O and δ
18O-N
2O produced by nitrification versus denitrification (
Thuss et al. 2014). Similarly, ecosystem metabolism techniques have recently been improved by the use of diel δ
18O-O
2 and δ
13C-DIC modelling (
Tobias et al. 2007;
Venkiteswaran et al. 2007;
Holtgrieve et al. 2010;
Parker et al. 2010) leading to a greater understanding of nutrient dynamics (
Fourqurean et al. 1997;
Fry et al. 2000;
Savage and Elmgren 2004;
Murray 2008). The isotopic labelling of benthic biofilm by differing
and
sources has recently been described (
Hood et al. 2014;
Loomer et al. 2014;
Peipoch et al. 2014). Here, we build on these studies by developing and testing a model that uses changes in concentrations and natural abundance (that is, not isotopically labelled) stable isotopic ratios to quantify the contributions of the various nitrogen-removal pathways in nutrient-impacted rivers. We applied this model to quantify the fate of N from the WWTP effluent discharges in a river highly impacted by both agricultural and WWTP nutrients.
The objectives of this research are to (i) quantify changes in concentrations and δ15N values of TAN and with distance downstream from WWTPs; (ii) develop a parsimonious process-based model for N cycling and the fate of WWTP N in rivers, and assess model performance with field measurements; and (iii) provide model-based estimates of the rates of nitrification, denitrification, NH3 volatilization, and N assimilation in WWTP plumes in a river impacted by both WWTP and agricultural nutrient inputs.
Discussion
The process-based NANNO model was able to reproduce the observed dynamics in concentrations and the δ
15N values of TAN and
(
Tables S3–
S6). Results from two seasons, with different proportional fates of N processing, at two different WWTPs with different TAN:
ratios in their effluent indicate a good degree of coherence between model results and field data (
Figs. 3–
6 and
Tables S3–
S6). Additionally, the shapes of the curves (increases, decreases, and plateaus) were all generally reproducible by the model. The model was least successful in reproducing behaviour when there were increases in
concentration without a change in δ
15N-
. This scenario suggests nitrification where the new
has the same δ
15N-
as the extant
.
Although process-based models provide several advantages over purely empirical ones, the goodness of fit and performance of such models depends on their particular application (
Saloranta et al. 2003). The variables in this underdetermined model system were constrained in several ways. The variables in this underdetermined model system were constrained in several ways, but the available data were insufficient to identify a unique best-fit model parameterization within the constrained parameter space. The local sensitivity of specific behaviours on the resulting modelling parameters can only be confidently assessed once best-fit model solutions have been developed (
Moore and Doherty 2005). This issue is typical of ecosystem models, which nevertheless are often used to predict the fate of nitrogen, oxygen, and biological oxygen demand in rivers (e.g., the Grand River Simulation Model that is used by the local Grand River Conservation Authority when making decisions about nutrient loads and river health).
In both Waterloo cases, denitrification played a modest role in reducing N concentrations (
Table 2). N
2O concentrations in and fluxes from the Grand River are high downstream of these WWTPs (
Rosamond et al. 2011,
2012;
Venkiteswaran et al. 2014). More detailed sampling of N
2O and its δ
15N values may provide additional constraints to improve the model fit.
Nitrification played a moderate role in N cycling in all four cases. There were no clear correlations between nitrification rates and rates of other N processes suggesting that predictions about the fate of N in the Grand River cannot be simply derived from other components of ecosystem metabolism. Where measurable, concentrations and δ15N values may provide additional information to the model by constraining nitrification.
The δ
15N values of benthic periphyton and invertebrates (
Loomer 2008;
Loomer et al. 2014) as well as macrophtes (
Hood 2012;
Hood et al. 2014) are often used as indicators of different N sources and N pollution because they form the base of the food web. Interpreting these data requires an ability to understand and predict the fate of large isotopically distinct N sources like WWTP effluent since the δ
15N values measured in biota ultimately depend on the source of N and isotopic fractionation during biological assimilation. Moreover, macrophytes integrate N over a much longer time scale than the effluent-plume travel time or diel variability (
Hood et al. 2014;
Loomer et al. 2014).
There are several key model parameters that are insufficiently characterized, such as isotopic fractionation during biological assimilation of both TAN and
, preferential use of different N species, and release of TAN and
. The variability in isotopic fractionation during biological
assimilation is large and varies nonlinearly with concentration (
Hoch et al. 1992;
Pennock et al. 1996;
Yoneyama et al. 2001). This poses a vexing problem at the ecosystem scale since the isotopic enrichment—concentration relationship varies between species and both concentrations and species vary within ecosystems.
The mass and δ
15N of river biomass are difficult to capture in the parsimonious NANNO model structure; model fitting may be improved if the release of TAN and
by biomass contributes significantly to river N relative to WWTP effluent (
Loomer et al. 2014). Nitrogen assimilation and release rates can be estimated with nutrient spiralling techniques but this analysis often conflates TAN and
. It is therefore difficult to discern which N form is used, which is released, and how these results apply to a river with more than 100 km of upstream nutrient inputs. The degree of importance, if any, to dissolved organic N mineralization or N release from microbes and macrophytes in the nutrient-replete WWTP plumes is unknown.
Understanding the ecosystem effects of changes in nitrogen sources, such as altering WWTPs to produce only
instead of
to improve river O
2 concentrations requires knowledge about which N enters the base of the foodweb via primary producers and consumers. In cases where δ
15N-
and δ
15N-TAN values are far enough apart, or one is changing while the other is constant, the use of each by primary producers and consumers may be teased apart. Biological
assimilation is associated with little to no isotopic fractionation (
Mariotti et al. 1981;
Yoneyama et al. 1998,
2001) and in the WWTPs’ effluent plumes in the Grand River δ
15N-
values do not vary as much as δ
15N-TAN values. In such scenarios, response to increasing δ
15N-TAN may be observable as a concomitant increase in the δ
15N of primary producers and consumers (
Hood et al. 2014;
Loomer et al. 2014).
Since O
2, N, and phosphorus cycles are strongly linked, improving the understanding of nitrogen processes allows previous work on O
2 and phosphorus cycling in the Grand River (
Barlow-Busch et al. 2006;
Venkiteswaran et al. 2014,
2015) to be extended to process-based biogeochemical models that incorporate multiple elements and their isotopes. Components that may be added to NANNO to improve constraints on nitrogen processes include δ
18O-
values. However, recent work has demonstrated that predicting the δ
18O values of nitrogenous species is more complicated than originally thought because there are poorly understood abiotic factors that alter the δ
18O value of
and
as well as multiple pathways to produce N
2O (
Buchwald and Casciotti 2010;
Casciotti et al. 2010;
Snider et al. 2010,
2012,
2013,
2015;
Buchwald et al. 2012). Nevertheless, there are opportunities to produce a more constrainable model.