Introduction

Direct N2O emission from agricultural land is a major anthropogenic source of greenhouse gas affecting climate change (Myhre et al. 2013), as well as posing a threat to the stratospheric ozone layer (Ravishankara et al. 2009). The Intergovernmental Panel on Climate Change (IPCC 2006) guidelines based on data compiled by Bouwman (1996) and Bouwman et al. (2002) recommend estimating these direct emissions at approximately 1.25% of applied synthetic fertilizers and manure. However, lower emission factors (EF) have been reported under various regional contexts from temperate (Gabrielle et al. 2006) to Mediterranean climates (Cayuela et al. 2016). On the other hand, there is evidence that emission response to increasing N inputs is exponential rather than linear (Shcherbak et al. 2014). Thus, quantifying direct N2O emission from agriculture remains uncertain at both regional and global scales.

In addition to these direct agricultural emissions, “indirect” processes, occurring far away from agricultural fields, also contribute, but their estimation [roughly 0.75% of leached nitrogen according to IPCC (2006) guidelines, revised to 1.1% in the IPCC (2019) report] is still a subject of controversy. Using a top–down approach based on the observed N2O increase in atmospheric concentration (rather than the bottom–up approach extrapolating from land-based data as recommended by the IPCC), Crutzen et al. (2008) estimated that 3–5% of anthropogenically produced reactive nitrogen at the global scale is actually emitted as N2O, implying a much higher overall emission factor compared with the IPCC factor.

Other authors using top–down approaches at regional scales based on tall tower trace gas measurements (Corazza et al. 2011; Thomson et al. 2011; Griffis et al. 2013) came to the same conclusion. Later on, Davidson (2009) and Reay et al. (2012) tried to reconcile the top–down and bottom–up methodologies by considering that N2O emitted by recycled nitrogen, such as manure nitrogen, as well as emissions from nitrate-contaminated semi-natural environments, contributes more to global N2O emissions than agricultural soils themselves. Syakila and Kroeze (2011), revisiting the global nitrous oxide budget, confirmed that the observed increase in atmospheric N2O cannot be explained by inventories based on the current IPCC guidelines, highlighting the significance of poorly constrained indirect anthropogenic emissions linked to the nitrogen cascade.

There is in fact no comprehensive and formalized vision of the cascade of processes responsible for indirect N2O emissions. Systematic measurements at the landscape scale show the existence of N2O emission hot spots in wetlands and riparian areas, in places where nitrate-enriched groundwater comes into contact with the biogeochemically active upper soil layers (Clement et al. 2002; Oehler et al. 2007; Vilain et al. 2010, 2012; Anderson et al. 2014). However, given the limited areal extent of riparian wetlands, their contribution to agricultural N2O emission, though significant, was lower than direct emissions.

Rivers and streams appear as other hotspots of N2O emissions to the atmosphere. Garnier et al. (2009) and Outram et al. (2012) showed that surface waters in agricultural drainage networks are systematically oversaturated with respect to atmospheric concentration, particularly in small stream orders, showing that rivers in general and headwaters in particular are significant N2O emitters. Interestingly, domestic wastewater-impacted rivers are also major emitters (Garnier et al. 2007, 2009; Tallec et al. 2006). Beaulieu et al (2011) directly measured N2 and N2O emission using whole-stream 15NO3 -tracer experiments in 72 headwater streams draining different land-use types across the United States. Extrapolating their results, they concluded that river network emissions can represent about 10% of global anthropogenic N2O emissions. Turner et al. (2015) reported measurements of N2O fluxes across the water–atmosphere interface of streams of different Strahler order in the Mississippi basin and suggested that taking into account headwater streams (including episodic zero-order microflow stream channels or tile drain collectors) can more than double the agricultural N2O budget compared with bottom–up estimates based on IPCC guidelines.

River emission results from two processes: (i) degassing of N2O produced in soils and aquifers once groundwater reaches surface water, and (ii) in-stream production of N2O in oxygen-depleted sectors either in the water column or in the benthos, through both nitrification and denitrification (Mulholland et al. 2008). In contrast to the former, the latter processes are particularly significant in large river stretches and estuaries, as reported by Barnes and Owens (1999), Kroeze et al. (2005), Garnier et al. (2006, 2007, 2009), and Rodrigues et al. (2007).

Here we develop a comprehensive modeling framework embracing the whole N cascade from cropland to river waters in order to relate direct N2O emission from agricultural soils and the indirect processes listed above. Several mechanistic models have been developed to calculate N2O emissions from agricultural soil profiles, such as DNDC (Li et al. 1992; Beheydt et al. 2007), DayCENT (Parton et al. 1996, 2001; Del Grosso et al. 2005), NOE (Hénault et al. 2005), and STICS (Bessou et al. 2010). These models perform well for estimating emissions linked to given pedo-climatic context and agricultural practices, but their one-dimensional (vertical) character prevents them from assessing the landscape interactions responsible for indirect emissions. A number of models were developed to simulate nitrogen cycling at the landscape scale (see the review by Cellier et al. 2011): these tools are able to calculate nitrate transfers and associated N2O emissions through the different compartments of a small terrestrial landscape (typically 1–10 km2) at very fine spatial and temporal resolution, but they require a comprehensive and detailed data set on topography, hydrology, soil conditions, and farming practices.

Our approach requires far fewer data; it is closely connected to our modeling approach of nitrate elimination at the river–watershed interface (Billen et al. 2018). It is applied here to three temperate watersheds of increasing size, characterized by intensive agriculture: the Orgeval, the Haut Loir, and the Seine. The purpose of this work is to assess, using a coherent modeling procedure, the share and location of the different processes responsible for indirect N2O emissions in the land- and waterscape of agricultural basins of different sizes.

Modeled watersheds

The Orgeval is a small catchment (106 km2) in the Brie Laitière agricultural region located approximately 70 km east of Paris, nested in the Seine basin. Agricultural land, mainly devoted to cereal cropping, occupies 82% of the total area, while forests represent 17%, hence with very few grasslands and livestock left. Active riparian wetlands, defined as potential wetlands covered by grasslands or forests, represent 3.2% of the watershed. In all, 90% of the agricultural area is equipped with tile drains. The geological, hydrological, and land use features are described in detail by Garnier et al. (2014, 2016). Direct N2O emissions from the main cropping systems have been measured with manual and automatic chambers by Vilain et al. (2010, 2011, 2012) and Benoit et al. (2014, 2015a, b) for 10 years. Annual mean values integrated over crop rotation cycles ranged between 0.65 and 0.9 kgN/ha/yr. N2O concentrations in groundwater (piezometers) and surface water have been regularly measured in parallel (Vilain et al. 2011; Garnier et al. 2007, 2009). The basin is part of the OZCAR network of experimental watersheds (https://ozcar-ri.prod.lamp.cnrs.fr/oracle/).

The Haut Loir basin (3600 km2, upstream from Chateaudun) drains two distinct agricultural regions with contrasting geological and pedological substrates. The western part of the basin is representative of the Perche region, with clayey soils, heavily tile-drained, and a dense tributary network, while the eastern part is representative of the Beauce region, with a very low river drainage density mainly fed by a deep groundwater system. Direct N2O emissions were regularly measured with manual chambers at two locations. Annual mean emission was approximately 1–2 kgN/ha/yr. The main crop rotations in the whole watershed and their N budget were recorded (see Billen et al. 2018).

N2O concentration in surface water, as well as emission fluxes at the water–atmosphere interface, were measured throughout the 2017–2018 period by Grossel et al. (in preparation).

The Seine watershed, upstream from Caudebec (73,500 km2, including the freshwater estuary), has been described in detail by Meybeck et al. (1998) and Flipo et al. (2019). Its land use comprises 57% cropland, 26% forest, 10% grassland, and 7% urban or artificialized areas. The current structure of its agricultural system is described in detail by Le Noë et al. (2017) and Billen et al. (2019). Direct and indirect N2O emissions were estimated by Garnier et al. (2009, 2019) using bottom–up methodologies. Numerous measurements of N2O concentrations in the drainage network are available. They show a clear stream order structure, with higher concentrations in headwaters, decreasing to near-saturation values in third- to fifth-order streams, then increasing again in higher-order rivers (Garnier et al. 2009).

The model described below was run on these three watersheds for three contrasting hydrological years: 2016, 2017, and 2018. The first is characterized by a relatively atypical flood in June, while the other two years present a more usual temperate regime, with high flow in winter and low flow from June to October. The average specific runoff over the 3 years was 175, 112, and 145 mm/yr for the Orgeval, Haut Loir, and Seine watersheds, respectively.

Model description

The model of N2O fluxes developed in this paper closely relies on the description of nitrate flow from cropland subroot zone to river outlet, including transfer through riparian wetlands, recently published by Billen et al. (2018). In brief, this model first assesses a fixed nitrate concentration to subroot water from each land use class of the basin, based on the N balance of typical crop rotations. This concentration is assigned to pre-calculated subsurface and base flow runoff. Unless tile drainage allows some fraction of the subsurface runoff to bypass riparian zones, the nitrate flow issued from the watershed has to cross an active denitrification area before reaching surface water. The resulting nitrate reduction that occurs there is calculated taking into account the extent of riparian wetlands in each elementary watershed, their potential denitrification rate, and the seasonal variations in temperature. Other processes affecting nitrate in the river drainage include algal uptake, water column, and benthic nitrification and denitrification as described in the RIVE model (www.fire.upmc.fr/rive), as well as inputs from urban point sources.

The same approach is considered here to model the corresponding flows of N2O. The model: (1) assigns a N2O concentration to subroot water for each land use class in the watershed; (2) assesses N2O production and emission in the riparian wetlands; and (3) calculates in-stream N2O production and the resulting emissions across the water–atmosphere interface (Fig. 1).

Fig. 1
figure 1

Principles of the model of N2O production and transfer through the land–river continuum. (Cbas, Crip and Criv stand for N2O concentration in the water phase of basin soil, riparian soil and river respectively; Sbas, Srip, Sriv are the suface area of watershed, riparian wetlands and river network.; Ebas, Erip, Eriv stand for the N2O emission from watershed soils, riparian wetlands and river surface, respectively. kvs is the gas transfer coefficient from soil to atmosphere; kr is the gas transfer coefficient from river water to the atmosphere)

Full size image

N 2 O concentration in subroot water

Initialization of the model of N2O transfer requires an assessment of a representative N2O concentration to subroot water of each land use class considered in the watershed. Direct measurement of this concentration is technically extremely difficult. However, N2O concentration could be measured in groundwater collected from piezometer or groundwater resurgences, and showed large over-saturation with respect to atmospheric partial pressure. In fact, this level of over-saturation is the very driver of direct emission of N2O from soils to the atmosphere. As mentioned above, this latter flux is commonly measured using manual or automatic chambers, or can be estimated at an annual time scale through empirical relationships involving agricultural practices and climatic features, such as the one proposed by Garnier et al. (2019). An average N2O concentration of soil water (Cbas, in µgN-N2O/l) can therefore be associated with the average annual N2O direct emission from watershed soil to the atmosphere (Ebas, in mgN-N2O/m2/h), considering a simple soil water–atmosphere transfer coefficient (kvs):

$${\text{Ebas }}\left( {{\text{mgN}} - {\text{N}}_{2} {\text{O}}/{\text{m}}^{2} /{\text{h}}} \right) \, = {\text{ kvs }} \times \, \left( {{\text{Cbas }}{-}{\text{ C}}_{{{\text{eq}}}} } \right)$$
(1)

where Ceq (mgN-N2O/m3) is the N2O concentration at equilibrium with the atmospheric partial pressure at the mean annual temperature, i.e. 0.33 mgN-N2O/m3 à 12 °C (Weiss and Price 1980). The parameter kvs (m/h) is defined as an average apparent ventilation coefficient between the soil water phase and the atmosphere. It can be calibrated on situations where both the flux of soil nitrous oxide emission and the soil water concentration are known. This is the case in a few spots in the investigated basins (Table 1), leading to a kvs value of approximately 0.0025 m/h.

Table 1 Calibration of the apparent ventilation coefficient between the soil water phase and the atmosphere. kvs is calculated as Ebas/(Cbas-Ceq)
Full size table

N2O production and emission in riparian wetlands

During the water transfer from the watershed to the river across the riparian zone, the process of denitrification, while reducing nitrate concentration of the inflowing water discharge (Q, m3/s), also tends to increase its N2O concentration. The nitrate module of the Riverstrahler model (Billen et al. 2018) explicitly calculates the reduction of nitrate concentration (ΔCNO3denrip) due to denitrification. During this process, nitrate is partly converted into N2, and partly into N2O as an intermediate of denitrification. Although the share of denitrification flux converted to N2O can vary widely according to several factors, including availability of nitrate, organic matter, as well as pH (Weier et al. 1993; Saggar et al. 2013), the assumption is made of a constant fraction (pNden) expressing the N2O/(N2O + N2) ratio. Because the value of this ratio is not very well constrained by empirical data, it is considered in the model as an adjustment parameter and a sensitivity analysis toward the chosen value is carried out (see”Discussion” section).

The rate of production of N2O (PN2Odenrip, mg/m3/h) through riparian denitrification can thus be estimated from the corresponding nitrate concentration reduction (ΔCNO3denrip, mg/m3) and the discharge Q (m3/h) flowing through the riparian zone:

$${\text{P}}_{{{\text{N}}2{\text{O}}}} {\text{denrip }}\left( {{\text{mg}}/{\text{m}}^{3} /{\text{h}}} \right) \, = \, - {\text{pNden }} \times \, \Delta {\text{C}}_{{{\text{NO}}3}} {\text{denrip }} \times {\text{ Q}}$$
(2)

As a result of this additional N2O production, a new equilibrium is established in the riparian zone, with increased emission across the riparian soil–atmosphere interface (Erip, mgN-N2O/m2/h), at a rate that can be expressed as:

$${\text{Erip }} = {\text{ kvs }} \times \, \left( {{\text{Crip }}{-}{\text{ C}}_{{{\text{eq}}}} } \right)$$
(3)

where kvs is the soil–atmosphere transfer coefficient defined above, and Crip is the average N2O concentration in the riparian zone (here considered as a perfectly mixed reactor).

This new equilibrium concentration (Crip) in the riparian zone can be calculated by considering the material balance equation between input and output of the water flux Q through the riparian reactor, considering production of N2O by specifically riparian processes (PN2Odenrip) as well as by nonspecifically riparian processes, the latter corresponding to the average N2O emission in the watershed soils and thus equal to kvs×Srip×(Cbas-Ceq)) (see Eq. 1):

$${\text{Q}} \times {\text{Crip }} = {\text{ Q}} \times {\text{Cbas }} + {\text{ P}}_{{{\text{N}}2{\text{O}}}} {\text{denrip }} + {\text{ kvs}} \times {\text{Srip}} \times \left( {{\text{Cbas}} - {\text{C}}_{{{\text{eq}}}} } \right) \, - {\text{ kvs}} \times {\text{Srip }} \times \, \left( {{\text{Crip }}{-}{\text{ C}}_{{{\text{eq}}}} } \right)$$
(4)

where Srip (m2) is the area of the riparian zone in the element of watershed considered, and Cbas here is the average N2O concentration in soil water from the watershed resulting from non-riparian processes.

From Eq. 4 it follows that:

$${\text{Crip }} = {\text{ Cbas }} + {\text{ P}}_{{{\text{N}}2{\text{O}}}} {\text{denrip}} / \, \left( {{\text{Q}} + {\text{kvs}} \times {\text{Srip}}} \right)$$
(5)

Or

$${\text{Crip }} = {\text{ Cbas}} - {\text{pNden }} \times \, \Delta {\text{C}}_{{{\text{NO}}3}} {\text{denrip }}/ \, \left( {1 + {\text{kvs}} \times {\text{Srip}}/{\text{Q}}} \right)$$
(6)

This procedure also allows calculation of the emission flux from riparian soil to the atmosphere (Eq. 3), thus providing a complete budget of N2O in the riparian zone.

In-stream processes

Point sources of N2O linked to urban wastewater inputs must be considered in modeling N2O in the river network. Based on a number of surveys of N2O concentration in effluents from different types of wastewater treatment plants, as well as on detailed specific studies (Tallec et al. 2006, 2008), the model considers that wastewater inputs are at near saturation N2O concentration, except when a tertiary treatment of nitrification and/or denitrification treatment is involved, in which case the N2O concentration is 10 times the saturation level.

In-stream production of N2O can also occur by nitrification and denitrification in the river water column and in the benthic phase. These are part of the RIVE module of the Riverstrahler model, described in detail at www.fire.upmc.fr/rive, as well as by Garnier et al. (2006, 2007) for the special case of the lower Seine River. The same parameterization has been used in the present study. Regarding benthic processes, the calculation procedure used is the one resulting from the meta-model established by Billen et al. (2015), which allows calculation of the integrated benthic denitrification flux, from the calculation of the organic carbon deposition, together with the oxygen and nitrate concentrations in the water column. The resulting N2O emission flux from the benthos is assumed to be a fraction pNden of the integrated benthic denitrification.

Transfer of N2O through the water–atmosphere interface is classically represented by a relation linking the flux Eriv (mgN-N2O/m2/h) to the oversaturation of the water phase, with a transfer coefficient krN2O (m/h) depending on the flow velocity (v, m/h), the water column depth (dpth, m), and the slope (dimensionless) (Raymond et al. 2012), as done by Marescaux et al. (2018) (Weiss and Price 1980):

$${\text{Eriv }} = {\text{ kr}}_{{{\text{N}}2{\text{O}}}} \times \, \left( {{\text{C}}_{{{\text{N}}2{\text{O}}}} {-}{\text{ Ceq}}} \right)$$
(7)
$${\text{with Ceq }} = \, 0.0002 \times {\text{T}}^{{2}} \, {-}0.0167 \times {\text{T }} + \, 0.5038$$
(8)
$${\text{kr}}_{{{\text{N}}2{\text{O}}}} \left( {{\text{m}}/{\text{h}}} \right) \, = {\text{ k}}_{600} \times {\text{SQR}}\left( {600/{\text{Schmidt}}_{{{\text{N}}2{\text{O}}}} } \right)$$
(9)
$${\text{Schmidt}}_{{{\text{N}}2{\text{O}}}} = \, 2056 \, {-} \, 137 \times {\text{T }} + \, 4.317 \, \times {\text{ T}}^{2} \, {-}0.05435 \, \times {\text{ T}}^{3}$$
(10)
$${\text{k}}_{600} \left( {{\text{m}}/{\text{h}}} \right) \, = \, 5037/24 \, \times \, \left( {{\text{v}} \times {\text{slope}}} \right)^{0.89} \times {\text{ dpth}}^{0.54}$$
(11)

Results

The model was first implemented at the scale of the small watersheds Orgeval and Haut Loir, for which the direct emissions estimated from measurements with static chambers on the main land use and crop rotations were available, thus allowing an estimation of an average subroot N2O concentration used as the main constraint variable of the model. The simulated N2O concentrations in the river network are calculated for different values (from 0 to 0.07) of the adjustable parameter pNden and compared with available observations at a number of stations (Fig. 2). A satisfactory agreement is obtained for a pNden fraction of approximately 0.015. Higher values lead to overestimation of N2O concentrations in all river stations. Lower values strongly underestimate the observations (Fig. 2).

Fig. 2
figure 2

Model-calculated time and space variations of N2O concentration in river water of a the Haut Loir basin and, b the Orgeval basin, for three values of the pNden parameter (representing the fraction of total denitrification emitted as N2O in riparian and stream processes). c Scatter plot of calculated versus observed river N2O concentrations for 3 values of pNden. The r2 of the correlation and the simulation RMSE (normalized to the mean) are indicated for each group of simulations. In all figures, red curves and points correspond to pNden = 0, blue is for pNden = 0.015, green for pNden = 0.07

Full size image

The model also provides a complete budget of N2O emissions along the continuum from watershed soils to the outlet of the river network (Table 2) and can be mapped at the resolution used for the calculations, i.e. the detailed land-use and potential wetland data within each of the sub-watershed areas and 1-km stretches for the river network (Fig. 3).

Table 2 Model-calculated budget of N2O emissions along the continuum from watershed soils to the river outlet in the three watersheds investigated
Full size table
Fig. 3
figure 3

Map of calculated N2O fluxes from watershed soils, riparian wetlands, and river surfaces. a. Haut Loir basin. b. Orgeval basin. For the Haut Loir basin, a magnification of two selected areas shows the riparian wetlands along the river, with their higher emissions than the surrounding fields

Full size image

In the case of the Haut Loir basin, a number of direct measurements of the emission flux of N2O across the river–atmosphere interface are also available (Grossel et al., in preparation). They range between 0.002 and 0.28 mgN/m2/h (mean 0.11 mgN/m2/h), which fits reasonably with the model estimates of 0.08 mgN/m2/h on average at the same stations for pNden = 0.015.

These results confirm that riparian wetlands as well as river surfaces can indeed be seen as hotspots of N2O emission, with per-area rates of emission sometimes more than twice the average direct emissions from the dominant agricultural land areas in the watershed, depending on the particular structure of the land- and waterscape. However, owing to their somewhat limited share in the total watershed area, these areas contribute to only 2–6% and 0.8–-2%, respectively, of the total emissions. In other words, indirect emissions from riparian and river surfaces together in fact only represent 3–7% of the direct emissions from the soils. The Haut Loir watershed, because of its lower relative areas of wetlands, has the lowest indirect emissions.

The model was also implemented at the scale of the much larger Seine watershed, using the estimates of N2O emission from cropland, grassland, and forest in the different agricultural regions made by Garnier et al. (2019) based on the record of fertilization practices (Le Noë et al. 2017) and regional climate data. Land use was geographically distributed based on Corine Land Cover, and the riparian wetlands were defined from the map of potential wetlands established on the basis of a topographical index (Curie et al. 2007; Berthier et al. 2014; https://geowww.agrocampus-ouest.fr/web/?p=1538), as described by Billen et al. (2018).

The calculated average N2O concentrations for different values of the pNden parameter were compared with those measured since 2012 in water from different Strahler orders (Fig. 4). Some unpublished measurements available in springs or in piezometers were also used for this comparison, the latter being referred to as “zero order,” the former as “subroot.” The measurements collectively show a U-shaped distribution of concentrations versus stream orders (see also Garnier et al. 2009 and Marescaux et al. 2018). The high concentrations at zero order (or subroot water) rapidly decrease with increasing stream orders, down to near saturation levels at order 3–5. N2O concentration increases again at higher orders. This U-shaped distribution is well reproduced by the model, which, however, predicts, on average, larger concentrations in riparian waters and headwaters; this quickly tends toward saturation once exposed to the rapid evasion processes across the water–atmosphere interface in surface river bodies. In-stream processes, fostered by increased point sources of organic pollution, are responsible for the further increase in N2O concentration at higher stream order. The results obtained by simulating N2O distribution with different values for the N2O/(N2O + N2) denitrification ratio (pNden) show again that a value of about 0.015 is the maximum compatible with available observations at the scale of the whole watershed (Fig. 4). A full map of N2O emissions from soils, riparian wetlands, and rivers in the Seine basin, similar to that shown in Fig. 3, can be viewed at address: https://doi.org/10.26047/PIREN.data.N2O.2019.

Fig. 4
figure 4

Model-calculated average N2O concentration in watershed and riparian water, as well as in surface water of different Strahler streams. Comparison with available measurements carried out in the Seine basin since 2012 (Garnier et al. 2009; Marescaux et al. 2018 and unpublished data). The level of N2O saturation with respect to atmospheric concentration for two temperatures is also indicated (blue dotted lines: 0.36 µgN-N2O at 10 °C, 0.21 µgN-N2O at 25 °C)). b Scatter plot of calculated versus observed N2O concentrations for 3 values of the N2O/(N2O + N2) ratio (pNden). The r2 of the correlation and the simulation RMSE (normalized to the mean) are indicated for each group of simulations. In both figures, red curves and points correspond to pNden = 0, blue is for pNden = 0.015, green for pNden = 0.07

Full size image

The complete budget of direct and indirect N2O emission at the scale of the Seine watershed shows that riparian emissions represent 17% of direct watershed emissions, while stream surface emission accounts for only 4%. Thus, indirect emissions represent almost 20% of total watershed emissions (Table 2).

Discussion

The approach developed in this paper pursues the attempt initiated by the introduction of a module of riparian denitrification in the Riverstrahler model (Billen et al. 2018) for describing the N cascade along the continuum from watershed soils to the river outlet within a coherent modeling framework. The model proposed is constrained by the direct N2O emissions from the different land use classes of watershed soils, assumed to be known on a yearly average basis, from either empirical measurements, using e.g. static chambers, or another independent modeling approach such as, for instance, the empirical relationship linking annual N2O emission to fertilizer inputs and climate data established by Garnier et al. (2019) based on data from the literature. From this mean annual direct N2O emission value, the average N2O concentration in groundwater and subsurface runoff is calculated using a simple empirically calibrated transfer coefficient. Riparian denitrification is first calculated by the model in terms of nitrate reduction, based on the hydrological behavior of each sub-basin and on a temperature-dependent denitrification potential of riparian wetlands, the extension of which is defined by topographic information. A fraction pNden of this denitrification is assumed to be converted into N2O, which allows the model to calculate both the enhanced emission from riparian areas and the N2O concentration flowing into surface water. Within the drainage network, the model calculates the in-stream N2O production by benthic and water column processes, including nitrification and denitrification, as well as transfer across the water–atmosphere interface. The model is validated by its capability to predict the distribution of N2O concentration in river water at different places in the river network, as well as the order of magnitude of N2O emission fluxes at the river–atmosphere interface.

The model involves the calibration of two parameters, the soil–atmosphere transfer coefficient (kvs) and the fraction of denitrification emitted as N2O versus N2 (pNden). The meaning and value chosen for these two parameters requires some discussion.

N2O emission by soils is the result of complex and episodic processes, depending on highly transient conditions in a relatively heterogeneous medium (Saggar et al. 2013). This explains the difficulties in measuring and modeling them. The model approach developed here does not intend to describe these processes; it starts from the knowledge of annual emissions, which can be assessed either by in situ measurements or by using an empirical modeling approach (e.g., Garnier et al. 2019) and considers that the soil water pool that generates runoff is buffered to the extent that its average N2O concentration reflects the integrated soil N2O production processes. A soil–atmosphere transfer coefficient (kvs, in m/h) is calibrated on local situations where both the annual N2O emission flux from the soil and the average N2O concentration in the soil water pool are known. The emission values found in two such situations, in the Haut Loir and Orgeval basins, are consistently about 0.0025 m/h (Table 1). Combined with the mean runoff rate (100–270 mm/yr), these values imply that no more than 1–2% of the N2O flux produced in the soil is transferred to the hydrosphere rather than to the atmosphere. This contribution to the indirect emissions is thus somewhat limited.

The larger contributor to indirect emissions, therefore, consists of denitrification processes during the transfer of water from soils to the outlet of the river system, including riparian zones and the rivers themselves. The model calculates these processes by mechanistic modules, well established and validated by their capacity to correctly simulate nitrate and ammonium concentrations (Billen et al. 2018) with, however, one additional calibrated parameter (pNden) representing the fraction of N2O emitted during nitrate reduction through denitrification. This fraction is assumed to be constant in space and time for rivers and riparian wetlands.

N2O emission from denitrification results from an imbalance between the rate of NO3− reduction to N2O and its further reduction into N2. Decreased pH (Rochester 2003; Simek et al. 2002; Simek and Cooper 2002) or increased temperature (Benoit et al. 2015a, b) has been shown to result in higher N2O/(N2 + N2O) ratios in agricultural soils. The review by Rochester (2003) of several studies, as well as the data of Hénault et al. (2005), suggests a consistent relationship between this ratio and the pH, starting at values close to 0.01 at pH 7–8 and increasing up to 0.75 at pH 4–5. The problem with this relationship is that it was established based on laboratory studies in which the pH was manipulated. Only two results included in the relationship of Rochester (2003) come from direct field measurements: they both show a N2O/(N2 + N2O) ratio of 0.01 at a pH of 5.8 and 8.2, respectively. Simek et al. (2002) emphasized the point that the imbalance of microbial or enzymatic activities leading to N2O emission during denitrification can be transient until the readjustment of the rates of the two enzymatic processes (nitrate and N2O reductions) or the adaptation of the microbial consortium. N2O emission can therefore be a transient process linked to changes in pH or temperature rather than controlled by constant conditions of pH and temperature. We therefore limited our literature search to studies reporting pNden values based on measurements carried out under field conditions, i.e., without manipulation of the soil samples, either using 15 N-labelled fertilizer tracing, acetylene block technique on entire core samples or the results of emission chambers combined with estimation of denitrification. Even so, the range of values reported for the N2O/(N2 + N2O) ratio is considerable, varying from less than 0.01 to 0.6, without a clear relationship to pH or any other environmental variables (Denmead et al. 1979; Rudaz et al. 1999; Rochester 2003; Clement et al. 2002; Oehler et al. 2007; Scheer et al. 2009; Autret et al. 2019).

The review by Schlesinger (2009), distinguishing records of N2O/(N2 + N2O) for different types of ecosystems from agricultural soils, soils under natural vegetation or freshwater wetlands shows that the latter (0.08 ± 0.02) are typically much lower than the formers (0.37 ± 0.03 and 0.49 ± 0.07), as also shown by Beaulieu et al (2011). In a more recent paper, Butterbach-Bahl et al (2013) brought arguments showing that part of the data gathered by Schlesinger (2009), based on the use of the acetylene block technique, are subject to a systematic underestimation of N2 production by denitrification, leading to overestimation of the nitrous oxide ratio. Discarding these data leads to an average N2O/(N2 + N2O) ratio of 0.016 ± 0.04 for wetlands and rivers (based however on 5 measurements only). The measurements by Beaulieu et al (2011) obtained by in situ 15 N tracing across 53 US rivers and streams, ranged from 0.0004 to 0.056, with an interquartile range between 0.003 and 0.01.

The pNden value selected in our model calculations (0.015) is in exactly the same range. Figures 2 and 4 show the sensitivity of the model results to the value of this parameter. A more detailed analysis of the model response in terms of riparian and riverine emissions is shown for the case of the Seine watershed in Fig. 5. For pNden set to zero, the emission of riparian wetlands corresponds to the background soil processes, similar to the average of those direct emissions occurring in the watershed, while river indirect emissions consist of the evasion of N2O transferred from watershed soils to surface waters. At increasing pNden values, additional emission occurs from both riparian wetlands and river surfaces precisely due to denitrification processes occurring there. For a pNden value of approximately 0.015, the calculated emissions fit the observed distribution of river concentrations and emission fluxes. Higher values of pNden, e.g., 5 to 10 times higher (0.07—0.15), would clearly result in much higher estimates of indirect N2O emissions by riparian zones and river networks, which would reach values close to the watershed direct emissions themselves. However, such values are not supported by empirical observations of N2O concentration and emission rates available in the investigated watersheds.

Fig. 5
figure 5

Sensitivity analysis of the effect of pNden value on the model estimate of N2O emissions by riparian and riverine processes in the Seine watershed

Full size image

The estimate of indirect N2O emissions obtained with the model for the Seine watershed is consistent with previous figures published by Garnier et al. (2009) Marescaux et al. (2018) based on empirical budgeting approaches (Table 3). They do not contradict the default parameters suggested in the IPCC guidelines. However, our estimate of indirect emissions in the Seine watershed does not take into account the possible contribution of ponds, lakes, and other stagnant water bodies.

Table 3 Budget of direct and indirect N2O emissions from the Seine watershed, published by Garnier et al. (2009) and Marescaux et al. (2018), compared with the results of the present study (with an estimated 25% range based on the uncertainty of pNden calibration)
Full size table

Our estimated contribution of indirect emissions to total emissions from the Seine watershed is much lower compared with those reported by Turner et al. (2015) for the Mississippi. Our work is based on two different watersheds (Orgeval nested in the Seine basin) sharing the same geological substrate, both characterized by carbonate-rich sedimentary rocks. We cannot exclude that under other contexts and, in particular, in more acidic soils and waters, different values of pNden would apply and lead to much higher N2O emission in both river and riparian wetlands.

Conclusion

The modeling approach developed in this paper offers a coherent framework to model direct and indirect N2O emissions at the scale of the whole continuum from agricultural soils to river system. Its application to three watersheds of increasing size, from a few hundred to nearly 100,000 km2 in a temperate climate and sedimentary rock context, confirms the view that, in spite of high N2O emission rates per surface area in riparian wetlands and headwater streams, indirect emissions represent barely 20% of direct soil emissions when expressed at the scale of the whole watershed area. The largest part of these indirect emissions occurs in riparian areas where nitrate leached from the watershed soils undergoes partial denitrification. Headwaters receiving oversaturated groundwater from the watershed are also significant emitters. Downstream in the large river network, in-stream processes, possibly amplified by point source inputs, may contribute to enhanced emissions through the river–atmosphere interface.