# Disparate Seasonal Nitrate Export From Nested Heterogeneous Subcatchments Revealed With StorAge Selection Functions

Tam V. Nguyen, Rohini Kumar, Andreas Musolff, Stefanie R. Lutz, Fanny J. Sarrazin, Sabine Attinger, Jan H. Fleckenstein

<figure><img src="image_1.png" /></figure>

► To cite this version:

Tam V. Nguyen, Rohini Kumar, Andreas Musolff, Stefanie R. Lutz, Fanny J. Sarrazin, et al.. Disparate Seasonal Nitrate Export From Nested Heterogeneous Subcatchments Revealed With StorAge Selection Functions. Water Resources Research, 2022, 58 (3), 10.1029/2021WR030797. http://dx.doi.org/10.1029/2021WR030797

HAL Id: hal-04760727

https://hal.science/hal-04760727v1

Submitted on 31 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Water Resources Research®

RESEARCH ARTICLE

10.1029/2021WR030797

**Key Points:**

* We introduce a novel spatially varying StorAge Selection (SAS)-based model to explore nitrate export from nested subcatchments with heterogeneous settings

* Age selection preference for discharge, transit times of discharge, and nitrate export dynamics varied spatially among subcatchments

* Contrasting SAS functions between subcatchments seasonally shift the dominant source contributions to overall catchment nitrate export

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

T. V. Nguyen,
tam.nguyen@ufz.de

**Citation:**

Nguyen, T. V., Kumar, R., Musolff, A., Lutz, S. R., Sarrazin, F., Attinger, S., & Fleckenstein, J. H. (2022). Disparate seasonal nitrate export from nested heterogeneous subcatchments revealed with StorAge Selection functions. *Water Resources Research*, *58*, e2021WR030797. https://doi.org/10.1029/2021WR030797

Received 9 JUL 2021
Accepted 3 MAR 2022

**Author Contributions:**

**Conceptualization:** Tam V. Nguyen, Rohini Kumar, Andreas Musolff, Jan H. Fleckenstein

**Formal analysis:** Fanny Sarrazin, Sabine Attinger, Jan H. Fleckenstein

**Funding acquisition:** Sabine Attinger, Jan H. Fleckenstein

**Investigation:** Rohini Kumar, Andreas Musolff, Jan H. Fleckenstein

**Methodology:** Tam V. Nguyen, Rohini Kumar, Andreas Musolff, Stefanie R. Lutz, Jan H. Fleckenstein

© 2022. The Authors.

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

<figure><img src="image_3.png" /></figure>

Disparate Seasonal Nitrate Export From Nested Heterogeneous Subcatchments Revealed With StorAge Selection Functions

Tam V. Nguyen¹ ID, Rohini Kumar² ID, Andreas Musolff¹ ID, Stefanie R. Lutz¹,³ ID, Fanny Sarrazin² ID, Sabine Attinger²,⁴ ID, and Jan H. Fleckenstein¹,⁵ ID

¹Department of Hydrogeology, Helmholtz Centre for Environmental Research – UFZ, Leipzig, Germany, ²Department of Computational Hydrosystems, Helmholtz Centre for Environmental Research – UFZ, Leipzig, Germany, ³Copernicus Institute of Sustainable Development, Utrecht University, Utrecht, The Netherlands, ⁴Institute of Earth and Environmental Sciences, University of Potsdam, Potsdam, Germany, ⁵Hydrologic Modelling Unit, Bayreuth Center of Ecology and Environmental Research (BayCEER), University of Bayreuth, Bayreuth, Germany

**Abstract** Understanding catchment controls on catchment solute export is a prerequisite for water quality management. StorAge Selection (SAS) functions encapsulate essential information about catchment functioning in terms of discharge selection preference and solute export dynamics. However, they lack information on the spatial origin of solutes when applied at the catchment scale, thereby limiting our understanding of the internal (subcatchment) functioning. Here, we parameterized SAS functions in a spatially explicit way to understand the internal catchment responses and transport dynamics of reactive dissolved nitrate (N-NO₃). The model was applied in a nested mesoscale catchment (457 km²), consisting of a mountainous partly forested, partly agricultural subcatchment, a middle-reach forested subcatchment, and a lowland agricultural subcatchment. The model captured flow and nitrate concentration dynamics not only at the catchment outlet but also at internal gauging stations. Results reveal disparate subsurface mixing dynamics and nitrate export among headwater and lowland subcatchments. The headwater subcatchment has high seasonal variation in subsurface mixing schemes and younger water in discharge, while the lowland subcatchment has less pronounced seasonality in subsurface mixing and much older water in discharge. Consequently, nitrate concentration in discharge from the headwater subcatchment shows strong seasonality, whereas that from the lowland subcatchment is stable in time. The temporally varying responses of headwater and lowland subcatchments alternate the dominant contribution to nitrate export in high and low-flow periods between subcatchments. Overall, our results demonstrate that the spatially explicit SAS modeling provides useful information about internal catchment functioning, helping to develop or evaluate spatial management practices.

1. Introduction

Agricultural practices have been identified as the main cause of poor water quality in many areas worldwide. High nitrate (N-NO₃) concentrations are commonly found in groundwater and surface water in areas with intensive agriculture (Randall & Mulla, 2001; Thorburn et al., 2003). Groundwater and surface water with high nitrate concentrations can negatively affect human health and the ecosystem (Boeykens et al., 2017; Knobeloch et al., 2000). In Europe, despite implemented regulations on agricultural practices (e.g., Council Directive 91/676/EEC), high nitrate concentrations in groundwater and surface water in many areas have persisted for several decades (European Commission, 2018; Knoll et al., 2019). To further develop and evaluate such regulations, understanding how catchments retain and release water and solutes (e.g., nitrate) plays an important role, especially for mesoscale catchments (10¹⁻¹⁰ km², Breuer et al., 2008) since management is often implemented at this scale (European Environment Agency, 2012).

At the mesoscale, catchments characteristics (e.g., land use, management practices, soil, topography, geological settings, and climatic conditions) are often heterogeneous (Dupas et al., 2020; Ebeling et al., 2021; Wollschläger et al., 2017). These characteristics were found to be linked to archetypal catchment solute export regimes (Ebeling et al., 2021; Musolff et al., 2015, 2017). However, in highly heterogeneous catchments, the internal (subcatchment) responses could be significantly different from the integrated catchment response, such that the integrated catchment response cannot be used to infer subcatchment behavior (Ehrhardt et al., 2019; Lassaletta et al., 2009; Scanlon et al., 2010; Winter et al., 2021). Therefore, effective spatial management of nutrient export in mesoscale catchments calls for an understanding of subcatchment functioning and its spatial integration.
**Project Administration:** Sabine Attinger, Jan H. Fleckenstein

**Resources:** Jan H. Fleckenstein

**Software:** Tam V. Nguyen

**Supervision:** Jan H. Fleckenstein

**Validation:** Tam V. Nguyen

**Visualization:** Tam V. Nguyen

**Writing – original draft:** Tam V. Nguyen

**Writing – review & editing:** Tam V.

Nguyen, Rohini Kumar, Andreas Musolff,

Stefanie R. Lutz, Fanny Sarrazin, Sabine

Attinger, Jan H. Fleckenstein

In recent years, the StorAge Selection (SAS) functions concept has emerged as a useful tool to improve our mechanistic understanding of catchment functioning (Botter et al., 2011; Harman, 2019; Hrachowitz et al., 2016; J. Yang, Heidbüchel, et al., 2018; Nguyen et al., 2021; Rinaldo et al., 2015). SAS functions describe catchment mixing and release of water and dissolved solutes of different ages, thus regulating the transit time distributions (TTDs) and solute composition of outflows (Botter et al., 2011; Harman, 2015; van der Velde et al., 2012). It is noted that the term “catchment mixing” (hereafter also called subsurface mixing) within the SAS function concept refers to the mixing at the catchment outlet, where water and solutes from different flow paths/ages eventually exit the catchment. SAS functions are typically incorporated into catchment-scale transport models in a lumped approach and rarely used in a distributed approach. The lumped approach (catchment-scale SAS functions) represents the integrated response of the catchment (Benettin et al., 2013; Nguyen et al., 2021), tracing the temporal dynamics of dissolved solutes in discharge at the catchment outlet, but not their explicit spatial origin. It remains unclear, to what extent parameters obtained from the spatially lumped approach are transferable to the subcatchment scale given potentially different subcatchment responses, as previously mentioned.

A spatially distributed SAS approach accounts for spatial heterogeneity in mesoscale catchments and can thus provide insights into subcatchment functioning and the spatiotemporal origin of solutes in outflows. In the distributed approach, SAS functions are applied for each model grid cell. Different implementations of the distributed SAS approach have been proposed. For example, Nguyen et al. (2021) used the non-well mixed SAS functions for each individual grid cell. Remondi et al. (2018) used several well-mixed SAS functions for different vertical storage compartments within a grid cell. Although the well-mixed assumption is applied for each vertical storage compartment, the overall response of the grid cell could be far from well-mixed (Benettin et al., 2017; Remondi et al., 2018). These approaches could reasonably represent solute export at the catchment outlet (Nguyen et al., 2021) as well as the internal gauging stations (Remondi et al., 2018). The aforementioned applications of the distributed approach are limited to either catchment with homogeneous geological settings (Nguyen et al., 2021) or to transport of conservative solutes (Remondi et al., 2018), while applications of these approaches for catchments with heterogeneous geological settings and non-conservative solutes (e.g., nitrate) are still lacking. While numerical studies have been able to provide insights into the functional forms of SAS functions (which represent subsurface mixing dynamics) at the catchment scale (e.g., J. Yang, Heidbüchel, et al., 2018), the functional forms and spatial variability of SAS functions at the grid-scale largely remain unknown. Furthermore, direct verification of the functional forms of SAS functions for each grid cell (e.g., using numerical groundwater models with particle tracking) would be technically/computationally very demanding if at all feasible. Therefore, a semi-distributed SAS approach, in which a few SAS compartments represent distinct subcatchments, may represent a reasonably sized modeling unit for which we can establish sufficient process understanding to verify SAS functions and solute concentrations.

We hypothesize that a semi-distributed SAS approach can capture the spatial heterogeneity of the catchment at an intermediate level and provide an understanding of subcatchment functioning. With the semi-distributed SAS approach, SAS functions at the subcatchment level can be validated (a) indirectly using instream solute/tracer concentrations at the internal gauging stations or (b) directly using numerical groundwater models (if necessary). Despite the potential benefits of the semi-distributed approach as mentioned above or elsewhere (Hrachowitz et al., 2016; Nguyen et al., 2021), an application or implementation of this concept has not yet been attempted. In addition, the temporal dynamics of SAS functions in large catchments have not been given enough attention with SAS-based models. Previous studies have restricted the temporal changes of SAS functions between (a) young (b) old, or (c) both young and old water selection preference schemes (Nguyen et al., 2021; van der Velde et al., 2015), while more selection preference schemes could exist (J. Yang, Heidbüchel, et al., 2018).

Considering the aforementioned issues, the main objective of this research is to provide a mechanistic understanding of nitrate export dynamics from a nested mesoscale catchment using the SAS approach. For this purpose, we modified the mesoscale Hydrologic Model (mHM)-SAS model (Nguyen et al., 2021) to enable its application in a semi-distributed manner and to improve the representation of the temporal dynamics of SAS functions. The modified model is used to explore subcatchment functioning in terms of nitrate export dynamics in a mesoscale catchment with three nested subcatchments located in Central Europe with a total area of 457 km². We also evaluate if a spatially lumped SAS approach could be used for understanding subcatchment functioning, especially in terms of nitrate export. Through this study, we aim at advancing the application of spatially explicit SAS-based
<figure><img src="image_4.png" /><figcaption>Figure 1. The modified mesoscale Hydrologic Model-StorAge Selection model (Nguyen et al., 2021) with added instream processes.</figcaption></figure>

models for mesoscale heterogeneous catchments, thereby informing the design of management strategies that tackle nitrate-related issues at both local and regional scales.

## 2. Methodology

### 2.1. The mHM-SAS Model

The mHM-SAS model (Nguyen et al., 2021) consists of a spatially distributed soil nitrogen model and a spatially lumped or distributed nitrate transport model for the subsurface below the soil/root zone (Figure 1). The mHM-SAS model uses the hydrological model of the mHM (Kumar et al., 2013; Samaniego et al., 2010), the soil nitrogen model of the HYdrological Predictions for the Environment model (Lindström et al., 2010; X. Yang, Jomaa, et al., 2018), and the subsurface transport model with SAS functions (van der Velde et al., 2012). The mHM-SAS model allows applying SAS functions for (a) the subsurface (representing the saturated and unsaturated zones below the soil/root zone) over the entire catchment (lumped SAS approach) or (b) the subsurface of each model grid cell (distributed SAS approach). Below are the brief technical descriptions of each compartment of the mHM-SAS. A more detailed description of the mHM-SAS model can be found elsewhere (Nguyen et al., 2021; Samaniego et al., 2010; X. Yang, Jomaa, et al., 2018).

The water balance for each soil layer in the mHM-SAS was calculated as follows (Samaniego et al., 2010):

$$ \frac{dx_i^j(t)}{dt} = (1 - \rho^j) \cdot I_i^{j-1}(t) - E_i^j(t) - I_i^j(t) \quad (1) $$

where i and j denote the indices of the grid cell and soil layer, respectively, $x_i^j(t)$ [L³] is the depth of the soil moisture content, $E_i^j(t)$ [L³T⁻¹] is the actual evapotranspiration, $\rho^j$ [-] is the impervious surface fractions, and $I_i^{j-1}(t)$ [L³T⁻¹] and $I_i^j(t)$ [L³T⁻¹] are the infiltrated water from the soil layer j - 1 to the soil layer j and the exfiltrated water from the soil layer j to the lower layer, respectively. Exfiltrated water from the last soil layer enters the subsurface. The water and water-age balance equations for the subsurface are expressed as follows (Benettin & Bertuzzo, 2018; Nguyen et al., 2021):

$$ S(t) = S_0 + V(t) \quad (2) $$

$$ \frac{dV(t)}{dt} = I(t) - Q(t) \quad (3) $$
$$ \frac{\partial S_T(t, T)}{\partial t} = I(t) - Q(t) \cdot P_Q(t, T) - \frac{\partial S_T(t, T)}{\partial T} \quad (4) $$

where $S_0$ [L³] and $S(t)$ [L³] are the initial storage and the storage at time $t$, respectively, $V(t)$ [L³] is the variation of storage based on the total inflow $I(t)$ [L³T⁻¹] to the SAS compartment and the total discharge $Q(t)$ [L³T⁻¹] out of the SAS compartment, $P_Q(t, T) = \int_0^\infty p_Q(t, T) \cdot dT$ [−] is the cumulative TTD of discharge, $p_Q(t, T)$ [T⁻¹] is the TTD at time $t$, and $T$ [T] is the age of water since its entry to the SAS compartment. The TTD, $p_Q(t, T)$ [T⁻¹], is linked with the residence time distribution (in form of the normalized age-ranked storage, $P_S$ [−]; Benettin & Bertuzzo, 2018) via the SAS function $\omega_Q(P_S, t)$ [−]:

$$ p_Q(t, T) = \omega_Q(P_S, t) \cdot \frac{\partial P_S}{\partial T} \quad (5) $$

Within the soil zone, the mHM-SAS model considers the transformation of nitrogen (N) between different N pools (dissolved inorganic nitrogen—DIN, dissolved organic nitrogen, active organic nitrogen, and inactive organic nitrogen) via mineralization, dissolution, and degradation. The mineralization, dissolution, and degradation rates depend on the first-order rate constants and environmental factors (Lindström et al., 2010; X. Yang, Jomaa, et al., 2018). DIN is assumed to be exclusively composed of nitrate (N-NO₃) (X. Yang, Jomaa, et al., 2018). Nitrate is transported with water from the soil zone to the subsurface (below the soil zone) and eventually to the stream. In this study, we focus on the transport of nitrate in the subsurface. Using a first-order reaction for subsurface denitrification, the nitrate concentration in discharge is calculated as follows:

$$ C_Q(t) = \int_0^{+\infty} C_S(T, t) \cdot \exp(-k \cdot T) \cdot p_Q(t, T) \cdot dT \quad (6) $$

where $C_S(T, t)$ [ML⁻³] is the nitrate concentration in the SAS compartment at time $t$ [T], the initial concentration in storage is $C_0 = C(T, t = 0)$, and $k$ [T⁻¹] is the first-order denitrification rate constant.

## 2.2. The Modified mHM-SAS Model

In this study, we modified the mHM-SAS model to enable a semi-distributed SAS approach. The subcatchments were used as the spatial units for which the SAS functions are applied. This is in line with the common idea that SAS functions are catchment-scale descriptors (Botter et al., 2011; van der Velde et al., 2012). Subcatchment delineation should not only be based on the surface or subsurface drainage area but also ensure a certain uniformity in topography, land use, and geological settings. Therefore, there is no unique way to define the subcatchment size, which is further discussed in detail in the case study (Section 2.3).

In the semi-distributed SAS approach, incoming fluxes from the fully distributed soil zone (i.e., the grid cells) to the SAS compartment need to be aggregated following the subcatchment delineation. Discharge and nitrate export from the SAS compartment enter the stream network and are routed from the upstream reach to the downstream reach (Figure 1). For streamflow routing, we adopted the Muskingum-Cunge method (Cunge, 1969), which was implemented in the Soil and Water Assessment Tool (Neitsch et al., 2011). Further modifications were added to account for instream processes as instream nitrate removal can be significant (e.g., Alexander et al., 2000) and instream nitrate dynamics are different from subsurface solute dynamics. In this study, instream nitrate removal via denitrification and uptake are lumped into net instream denitrification as follows (Lindström et al., 2010; X. Yang, Jomaa, et al., 2018):

$$ L_N^{\text{removal}} = k_{\text{instr}} \cdot f_{\text{temp}} \cdot f_{\text{conc}} \cdot c_{\text{area}} \quad (7) $$

where $L_N^{\text{removal}}$ [MT⁻¹] is the net instream nitrate removal, $k_{\text{instr}}$ [ML⁻²T⁻¹] is the instream N removal flux per unit area, $f_{\text{temp}}$ [−] and $f_{\text{conc}}$ [−] are the reduction factors taking into account temporal varying instream water temperature and nitrate concentration, respectively, and $c_{\text{area}}$ [L²] is the streambed area. A more detailed description of the variables $f_{\text{temp}}$, $f_{\text{conc}}$, and $c_{\text{area}}$ can be found elsewhere (Lindström et al., 2010; X. Yang, Jomaa, et al., 2018).

Furthermore, we modified the parameterization of SAS functions in the mHM-SAS model. In this study, we focus on the two-parameter beta function (Equation 8) because of its flexibility in representing different types of
selection preferences for outflows and its practical use (Buzacott et al., 2020; J. Yang, Heidbüchel, et al., 2018; Nguyen et al., 2021; van der Velde et al., 2015). In previous studies, the temporal variability of the beta function parameters was restricted to certain limited types of selection preferences. For example, van der Velde et al. (2012) fixed one parameter of the beta function as a constant, limiting the selection preference to either (a) young or (b) old water preferences according to catchment storage. Nguyen et al. (2021) used a step function to represent the temporal changes of the selection preference scheme (the beta function) for young or old (or both young and old) water in storage based on changes in the antecedent hydrologic conditions (the ratio between the accumulated inflow and outflow over previous time steps). In this study, we generalized the concept proposed by Nguyen et al. (2021) by allowing the selection preference scheme to change continuously based on antecedent hydrologic conditions (Equations 9–11). The temporal changes in the parameters of the beta function are expressed as follows:

$$\text{beta}(P_S, a, b) = \frac{\Gamma(a + b)}{\Gamma(a) \cdot \Gamma(b)} \cdot P_S^{a-1} \cdot (1 - P_S)^{b-1} \quad (8)$$

$$r(t) = \frac{\int_{t-n}^{t} I(t) \cdot dt}{\int_{t-n}^{t} Q(t) \cdot dt} \quad (9)$$

$$a = \frac{\alpha}{r(t)} \quad (10)$$

$$b = \beta \cdot r(t) \quad (11)$$

where beta($P_S$, $a$, $b$) is the beta function with two positive shape parameters $a$ [–] and $b$ [–], $\Gamma$ is the gamma function, $r(t)$ [–] is the ratio between inflow and outflow to the SAS compartment during the time [$t-n, t$], $n$ [T] is the time window to account for antecedent hydrologic conditions, $Q(t)$ [L$^3$T$^{-1}$] is the outflow from the SAS compartment at time $t$, and $\alpha$ [–] and $\beta$ [–] are time-invariant parameters that control the rate of change of $a$ and $b$ with $r(t)$. In this approach, $\alpha$, $\beta$, and $n$ are model parameters ($\alpha$, $\beta$, $n > 0$). Equations 8–11 show that an increase in $r(t)$ will result in a decrease in $a$ and an increase in $b$, indicating a stronger preference for younger water. This reflects that an increase in $r(t)$ represents an increase in catchment storage (or wetness), leading to the selection of younger water from storage. Compared to previous approaches (Nguyen et al., 2021; van der Velde et al., 2015), this approach does not restrict the parameter range of the beta function, allowing for all selection preference schemes that the beta function could represent (Figure S1 in Supporting Information S1).

## 2.3. Study Area and Data

The study area is the Selke catchment located in the northeastern Harz Mountains, Germany. The Selke catchment has an area of about 457 km$^2$ with diverse landscapes and hydrogeological settings (Figures 2a–2d). The catchment consists of both lowland and mountainous areas with elevation ranging between 106 and 592 m above mean sea level (a.m.s.l) (Figure 2a). In the mountainous part, agricultural lands are patchy. The lowland areas are characterized by extensive agricultural land use (Figure 2b). Both soil and geological maps show that the mountainous areas are less heterogeneous than the lowland areas (Figures 2c and 2d). In the mountainous areas (steeper slope), cambisols with high permeability overlaying low permeable schist and claystone layers result in predominantly shallow flow paths (Jiang et al., 2014). In the lowland areas (mild slope), chernozems with low permeability overlaying sedimentary deposit layers allow for the development of deeper flow paths. The mountainous areas have shallower aquifers compared to the lowland areas with deeper aquifers.

Based on the distinct catchment characteristics and to make use of the observed data from the three gauging stations (Silberhütte, Meisdorf, and Hausneindorf) for model evaluation, we delineated the Selke catchment into three subcatchments, namely the upper, middle, and lower Selke (Figure 2). The upper Selke is a mixed agriculture-forest subcatchment with high altitude, high average annual rainfall, steep slope, shallow aquifer, and shallow flow paths. The middle Selke is a forest-dominated subcatchment with hydrogeological settings similar to the upper Selke. The lower Selke is an agriculturally dominated subcatchment with gentle topography (mild slopes), deeper aquifers, and deep subsurface flow paths. Detailed information about these subcatchments is presented in Table 1 (see also Figure 2 for the spatial arrangement of different landscape attributes).
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Location of the study area and subcatchment delineation with (a) elevation, (b) land use, (c) soil types, and (d) geological units.</figcaption></figure>

In this study, model input and evaluation data were combined from different sources. Daily precipitation, temperature, and potential evapotranspiration were provided by the German Weather Service (DWD). Daily streamflow and instream nitrate concentration were obtained from the State Office of Flood Protection and Water Management of Saxony-Anhalt (LHW) and Helmholtz Center for Environmental Research (UFZ), respectively.

<caption>**Table 1**</caption>

Information About the Upper, Middle, and Lower Selke Subcatchments

<table><thead><tr><th>Subcatchment</th><th>Upper Selke</th><th>Middle Selke</th><th>Lower Selke</th></tr></thead><tbody><tr><td>Outlet gauge</td><td>Silberhütte</td><td>Meisdorf</td><td>Hausneindorf</td></tr><tr><td>Area (km² and % catchment area)</td><td>100.9 (22.1%)</td><td>78.9 (17.2%)</td><td>277.6 (60.7%)</td></tr><tr><td>Forest (% subcatchment area)</td><td>61.5</td><td>87.5</td><td>12.1</td></tr><tr><td>Agriculture (% subcatchment area)</td><td>36.0</td><td>10.2</td><td>75.8</td></tr><tr><td>Average elevation (m a.m.s.l)</td><td>448.9</td><td>370.0</td><td>164.8</td></tr><tr><td>Average slope (%)</td><td>6.8</td><td>11.5</td><td>2.6</td></tr><tr><td>Dominant soil types</td><td>Dystric/spodic cambisols</td><td>Haplic chernozems</td><td></td></tr><tr><td>Dominant geological units</td><td>Mississippian wacke/shale</td><td>Sedimentary material</td><td></td></tr><tr><td>Annual average precipitation (mm/year) (data from 2012 to 2019)</td><td>515.6</td><td>457.5</td><td>432.5</td></tr><tr><td>Average annual contribution to total catchment discharge (%)</td><td>50</td><td>25</td><td>25</td></tr></tbody></table>
Estimated nitrate load from wastewater treatment plants (WWTP) as well as their locations were taken from X. Yang, Jomaa, et al. (2018). Land use management practices (fertilizer, manure application, and crop rotation) are based on field surveys and interviews (J. Yang, Heidbüchel, et al., 2018; X. Yang, Jomaa, et al., 2018). Other data (digital elevation model, land use, soil, and geological map) were provided by the Federal Institute for Geosciences and Natural Resources, Germany. Meteorological forcing constitutes of daily total precipitation and average air temperature were acquired from the German Weather Service (DWD). The point station data were gridded at a spatial resolution of 1 × 1 km² using the external drift kriging interpolation approach with terrain elevation as an external variable (X. Yang, Jomaa, et al., 2018; Zink et al., 2017). The potential evapotranspiration was estimated with the Hargreaves and Samani (1985) method.

## 2.4. Parameter Sensitivity Analysis

The objective of parameter sensitivity analysis is to identify the parameters (or processes) that contribute most to the variability of streamflow and instream nitrate concentrations. This information is further used to select parameters for optimization. The Elementary Effect Test (EET; Campolongo et al., 2007; Morris, 1991) implemented in the Sensitivity Analysis For Everybody (Pianosi et al., 2015) toolbox was used for parameter sensitivity analysis. The EET is an effective tool for screening non-influential parameters for models with a high number of parameters (Campolongo et al., 2007; Pianosi et al., 2016). A further description of the EET is presented in Text S1 in Supporting Information S1.

In this study, all global (catchment) and local (subcatchment-specific) parameters (M = 75 parameters) were selected for sensitivity analysis (Table S1). Global parameters are catchment-scale parameters, while local parameters are SAS-related parameters that are defined for each subcatchment. The parameter ranges were selected based on previous studies (J. Yang, Heidbüchel, et al., 2018; Neitsch et al., 2011; Nguyen et al., 2021; X. Yang, Jomaa, et al., 2018) and parameter distributions were assumed to be uniform. Parameter sensitivity analysis was carried out for the period 2012–2019. All model runs were performed at a daily time step with a spatial resolution of 1 km². Detailed results of the parameter sensitivity analysis are shown in Text S2 and Figures S2 and S3 in Supporting Information S1.

## 2.5. Parameter Estimation and Uncertainty Analysis

In this study, parameters were optimized for the period 2012–2015 and validated for the period 2016–2019 using observed streamflow and instream nitrate concentrations at the Silberhütte, Meisdorf, and Hausneindorf gauging stations (Figure 2). Based on the results of parameter sensitivity analysis, we selected the 21 most sensitive parameters (8 hydrological parameters and 13 nitrate parameters) and the initial nitrate concentration in the subsurface for optimization (Table 3, Text S2 and Figure S2 in Supporting Information S1). These selected parameters include the different SAS-related parameters of all subcatchments, allowing for the quantification of the uncertainty in the subsurface mixing and TTs. The insensitive (uncalibrated) parameters were fixed to values that were either defined from the previous study (Nguyen et al., 2021) or randomly selected within their indicated ranges (Table S1).

For parameter optimization, we generated 400,000 parameter sets using the uniform Latin Hypercube Sampling (LHS) technique. LHS is an efficient approach for searching an ensemble of optimal solutions, accounting for parameter uncertainties (Abbaspour et al., 2004; Sarrazin et al., 2018). The same initial ranges of subsurface transport parameters (α, β, and k) in the three subcatchments were used (Section 3.2). This means that we did not impose any prior knowledge on subsurface mixing, water age, and denitrification conditions in these subcatchments. The model prediction uncertainty was characterized by the 95 percent prediction uncertainty (95PPU) band of behavioral simulations (Abbaspour et al., 2004). The lower and upper limits of the 95PPU band correspond to the 2.5% and 97.5% percentiles of the output variable at the respective time step. The 95PPU band was evaluated by the *p* factor [0, 1] (the percentage of measured data bracketed by the 95PPU band) and *r* factor [0, ∞) (the average thickness of the 95PPU band divided by the standard deviation of the measured data) (Abbaspour et al., 2004). In general, higher *p* and lower *r* factors indicate lower prediction uncertainty.
<table><caption>**Table 2**
List of Simulation Scenarios</caption><thead><tr><td>Simulation scenario (SC)</td><td>SAS approach (number of subcatchments)</td><td>Calibrated gauging station</td></tr></thead><tbody><tr><td>SC1: base case</td><td>Semi-distributed (3)</td><td>Silberhütte, Meisdorf, and Hausneindorf</td></tr><tr><td>SC2: lumped SAS</td><td>Lumped (1)</td><td>Hausneindorf</td></tr><tr><td>SC3: semi-distributed SAS</td><td>Semi-distributed (3)</td><td>Using calibrated parameters from SC 2</td></tr></tbody></table>

*Note.* All simulation scenarios use the conceptual model as shown in Figure 1 with the number of subcatchments varies from 1 to 3, depending on the simulation scenario.

The model performance was evaluated using the Nash-Sutcliffe Efficiency (NSE; Nash & Sutcliffe, 1970), its logarithmic transformation (lnNSE), and the bias (BIAS) (Text S3 in Supporting Information S1). Behavioral simulations were selected using “soft rules” (e.g., Choi & Beven, 2007; Hartmann et al., 2017; Sarrazin et al., 2018) by defining different threshold values for NSE, lnNSE, and BIAS for streamflow (Q) and instream nitrate concentrations (C). This ensures that the simulated results for both Q and C at all gauging stations meet a certain quality. The threshold values for NSE, lnNSE, and BIAS were defined based on the simulated results, in a way that allows uncertainty to be quantified (Hartmann et al., 2017), and are presented in Text S3 in Supporting Information S1.

## 2.6. Evaluating the Spatial Model Structure

Besides the aforementioned simulations (hereinafter referred to as simulation scenario 1: base case), we performed additional simulation scenarios (SC 2 and SC 3; Table 2) to evaluate the spatial model structure (semi-distributed and lumped SAS). Specifically, we determined whether the subsurface transport parameters obtained from the lumped SAS approach (catchment-scale SAS functions, SC 2) are applicable for the subcatchments, providing a similar understanding of the subcatchment functioning as the semi-distributed SAS-based approach (SC 3). In the lumped SAS approach (SC 2), we conceptualized the entire subsurface of the Selke catchment as a single storage compartment and applied the SAS concept to model nitrate export from this compartment. For this evaluation, the lumped SAS model was calibrated at the catchment outlet (Hausneindorf gauging station). The parameters selected for optimization were based on the result of sensitivity analysis from the semi-distributed SAS model (SC 1, Table S2). Then, the calibrated model parameters from the lumped approach (SC 2) were used for the subcatchments (SC 3) to validate their applicability. In both semi-distributed (SC 1) and lumped (SC 2) SAS approaches, the same criteria were applied to select behavioral simulations.

# 3. Results and Discussion

## 3.1. Streamflow, Instream Nitrate Concentrations, and Instream Nitrate Removal

Figure 3 shows the simulated streamflow and instream nitrate concentrations at the three gauging stations from the base case scenario SC1 (Table 3). The model could well capture the seasonality of streamflow and instream nitrate concentrations at the internal gauging stations (Silberhütte and Meisdorf) as well as at the catchment outlet (Hausneindorf). The model could represent high instream nitrate concentrations during the exceptional drought years 2018 and 2019 (Hari et al., 2020), which were not part of the model calibration. However, high flows are consistently underestimated by the model, which is a common issue with hydrological models driven by daily meteorological forcing (e.g., Mizukami et al., 2019).

Statistical indices (the median NSE, lnNSE, and BIAS) show that the model performance was satisfactory (Figure S4a in Supporting Information S1). In general, the model performance for the validation period is slightly better than for the calibration period (except for instream nitrate concentrations at the catchment outlet), indicating a slight underfitting in the calibration period. Considering differences in hydrological conditions between the calibration and validation periods, in which the validation period is drier with a multi-year drought period, the slight underfitting in the calibration period is acceptable. The NSE for instream nitrate concentrations at the catchment outlet during the calibration period is low due to the low seasonality of the observed data (Figure 3f or Figure S3a in Supporting Information S1). In this case, the NSE is high only if it can explain the short time-scale (e.g.,
<table><caption>**Table 3**
List of the Selected Parameters for Optimization and the Statistical Characteristics of Behavioral Parameter Sets of the Base Case Scenario SC1</caption><thead><tr><th rowspan="2">Parameter</th><th rowspan="2">Description</th><th colspan="2">Initial range</th><th>Calibrated</th></tr><tr><th>Min</th><th>Max</th><th>Median [min, max]</th></tr></thead><tbody><tr><td colspan="5"><strong>Global (catchment-scale) parameter</strong></td></tr><tr><td>soil<sub>4</sub></td><td>Pedotransfer function parameters for soil hydrology routines of mHM</td><td>0.65</td><td>0.95</td><td>0.78 [0.65, 0.95]</td></tr><tr><td>soil<sub>6</sub></td><td></td><td>-0.37</td><td>-0.18</td><td>-0.32 [-0.37, -0.26]</td></tr><tr><td>soil<sub>7</sub></td><td></td><td>0.54</td><td>1.12</td><td>0.81 [0.57, 1.08]</td></tr><tr><td>soil<sub>9</sub></td><td></td><td>-0.55</td><td>-0.09</td><td>-0.27 [-0.55, -0.11]</td></tr><tr><td>soil<sub>14</sub></td><td>Fraction of roots in forest areas</td><td>0.90</td><td>0.99</td><td>0.98 [0.96, 0.99]</td></tr><tr><td>soil<sub>17</sub></td><td>Shape factor for calculating infiltration</td><td>1.00</td><td>4.00</td><td>2.49 [1.69, 3.23]</td></tr><tr><td>runoff</td><td>Direct surface runoff parameter</td><td>0.00</td><td>5.00</td><td>3.42 [0.08, 5.00]</td></tr><tr><td>pet<sub>1</sub></td><td>Correction factor for potential evapotranspiration</td><td>0.70</td><td>1.30</td><td>0.96 [0.92, 1.00]</td></tr><tr><td>k<sub>na</sub></td><td>Denitrification rate in nonagricultural soil (day<sup>-1</sup>)</td><td>1.00e-8</td><td>1.00e-1</td><td>1.00e-2 [4.13e-3, 2.42e-2]</td></tr><tr><td>k<sub>a</sub></td><td>Denitrification rate in agricultural soil (day<sup>-1</sup>)</td><td>1.00e-8</td><td>1.00e-1</td><td>2.18e-2 [5.68e-3, 4.11e-2]</td></tr><tr><td>k<sub>str</sub></td><td>Instream N removal flux per unit area (kg m<sup>-2</sup> day<sup>-1</sup>)</td><td>1.00e-8</td><td>1.00e-3</td><td>2.72e-6 [3.07e-8, 2.76e-4]</td></tr><tr><td>C<sub>0</sub></td><td>Initial nitrate concentration in the subsurface (mg/L)</td><td>0.5</td><td>10.0</td><td>7.86 [4.43, 8.85]</td></tr><tr><td colspan="5"><strong>Local (subcatchment-specific) parameter</strong></td></tr><tr><td>α<sub>up</sub></td><td>Parameters of the SAS function (upper Selke)</td><td>0.01</td><td>5.00</td><td>0.36 [0.10, 0.97]</td></tr><tr><td>β<sub>up</sub></td><td></td><td>0.01</td><td>5.00</td><td>4.28 [0.77, 4.84]</td></tr><tr><td>S<sub>0_up</sub></td><td>Initial subsurface storage of the upper Selke (mm)</td><td>500.00</td><td>5,000.00</td><td>798.0 [565.0, 4,959.9]</td></tr><tr><td>k<sub>up</sub></td><td>Subsurface denitrification rate in the upper Selke (day<sup>-1</sup>)</td><td>1.00e-8</td><td>1.00e-2</td><td>9.06e-3 [3.42e-3, 9.60e-3]</td></tr><tr><td>α<sub>mid</sub></td><td>Parameters of the SAS function (middle Selke)</td><td>0.01</td><td>5.00</td><td>0.44 [0.10, 1.29]</td></tr><tr><td>β<sub>mid</sub></td><td></td><td>0.01</td><td>5.00</td><td>3.29 [1.06, 4.00]</td></tr><tr><td>k<sub>mid</sub></td><td>Subsurface denitrification rate in the middle Selke (day<sup>-1</sup>)</td><td>1.00e-8</td><td>1.00e-2</td><td>1.05e-3 [3.01e-4, 7.87e-3]</td></tr><tr><td>α<sub>low</sub></td><td>Parameters of the SAS function (lower Selke)</td><td>0.01</td><td>5.00</td><td>1.84 [0.22, 4.78]</td></tr><tr><td>β<sub>low</sub></td><td></td><td>0.01</td><td>5.00</td><td>1.95 [0.12, 4.71]</td></tr><tr><td>k<sub>low</sub></td><td>Subsurface denitrification rate in the lower Selke (day<sup>-1</sup>)</td><td>1.00e-8</td><td>1.00e-2</td><td>4.96e-6 [3.43e-7, 8.40e-5]</td></tr></tbody></table>

daily) fluctuations in the observed data (Schaefli & Gupta, 2007). Such short time-scale fluctuations may be interpreted as noise in the data due to measurement/observational errors. Nevertheless, other statistical indices, for example, the Kling-Gupta efficiency (Gupta et al., 2009) and the Pearson correlation coefficient, indicate good model performance for instream nitrate concentrations at the catchment outlet (Figure S4a in Supporting Information S1). The *r* factors for instream nitrate concentrations (*C*) tend to be higher than the *r* factors for streamflow (*Q*), indicating higher uncertainty for modeling instream nitrate concentrations (Figure S4b in Supporting Information S1). This is expected because the nitrate submodel is affected by additional uncertainties in model structure and input data related to the agricultural management practices. The *p* factors for both *C* and *Q* show that less than 60% of the observed values are inside the 95PPU band. This is acceptable considering the narrow width of the 95PPU band (reflected in small *r* factors) and strict criteria for NSE, lnNSE, and BIAS for behavioral solutions (Text S3 in Supporting Information S1).

The results show that the fraction of instream nitrate removal (the instream nitrate removal mass divided by the total instream nitrate loading) is highly seasonal, namely high during summer and low during winter (Figure S5 in Supporting Information S1). This is consistent with findings from previous studies in the area (X. Yang, Jomaa, et al., 2018). The fraction of instream nitrate removal during the drought periods in 2018 and 2019 was high due to unusually high water temperature and low-flow conditions in these periods. Although the fraction of instream nitrate removal could be up to about 50% during dry periods, the maximum cumulative instream nitrate removal among all behavioral simulations for the entire simulation period (2012–2019) accounts for a maximum of 3% of
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Simulated streamflow and instream nitrate concentrations at (a and d) the Silberhütte, (b and e) the Meisdorf, and (c and f) the Hausneindorf gauging stations in the base case scenario SC1. Solid lines indicate the median values, while bands indicate the 95 percent prediction uncertainty bands.</figcaption></figure>

the total nitrate export. The overall instream nitrate removal, however, could be significant for other areas (e.g., Alexander et al., 2000).

## 3.2. Behavioral Parameter Ranges

Statistical information about the behavioral parameter sets are shown in Table 3. Among the calibrated parameters, only local parameters provide information about subcatchment functioning. The calibrated subsurface
mixing parameters for upper ($α_{up}$, $β_{up}$) and middle ($α_{mid}$, $β_{mid}$) Selke are in a similar and much narrower range, while those for the lower Selke cover a wider ranges (Table 3). Our model results indicate similar subsurface mixing dynamics and reaction rates between the upper and middle Selke, but different values for the lower Selke. This reflects similar hydrogeological settings in the upper and middle Selke and the distinct hydrogeological setup of the lower Selke (Figure 2). The behavioral ranges of $α_{low}$ and $β_{low}$ parameters in the lower Selke are not significantly reduced compared to their initial ranges, indicating a relatively high uncertainty of these parameters. Similar to previous works in the study area (Nguyen et al., 2021) and elsewhere (Benettin et al., 2015, 2017), we found that using the observed streamflow and instream solute concentrations is not sufficient to constrain the initial subsurface storage (Table 3).

The subsurface denitrification rates in the upper and middle Selke ($k_{up}$ and $k_{mid}$) are at least an order of magnitude higher compared to those in the lower Selke ($k_{low}$, very low or negligible). These findings are in line with previous studies in the area (Hannappel et al., 2018; J. Yang et al., 2021; Winter et al., 2021) and nearby areas (Ehrhardt et al., 2019). Hannappel et al. (2018) evaluated the potential for subsurface denitrification based on nitrogen input, and on observed redox potential, oxygen, iron, as well as nitrate concentrations in several wells in the same area - four wells in the upper and middle Selke, and three wells in the lower Selke. Out of four wells in the upper and middle Selke, one well shows a clear sign of denitrification, two wells do not show a clear sign of denitrification, and one well shows no sign of denitrification. However, all three wells in the lower Selke show no sign of denitrification. Therefore, hydrochemical subsurface conditions in the lower Selke seem to prevent effective denitrification, while nitrate removal by denitrification seems more likely in the upper and middle Selke. This is further confirmed by J. Yang et al. (2021), using the SAS-based model for simulating nitrate transport in a small catchment located in the upper Selke. Their estimated subsurface denitrification rate is within the behavioral range reported in our study. Data-driven analyses also showed that subsurface denitrification in the lower Selke (Winter et al., 2020) and its nearby catchments (Ehrhardt et al., 2019) is unlikely to be an important process.

### 3.3. Subcatchment Discharge and Nitrate Export

Figure 4a shows the contribution of discharge from each subcatchment to the total catchment discharge. Overall, the simulated results show that a dominant fraction of catchment discharge (about 48%–51% considering the 95PPU band) originates from the upper Selke although it only accounts for 22.1% of the catchment area. The middle and lower Selke contribute a comparable amount of discharge (23%–25% and 24%–29% of catchment discharge, respectively) despite having significantly different areal percentages (17.2% and 61.7%, respectively). These results are comparable with those obtained from observed data (Table 1). Although the fraction of total discharge from the upper and middle Selke varies seasonally in a wide range, it remains mostly above 50%, and thus constitutes a dominant source of catchment discharge even during low-flow periods (Figure 4a).

In terms of exported nitrate load, the lower Selke contributes a substantial portion of nitrate load (about 44%–55%) despite its relatively low discharge contribution (Figures 4a and 4b). The exported nitrate loads from the upper and middle Selke account for 31%–38% and 13%–18% of the catchment nitrate export, respectively. During high-flow periods, the exported nitrate load from the upper and middle Selke (predominantly the upper Selke) is much higher than that from the lower Selke (Figure 4b). During low-flow periods, however, the lower Selke contributes the major fraction of the catchment nitrate export. This is because during low-flow periods (a) instream nitrate concentrations in discharge from the lower Selke are much higher than that from the upper and middle Selke (Figure 4c), and (b) discharge contribution from the lower Selke could increase up to 50%. The results also show that instream nitrate concentrations from the upper and middle Selke have a clear seasonal pattern (high during high-flow and low during low-flow periods), while that from the lower Selke is relatively stable (Figure 4c). This is related to the differences in the subsurface mixing, transport time, and denitrification timescale (Section 3.4). The dynamics of instream concentrations in discharge from individual subcatchments found in this study are in line with results from a data-driven analysis (Winter et al., 2020). The uncertainty in the simulated nitrate concentrations in discharge from the lower Selke is relatively large during low-flow periods in 2012, 2016, and 2018 compared to other periods (Figure 4c). This is due to the uncertainty in the estimated nitrate concentrations in the oldest water pool (or initial nitrate concentration $C_0$) and the interplay between denitrification and transport timescales (Section 3.4).

The reliability of the simulated nitrate concentrations in discharge from the middle Selke (Figure 4c) is lower compared to that from the upper Selke. This is because most of the nitrate at the Meisdorf gauging station
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Contribution of (a) discharge *Q*, (b) exported nitrate load $L_{N-NO_3}$, and (c) nitrate concentrations in discharge $C_Q$ from individual subcatchments (scenario SC1). The superscripts “up,” “mid,” and “low” indicate the upper, middle, and lower Selke, respectively. Discharge and exported nitrate load were aggregated from daily to monthly for better visualization. Light (blue, red, and green) color bands in panels (a–c) and gray band in panels (a and b) indicate the 95 percent prediction uncertainty (95PPU) bands from behavioral simulations, darker (blue, red, and green) color bands in panels (a–c) indicate the area (volume of discharge, and mass of nitrate) under the 95PPU bands, and solid lines in panels (a–c) indicate the median values.</figcaption></figure>

originates from the upper Selke, so the nitrate concentration data at the Meisdorf gauge do not contain much additional information for the model calibration. Considering the aforementioned reason and the comparable instream dynamics (Figure 4c) as well as the behavioral subsurface parameter ranges (Table S2) between the upper and middle Selke, those subcatchments could be considered as one subcatchment for future studies.

### 3.4. Linking Subsurface Mixing With Water Age and Nitrate Export

Figure 5 shows the relation between subsurface mixing and the TTs of discharge as well as nitrate export from the three subcatchments. The upper Selke tends to select young water ($a/b < 1$ in Equations 8–11) for discharge, apart from low-flow periods ($a/b > 1$; Figure 5a). The catchment progressively shifts from young (or old) to old (or young) water selection preference with decreasing subsurface storage (catchment wetness, Figure 5a, Figure S6 in Supporting Information S1). The Spearman correlation coefficient between *r*, which represents subsurface mixing dynamics, and subsurface storage in the upper Selke is relatively high with the median value of 0.71 (Figure S6 in Supporting Information S1). The simulated instream nitrate concentrations at Silberhütte are not notably better than those obtained after recalibrating the model using a time-invariant SAS approach (Figure S7
<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 5.** Relation between subsurface mixing dynamics (characterized by the *a/b* ratio of the beta function), TTs (characterized by the median transit time *TT*<sub>50</sub>), and nitrate concentration dynamics in discharge (*C*<sub>Q</sub>) from (a and b) the upper Selke and (c and d) the lower Selke, and (e) the interplay between transport time and denitrification timescale (characterized by the Damköhler number, *D*<sub>a</sub>) in the upper and lower Selke. Solid lines indicate the median values, while bands indicate the 95 percent prediction uncertainty bands. The superscripts “*up*” and “*low*” mean the upper and lower Selke, respectively. *C*<sub>percol</sub> is the nitrate concentration in percolated water. The Damköher number is the ratio between *TT*<sub>50</sub> of discharge and the reaction time (1/*k*). Results are from individual subcatchments.</figcaption></figure>

in Supporting Information S1). This means that, in the time-invariant approach, the calibrated parameter values may compensate for the model structural deficiencies. The time-variant SAS approach is more consistent with mechanistic understanding of subsurface mixing dynamics obtained from a detailed physically based model in the study area (J. Yang, Heidbüchel, et al., 2018; X. Yang, Jomaa, et al., 2018) or a SAS-based model in another area (Benettin et al., 2015). For example, in the time-variant approach, the selection preference for young water during high-flow periods is consistent with our understanding that fast shallow flow paths dominate under these flow conditions. These flow paths can be activated due to a combination of (a) high precipitation, (b) high permeability of the uppermost cambisols layer, and (c) low percolation rate of the lower schist and claystone layers (Section 2.3 and Figure 2c). During low-flow periods, we found a dominance of old water in the simulated discharge, causing a strong difference in the TTs of discharge between high and low-flow periods (Figure 5b). This mirrors that the relative contribution of discharge from deeper and longer flow paths to streamflow becomes more pronounced
in low-flow periods, because less flow from the shallower zone with shorter flow paths is generated when those shallow flow paths increasingly cease. It is noteworthy that the maximum TT is restricted by the time frame of the simulation rather than the actual age of the oldest water, which is unknown (Figure 5b). Discharge with older water has less nitrate compared to discharge with younger age due to longer time for denitrification, creating a pronounced seasonality in instream nitrate concentrations (Figure 5b). In addition, the seasonality in instream nitrate concentrations is also due to the seasonality of nitrate concentrations in the percolation water (Figure 5b). However, due to denitrification and subsurface mixing, the range of nitrate concentrations in discharge is buffered compared to that in the percolating water.

The middle Selke shows a similar behavior to the upper Selke in terms of subsurface mixing and nitrate export dynamics (Figure S8 in Supporting Information S1). In addition, the subsurface denitrification rates in the upper and middle Selke are comparable (Table 3). This is because the upper and middle Selke have similar hydrogeological settings (Table 1 and Figure 2). Visual assessment shows that the model prediction uncertainties (the 95PPU of *a/b*, $TT_{50}$, and $C_Q$) for the middle Selke tend to be higher than that for the upper Selke (Figure S4 in Supporting Information S1) for the reason mentioned earlier (Section 3.3).

Compared to the upper and middle Selke, the selection preference for discharge (*a/b* ratio) in the lower Selke varies over a smaller range (Figures 5a and 5c, Figure S5b in Supporting Information S1). This is not surprising considering that the lower Selke has smaller topographic gradients (flatter terrain) and a deeper aquifer system with more steady, less dynamic subsurface flow field (J. Yang, Heidbüchel, et al., 2018; Nixdorf & Trauth, 2018). The median Spearman correlation coefficient between *r* and subsurface storage in the lower Selke is 0.53. The median *a/b* ratio shows that subsurface mixing in the lower Selke varies around the symmetrical mixing schemes (*a/b* = 1, Figure S1 in Supporting Information S1) except during the very dry periods in which the system discharges only old water. As a result, the TTs of discharge from the lower Selke are much higher than those from the upper and middle Selke, which preferably discharge young water most of the time (Figures 5a and 5c and Figure S5b in Supporting Information S1). Data-driven analyses in the area also indicated higher TTs in the lower Selke compared to the upper and middle Selke (Winter et al., 2021). The relation between nitrate concentrations in discharge and TTs of discharge from the lower Selke is unclear, as nitrate concentrations in discharge from the lower Selke seem to be relatively steady throughout the years (Figure 5d). This is because subsurface mixing in the lower Selke (Figure 5c) is relatively stable around the *a/b* ratio of 1, which describes the symmetrical mixing behavior. Interestingly, the median subsurface transport time in the lower Selke subcatchment is faster compared to the denitrification timescale defined by the very low denitrification rate (Figure 5e and Table 3). During low-flow periods, the initial nitrate concentration in the oldest water pool has negligible impacts on the nitrate export from the upper Selke compared to that from the lower Selke (Figures 5b and 5d). This is because the upper Selke during those periods is characterized by relatively long subsurface transport times compared to the denitrification timescale so that denitrification is controlled by the high denitrification rates and most nitrate is removed along the deeper flow paths. In contrast, in the lower Selke subsurface transport times, although generally longer than in the upper Selke, are short relative to the very long reaction time scales caused by the very low denitrification rates, making the system transport-controlled as indicated by the $D_a$ numbers during low-flow periods (Figure 5e). In general, subsurface transport in the upper Selke is characterized by a strong variability of transport time-scales over the denitrification timescale (shown by the $D_a$ numbers, Figure 5e), while subsurface transport in the lower Selke is more steady and characterized by transport time scales that are shorter than the respective reaction time scales.

The TTs of discharge at the internal gauging stations (Silberhütte and Meisdorf) and at the catchment outlet (Hausneindorf) are not substantially different (Figure S9 in Supporting Information S1). This results from (a) the similarity in the TTs of discharge from the upper and middle Selke (as discussed in Section 3.4) and from the fact that (b) discharge at the catchment outlet is dominated by discharge from the upper and middle Selke (Figure 4a). This highlights difficulties in understanding subcatchment TTs dynamics at downstream gauging stations where discharge from different subcatchments is mixed.

### 3.5. Semi-Distributed Versus Lumped SAS Approach

Comparing results from the semi-distributed (SC 1) and lumped (SC 2) SAS approaches show that model performances for instream nitrate concentrations at the catchment outlet are somewhat different (Figure 6a). The median statistical indices (NSE and lnNSE) indicate that the lumped approach calibrated with data from the
<figure><img src="image_10.png" /><figcaption>(b) Instream nitrate concentration ($C_Q$) at the Silberhütte gauging station (Scenarios SC 1 and SC 3)</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>(c) Instream nitrate concentration ($C_Q$) at the catchment outlet (Scenarios SC 1 and SC 3)</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Figure 6.** Observed and simulated (a and c) instream nitrate concentrations at the catchment outlet (Hausneindorf) and (b) the internal gauging station (Silberhütte) from different simulation scenarios (Table 2). Solid lines indicate the median values, while bands indicate the 95 percent prediction uncertainty bands.</figcaption></figure>

catchment outlet only has a better model performance than the semi-distributed approach in the calibration period (Figure 6a). The slightly poorer model performance of the semi-distributed model, despite having higher degrees of freedom, is because the semi-distributed SAS model is constrained with streamflow and instream nitrate concentration data not only from the catchment outlet, but also from the internal gauging stations (Table 2). In the validation period, however, it can be seen clearly that the semi-distributed approach performs significantly better than the lumped approach. This suggests that the dynamics of age selection for discharge and the associated turnover of nitrate are indeed distinctly different between the three sub-catchments and cannot be adequately represented with the simpler lumped approach. A better model performance in the calibration period with the lumped approach could be an artifact of the optimization. In fact, a misrepresentation of three SAS functions by one SAS function can be compensated by other non-SAS-related parameters in the calibration period but not in the validation period. The results further suggest that for model applications beyond the calibration period (e.g., climate change and land-use change impact studies), the semi-distributed approach should be preferred over the lumped approach. For streamflow simulation, the two approaches have comparable results with median NSE values from both approaches being within the range [0.73, 0.89] for both calibration and validation periods (not shown). This indicates that subsurface discretization is more relevant for water quality modeling than hydrological modeling, which is in line with other studies (e.g., Jha et al., 2004).
Next, we compared the simulations of nitrate dynamics based on the semi-distributed and lumped SAS approaches to understand how well the internal catchment functioning can be represented with spatially lumped SAS functions. Figures 6b and 6c shows the simulated instream nitrate concentrations at the internal gauging station (e.g., Silberhütte) and the catchment outlet (Hausneindorf) using the calibrated catchment-scale subsurface parameters (α, β, k, S₀; Table S2) obtained from the lumped approach for all subcatchments (SC 3). It is clearly visible that these parameters cannot be used for the subcatchments as they provide a false understanding of the subcatchment functioning (Figure 6b). For example, the simulated nitrate concentrations in discharge from the upper Selke (SC 3) are relatively high and steady throughout the years, while those from the observed data and the semi-distributed SAS approach (SC 1) show a strong seasonality (Figure 6b). The relatively steady simulated nitrate concentrations in discharge (SC 3) are due to (a) faster TTs and higher nitrate concentrations in young water (percolated water) during high-flow periods and (b) longer TTs and high nitrate concentration in the old water pool due to the very low denitrification rate during low-flow periods (Figure S10 in Supporting Information S1, Table S2). In the dry periods, the simulated nitrate concentrations (SC 3) are even slightly higher than those in the high-flow periods (e.g., during 2018–2019; Figure 6b), suggesting contrasting catchment functioning compared to observed data and results from the semi-distributed approach (SC 1). The simulated subsurface mixing and TT dynamics from the upper Selke (SC 3) indicate that using parameters from the lumped approach (SC 2) will shift the selection preference for discharge to much older water compared to the semi-distributed approach (SC 1) (Figure S10 in Supporting Information S1).

Despite a clear mismatch at the internal gauging station, the simulated instream nitrate concentrations at the catchment outlet match quite well the observations (Figure 6c). This indicates that taking the same subsurface transport parameters for all subcatchments in a highly heterogeneous catchment could provide the right results at the catchment outlet for the wrong reasons. For spatially explicit SAS models, this means that a parameter regionalization technique could be needed to parameterize the subsurface transport parameters of each spatial modeling unit (e.g., sub-catchments or HRUs or grid-cells) to be applicable in heterogeneous catchments, thus assisting (land use) management decisions.

Results from the lumped and semi-distributed models also imply that if the lower Selke is further divided into smaller modeling units, individual responses from these modeling units can be different from the integrated response of the lower Selke. In addition, the subsurface denitrification rates are expected to vary in a wider range when smaller modeling units are used (as demonstrated from two scenarios SC 1 and SC3). This is because the geological setting of the lower Selke is highly heterogeneous (Figure 2d). In this case, additional data (internal gauging stations) are required for further understanding the internal functioning of different modeling units within the lower Selke. In contrast, further discretization of the upper and middle Selke into smaller modeling units might not change our understanding of the internal subcatchment functioning, as the soil and geological conditions in these areas are rather homogeneous. Therefore, the responses of these smaller modeling units are expected to be similar (as shown by the similar responses of the upper and lower Selke; Section 3.4).

## 4. Model Capabilities, Implications for Management Practices, and Limitations

This study demonstrated that the spatially explicit (e.g., semi-distributed) SAS approach can provide valuable additional insights into the functioning of each subcatchment with internally consistent process descriptions, while at the same time it does not compromise the quality of the model fit at the integral point of the main catchment outlet. In contrast, the lumped SAS approach could only yield robust results at the main catchment outlet and yielded inadequate results at internal points in the model domain. Our application of the semi-distributed SAS model in a nested mesoscale heterogeneous catchment has demonstrated the model's ability to capture nitrate dynamics at internal gauging stations as well as at the main catchment outlet. Applying SAS functions in a semi-distributed framework as presented here, helps to overcome some of the limitations of the spatially lumped characteristics of the general SAS concept. Results from a semi-distributed model can provide not only additional spatial information, such as subcatchment nitrate export, but also temporal information on the age of water and potentially nitrate, which is) related to the source and origin of the exported nitrate.

The spatially explicit SAS approach is especially relevant for planning and evaluating spatial management practices as (a) parameters infer from the lumped approach could fail to represent the subcatchment functioning, (b) the lumped approach is less robust than the semi-distributed approach, and (c) the lumped approach does not
provide information about both spatial and temporal origins of nitrate in discharge for effective management. Results from the Selke with the semi-distributed SAS approach show that the lowland catchments (lower Selke) should have different management practices compared to the mountainous headwater catchments (middle and upper Selke). Agricultural management practices that aim to quickly reduce nitrate export during high-flow periods should be implemented in the mountainous headwater catchments rather than in the lowland catchment. This is because of the short TTs and transport-limited characteristics of these catchments during high-flow periods. However, management practices that aim to reduce exported nitrate loads (a) during low-flow periods or (b) in the coming decade(s) should be implemented in the lowland catchments with longer TTs and transport-limited characteristics. Due to the short median TTs in the mountainous catchments (~1 year), the effectiveness of management practices in these catchments can be evaluated in the following years. In contrast, long median TTs in the lowland catchment would require decades for the effects of a certain management practice to become effective and visible.

Despite the advantages of the semi-distributed SAS approach, the application of this approach in larger catchments with more diverse hydrogeological settings could face several challenges. In such catchments, the number of subsurface parameters could be high due to a high number of subcatchments. In this case, understanding the linkage between key catchment characteristics (e.g., topography, geology, land use, and meteorological conditions) with subcatchment functioning (parameters of the SAS function) could avoid unnecessary small spatial resolution and model overparameterization. This can provide useful insights into the optimal spatial modeling resolution, in which the number of modeling units is at a minimum while the spatial heterogeneity of subcatchment responses is adequately captured. For such understanding, applications of the semi-distributed SAS approach in much larger catchments with diverse settings are required. In addition, the calculation of instream nitrate removal in the current version of the model is relatively simple using a removal rate constant and other reduction factors to account for instream temperature and nitrate concentrations, more elaborated approaches would be useful to disentangle different nitrate removal and transformation processes in the stream network (Mulholland et al., 2008; X. Yang et al., 2019).

## 5. Conclusions

In this study, we developed a semi-distributed SAS-based model, in which SAS functions are applied at the subcatchment level. The proposed model was applied in a mesoscale nested catchment, namely the Selke catchment located in Germany. The catchment was delineated into three subcatchments for application of SAS functions, consisting of (a) a upper mountainous headwater subcatchment (upper Selke) with a mixture of forest and agricultural land, (b) a middle mountainous subcatchment (middle Selke) dominated by forest land, and (c) a lowland subcatchment (lower Selke) dominated by agricultural land. The main results from this study are as follows:

1. The semi-distributed SAS approach could represent instream nitrate concentration dynamics not only at the catchment outlet but also at the internal gauging stations.

2. The headwater subcatchment has high seasonal variations in the subsurface mixing schemes, while that in the lowland catchment is less pronounced. Nitrate concentrations in discharge from the headwater subcatchment show a strong seasonality, while those from the lowland subcatchment are relatively steady over different seasons.

3. Instream denitrification only removes a minor part of the exported nitrate loads.

4. The median age of water in discharge ($TT_{50}$) from the headwater subcatchment is much younger than that from the lowland subcatchment.

5. The headwater and lowland subcatchments take turns at dominating catchment nitrate export in high and low-flow periods.

6. Parameters infer from the lumped approach fail to represent the subcatchment functioning and the lumped approach is less robust than the semi-distributed approach.

Results from this study have demonstrated that the proposed model can provide useful insights into the functioning of each subcatchment, unlike the lumped SAS approach. The proposed model concept in combination with an appropriate regional parameterization approach could help to extend the application of the SAS concept in
larger catchments. Results from such model applications could help understand both spatial and temporal origins of nitrate in rivers, contributing toward efforts to reduce nitrate pollution.

## Data Availability Statement

The source code and input data to reproduce this work are available at https://git.ufz.de/nguyenta/mhm-sas and https://git.ufz.de/yangx/mHM-Nitrate. The raw meteorological data can be downloaded from the DWD (https://opendata.dwd.de/climate_environment/CDC/observations_germany/climate/daily/). Discharge data are available at https://gld-sa.dhi-wasy.de/GLD-Portal/.

## Acknowledgments

The authors would like to thank the Deutscher Wetterdienst (DWD), the Federal Institute for Geosciences and Natural Resources (BGR), and the State Agency for Flood Protection and Water Management of Saxony-Anhalt (LHW) for providing the required input data. The authors would also like to thank Xiaojiang Yang and Michael Rode for providing the model code, land use/land cover management practices data, and high-frequency nitrate data for model calibration. F. Sarrazin, R. Kumar, and S. Attinger acknowledge the Advanced Earth Modelling Capacity (ESM) project funded by the Helmholtz Association. The authors thank the Editor and two anonymous reviewers for their constructive comments that helped to considerably improve the quality of the manuscript. Open access funding enabled and organized by Projekt DEAL.

## References

Abbaspour, K. C., Johnson, C. A., & van Genuchten, M. T. (2004). Estimating uncertain flow and transport parameters using a sequential uncertainty fitting procedure. *Vadose Zone Journal*, 3(4), 1340–1352. https://doi.org/10.2136/vzj2004.1340

Alexander, R. B., Smith, R. A., & Schwarz, G. E. (2000). Effect of stream channel size on the delivery of nitrogen to the Gulf of Mexico. *Nature*, 403(6771), 758–761. https://doi.org/10.1038/35001562

Benettin, P., & Bertuzzo, E. (2018). Tran-SAS v1.0: A numerical model to compute catchment-scale hydrologic transport using StorAge Selection functions. *Geoscientific Model Development*, 11(4), 1627–1639. https://doi.org/10.5194/gmd-11-1627-2018

Benettin, P., Kirchner, J. W., Rinaldo, A., & Botter, G. (2015). Modeling chloride transport using travel time distributions at Plynlimon, Wales. *Water Resources Research*, 51(5), 3259–3276. https://doi.org/10.1002/2014WR016600

Benettin, P., Soulsby, C., Birkel, C., Tetzlaff, D., Botter, G., & Rinaldo, A. (2017). Using SAS functions and high-resolution isotope data to unravel travel time distributions in headwater catchments. *Water Resources Research*, 53(3), 1864–1878. https://doi.org/10.1002/2016WR020117

Benettin, P., Van Der Velde, Y., Van Der Zee, S. E. A. T. M., Rinaldo, A., & Botter, G. (2013). Chloride circulation in a lowland catchment and the formulation of transport by travel time distributions. *Water Resources Research*, 49(8), 4619–4632. https://doi.org/10.1002/wcr.20309

Boeykens, S. P., Piol, M. N., Samudio Legal, L., Saralegui, A. B., & Vázquez, C. (2017). Eutrophication decrease: Phosphate adsorption processes in presence of nitrates. *Journal of Environmental Management*, 203, 888–895. https://doi.org/10.1016/j.jenvman.2017.05.026

Botter, G., Bertuzzo, E., & Rinaldo, A. (2011). Catchment residence and travel time distributions: The master equation. *Geophysical Research Letters*, 38(11), 1–6. https://doi.org/10.1029/2011GL047666

Breuer, L., Vaché, K. B., Julich, S., & Frede, H. G. (2008). Current concepts in nitrogen dynamics for mesoscale catchments. *Hydrological Sciences Journal*, 53(5), 1059–1074. https://doi.org/10.1623/hysj.53.5.1059

Buzacott, A. J. V., van der Velde, Y., Keitel, C., & Vervoort, R. W. (2020). Constraining water age dynamics in a south-eastern Australian catchment using an age-ranked storage and stable isotope approach. *Hydrological Processes*, 34(23), 4384–4403. https://doi.org/10.1002/hyp.13880

Campolongo, F., Cariboni, J., & Saltelli, A. (2007). An effective screening design for sensitivity analysis of large models. *Environmental Modeling & Software*, 22(10), 1509–1518. https://doi.org/10.1016/j.envsoft.2006.10.004

Choi, H. T., & Beven, K. (2007). Multi-period and multi-criteria model conditioning to reduce prediction uncertainty in an application of TOPMODEL within the GLUE framework. *Journal of Hydrology*, 332(3–4), 316–336. https://doi.org/10.1016/j.jhydrol.2006.07.012

Cunge, J. A. (1969). On the subject of a flood propagation computation method (Muskingum method). *Journal of Hydraulic Research*, 7(2), 205–230. https://doi.org/10.1080/00221686909500264

Dupas, R., Ehrhardt, S., Musolff, A., Fovet, O., & Durand, P. (2020). Long-term nitrogen retention and transit time distribution in agricultural catchments in western France. *Environmental Research Letters*, 15(11), 115011. https://doi.org/10.1088/1748-9326/abe47

Ebeling, P., Kumar, R., Weber, M., Knoll, L., Fleckenstein, J. H., & Musolff, A. (2021). Archetypes and controls of riverine nutrient export across German catchments. *Water Resources Research*, 3(4), 1–29. https://doi.org/10.1029/2020wr021834

Ehrhardt, S., Kumar, R., Fleckenstein, J. H., Attinger, S., & Musolff, A. (2019). Trajectories of nitrate input and output in three nested catchments along a land use gradient. *Hydrology and Earth System Sciences*, 23(9), 3503–3524. https://doi.org/10.5194/hess-23-3503-2019

European Commission. (2018). *Report from the commission to the council and the European Parliament on the implementation of Council Directive 91/676/EEC concerning the protection of waters against pollution caused by nitrates from agricultural sources based on Member State reports fo.* European Commission.

European Environment Agency. (2012). EEA catchments and rivers network system ECRINS v1.1 – Rationales, building and improving for widening uses to water accounts and WISE applications (Issue 7).

Gupta, H. V., Kling, H., Yilmaz, K. K., & Martinez, G. F. (2009). Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling. *Journal of Hydrology*, 377(1–2), 80–91. https://doi.org/10.1016/j.jhydrol.2009.08.003

Hannappel, S., Köpp, C., & Bach, T. (2018). Characterization of the denitrification potential of aquifers in Saxony-Anhalt. *Grundwasser*, 23(4), 311–321. https://doi.org/10.1007/s00767-018-0402-7

Hargreaves, G. H., & Samani, Z. A. (1985). Reference crop evapotranspiration from temperature. *Applied Engineering in Agriculture*, 1(2), 96–99. https://doi.org/10.13031/2013.26773

Hari, V., Rakovec, O., Markonis, Y., Hanel, M., & Kumar, R. (2020). Increased future occurrences of the exceptional 2018–2019 Central European drought under global warming. *Scientific Reports*, 10(1), 12207. https://doi.org/10.1038/s41598-020-68872-9

Harman, C. J. (2015). Time-variable transit time distributions and transport: Theory and application to storage-dependent transport of chloride in a watershed. *Water Resources Research*, 51(1), 1–30. https://doi.org/10.1002/2014WR015707

Harman, C. J. (2019). Age-ranked Storage-discharge relations: A unified description of spatially lumped flow and water age in hydrologic systems. *Water Resources Research*, 55(8), 7143–7165. https://doi.org/10.1029/2017WR022304

Hartmann, A., Antonio Barberá, J., & Andreo, B. (2017). On the value of water quality data and informative flow states in karst modelling. *Hydrology and Earth System Sciences*, 21(12), 5971–5985. https://doi.org/10.5194/hess-21-5971-2017

Hrachowitz, M., Benettin, P., van Breukelen, B. M., Fovet, O., Howden, N. J. K., Ruiz, L., et al. (2016). Transit times-the link between hydrology and water quality at the catchment scale. *Wiley Interdisciplinary Reviews: Water*, 3(5), 629–657. https://doi.org/10.1002/wat2.1155

Jha, M., Gassman, P. W., Secchi, S., Gu, R., & Arnold, J. (2004). Effect of watershed subdivision on swat flow, sediment, and nutrient predictions. *Journal of the American Water Resources Association*, 40(3), 811–825. https://doi.org/10.1111/j.1752-1688.2004.tb04460.x
Jiang, S., Jomaa, S., & Rode, M. (2014). Modelling inorganic nitrogen leaching in nested mesoscale catchments in central Germany. *Ecohydrology*, 7(5), 1345–1362. https://doi.org/10.1002/eco.1462

Knobeloch, L., Salna, B., Hogan, A., Postle, J., & Anderson, H. (2000). Blue babies and nitrate-contaminated well water. *Environmental Health Perspectives*, 108(7), 675–678. https://doi.org/10.1289/ehp.00108675

Knoll, L., Breuer, L., & Bach, M. (2019). Large scale prediction of groundwater nitrate concentrations from spatial data using machine learning. *The Science of the Total Environment*, 668, 1317–1327. https://doi.org/10.1016/j.scitotenv.2019.03.045

Kumar, R., Samaniego, L., & Attinger, S. (2013). Implications of distributed hydrologic model parameterization on water fluxes at multiple scales and locations. *Water Resources Research*, 49(1), 360–379. https://doi.org/10.1029/2012WR012195

Lassaletta, L., García-Gómez, H., Gimeno, B. S., & Rovira, J. V. (2009). Agriculture-induced increase in nitrate concentrations in stream waters of a large Mediterranean catchment over 25 years (1981–2005). *The Science of the Total Environment*, 407(23), 6034–6043. https://doi.org/10.1016/j.scitotenv.2009.08.002

Lindström, G., Pers, C., Rosberg, J., Strömqvist, J., & Arheimer, B. (2010). Development and testing of the HYPE (Hydrological Predictions for the Environment) water quality model for different spatial scales. *Hydrology Research*, 41(3–4), 295–319. https://doi.org/10.2166/nh.2010.007

Mizukami, N., Rakovec, O., Newman, A. J., Clark, M. P., Wood, A. W., Gupta, H. V., & Kumar, R. (2019). On the choice of calibration metrics for “high-flow” estimation using hydrologic models. *Hydrology and Earth System Sciences*, 23(6), 2601–2614. https://doi.org/10.5194/hess-23-2601-2019

Morris, M. D. (1991). Factorial sampling plans for preliminary computational experiments. *Technometrics*, 33(2), 161–174. https://doi.org/10.2307/1269043

Mulholland, P. J., Helton, A. M., Poole, G. C., Hall, R. O., Hamilton, S. K., Peterson, B. J., et al. (2008). Stream denitrification across biomes and its response to anthropogenic nitrate loading. *Nature*, 452(7184), 202–205. https://doi.org/10.1038/nature06686

Musolff, A., Fleckenstein, J. H., Rao, P. S. C., & Jawitz, J. W. (2017). Emergent phenotype patterns of coupled hydrologic and biogeochemical responses in catchments. *Geophysical Research Letters*, 44(9), 4143–4151. https://doi.org/10.1002/2017GL072630

Musolff, A., Schmidt, C., Selle, B., & Fleckenstein, J. H. (2015). Catchment controls on solute export. *Advances in Water Resources*, 86, 133–146. https://doi.org/10.1016/j.adwatres.2015.09.026

Nash, J. E., & Sutcliffe, J. V. (1970). River flow forecasting through conceptual models part I – A discussion of principles. *Journal of Hydrology*, 10(3), 282–290. https://doi.org/10.1016/0022-1694(70)90255-6

Neitsch, S., Arnold, J., Kiniry, J., & Williams, J. (2011). *Soil & water assessment tool theoretical documentation version 2009*. Texas Water Resources Institute. https://doi.org/10.1016/j.scitotenv.2015.11.063

Nguyen, T. V., Kumar, R., Lutz, S. R., Musolff, A., Yang, J., & Fleckenstein, J. H. (2021). Modeling nitrate export from a mesoscale catchment using StorAge Selection functions. *Water Resources Research*, 57(2). https://doi.org/10.1029/2020wr028490

Nixdorf, E., & Trauth, N. (2018). Evaluating the reliability of time series analysis to estimate variable riparian travel times by numerical ground-water modelling. *Hydrological Processes*, 32(3), 408–420. https://doi.org/10.1002/hyp.11428

Pianosi, F., Beven, K., Freer, J., Hall, J. W., Rougier, J., Stephenson, D. B., & Wagener, T. (2016). Sensitivity analysis of environmental models: A systematic review with practical workflow. *Environmental Modelling & Software*, 79, 214–232. https://doi.org/10.1016/j.envsoft.2016.02.008

Pianosi, F., Sarrazin, F., & Wagener, T. (2015). A MATLAB toolbox for global sensitivity analysis. *Environmental Modelling & Software*, 70, 80–85. https://doi.org/10.1016/j.envsoft.2015.04.009

Randall, G. W., & Mulla, D. J. (2001). Nitrate nitrogen in Surface waters as influenced by climatic conditions and agricultural practices. *Journal of Environmental Quality*, 30(2), 337–344. https://doi.org/10.2134/jeq2001.302337x

Remondi, F., Kirchner, J. W., Burlando, P., & Fatichi, S. (2018). Water flux tracking with a distributed hydrological model to quantify controls on the Spatiotemporal variability of transit time distributions. *Water Resources Research*, 54(4), 3081–3099. https://doi.org/10.1002/2017WR021689

Rinaldo, A., Benettin, P., Harman, C. J., Hrachowitz, M., McGuire, K. J., Van Der Velde, Y., et al. (2015). StorAge Selection functions: A coherent framework for quantifying how catchments store and release water and solutes. *Water Resources Research*, 51(6), 4840–4847. https://doi.org/10.1002/2015WR017273

Samaniego, L., Kumar, R., & Attinger, S. (2010). Multiscale parameter regionalization of a grid-based hydrologic model at the mesoscale. *Water Resources Research*, 46(5), 1–25. https://doi.org/10.1029/2008WR007327

Sarrazin, F., Hartmann, A., Pianosi, F., Rosolem, R., & Wagener, T. (2018). V2Karst V1.1: A parsimonious large-scale integrated vegetation-recharge model to simulate the impact of climate and land cover change in karst regions. *Geoscientific Model Development*, 11(12), 4933–4964. https://doi.org/10.5194/gmd-11-4933-2018

Scanlon, T. M., Ingram, S. M., & Riscassi, A. L. (2010). Terrestrial and in-stream influences on the spatial variability of nitrate in a forested headwater catchment. *Journal of Geophysical Research*, 115(G2). https://doi.org/10.1029/2009jg001091

Schaefl, B., & Gupta, H. V. (2007). Do Nash values have value? *Hydrological Processes*, 21(15), 2075–2080. https://doi.org/10.1002/hyp.6825

Thorburn, P. J., Biggs, J. S., Weier, K. L., & Keating, B. A. (2003). Nitrate in groundwaters of intensive agricultural areas in coastal Northeastern Australia. *Agriculture, Ecosystems & Environment*, 94(1), 49–58. https://doi.org/10.1016/S0167-8809(02)00018-X

van der Velde, Y., Heidbüchel, I., Lyon, S. W., Nyberg, L., Rodhe, A., Bishop, K., & Troch, P. A. (2015). Consequences of mixing assumptions for time-variable travel time distributions. *Hydrological Processes*, 29(16), 3460–3474. https://doi.org/10.1002/hyp.10372

van der Velde, Y., Torfs, P. J. J. F., Van Der Zee, S. E. A. T. M., & Uijlenhoet, R. (2012). Quantifying catchment-scale mixing and its effect on time-varying travel time distributions. *Water Resources Research*, 48(6), 1–13. https://doi.org/10.1029/2011WR011310

Winter, C., Lutz, S. R., Musolff, A., Kumar, R., & Weber, M. (2020). Disentangling the impact of catchment heterogeneity on nitrate export dynamics from event to long-term time scales. *Water Resources Research*.

Winter, C., Lutz, S. R., Musolff, A., Kumar, R., Weber, M., & Fleckenstein, J. H. (2021). Disentangling the impact of catchment heterogeneity on nitrate export dynamics from event to long-term time Scales. *Water Resources Research*, 57(1). https://doi.org/10.1029/2020WR027992

Wollschläger, U., Attinger, S., Borchardt, D., Brauns, M., Cuntz, M., Dietrich, P., et al. (2017). The Bode hydrological observatory: A platform for integrated, interdisciplinary hydro-ecological research within the TERENO Harz/Central German Lowland Observatory. *Environmental Earth Sciences*, 76(1), 29. https://doi.org/10.1007/s12665-016-6327-5

Yang, J., Heidbüchel, I., Musolff, A., Reinstorf, F., & Fleckenstein, J. H. (2018). Exploring the dynamics of transit times and subsurface mixing in a small agricultural catchment. *Water Resources Research*, 54(3), 2317–2335. https://doi.org/10.1002/2017WR021896

Yang, J., Heidbüchel, I., Musolff, A., Xie, Y., Lu, C., & Fleckenstein, J. H. (2021). Using nitrate as a tracer to constrain age selection preferences in catchments with strong seasonality. *Journal of Hydrology*, 603, 126889. https://doi.org/10.1016/j.jhydrol.2021.126889
AGU

ADVANCING EARTH
AND SPACE SCIENCE

Yang, X., Jomaa, S., Büttner, O., & Rode, M. (2019). Autotrophic nitrate uptake in river networks: A modeling approach using continuous high-frequency data. *Water Research*, **157**, 258–268. https://doi.org/10.1016/j.watres.2019.02.059

Yang, X., Jomaa, S., Zink, M., Fleckenstein, J. H., Borchardt, D., & Rode, M. (2018). A new fully distributed model of nitrate transport and removal at catchment scale. *Water Resources Research*, **54**(8), 5856–5877. https://doi.org/10.1029/2017WR022380

Zink, M., Kumar, R., Cuntz, M., & Samaniego, L. (2017). A high-resolution dataset of water fluxes and states for Germany accounting for para-metric uncertainty. *Hydrology and Earth System Sciences*, **21**(3). https://doi.org/10.5194/hess-21-1769-2017