Global agricultural ammonia emissions simulated with the ORCHIDEE land surface model

Ammonia (NH₃) is an important atmospheric constituent. It plays a role in air quality and climate through the formation of ammonium sulfate and ammonium nitrate particles. It has also an impact on ecosystems through deposition processes. About 85 % of NH₃ global anthropogenic emissions are related to food and feed production and, in particular, to the use of mineral fertilizers and manure management. Most global chemistry transport models (CTMs) rely on bottom-up emission inventories, which are subject to significant uncertainties. In this study, we estimate emissions from livestock by developing a new module to calculate ammonia emissions from the whole agricultural sector (from housing and storage to grazing and fertilizer application) within the ORCHIDEE (Organising Carbon and Hydrology In Dynamic Ecosystems) global land surface model. We detail the approach used for quantifying livestock feed management, manure application, and indoor and soil emissions and subsequently evaluate the model performance. Our results reflect China, India, Africa, Latin America, the USA, and Europe as the main contributors to global NH₃ emissions, accounting for 80 % of the total budget. The global calculated emissions reach 44 Tg N yr⁻¹ over the 2005–2015 period, which is within the range estimated by previous work. Key parameters (e.g., the pH of the manure, timing of N application, and atmospheric NH₃ surface concentration) that drive the soil emissions have also been tested in order to assess the sensitivity of our model. Manure pH is the parameter to which modeled emissions are the most sensitive, with a 10 % change in emissions per percent change in pH. Even though we found an underestimation in our emissions over Europe (−26 %) and an overestimation in the USA (+56 %) compared with previous work, other hot spot regions are consistent. The calculated emission seasonality is in very good agreement with satellite-based emissions. These encouraging results prove the potential of coupling ORCHIDEE land-based emissions to CTMs, which are currently forced by bottom-up anthropogenic-centered inventories such as the CEDS (Community Emissions Data System).

1 Introduction

Ammonia (NH₃) is a crucial species in the atmosphere, playing a role in the alteration of air quality and climate through its implication in airborne particle matter formation (PM or aerosols) (Anderson et al., 2003; Bauer et al., 2007). The NH₃ lifetime is short and has been reported to range from a few hours to a few days (Pinder et al., 2008; Behera et al., 2013), as ammonia mostly originates from surface emissions and its deposition velocity is high over most surfaces (Hov et al., 1994; Evangeliou et al., 2020). Due to this characteristic, NH₃ is transported over relatively short distances and readily reacts with abundant gases such as nitric and sulfuric acids to form secondary aerosols (Malm et al., 2004). The resulting aerosols, such as ammonium nitrates or ammonium sulfates, have important impacts on the Earth's radiative budget due to their ability to scatter incoming radiation, act as cloud condensation nuclei, and indirectly increase cloud lifetime (Abbatt et al., 2006; Henze et al., 2012; Behera et al., 2013; Evangeliou et al., 2020). The impact of NH₃ on the total radiative forcing is estimated to be −0.06 W m⁻², with respective contributions of about −0.07 and 0.01 W m⁻² to the radiative forcing of nitrate and sulfate (Myhre et al., 2013). By analyzing different Representative Concentration Pathway (RCP) scenarios, Hauglustaine et al. (2014) showed the importance of ammonia with respect to the direct aerosol forcing in the future due to the potentially significant increase in agricultural emissions. In the most extreme scenario, emissions could increase by 50 % by 2100 compared with their present-day level.

In addition to its contribution to the radiative budget, the balance between NH₃, SO₂, and NO_x emissions controls the formation of secondary inorganic aerosols (SIAs), which are important components of fine particles (PM_2.5) (Paulot et al., 2016; Fu et al., 2017; Sutton et al., 2020). Quantifying ammonia emissions is of high interest for air quality policies, as it appears that NH₃ emission reductions would also be efficient to reduce inorganic aerosol formation (Lachatre et al., 2019).

There are many issues in the development of reliable NH₃ emission inventories, as analyzed by Nair and Yu (2020), such as the lack of emission measurements, the difficulties in validating the NH₃ concentration with measurements, and the critical assumptions behind the modeling approaches in terms of emission factors and activity rates. Even though ammonia emissions are challenging to estimate, several studies have aimed at quantifying global emissions and their associated uncertainty. For example, Dentener and Crutzen (1994) estimated a global NH₃ emission of 45 Tg N yr⁻¹ (±50 %), which is a low estimate compared with the 54 Tg N yr⁻¹ (±25 %) of Bouwman et al. (2005) and the 75 Tg N yr⁻¹ (±50 %) of Schlesinger and Hartley (1992) (Zhu et al., 2015). Agricultural activities are among the significant sources of ammonia in the world, accounting for about 85 % of the global anthropogenic NH₃ emissions (Behera et al., 2013). Agricultural emissions originate from fertilizer application and livestock management, with the latter including livestock housing, manure storage, and manure application. Globally, recent studies have developed methodologies in order to quantify emissions from this sector. For example, Beusen et al. (2008), Paulot et al. (2014), and McDuffie et al. (2020) estimated similar emissions of about 32–35 Tg N yr⁻¹, which is less than the 41–47 Tg N yr⁻¹ estimates of Crippa et al. (2018) and Vira et al. (2020).

Modeling NH₃ sources from agriculture is especially difficult because it depends on several factors related to the environment (e.g., atmospheric conditions and soil properties) and to agricultural practices, which are also crucial to capture the temporal and spatial variability in emissions correctly. Emissions from manure management are driven by the amount of N in the feed, animal body characteristics, animal housing conditions, temperature, and animal waste handling practices (Anderson et al., 2003). The soil NH₃ emissions originate from N application, either from fertilizer or manure, and are controlled by the soil pH, temperature, water content, surface wind speed, and the atmospheric NH₃ concentration (Kirk and Nye, 1991; Hertel et al., 2011; Behera et al., 2013). Other factors such as the ammonium content of the fertilizer and the timing of N application are also crucial for emission estimates (Riddick et al., 2016; Vira et al., 2020).

The first type of approach to the quantification of agricultural ammonia emissions is bottom-up inventories. Most global inventories, such as the CEDS (Community Emissions Data System; McDuffie et al., 2020), EDGAR (Emissions Database for Global Atmospheric Research; Crippa et al., 2018), and HTAP (Hemispheric Transport of Air Pollution; Janssens-Maenhout et al., 2015), are based on activity data associated with corresponding emission factors (EFs). Chemistry transport models (CTMs) are usually forced with these global emission inventories. As examples, the inventory described by Bouwman et al. (1997) is prescribed in the study of Xu and Penner (2012), and the CEDS inventory (McDuffie et al., 2020) is used in Paulot et al. (2016) and Pai et al. (2021). Emission inventories do not account for environmental factors, such as the temperature or soil humidity, which is an important limitation for studying spatial–temporal variability in the atmospheric NH₃ and ${\mathrm{NH}}_{\mathrm{4}}^{+}$ concentrations. Most inventories rely on the fertilizer application period to represent the seasonality of emissions but are based on few studies and usually use the same temporal profile (most of the time reflecting European agricultural practices), which is extrapolated to the whole globe. More complex inventories exist, such as the updated version of the Global Livestock Environmental Assessment Model (GLEAM; Uwizeye et al., 2020) or the comprehensive food system developed by Conijn et al. (2018), and combine more detailed agricultural information (e.g., animal requirements, livestock system types, manure management, and the surface types receiving manure) with EFs but consider yearly emissions. Even though this type of approach is more accurate due to the detailed consideration of agricultural practices, it shows limitations for studying the temporal variability in emissions due to the static representation of the agricultural practices when using unique EFs or only one seasonal profile for the whole globe. Recently, more complex models based on an explicit description of processes that control the volatilization from soil have been developed. The FAN (Flow of Agricultural Nitrogen) model, initially developed by Riddick et al. (2016) and largely improved by Vira et al. (2020), combines information on agricultural practices, emission factors for manure management emissions, and physical processes for soil volatilization to compute NH₃ emissions from the different agricultural sources. When soil processes are tightly coupled to the main meteorological drivers, the related emissions respond to environmental changes, which is particularly interesting in the case of climate–surface interaction studies. Even if the FAN model is integrated into the Community Earth System Model (CESM), the manure produced by livestock is not directly linked to the biomass productivity, which can represent uncertainty in the N content of the manure and, therefore, in the resulting emissions.

In this study, in order to better account for the key parameters in the estimate of the NH₃ emissions, we implement a module representing the agricultural sector within the ORCHIDEE (Organising Carbon and Hydrology In Dynamic Ecosystems) land surface model (LSM). Our methodology is based on the integration of a complete dynamical agricultural module (CAMEO, Calculation of AMmonia Emissions in ORCHIDEE) within ORCHIDEE, which details a feed management module linked to the biomass productivity of the model and animal characteristic information, a manure management representation that combines regional agricultural handling practices, and a complex soil emission component based on key environmental parameters such as the vegetation growth, temperature, and soil humidity.

Section 2 describes the agricultural model within ORCHIDEE and the model setup of the 11 year control simulation (2005–2015) as well as with the sensitivity analysis simulation set. Global and regional results comparing previous works (CEDS and the FAN model from Vira et al., 2020) and seasonal analysis using airborne measurements (Infrared Atmospheric Sounding Interferometer-derived emissions) are presented and discussed in Sect. 3. The conclusions are provided in Sect. 4.

2 Methods

This section describes the process-based model for the N flow from agriculture within the ORCHIDEE LSM. The new module implemented here aims at calculating two types of emissions from agriculture: the manure management chain emissions (livestock housing and yard emissions as well as manure storage emissions) and soil emissions (accounting for fertilizer and manure application). The ORCHIDEE model framework is described in Sect. 2.1.1 followed by the different interactive components (shown in Fig. 1): the feeding of livestock, the whole manure management chain, fertilizer surface application, and the soil–plant–atmosphere continuum processes leading to soil emissions. Section 2.2 describes the setup of the simulations and the model evaluation protocol.

2.1 The ORCHIDEE LSM

2.1.1 General description

ORCHIDEE is a global-scale terrestrial ecosystem model coupling energy, water, and both the carbon and nitrogen cycles (Ciais et al., 2005; Krinner et al., 2005; Piao et al., 2007). Vegetation comprises 15 plant functional types (PFTs), among which two crop types (C₃ and C₄) and four grass types (temperate, boreal, and tropical C₃ grasses as well as a single C₄ class) are represented. The initial version used in this study includes a simple management of the crop biomass (which assumes that 45 % of the net primary productivity, NPP, is harvested) but no grassland management.

The main N processes within the soil–plant–atmosphere continuum are based on the OCN model (Zaehle and Friend, 2010; Zaehle et al., 2010). The representation of nitrification and denitrification processes are based on the DNDC (DeNitrification–DeComposition) model (Li et al., 1992; Li, 2000; Zhang et al., 2002). It accounts for ammonia/ammonium (${\mathrm{NH}}_{\mathrm{3}}/{\mathrm{NH}}_{\mathrm{4}}^{+}$), nitrate (${\mathrm{NO}}_{\mathrm{3}}^{-}$), nitrogen oxides (NO_x), and nitrous oxide (N₂O) pools as well as the related emissions. In addition to NH₃, NO_x, N₂O, and N₂ emissions, N is lost through runoff and leaching processes. The N inputs to soil mineral pools include atmospheric NO_y and NH_x deposition, biological nitrogen fixation (BNF), and the application of synthetic and organic fertilizers over agricultural lands. The total ammonia nitrogen (TAN) pool is also updated according to plant uptake, as described in Zaehle and Friend (2010). The version of ORCHIDEE used for this study is ORCHIDEEv3, revision 6863. It was part of the ensemble of terrestrial ecosystem models used for the 2019 Global Carbon Budget (Friedlingstein et al., 2019) and was recently evaluated by Seiler et al. (2022). Overall, ORCHIDEEv3 shows good agreement with observation-based data for carbon fluxes and vegetation state. Former revisions of ORCHIDEEv3 have also been used to quantify the global gross primary production (GPP) flux (Vuichard et al., 2019) and soil N₂O emissions (Tian et al., 2018). In the initial version, the organic fertilizer (i.e., manure) amount was prescribed annually (Zhang et al., 2017a), and the corresponding quantity of N was applied at a constant rate daily over the whole year. In addition, the emissions from the whole manure management were missing, and only soil emissions were taken into account. A description of the ORCHIDEE model, including the nitrogen cycle and its interaction with the carbon cycle, is detailed in Vuichard et al. (2019). At the global scale, the model evaluation shows good agreement between the gross primary production simulated with the carbon–nitrogen interaction version and the observational validation set.

In this paper, we integrate the following new developments within a new module called CAMEO for the Calculation of AMmonia Emissions in ORCHIDEE:

• a new grassland and cropland management module dedicated to livestock feeding (Sect. 2.1.3);

• a module computing manure production and the associated emissions from indoor farming livestock activities (housing, yard, manure storage) (Sect. 2.1.4);

• a new parametrization for agricultural N application onto croplands and grasslands (Sect. 2.1.5);

• an improved soil emission scheme based on a more realistic representation of the soil–plant–atmosphere continuum (Sect. 2.1.6).

Figure 1 Scheme of the agroecosystem representation developed within ORCHIDEE. It is composed of four main components describing feed management, manure management, surface N application, and the soil–plant–atmosphere continuum processes leading to the soil emissions. The livestock distributions and synthetic fertilizers are the main forcing files in the system and are represented in the dashed frames. The time steps (1 year, 1 d, and 30 min) of the processes are indicated using the arrows.

2.1.2 The agricultural N-flow module within ORCHIDEE: CAMEO

Livestock feed management

Both the livestock feeding (BM_ing, $\mathrm{kg}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) and bedding (BM_bedding, $\mathrm{kg}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) needs are calculated within each grid cell from livestock density distribution maps, for different livestock categories. The livestock types considered in our study are non-dairy cattle, dairy cattle, pigs, small ruminants, and chickens, which are the main contributors to global livestock NH₃ emissions. D_a, the distribution of each livestock category a, is taken from the Gridded Livestock of the World (GLW 2; Robinson et al., 2014) for the year 2006.

BM_ing for livestock category a is calculated as follows:

$$\begin{array}{}\text{(1)}& {\mathrm{BM}}_{\mathrm{ing}\left(\mathrm{a}\right)}=D_a\times \mathrm{SI}\times W_{\mathrm{a}},\end{array}$$

where W_a is the animal weight (kg), and SI is the specific intake (the intake per animal weight unit, $\mathrm{kg}\mathrm{C}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{yr}}^{-\mathrm{1}}$).

A daily dry matter intake equal to 2.5 % of the livestock weight (Paustian et al., 2006) is considered for every livestock category, and a factor of 0.45 is used to convert the dry matter into carbon matter (Paustian et al., 2006), leading to an SI value of 0.01 $\mathrm{kg}\mathrm{C}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{yr}}^{-\mathrm{1}}$. Regarding livestock weights, we use regional values adapted from the Supplement of FAO (2018), as listed in Table 1.

Table 1 Regional weights of the animal (kg) used in the intake needs calculation. Data have been adapted from FAO (2018). The regions listed in the table are as follows: NA (North America), RUS (Russian Federation), WE (western Europe), EE (eastern Europe), NENA (Near East and North Africa), ESEA (East and Southeast Asia), OCE (Oceania), SA (South Asia), LAC (Latin America and the Caribbean), and SSA (sub-Saharan Africa).

The livestock feeding and bedding needs are provided by a fraction of the crop and grass NPP which is harvested (called ${\mathrm{BM}}_{\mathrm{harv}/\mathrm{graz}}$). In order to quantify the amount of grassland and cropland biomass needed to feed each livestock category (BM_{ing,grass (a)} and BM_{ing,crop (a)}, respectively), we use the fractions of grass and crop that constitute the diet composition of each animal (d_grass (a) and d_crop (a), respectively). The ruminant animals (i.e., cattle and small ruminants) have a diet composed of a portion of grass and crop, whereas the non-ruminant animals (i.e., pigs and chickens) have a crop-only diet. The bedding needs are taken from crop residues only.

$$\begin{array}{}\text{(2)}& \begin{array}{rl}& {\mathrm{BM}}_{\mathrm{ing},\mathrm{grass}\left(\mathrm{a}\right)}={\mathrm{BM}}_{\mathrm{ing}\left(\mathrm{a}\right)}\times d_{\mathrm{grass}\left(\mathrm{a}\right)}\\ & {\mathrm{BM}}_{\mathrm{ing},\mathrm{crop}\left(\mathrm{a}\right)}={\mathrm{BM}}_{\mathrm{ing}\left(\mathrm{a}\right)}\times d_{\mathrm{crop}\left(\mathrm{a}\right)}\end{array}\end{array}$$

The diet composition $d_{\mathrm{crop}/\mathrm{grass}\left(\mathrm{a}\right)}$ is calculated from regional feeding information detailed in the Global Livestock Environmental Assessment Model (FAO, 2018) and is described in the Supplement. The bedding is estimated as follows EMEP/EEA (2019):

$$\begin{array}{}\text{(3)}& {\mathrm{BM}}_{\mathrm{bedding}\left(\mathrm{a}\right)}=D_a\times \mathrm{0.32}\times {\mathrm{Straw}}_{\mathrm{a}}.\end{array}$$

Here, Straw_a corresponds to the amount of straw ($\mathrm{kg}{\mathrm{head}}^{-\mathrm{1}}{\mathrm{yr}}^{-\mathrm{1}}$) used as bedding for each livestock type (Table 2), and the 0.32 factor corresponds to the C content of straw, assuming a C:N ratio of 80 (USDA, 2022) for the straw material and an N content of 4 g N kg⁻¹ (EMEP/EEA, 2019). This value is consistent with recent experimental studies by Su et al. (2020), who found 0.35 kg kg⁻¹ of total carbon in wheat straw.

BM_{bedding (a)} and BM_ing,crop(a) constitute the total demand for crops. We assume that the demand for crop biomass in each grid cell is satisfied by the amount of crop biomass harvested globally (global market). In contrast, the grass biomass needs are satisfied locally. Indeed, the grass biomass needs define the grassland management intensity through a grazing indicator (GI, unitless) that corresponds to the fraction of grass NPP for the year y that is harvested. Hence, the GI is defined as follows:

$$\begin{array}{}\text{(4)}& {\mathrm{GI}}_{\left(y\right)}=min\left({\displaystyle \frac{{\mathrm{BM}}_{\mathrm{ing}\mathrm{grass}}}{{\mathrm{NPP}}_{\mathrm{grass}\mathrm{above}\left(y\right)}}};{\max }_{\mathrm{above}}\right),\end{array}$$

where NPP_{grass above(y)} corresponds to the above NPP of grasslands at the grid cell level (kg C m⁻² per grid cell per year) and max_above, which is a parameter equal to 0.7 and defined as the maximum of the above biomass available for grazing/cutting.

NPP_{grass above(y)} is a function of the grassland NPP (kg C m⁻² per grassland per year) but also of the grassland area defined in each grid cell. Due to an inconsistency between the land use map and livestock density map, the targeted BM_ing grass value may not be reached by the use of the GI. To ensure that BM_ing grass demand is always satisfied, we adjust the diet composition of ruminants in some grid cells by increasing d_crop (a) as much as needed (and by reducing d_grass (a) by the same factor). The adjusted value of d_grass (a) is named d_{grass, adjusted (a)} and is depicted in the Supplement (Fig. S1). The GI is then applied to NPP_{grass above(y)} on a daily basis in order to obtain the total effective grazed biomass. For each animal, d_{grass, adjusted (a)} is used to deduce the effective crop biomass from the effective grazed biomass. Finally, each grid cell's effective crop biomass is constrained by the global crop harvested NPP. To do so, we compute the ratio between the global effective crop biomass and the global crop harvested NPP (HI) at a yearly time step. When HI >1, we impose the same constraint locally by dividing the effective crop biomass by HI.

As our methodology is based on an N-flow scheme, the C:N ratio imposed by the model for the crop and grass products is used to convert the carbon into the N biomass ingested (in $\mathrm{kg}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$). The grassland C:N ratio is unique for each grid cell and varies spatially from 23 to 62, whereas the cropland C:N ratio is fixed for the whole globe and is estimated to be ∼38. BM_{ing tot(y,a)} represents the total (including crop and grass products) N biomass ingested, which is used to compute the resulting manure emissions (described in the next section). Concerning the crop used as straw, a fixed C:N ratio of 80 is chosen (EMEP/EEA, 2019).

Indoor N flows and ammonia emissions

We adapt the scheme developed by Dämmgen and Hutchings (2008) which defines indoor ammonia emissions for each animal category. These pathways have also been used in the Tier 2 methodology of the manure management part of the EMEP/EEA Air Pollutant Emission Inventory Guidebook EMEP/EEA (2019). It is based on an N-flow model with mass transfers and emissions proportional to the TAN.

The main output of this module is the N emissions that occur during housing, yarding, and storage of animals, along with the resulting manure produced. The seasonal variability in indoor N emissions is neglected, and the emissions and the manure flow are calculated yearly.

Firstly, we compute the N biomass excreted by each animal category (m_excreted) based on the excretion rates estimated by Paustian et al. (2006) for the Intergovernmental Panel on Climate Change (IPCC) Tier 2 recommendations (see Table 2).

$$\begin{array}{}\text{(5)}& m_{\mathrm{excreted}\left(\mathrm{a}\right)}={\mathrm{BM}}_{\mathrm{ing}\mathrm{tot}(y,\mathrm{a})}\times (\mathrm{1}-N_{\mathrm{retained}\left(\mathrm{a}\right)})\end{array}$$

Secondly, we compute the manure excreted during the different livestock activities as a proportion of the year spent in housing, in the yard, and grazing, based on EMEP/EEA (2019). The fraction of time spent in the yard (x_yard, Table 2) is prescribed. The remaining time fraction is split into grazing (x_graz, Table 2) and housing periods.

$$\begin{array}{}\text{(6)}& \begin{array}{rl}& m_{\mathrm{N}(\mathrm{yard},\mathrm{a})}=x_{\mathrm{yard},\mathrm{a}}\times m_{\mathrm{excreted}\left(\mathrm{a}\right)}\\ & m_{\mathrm{N}(\mathrm{graz},\mathrm{a})}=x_{\mathrm{graz},\mathrm{a}}\times (\mathrm{1}-x_{\mathrm{yard},\mathrm{a}})\times m_{\mathrm{excreted}\left(\mathrm{a}\right)}\\ & m_{\mathrm{N}(\mathrm{house},\mathrm{a})}=(\mathrm{1}-x_{\mathrm{graz},\mathrm{a}})\times (\mathrm{1}-x_{\mathrm{yard},\mathrm{a}})\times m_{\mathrm{excreted}\left(\mathrm{a}\right)}\end{array}\end{array}$$

Default values of the TAN fraction contained in the excretal N (x_TAN,a) from the manure management part of the EMEP/EEA Air Pollutant Emission Inventory Guidebook 2019 (EMEP/EEA, 2019) (see Table 2) are used to calculate the amount of TAN produced during each activity, i (housing, yarding, and grazing):

$$\begin{array}{}\text{(7)}& m_{\mathrm{TAN}(i,\mathrm{a})}=x_{\mathrm{TAN},\mathrm{a}}\times m_{\mathrm{N}(i,\mathrm{a})}.\end{array}$$

Table 2 Default values for the fractions of the year spent grazing and in the yard, the proportion of TAN in the N mass retained, and the straw used as bedding. The information has been taken and adapted from EMEP/EEA (2019). The straw used as bedding for pigs is the average between different livestock types. The N retained is taken from Paustian et al. (2006).

m_{TAN (graz,a)} and m_N (graz,a) are used in Sect. 2.1.5 for N application on cultivated areas.

The ammonia emissions from housing $E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{a})}$ (in $\mathrm{kg}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) combine the volatilization from liquid and solid TAN masses, with specific emission factors: ${\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{liq},\mathrm{a})}$ (kilograms of NH₃-N per kilogram of TAN) and ${\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{sol},\mathrm{a})}$ (kilograms of NH₃-N per kilogram of TAN). ${\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{liq},\mathrm{a})}$ and ${\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{sol},\mathrm{a})}$ values are taken from Sommer et al. (2019) for each animal a except for the small ruminants category, which is assigned a value from EMEP/EEA (2019) (see Table 3). $E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{a})}$ is written as follows:

$$\begin{array}{}\text{(8)}& \begin{array}{rl}{E}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{a})}& =(x_{\mathrm{liq},\mathrm{a}}\times {\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{liq},\mathrm{a})}+(\mathrm{1}-x_{\mathrm{liq},\mathrm{a}})\\ & \times {\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{sol},\mathrm{a})})\times {\mathrm{m}}_{\mathrm{TAN}(\mathrm{house},\mathrm{a})},\end{array}\end{array}$$

where x_liq,a (unitless) is the proportion of manure handled as liquid for livestock type a, adapted from the Global Livestock Environmental Assessment Model (FAO, 2018) (see Supplement).

Emissions from yarding ($E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{yard},\mathrm{a})}$) are calculated from the mass excreted in the yard, and there is no distinction between liquid and solid handling.

$$\begin{array}{}\text{(9)}& E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{yard},\mathrm{a})}={\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{yard},\mathrm{a})}\times {\mathrm{m}}_{\mathrm{TAN}(\mathrm{yard},\mathrm{a})}\end{array}$$

Table 3 Emission factors (EFs) given as a percentage of the TAN content in the manure. EFs for the yard and the other N species come from EMEP/EEA (2019). The other EFs are taken from Sommer et al. (2019). There is no distinction between liquid and solid manure for yard EFs. Numbers in parentheses are the standard deviations given in Sommer et al. (2019) and used for the sensitivity analysis.

We compute the amounts of N and TAN that are stored as liquid and solid before application ($m_{\mathrm{N}(\mathrm{stor},\mathrm{type},\mathrm{a})}$ and $m_{\mathrm{TAN}(\mathrm{stor},\mathrm{type},\mathrm{a})}$ for type=liq,sol, respectively; ${\mathrm{kgNm}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$; Eq. 10). For storage, we assume that all of the manure from housing and yarding is stored, except the N lost by ammonia emissions in the house and yard ($E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{a})}$ and $E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{yard},\mathrm{a})}$, respectively). Manure from yarding is considered liquid and goes in the liquid manure storage. Concerning the liquid storage (i.e., slurries), a fraction f_min of the organic N (N-TAN) is converted into TAN through mineralization. A value of 0.1 is used for f_min (Dämmgen and Hutchings, 2008; EMEP/EEA, 2019). For solid storage, we account for an additional N source from bedding ($m_{\mathrm{bed},\mathrm{N},\mathrm{a}}$). Incorporation of bedding in the manure storage induces an immobilization of TAN in the organic matter when manure is handled as straw-based solid manure, at a rate f_imm proportional to $m_{\mathrm{bed},\mathrm{N},\mathrm{a}}$. An f_imm value of 0.0067 kg kg⁻¹ is used (Kirchmann and Witter, 1989; Webb and Misselbrook, 2004; EMEP/EEA, 2019). This immobilization greatly reduces the resulting NH₃ emissions.

$$\begin{array}{}\text{(10)}& \begin{array}{rl}& m_{\mathrm{N}(\mathrm{stor},\mathrm{liq},\mathrm{a})}=(m_{\mathrm{N}(\mathrm{house},\mathrm{liq},\mathrm{a})}-E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{liq},\mathrm{a})})\\ & +(m_{\mathrm{N}(\mathrm{yard},\mathrm{a})}-E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{yard},\mathrm{a})});\\ & m_{\mathrm{TAN}(\mathrm{stor},\mathrm{liq},\mathrm{a})}=(m_{\mathrm{TAN}(\mathrm{house},\mathrm{liq},\mathrm{a})}-E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{liq},\mathrm{a})})\\ & \times (\mathrm{1}-f_{\min })+\\ & m_{\mathrm{TAN}(\mathrm{stor},\mathrm{liq},\mathrm{a})}\times f_{\min }+(m_{\mathrm{TAN}(\mathrm{yard},\mathrm{a})}-E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{yard},\mathrm{a})});\\ & m_{\mathrm{N}(\mathrm{stor},\mathrm{sol},\mathrm{a})}=m_{\mathrm{N}(\mathrm{house},\mathrm{sol},\mathrm{a})}-E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{sol},\mathrm{a})}\\ & +m_{\mathrm{bed},\mathrm{N},\mathrm{a}};\\ & m_{\mathrm{TAN}(\mathrm{stor},\mathrm{sol},\mathrm{a})}=m_{\mathrm{TAN}(\mathrm{house},\mathrm{sol},\mathrm{a})}-E_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{house},\mathrm{sol},\mathrm{a})}\\ & -m_{\mathrm{bed},\mathrm{a}}\times f_{\mathrm{imm}}.\end{array}\end{array}$$

Here, m_bed,N is the N mass of bedding ($\mathrm{kg}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$), and m_bed is the dry matter mass of bedding.

Manure from storage is supposed to be entirely used as fertilizer. The quantities m_{TAN (applic)} and m_N (applic) are the respective TAN and N manures that will be applied to the surface, as described in Sect. 2.1.5. They are obtained by removing the total N emissions from the stored manure (Eq. 11).

$$\begin{array}{}\text{(11)}& \begin{array}{rl}\text{liquid}& \left\{\begin{array}{l}{m}_{\mathrm{TAN}(\mathrm{applic},\mathrm{liq},\mathrm{a})}=m_{\mathrm{TAN}(\mathrm{stor},\mathrm{liq},\mathrm{a})}-E_{\mathrm{stor}(\mathrm{liq},\mathrm{a})}\\ m_{\mathrm{N}(\mathrm{applic},\mathrm{liq},\mathrm{a})}=m_{\mathrm{N}(\mathrm{stor},\mathrm{liq},\mathrm{a})}-E_{\mathrm{stor}(\mathrm{liq},\mathrm{a})}\end{array}\right.\\ \text{solid}& \left\{\begin{array}{l}{m}_{\mathrm{TAN}(\mathrm{applic},\mathrm{sol},\mathrm{a})}=m_{\mathrm{TAN}(\mathrm{stor},\mathrm{sol},\mathrm{a})}-E_{\mathrm{stor}(\mathrm{sol},\mathrm{a})}\\ m_{\mathrm{N}(\mathrm{applic},\mathrm{sol},\mathrm{a})}=m_{\mathrm{N}(\mathrm{stor},\mathrm{sol},\mathrm{a})}-E_{\mathrm{stor}(\mathrm{sol},\mathrm{a})}\end{array}\right.\end{array}\end{array}$$

In addition to emissions of NH₃, other N species (N₂O, NO, and N₂) emissions can occur from storage and, thus, are required to calculate the final manure mass from storage. These emissions are obtained using the EFs listed in Table 3 as follows:

$$\begin{array}{}\text{(12)}& \begin{array}{rl}{E}_{\mathrm{stor},\mathrm{liq},\mathrm{a}}& =m_{\mathrm{TAN}(\mathrm{stor},\mathrm{liq},\mathrm{a})}\times \\ & ({\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{stor},\mathrm{liq},\mathrm{a})}+{\mathrm{EF}}_{{\mathrm{N}}_{\mathrm{2}}\mathrm{O}(\mathrm{store},\mathrm{liq},\mathrm{a})}\\ & +{\mathrm{EF}}_{\mathrm{NO}(\mathrm{stor},\mathrm{liq},\mathrm{a})}+{\mathrm{EF}}_{{\mathrm{N}}_{\mathrm{2}}(\mathrm{stor},\mathrm{liq},\mathrm{a})});\\ E_{\mathrm{stor},\mathrm{sol},\mathrm{a}}& =m_{\mathrm{TAN}(\mathrm{stor},\mathrm{sol},\mathrm{a})}\times \\ & ({\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}(\mathrm{stor},\mathrm{sol},\mathrm{a})}+{\mathrm{EF}}_{{\mathrm{N}}_{\mathrm{2}}\mathrm{O}(\mathrm{stor},\mathrm{sol},\mathrm{a})}\\ & +{\mathrm{EF}}_{\mathrm{NO}(\mathrm{stor},\mathrm{sol},\mathrm{a})}+{\mathrm{EF}}_{{\mathrm{N}}_{\mathrm{2}}(\mathrm{stor},\mathrm{sol},\mathrm{a})}).\end{array}\end{array}$$

The remaining manure after storage m_(applic,a) and the manure produced during grazing m_(graz,a) are the main output of this specific module. Both quantities are the input for the surface application component of the model (described in the next section).

Robinson et al. (2014)FAOSTAT (2020)FAO (2018)FAO (2018)EMEP/EEA (2019)Paustian et al. (2006)EMEP/EEA (2019)EMEP/EEA (2019)FAO (2018)EMEP/EEA (2019)Sommer et al. (2019)

Table 4 Summary of the data sources used in Sects. 2.1.3 and 2.1.4 for the calculation of the indoor emissions. All of the data are used for each livestock type a except for the D_{dairy cattle,i} variable.

Organic application onto land

This section contains the description of the manure application to soil. m_(applic,a) and m_(graz,a) are the manure remaining after storage and the manure produced during grazing, respectively (a description is given in Sect. 2.1.4). Both are yearly stocks applied daily at a constant rate and during a specific period, driven mainly by environmental conditions, as described below. This assumption may neglect the actual seasonal patterns of N application usually defined by local governance in some regions. For instance, as discussed in Van Damme et al. (2022), the time of the year when fertilizers can be applied in Europe is strongly dependent on local regulations. Synthetic fertilizers are also considered in our representation and follow the same temporal distribution as the manure from storage. Note that ORCHIDEE does not differentiate between natural and managed grassland. Only grid cells with either the presence of livestock or fertilizer application are considered for the emission calculation, so we can assume that the pixel is managed.

As previously stated, the following two types of manure are considered:

• Manure from storage and applied to soil as fertilizer. The manure from storage is applied daily at a constant rate for 6 months from the beginning of vegetation growth, corresponding to the first leaf development depending on the PFT. The intermediate period of application (L_application=6 months) has been chosen in order to account for the heterogeneity in the agricultural practices, as the model only represents C₃ and C₄ crop types within the grid cell. Moreover, there is a lack of information in the literature about N application onto grassland at the global scale. We assume that cropland and grassland PFTs receive stored manure with a 2 times higher preference for cropland fractions.

• Manure deposited during grazing activity by the ruminants. The manure from grazing activity m_(graz,a) is calculated in Eq. (7) and is assumed to be only deposited on grassland PFTs by ruminants. The first day of manure deposition for grazing also corresponds to the beginning of vegetation growth. The amount and period of manure deposited during grazing are animal-specific and are determined by the fraction of time passed grazing (x_graz,a).

The soil–plant–atmosphere processes leading to soil emissions

In this section, we describe the physical processes in the soil that influence ammonia emissions. A single soil TAN pool (TAN_(soil), g N m⁻²) is considered. The soil TAN pool is dynamically updated depending on the processes implemented in the model. These processes are described in Zaehle and Friend (2010). The processes corresponding to the creation of ${\mathrm{NH}}_{\mathrm{4}}^{+}$ are related to mineralization, N application, and NH_x deposition, whereas the losses include nitrification, leaching, and volatilization.

TAN_(soil,aq) corresponds to the ammonium pool TAN_(soil), which is assumed to be diluted in the soil water at a different heights in the soil according to the z_activity parameter.

The z_activity parameter is regulated by all TAN sources, called “input” (“min” – mineralization, “dep” – deposition, BNF, “fert” – mineral fertilizer, and “manure” – applied manure and manure from grazing), in soil as follows:

$$\begin{array}{}\text{(13)}& \begin{array}{rl}{z}_{\mathrm{activity}}\left(t\right)=& (p_{\mathrm{zact}\mathrm{\_}\mathrm{deep}}\\ & \times {\mathrm{input}}_{\min }+p_{\mathrm{zact}\mathrm{\_}\mathrm{deep}}\times {\mathrm{input}}_{\mathrm{dep}}\\ & +p_{\mathrm{zact}\mathrm{\_}\mathrm{deep}}\times {\mathrm{input}}_{\mathrm{bnf}}+\\ & p_{\mathrm{zact}\mathrm{\_}\mathrm{surf}}\left(\mathrm{fert}\right)\times {\mathrm{input}}_{\mathrm{fert}}\\ & +p_{\mathrm{zact}\mathrm{\_}\mathrm{surf}}\left(\mathrm{manure}\right)\\ & \times {\mathrm{input}}_{\mathrm{manure}}+z_{\mathrm{activity}}(t-\mathrm{1})\\ & \times {\mathrm{TAN}}_{\left(\mathrm{soil}\right)})\\ & \times {\displaystyle \frac{\mathrm{1}}{{\mathrm{input}}_{\mathrm{tot}}+{\mathrm{TAN}}_{\left(\mathrm{soil}\right)}}},\end{array}\end{array}$$

where input_tot denotes the total TAN sources in soil. We assume that the fertilization and the application of manure are surface N additions to soil, whereas the other sources of TAN (mineralization, deposition, and BNF) are deeply added into soil (p_{zact_deep}=1.0 m). It is worth noting that input_fert is given as the ammonium content of the total N mineral fertilizer applied (the parameter ${\mathrm{Frac}}_{{\mathrm{NH}}_{\mathrm{4}}^{+},\mathrm{fert}}$ is the fraction of the ammonium content of the N fertilizer used to make the conversion, and this parameter is tested in the sensitivity analysis).

p_{zact_surf} is obtained as described in Riddick et al. (2016):

$$\begin{array}{}\text{(14)}& \begin{array}{rl}& p_{\mathrm{zact}\mathrm{\_}\mathrm{surf}\left(\mathrm{manure}\right)}=s_{\mathrm{W}}\left(m\right)\times {\mathrm{input}}_{\mathrm{manure}}/\mathrm{SWC},\\ & p_{\mathrm{zact}\mathrm{\_}\mathrm{surf}\left(\mathrm{fert}\right)}=s_{\mathrm{W}}\left(f\right)\times {\mathrm{input}}_{\mathrm{fert}}/\mathrm{SWC}.\end{array}\end{array}$$

Here, SWC is the soil water content computed by ORCHIDEE, s_W(m) is the specific water volume of manure ($\mathrm{5.67}\times {\mathrm{10}}^{-\mathrm{4}}$ m³ of water per gram of nitrogen; Sommer and Hutchings, 2001; Riddick et al., 2016), and s_W(f) is the specific water volume of synthetic fertilizers. s_W(f) depends on the soil temperature T_g and is given by United Nations Industrial Development Organization (UNIDO) (1988):

$$\begin{array}{}\text{(15)}& S_{\mathrm{W}}\left(f\right)={\displaystyle \frac{\mathrm{1}\times {\mathrm{10}}^{-\mathrm{6}}}{\mathrm{0.466}\times \mathrm{0.66}\times e^{\mathrm{0.0239}\times \left(T_{\mathrm{g}}-\mathrm{273}\right)}}}.\end{array}$$

The emissions of NH₃ ($E_{{\mathrm{NH}}_{\mathrm{3}}}$, $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$) are obtained following the resistive scheme used in the FAN model (Riddick et al., 2016; Vira et al., 2020):

$$\begin{array}{}\text{(16)}& E_{{\mathrm{NH}}_{\mathrm{3}}}={\displaystyle \frac{{\mathrm{NH}}_{\mathrm{3}}\left(g\right)-{\mathit{\chi }}_{\mathrm{a}}}{R_{\mathrm{a}}\left(z\right)+R_{\mathrm{b}}}},\end{array}$$

where NH₃(g) is the NH₃ concentration at the surface (g N m⁻³), χ_a is the free-atmosphere concentration (g N m⁻³), R_a(z) is the aerodynamical resistance (s m⁻¹), and R_b is the quasi-boundary-layer resistance (s m⁻¹).

χ_a is prescribed as a monthly field averaged over 11 years from a run of the global LMDZ-INCA (Laboratoire de Météorologie Dynamique–INteraction with Chemistry and Aerosols) model at a $\mathrm{2.5}{}^{\circ }\times \mathrm{1.3}{}^{\circ }$ resolution (39 vertical levels) over the 2005–2015 period (Hauglustaine et al., 2014). The spatial distribution of χ_a is presented in Fig. S4 (Supplement) for both May and December (2005–2015 climatology).

R_a(z) is computed interactively by the biophysics module of the ORCHIDEE model. R_b has been implemented according to Xu et al. (2019) as follows:

$$\begin{array}{}\text{(17)}& R_{\mathrm{b}}={\displaystyle \frac{v}{D_{{\mathrm{NH}}_{\mathrm{3}}}}}\times {\left[{\displaystyle \frac{c}{(\mathrm{LAI}{)}^{\mathrm{2}}}}\times \left({\displaystyle \frac{l\times {\mathit{\mu }}_{*}}{v}}\right)\right]}^{\mathrm{1}/\mathrm{3}},\end{array}$$

where $D_{{\mathrm{NH}}_{\mathrm{3}}}$ is the molecular diffusivity of NH₃ in air (m² s⁻¹; Massman, 1998), c is an empirical constant equal to 3, l is the leaf width (0.02 m; Massad et al., 2010), v is the kinematic viscosity of air (1.56×10⁻⁵ m² s⁻¹ at 25 ^∘C), T is the air temperature (in K), and LAI is the leaf area index (m² m⁻²), which is computed by the ORCHIDEE model. The resulting annual mean R_b ranges between 0 and 1.14 s m⁻¹ over the globe. $D_{{\mathrm{NH}}_{\mathrm{3}}}$ is a function of temperature and is written as follows:

$$\begin{array}{}\text{(18)}& D_{{\mathrm{NH}}_{\mathrm{3}}}=\mathrm{0.1978}\times {\left({\displaystyle \frac{T}{\mathrm{273.13}}}\right)}^{\mathrm{1.81}}\times {\mathrm{10}}^{-\mathrm{4}}.\end{array}$$

Henry's law coefficient (K_H) and the dissociation constant of ${\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)$ in water ($K_{{\mathrm{NH}}_{\mathrm{4}}}$) (Sutton et al., 1994) are used for the speciation between the different TAN species (${\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right),{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right),{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\left(\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{3}}\right)$).

$$\begin{array}{}\text{(19)}& {\displaystyle }{K}_{\mathrm{H}}={\displaystyle \frac{\left[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)\right]}{\left[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)\right]}}\text{(20)}& {\displaystyle }{K}_{{\mathrm{NH}}_{\mathrm{4}}}={\displaystyle \frac{\left[{\mathrm{H}}^{+}\right]\left[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)\right]}{\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]}}\end{array}$$

By combining Eqs. (19) and (20), we can compute the gaseous phase of ammonia NH₃(g), which is the fraction that will be volatilized. TAN_(soil,aq) corresponds to the aqueous phase of TAN in the soil, which is modulated by the height of the soil through the z_activity parameter.

$$\begin{array}{}\text{(21)}& {\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)={\displaystyle \frac{{\mathrm{TAN}}_{(\mathrm{soil},\mathrm{aq})}}{\frac{\mathit{\theta }}{{\mathrm{K}}_{\mathrm{fact}}}+\mathit{ϵ}}},\end{array}$$

where θ is the volumetric soil water content (in cubic meters of water per cubic meter of soil), and ϵ is the fraction of the air-filled soil volume computed by the ORCHIDEE model.

K_fact is calculated as

$$\begin{array}{}\text{(22)}& {\mathrm{K}}_{\mathrm{fact}}=\mathrm{1}/(\mathrm{1}+K_{\mathrm{H}}+K_{\mathrm{H}}\left[{\mathrm{H}}^{+}\right]/K_{{\mathrm{NH}}_{\mathrm{4}}}),\end{array}$$

and K_H, Henry's law constant for NH₍₃₎, depends on the surface temperature T_g as follows:

$$\begin{array}{}\text{(23)}& K_{\mathrm{H}}=\mathrm{H}\times T_{\mathrm{g}}\times e^{\mathrm{4092}\left(\mathrm{1}/T_{\mathrm{g}}-\mathrm{1}/T_{\mathrm{ref}}\right)},\end{array}$$

where T_ref is the reference temperature (298.15 K), and H is a conversion factor equal to 4.905. We use the value of 0.59 $\mathrm{mol}{\mathrm{m}}^{-\mathrm{3}}{\mathrm{Pa}}^{-\mathrm{1}}$ described in Sander (2015) by which the perfect gas constant has been multiplied in order to get a unitless constant. $K_{{\mathrm{NH}}_{\mathrm{4}}}$ is the dissociation equilibrium and also depends on the surface temperature T_g as follows:

$$\begin{array}{}\text{(24)}& K_{{\mathrm{NH}}_{\mathrm{4}}}=\mathrm{5.67}\times {\mathrm{10}}^{-\mathrm{10}}\times e^{-\mathrm{6286}\left(\mathrm{1}/T_{\mathrm{g}}-\mathrm{1}/T_{\mathrm{ref}}\right)}.\end{array}$$

The hydrogen ion concentration [H⁺] is assumed to be constant and equal to 10⁻⁷, which approximately corresponds to the pH given in Massad et al. (2010) for cattle manure, diammonium phosphate fertilizers in acidic soils, and ammonium nitrate fertilizers. A pH of 7 is also adopted in Riddick et al. (2016). In our simulations, the pH does not impact the surrounding soil pH in the model, in contrast to Vira et al. (2020), where the pH varies according to different TAN age classes.

2.2 Modeling setup

The ORCHIDEE model, including all of the developments described in Sect. 2, was run at a spatial resolution of 2^∘ (180×90). This spatial resolution is relatively low but enables one to perform an ensemble of sensitivity tests at a reasonable computing cost. We also performed a reference simulation at a 0.5^∘ resolution to ensure that the model resolution does not affect the results. We performed a 10-year reference simulation over the 2005–2015 period. This simulation starts in January 2005 from a simulation done with an ORCHIDEE version similar to the one presented in this paper but without the developments presented in Sect. 2.1.2. In the reference simulation, all annual forcing data are updated every year except those related to BNF and the livestock density, which are constant over time. A set of nine sensitivity test simulations characterized by specific changes in the parametrization were conducted to evaluate the impact of parameter uncertainty on agricultural ammonia emissions. The parameters that have been tested are the atmospheric ammonia concentration (χ_a), the pH of the manure (pH, default value of 7), the timing period of N application (L_application, default value of 183 d), the emission factor for housing and storage activities (${\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{indoor}\right)}$), the fraction of ammonium in fertilizer (${\mathrm{Frac}}_{{\mathrm{NH}}_{\mathrm{4}}^{+},\mathrm{fert}}$, default value of 0.6), and the N deep processes regulation parameter (p_{zact_deep}, default value of 1 m). Table 5 summarizes the set of simulations and the key parameters tested.

In Riddick et al. (2016), the value of χ_a was set to 0.3 µg N m⁻³, as it is representative of the concentration over low-activity agricultural sites (Zbieranowski and Aherne, 2012). Little sensitivity of the emissions to this parameter was found because χ_a is much smaller than NH₃(g). However, this parameter has been tested in our implementation through a sensitivity analysis.

Table 5 Summary of the simulations performed and the parameters tested. The ${\mathrm{EF}}_{{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{indoor}\right)}$ parameter is described in Sect. 2.1.4 for housing and storage emissions. The equations and sections related to the individual parameters are also indicated.

The ORCHIDEE model requires the following set of forcing data:

• meteorological data, including the near-surface air temperature and specific humidity, wind speed, pressure, short- and long-wave incoming radiation, rainfall, and snowfall, from the Climatic Research Unit (CRU) and Japanese reanalysis (JRA) dataset (CRU-JRA V2.1) (Harris et al., 2014) (preprocessed and adapted by Vladislav Bastrikov, LSCE, July 2020), provided at 6 h time steps;

• the global average annual atmospheric CO₂ concentration, which is provided by TRENDY (Le Quéré et al., 2018);

• the global annual land cover distribution, based on combined information from the Land-Use Harmonization 2 (LUH2v2) dataset at a 0.25^∘ resolution (Hurtt et al., 2020) and the ESA (2022) CCI Land Cover (see Lurton et al., 2020, for more details);

• atmospheric N deposition fluxes (NH_x and NO_y) from the IGAC/SPARC Chemistry-Climate Model Initiative (CCMI; Eyring et al., 2013), which have been used in the N₂O Model Intercomparison Project (NMIP) project (Tian et al., 2018), corresponding to information at a 0.5^∘ spatial resolution and a monthly temporal resolution;

• the annual mineral fertilizer rates over croplands and grasslands from an annual dataset developed by Lu and Tian (2017), which corresponds to a reconstruction from 1960 to 2014 for the global cropland matched with the HYDE 3.2 cropland distribution;

• the biological nitrogen fixation rate, which is provided as climatological data as a function of the evapotranspiration flux (see Vuichard et al., 2019, for more details);

• the distribution of each livestock category, which is taken from the Gridded Livestock of the World (GLW 2; Robinson et al., 2014) and represents the gridded animal densities (D_a) for the year 2006 at a 1 km resolution. Small-ruminant densities correspond to the sum of the sheep and goat densities. The dairy cattle distribution has been retrieved from the total cattle distribution combined with national dairy cattle densities given by FAOSTAT (2020). The calculation adopted is described in the Supplement.

2.3 Model evaluation dataset

Our integrated approach allows the computation of different variable levels before the final emission results, such as biomass productivity, animal excretion rate, and manure production. This set of variables offers the advantage of evaluating our emissions at different stages of the N flow against the previous works listed in Table 6.

Bouwman et al. (2013a)Billen et al. (2014); Bodirsky et al. (2014)Conijn et al. (2018); Uwizeye et al. (2020)Billen et al. (2014); Bodirsky et al. (2014)Conijn et al. (2018); Uwizeye et al. (2020)Crippa et al. (2018); Vira et al. (2020)Beusen et al. (2008); Potter et al. (2010)Bouwman et al. (2013b); Billen et al. (2014)Bodirsky et al. (2014); Zhang et al. (2017a)Conijn et al. (2018); Vira et al. (2020)Uwizeye et al. (2020)Bouwman et al. (2013b); Billen et al. (2014)Bodirsky et al. (2014); Zhang et al. (2017a)Conijn et al. (2018); Vira et al. (2020)Uwizeye et al. (2020)McDuffie et al. (2020); Vira et al. (2020)Evangeliou et al. (2020)

Table 6 Summary of the different simulated variables evaluated in this work by comparison with previous studies.

The decadal mean (2005–2015) values of global and regional calculated agricultural emissions, including indoor and soil emissions, are compared to the CEDS inventory (McDuffie et al., 2020) and the emissions simulated by FANv2 (Vira et al., 2020). From the LMDZ–ORCHIDEE–INCA coupling development perspective, it is interesting to compare our approach with CEDS (McDuffie et al., 2020), as it is a reference dataset offering a long period of data (1750–2019). FANv2 has been chosen for our evaluation because our work is based on a similar approach. The regional budget accounts for Africa; tropical southern Asia; Europe; China, Korea, and Japan (abbreviated as China–K–J in the figures); Oceania; India; the USA and Canada; and Latin America. The seasonal variations in ammonia emissions are also evaluated against satellite-derived emissions (Evangeliou et al., 2020). For that purpose, atmospheric NH₃ columns observed by the Infrared Atmospheric Sounding Interferometer (IASI) satellite have been combined with the NH₃ lifetime calculated by LMDZ-INCA in order to retrieve emissions. The NH₃ retrieval product used to derive emissions in our study comprises the 2011–2015 morning observations (Metop-A and -B) and follows a neural network retrieval approach (ANNI-NH3-v3R), as referred to in Van Damme et al. (2017, 2021). Both the lifetime and atmospheric columns are monthly products and share the exact same grid resolution (LMDZ-INCA grid resolution at $\mathrm{2.5}{}^{\circ }\times \mathrm{1.3}{}^{\circ }$). All three CAMEO, CEDS, and FANv2 seasonal variations are evaluated against IASI-derived emissions (defined as IASI^inv). In order to be consistent with IASI observations (where no source distinction is possible), the CAMEO, CEDS, and FANv2 agricultural emissions need to be complemented by fire emission data taken from van der Werf et al. (2017) (Global Fire Emissions Database, GFED s4) as well as by industrial and waste sources (McDuffie et al., 2020). The extended emissions are referred to as CAMEO₊, FAN₊, and CEDS₊, respectively, as described in Table 7. It is important to note that only the ORCHIDEE model can provide natural emissions; this source is not considered in the FAN₊ nor the CEDS₊ datasets.

(McDuffie et al., 2020)(van der Werf et al., 2017)(Vira et al., 2020)(McDuffie et al., 2020)(van der Werf et al., 2017)(McDuffie et al., 2020)(McDuffie et al., 2020)(van der Werf et al., 2017)

Table 7 Summary of the different datasets used in the comparison with the IASI-derived emissions. All of the emission sets (except FANv2 data, which are a 2010–2015 climatology) are taken from the 2011–2015 period and have been gridded onto the LMDZ-INCA default resolution of 144×142 pixels.

3 Results and evaluation

3.1 Evaluation of intermediate variables

Using the ORCHIDEE model, we estimate the total biomass produced, including grass and crop, to be about 103 Tg N yr⁻¹, which is slightly lower than previous estimates (110–152 Tg N yr⁻¹ estimated by Bouwman et al., 2013a; Billen et al., 2014; Bodirsky et al., 2014; Conijn et al., 2018; Uwizeye et al., 2020). The calculated global annual crop production (expressed in N) is about 74 Tg N yr⁻¹ and compares well with the 72 and 74 Tg N yr⁻¹ estimated by Billen et al. (2014) and Zhang et al. (2021), respectively. The calculated annual grass production (28.8 Tg N yr⁻¹) is more than 3 times lower than the 80.3 Tg N yr⁻¹ reported by Billen et al. (2014) and estimated from the difference between livestock ingestion and available feed resources. By doing so, uncertainties from several components (crop production, net import of vegetal proteins, and human consumption of vegetal proteins) are accumulated. Our resulting total biomass ingested by the livestock (88 Tg N yr⁻¹) is lower than the range found in the literature (122–167 Tg N yr⁻¹; Billen et al., 2014; Bodirsky et al., 2014; Conijn et al., 2018; Uwizeye et al., 2020) which can be attributed to the low grassland production calculated in our model.

However, if the grass N production is largely underestimated by ORCHIDEE, our grass C production estimate of 1.2 Pg C yr⁻¹ is close to the value of 1.95 Pg C yr⁻¹ reported in the Fifth Assessment Report (AR5) of the IPCC (Smith et al., 2014). In this respect, an overestimation of the C:N ratio may also explain part of the grass N production underestimation. It can be explained using a unique grass C:N ratio per pixel and a global C:N ratio for crop in the model. Because we assume that all of the manure stored is then applied to soil, for the evaluation phase, we only consider data from the literature that estimate the application rate of manure to soil. The global annual amount of manure application (66 Tg N yr⁻¹) is lower than the range of 99–129 Tg N yr⁻¹ estimated in recent studies (Beusen et al., 2008; Potter et al., 2010; Bouwman et al., 2013b; Billen et al., 2014; Zhang et al., 2017a; Conijn et al., 2018; Vira et al., 2020; Uwizeye et al., 2020).

Figure 2 compares the distribution of the N manure applied with the values retrieved by Zhang et al. (2017a) for the year 2006 in order to be consistent with the reference livestock distribution used in our approach. The spatial distribution of the manure application highlights the main livestock-raising regions, such as China, India, Europe, Latin America, and the USA. It shows good consistency with Zhang et al. (2017a), although it is higher in India, the USA, and Latin America and lower in China and Europe. These differences can be explained by the fact that we use different regional animal weights per livestock category (Table 1), instead of a fixed value as in Zhang et al. (2017a). This induces different N demands for similar livestock and, ultimately, different quantities of applied manure. Indeed, we use recently adapted data for animal weights from FAO (2018), whereas Zhang et al. (2017a) used IPCC guidelines (Tier1 IPCC, 2006; Paustian et al., 2006). For instance, for India, the non-dairy cattle weight is almost 4 times higher, which explains the differences observed in our calculation. Moreover, our study assumes a unique N excretion rate per livestock type and no livestock system distinctions as a simplification.

Figure 2 The manure application ($\mathrm{kg}\mathrm{N}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) (a) simulated by CAMEO and averaged over the 2005–2015 period and (b) calculated by Zhang et al. (2017a) for 2006.

Table 8 Global estimates of intermediate variables computed by the model in this study for the 2005–2015 period and the range of previous estimates.

3.2 Agricultural emissions at the global scale

We estimate global NH₃ agricultural emissions (averaged over 2005–2015) of about 44 Tg N yr⁻¹: 78 % from soil volatilization (driven by fertilizer and manure application) and the remainder from indoor emissions (from livestock housing, yarding, and storage). These global NH₃ emissions are within the range given by McDuffie et al. (2020) and Vira et al. (2020) of 39–47 Tg N yr⁻¹ (Fig. 3a).

Figure 3 Averaged global (a) and regional (b) NH₃ emissions (Tg N yr⁻¹) from indoor activities (“Indoor”, darker bars) and agricultural soil (“Agri. soils”, lighter bars) computed by CAMEO over the 2005–2015 period (“CAMEO”, blue bars) and FANv2 (Vira et al., 2020) over the 2010–2015 period (“FAN”, green bars). Total agricultural emissions (accounting for manure management and soil volatilization) estimated in the CEDS inventory (2005–2014 average) (McDuffie et al., 2020) are represented in purple (“CEDS-Agri. total”). China–K–J accounts for China, Korea, and Japan.

China, India, Africa, Latin America, the USA, and Europe appear as the main contributors to global NH₃ emissions, accounting for 80 % of the total budget (Fig. 3b). Most of these source areas, which have also been identified as agricultural regions by Van Damme et al. (2018), are regions with intensive crop cultivation (Fig. S2) and important livestock activities, inducing high N application rates (Fig. 2). The spatial distributions of the calculated agricultural NH₃ emissions show good agreement with the FANv2 and CEDS results (Fig. 4, b, c). In India and China, our emissions are slightly higher than FANv2 and CEDS estimates. Our values are lower than the FANv2 estimate in Latin America and Africa but high compared with CEDS emissions, which are particularly low in these two regions. In some parts of Africa and Latin America, where the use of synthetic fertilizer is low (never exceeding 2500 $\mathrm{kg}\mathrm{N}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$), livestock activity appears to be the main contributor to the emissions.

In intensive agricultural regions, data used for mineral fertilizer application rates can be a source of discrepancy between models. Vira et al. (2020) use the LUH2 dataset (Hurtt et al., 2020), which assumes that only croplands are fertilized. The amounts of fertilizer applied over croplands are comparable globally between Vira et al. (2020) and our study (respective minimum–maximum of 79–87 and 96–101 Tg N yr⁻¹ over 2010–2015) but differ in some regions (see Fig. S2a). In addition, in our study, grasslands are also fertilized with a global amount of 25.7 Tg Nyr⁻¹. This leads to differences in the simulated soil emissions, more specifically in India, the USA, and China, where grasslands are highly fertilized (Fig. S2b) and can be translated into high volatilization rates when compared with FANv2.

Figure 4 Simulated ammonia emissions ($\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) from total agricultural sources computed (a) by CAMEO (2005–2015 average), (b) by the FANv2 model (2010–2015 average) (Vira et al., 2020), and (c) from the CEDS inventory (2005–2015 average) (McDuffie et al., 2020).

NH₃ emissions peak in June–July–August for most regions (the USA, Europe, China, and Africa; Fig. 5) with maximum values reaching 16.4 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$ in eastern China. The peak in India appears earlier in spring, whereas two peaks occur in Latin America: one during December–January–February and one during September–October–November. Depending on the region, the seasonality of the emissions varies according to different factors, including environmental parameters and agricultural practices. This aspect will be analyzed in more detail in Sect. 3.4 and 3.5.

Figure 5 Regional seasonal agricultural ammonia emissions averaged over the 2005–2015 period ($\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) simulated by CAMEO. Boxes delimit the regions used in the analysis comparing CAMEO emissions with IASI^inv (Sect. 3.5). Note that the scales are different for each region.

The spatial pattern of the simulated indoor NH₃ emissions (Fig. 6b) is similar to that of manure application rates, with both being driven mainly by livestock density. Hot spot regions of indoor emissions are located in eastern China, eastern India, and northern Europe, with maximum values reaching up to 1.7 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$. The major sources of volatilization from soils are located in India, eastern China, and the USA, with a maximum value of 12 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$. The difference in spatial patterns between the two source categories is mainly due to the fact that soil emissions not only depend on livestock distribution (indoor emissions) but also on environmental conditions and mineral fertilizer application rates. The sensitivity of modeled emissions to some of these factors is presented in Sect. 3.4.

Figure 6 Simulated ammonia emissions averaged over the 2005–2015 period ($\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) from agricultural soil (fertilizer and manure application) (a) and manure management (b).

The agricultural NH₃ emissions from manure management (qualified as “indoor”) are poorly quantified at the global scale. While Vira et al. (2020) reported this emission source to be 18 Tg N yr⁻¹ for 2010, our estimate is half this value (9.6 Tg N yr⁻¹) but is in good agreement with the NH₃ emissions reported by Crippa et al. (2018) and Beusen et al. (2008) (9 Tg N yr⁻¹ for the year 2010 and 2000, respectively). Biomass excreted in our model is 40 % lower than what is produced in FANv2, which can partly explain the difference observed in the resulting indoor emissions. In addition, we use EFs from recent studies (Sommer et al., 2019; EMEP/EEA, 2019), whereas a parametrization relying on the temperature and the ventilation rate is used in FANv2 (Vira et al., 2020). However, the parametrization in FANv2 has been adjusted to reproduce default EFs for barns and stores from EMEP/EEA (2016) under European conditions. The use of updated EFs compared with Vira et al. (2020) largely explains the differences between the estimated indoor emission estimates. Moreover, in contrast to FANv2, our manure management module integrates a distinction between solid and liquid manure handling for each livestock type, with very different EF values. As discussed in Groen et al. (2016), Mu et al. (2017), and Uwizeye et al. (2017, 2020), uncertainties associated with EFs are large and can lead to over- or underestimates in indoor emissions and resulting soil emissions. The sensitivity of our calculated total emissions (manure management and soil) to this input parameter will be described in detail in Sect. 3.4, with a change of 14 % at the global scale demonstrating that the parameter has a significant impact on indoor emissions. Furthermore, we consider that each animal category has unique grazing, housing, and yarding periods, whereas Vira et al. (2020) consider regional livestock production systems.

3.3 Emissions at the regional scale

Good agreement is found for the NH₃ emissions in China between CAMEO, CEDS, and FAN estimates. In agreement with CEDS, India is the second biggest emitter region (Fig. 3b). However, FANv2 estimates much higher emissions in Africa and Latin America. As mentioned, there are important gaps between estimates given by CEDS and FANv2 for these two regions. In Africa and Latin America, our calculation leads to intermediate results between CEDS and FANv2 estimates.

Our estimates are quite similar to those given by FANv2 and CEDS, especially in Europe, China, and India. However, our indoor emissions are usually lower than what is computed by FANv2, except in the USA, where our emissions are slightly higher. In China, our total agricultural estimate of 9.3 Tg N yr⁻¹ is higher than those from CEDS and FANv2. However, our estimate is lower than those from several studies focusing on Chinese emissions: Kang et al. (2016), Li et al. (2021), Zhang et al. (2018), Crippa et al. (2018), and Zhang et al. (2017b) (9.7, 10.0, 10.3, 11.3, and 12.4 Tg N yr⁻¹, respectively). In India, the emissions that we compute (7.1 Tg N yr⁻¹) are closed to the FANv2 emissions (7.4 Tg N yr⁻¹). Our emissions in North America are 30 % higher than the FANv2 value and the other estimates. As shown in Fig. 6, the southern central part of the USA (mainly Texas) is a hot spot region, where the maximum value can reach 3 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$. The combination of low soil moisture and a high temperature simulated in this area can explain such high values of volatilization from the soil. In FANv2, emissions do not exceed 1.5 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$. Unlike the USA, our emissions in Europe are 30 % lower than FANv2 and EDGAR4.3 emissions, whereas they are 15 % lower than EMEP and CEDS emissions.

Unlike in FANv2, where three types of N fertilizers in the form of ammoniacal nitrogen, urea, and nitrates are considered, we assume a constant ammonium fraction of 0.6 in synthetic fertilizers. Even if the yearly fertilizer application is similar to the amount used in FANv2, the ammonium pool in soil from the mineral application can be different. This may imply differences in the emissions, especially in regions where the mineral application is intensive, such as Europe, China, and India (See Fig. S2). Concerning Africa and Latin America, our calculated emissions (∼5.4 Tg N yr⁻¹) are within the range of CEDS and FANv2 values. Africa and Latin America are characterized by specific environmental conditions along with different vegetation types, which may explain the uncertainties present in the estimates. There is also a lack of information regarding agricultural practices and resulting emissions in these regions. In Argentina, Castesana et al. (2018) estimated agricultural emissions of about 0.31 Tg N, whereas our emissions reach 0.91 Tg N and are closer to the estimates of Vira et al. (2020) (1.02 Tg N). The large differences mainly come from fertilizer use, reaching 1400 Gg N in their approach. The fertilizer use from the NMIP project (Tian et al., 2018) (752 Gg N) is in line with the reported values in Castesana et al. (2018) and consistent with the IFASTAT values for 2010–2015 (400–900 Gg N). We can not easily conclude whether the emissions differences come from EFs or manure production estimates. We can only compute a posteriori a single EF for soil emission from our process-based modeling, whereas no manure stock production is given in Castesana et al. (2018).

3.4 Sensitivity to model parameters

Among the parameters tested, the manure pH used in the calculation of the gaseous phase of ammonia is the strongest driver of NH₃ emissions (Fig. 7). At the global scale, the pH induces an increase of about 74 % when fixed at 7.5 compared with the reference value fixed at 7.

Figure 7 Global and regional differences between NH₃ emissions from the test simulation (TEST: CONC_0.3, CONC₃, pH_7.5, TIM₁₀, TIM₃₆₅, EF_max, $p_{{\mathrm{zact}}_{\mathrm{1.5}}}$, and FERT_0.75) and the reference simulation (REF) in teragrams of nitrogen (Tg N). Percentages indicate change relative to the REF value: $(\mathrm{TEST}-\mathrm{REF})/\mathrm{REF}\times \mathrm{100}$. China–K–J accounts for China, Korea, and Japan.

The impact of the pH is very variable from one region to another: it reaches up to 90 % in some regions such as Africa and the USA, whereas it is the lowest in India (49 %) (Fig. 7). In order to explain these regional differences, we explored the drivers of the spatial distribution of modeled NH₃ emission sensitivity to pH. The spatial distributions of the sensitivity of NH₃ emissions and the gaseous ammonia pool of the soil to pH are similar (Fig. 8a). In particular, the sensitivity is low in India compared with other regions, like Europe, for both variables. Figure 8b shows the spatial pattern of the dissociation constant of ammonium. The highest values are located in the warmest regions, such as India, where the temperature is one of the main drivers of the dissociation constant. As the dissociation reaction (Eq. 20) is favored in these regions (more NH₃ is available), it implies that volatilization is more likely to occur. Along with the high dissociation constant, India is characterized by an important soil ${\mathrm{NH}}_{\mathrm{4}}^{+}$ concentration (Fig. 8c) due to intensive agricultural input (mineral fertilizer and manure application), leading to a high quantity of TAN being available for emissions. In regions where conditions promote high NH₃ volatilization, pH is a weaker driver of emissions. Despite the regional differences in the pH sensitivity, it is an environmental parameter that is an essential driver in the emissions and can be a source of significant uncertainty in our model.

Riddick et al. (2016) also studied the emission sensitivity to pH change. They estimated an increase of 50 % and 70 % in the manure and fertilizer emissions, respectively, when changing the pH from 7 to 8. Even though the sensitivity that we describe seems higher than that in Riddick et al. (2016), no clear conclusion can be reached because we consider a unique pool of TAN. Indeed, we calculate the impact of a change in the total emissions, whereas Riddick et al. (2016) calculated changes in both the manure and fertilizer emissions by changing the pH of the two TAN pools (manure and the fertilizer) separately.

Figure 8 (a) The relative anomaly of gaseous ammonia in soil between the pH_7.5 simulation and the reference simulation as a percentage (%). (b) The ${\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)$/NH₃(aq) dissociation constant from the reference simulation (mol L⁻¹). (c) The soil NH₄ concentration modulated by the z_activity parameter (g N m⁻²).

Changing the duration of the N application (mineral fertilizer and stored manure) from 183 to 10 d induces a 30 % increase in global NH₃ emissions. The highest increase is calculated in China (86 %), whereas Europe's impact is slightly negative (−4 %). Reducing the duration of fertilization induces a significant change in the emission dynamic (Fig. 9), with emission peaks occurring right after the start of the vegetation growing season, considered in our model as the first day of the N application period (Fig. 10). However, the emission sensitivity to fertilization duration varies across regions, depending on the environmental conditions after the start of the vegetation growing season. In China, this signal is on average higher in April–May (Fig. 10). It is the period with the lowest soil moisture value and the highest soil temperature, which are conditions that maximize emissions. This could explain the high sensitivity that we observe in this region. On the contrary, in Europe, the growing season signal appears mainly in February and April, when the soil temperature is the lowest and the soil moisture the highest, indicating that these conditions are the least favorable for emissions, resulting in a negative sensitivity.

Figure 9 Regional NH₃ emissions ($\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) from different sets of simulations. Green lines are the emissions from the reference CAMEO simulations. Dark and light blue lines are emissions from the respective TIM₁₀ and TIM₃₆₅ simulations. Dark red and orange lines are emissions from the respective CONC_0.3 and CONC₃ simulations. China–K–J accounts for China, Korea, and Japan.

Figure 10 The regional monthly simulated soil moisture (m³ m⁻³, in blue) and soil temperature (K, in orange). The regional monthly signal for leaves to start to grow in ORCHIDEE averaged over the different PFTs is drawn using dotted gray lines (this metric is unitless and can be seen as a qualitative signal for the start of the vegetation growing season). Variables are for 2006. China–K–J accounts for China, Korea, and Japan.

When N is constantly applied throughout the year (365 d), the emissions are reduced by about 12 % globally, with India being the region with the strongest reduction (−28 %). The emissions are lower when N is applied throughout the year because it reduces the quantity of N emitted when conditions are the most favorable for volatilization. The variation in p_{zact_deep} from 1 to 1.5 m has a relatively constant impact on NH₃ emissions of about −20 % over every region. Increasing p_{zact_deep} increases the dilution of specific ammonium sources (that were assumed to be “deep sources” such as BNF, deposition, and mineralization) in soil, which in turn reduces emissions. Concerning the sensitivity to the content of ammonium in N fertilizers, when ${\mathrm{Frac}}_{{\mathrm{NH}}_{\mathrm{4}}^{+},\mathrm{fert}}$ is increased by 20 %, emissions increase by about 14 % on average. In China and India, where fertilizer use is the highest, the increase can reach +24 % (see Fig. S5). When fixing the atmospheric concentration at 0.3 and 3 µg N m⁻³, the global NH₃ emissions increase by 5 % and decrease by 11 %, respectively. In Africa, where the impact of using a concentration of 3 µg N m⁻³ is the highest, emissions are reduced by 26 %, whereas they slightly increase in India (3 %). Indeed, over India, atmospheric NH₃ concentrations from LMDZ-INCA are higher than 3 µg N m⁻³, in particular during May which experiences values up to 7 µg N m⁻³ (Fig. S4). Counterintuitively, using a fixed atmospheric NH₃ concentration (the CONC_0.3 and CONC₃ simulations) does not induce any important change in the seasonality of emissions (Fig. 9). The sensitivity to this parameter has also been tested in FANv1, and the same range of model response was found (Riddick et al., 2016).

Finally, higher emission factor values imply an increase in total NH₃ emissions of about 6 % globally. Although this change is not as significant as other factors, it is worth noting that this impact is only driven by indoor emissions, which account for 22 % of the total emissions globally. In Europe and Latin America, for instance, where the contribution of indoor emissions to the total emissions can reach more than 90 % (Fig. S3), the impact of using higher EF_max values calculated at the scale of these two regions (14 % and 12 %, respectively) is higher than the impact calculated at the global scale.

3.5 Emission seasonality

Seasonality patterns have been first explored by comparing our emissions against the CEDS inventory and the FANv2-simulated emissions. As shown in Fig. 11, we calculate maximum emissions during the spring and summer seasons, whereas the emissions peak almost everywhere only during spring in CEDS and FAN. In the three datasets, the lowest values are calculated during winter, when the meteorological conditions are not favorable for emissions and N application is lowest.

The summer peak observed in our emissions and the spring peak in CEDS and FAN in the USA, more specifically in the central and southern central parts of the USA, are also reported by the Magnitude And Seasonality of Agricultural Emissions (MASAGE) bottom-up inventory from Paulot et al. (2014). In MASAGE, two peaks are highlighted during the year: one in March and the other in June. Goebes et al. (2003) and Pinder et al. (2006) attributed these peaks to the timing of the mineral fertilizer and manure application. This is consistent with our approach, as indoor emissions do not vary over the year and only the N application is time dependent. In Europe, our emissions are higher in summer, whereas Paulot et al. (2014) estimated a clear peak in spring, like in FANv2 and CEDS. In addition, the analysis of Fortems-Cheiney et al. (2020), based on different emission inventories for France, showed the substantial contribution of mineral fertilizer application to emissions, leading to a peak in April. Paulot et al. (2014) demonstrated that April is when emissions reach a maximum over several European agricultural regions, such as Portugal and Spain or Benelux, Germany, and Denmark, due to local regulations preventing farmers from applying manure outside of the growing season. These regions are also characterized by large emissions in July, likely from livestock. This has recently been confirmed by top-down emissions based on Cross-track Infrared Sounder (CrIS) and IASI observation estimates in the UK (Marais et al., 2021). Our approach is constrained by the low level of detail in crop diversity, mainly due to the spatial resolution of the model. Thus, we choose a long enough N application period to catch the crop system diversity. In China, the highest emissions are calculated by our model in summer, which is supported by previous inventories (Streets et al., 2013; Kang et al., 2016; Xu et al., 2018) and satellite observations such as the Tropospheric Emissions Spectrometer (TES) instrument (Shephard et al., 2011) and the Atmospheric Infrared Sounder (AIRS) retrievals (Warner et al., 2017). Our emissions in India peak in spring but remain high during summer. This pattern is also highlighted in the HTAP emissions and IASI satellite data (Janssens-Maenhout et al., 2015; Van Damme et al., 2017), where no strong seasonality is shown over the Indo-Gangetic Plain, but higher emissions are noted from April to September.

Figure 11 Seasonal patterns of ammonia emissions ($\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) from total agricultural sources simulated by CAMEO (2005–2015 average, first column), estimated in the CEDS inventory (2005–2015 average, second column), and simulated by the FANv2 model (2010–2015 average, third column).

To complete our analysis, monthly emissions derived from the IASI satellite (IASI^inv) averaged over the 2011–2015 period are used as a comparison regarding different hot spot regions defined in Fig. 5. The same operation has been carried out for agricultural emissions from FANv2 and the CEDS inventory, and the details of the CAMEO₊, CEDS₊, and FAN₊ dataset constituents are listed in Table 6. First, it is worth noting that the seasonality of CAMEO₊, CEDS₊, and FAN₊ is primarily due to their respective agricultural emissions; the other sources have no important role in the seasonality. Indeed industrial and waste sources show very low annual standard deviations (in the ranges of 10⁻⁷ and 10⁻⁸ $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$, respectively) compared with the agricultural emissions estimated by CEDS, FANv2, or CAMEO (0–8 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$). Biomass burning emissions have a standard deviation reaching 0.6 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$ in areas characterized by high fire activity, whereas the deviation is intermediate in agricultural regions ($\sim {\mathrm{10}}^{-\mathrm{2}}$ $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$; Fig. S7). Overall, our emissions seem to be consistent with IASI^inv, with the general patterns and the absolute values being similar (Fig. 12), except for the equatorial Africa and central Brazil regions where CAMEO₊ largely underestimates IASI^inv. Seasonal emission patterns of CEDS₊ and FAN₊ are quite close to each other but very different from IASI^inv. FAN₊ and CEDS₊ usually depict a sharp peak in spring (1 month in advance in FAN₊) as well as another smaller peak only for CEDS₊ in September–October. While the seasonality in CEDS₊ is artificially retrieved from a specific profile, the temporal variability in FAN₊ is driven by the meteorological conditions and the crop types present within the pixel. Even though the representation of crops in version 5.0 of the Community Land Model (CLM5.0) used in FANv2 is more precise than in ORCHIDEE (eight crop types compared with two), we observe that the temporal variability is better represented in our approach when compared to IASI^inv. In FANv2, the fertilizer application is triggered 20 d after leaf emergence. This might explain the single sharp peak observed in spring, while the long period (6 months) that we have chosen seems to capture the general pattern of emissions better. More specifically, there is very good agreement between CAMEO₊ and IASI^inv in the central USA, where emissions peak in summer (>2 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$). In Europe, our calculated emissions show a clear and strong peak in July, whereas IASI^inv patterns are different. In northern Europe, IASI^inv presents two peaks, in March and in August, with very low emissions in winter, whereas emissions remain at around 0.5–0.8 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$ in CEDS₊, CAMEO₊ and FAN₊. These differences might be explained by high uncertainties associated with the IASI instrument during this period. Van Damme et al. (2014) highlighted limitations in the IASI (Metop-A) measurement availability over Europe in winter 2011. The number of cloud-free observations is low, especially in December and January, with only 4 % of the dataset being associated with an error lower than 50 % due to the thermal contrast (defined as the temperature difference between the Earth’s surface and the atmosphere at 1.5 km) and the amount of ammonia present in the atmosphere. In southern Europe (the Spanish region), IASI^inv emissions show a weaker seasonal cycle than CAMEO₊, peaking in summer. Spain is characterized by a diversity of agricultural production types (e.g., crops, fruits, and olives) which are not differentiated in our model. This highlights a limitation in our approach: the homogeneous representation of agricultural lands causes emissions to have the same seasonal pattern over this region. In the northern part of India, CAMEO₊ highlights a peak in May, which is 2 months earlier than the peak in the IASI^inv. In agreement with IASI^inv results, Tanvir et al. (2019) also depicted a clear peak in the NH₃ concentration in July using the TES data for the Indo-Gangetic Plain. It is worth noting that IASI observations in this specific region might be associated with a high level of uncertainties (Marais et al., 2021). Marais et al. (2021) indicated that, in addition to the intense biomass burning season and the relatively low abundance of acidic aerosols in northern India, warm temperatures may increase emissions and suppress the partitioning of NH₃ to aerosols, thereby inducing an enhancement in the spectral signal. In the Chinese hot spot, we observe the same peak in summer for both CAMEO₊ and IASI^inv, even though our maximum value is almost 2 times higher than that of IASI^inv (Fig. 12).

In Africa and Latin America, our emissions show a less pronounced variability over the year than IASI^inv. Results in Africa show three peaks (February, May, and October) in the IASI^inv, whereas our emissions highlight only one peak in June (Fig. 12). Equatorial Africa is a specific region with respect to emission seasonality that has recently been studied in Hickman et al. (2021). The aforementioned work revealed different ecoregion drivers of the atmospheric NH₃, explaining the seasonal patterns observed by the IASI satellite. The region chosen in our study is between the two northern ecoregions, wet and dry savannas, which are characterized by important livestock densities and intense biomass burning activities.

The June peak retrieved in CAMEO₊ and IASI^inv can be attributed to the seasonal pattern of the dry savanna. This time of the year corresponds to the rainy season in the dry ecoregions, and emissions from soil (from livestock excreta and natural processes) are expected to be stimulated through microbial activity enhanced by precipitation (Hickman et al., 2018). Hickman et al. (2021) demonstrated that precipitation and temperature are the most important predictors to explain the seasonality of the emissions in this region. The two other peaks in the IASI^inv during the dry season (February and October) could be a contribution from the wetter region located just below the dry savanna. Hickman et al. (2021) demonstrated that, in addition to the importance of the soil emissions in this region, vegetation fires during this period might explain additional emissions. This result is also supported by the long-term measurements from INDAAF (International Network to study Deposition and Atmospheric composition in Africa), which suggest that the seasonality in the wet savanna is the result of biomass burning, with a high increase in the concentrations during the dry season (Adon et al., 2010). The absence of these two peaks in our estimate can be explained by an underestimation of the biomass burning emissions from the GFEDs4 inventory used to complement our emissions (van der Werf et al., 2017). This inventory is based on MODIS biomass burning area, and recent analysis has suggested that MODIS underestimates fire emissions by a factor of 2–5 due to the non-detection of small fires (Roteta et al., 2019; Hickman et al., 2021; Ramo et al., 2021), mainly from agricultural practices.

The case of Latin America has been much less studied, and we observe strong seasonality in IASI^inv. The two regions in Latin America are characterized by croplands, intensive livestock farming, and biomass burning activities. Over central Brazil, IASI^inv reveals an important peak in September (>4 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) that is only represented in the CAMEO₊ and FAN₊ time series by a smooth increase. Andela et al. (2017) has shown a strong positive spatial correlation between burned area and cropland fractions in this ecoregion, probably suggesting significant agricultural waste burning. In addition, Castro Videla et al. (2013) have also shown a maximal biomass burning activity measured by MODIS via the monthly mean number of fires (MODIS fire dataset) in the Brazilian Caatinga shrublands and in the northeastern part of the Cerrado region in September. In the Pampas region, the general seasonality from IASI^inv and our emissions describes the highest emissions from September to March and low emissions during the rest of the year. However, our emissions do not highlight the clear peaks in March and September that are shown by IASI^inv. The work of Castro Videla et al. (2013) has pointed out sugarcane and soybean expansion in this region as the main drivers of biomass burning. This can lead to the same conclusion as for Africa with respect to the small fires that are often not detected by MODIS and, thus, might be not considered in the GEFDs4 inventory used in our work. It is also worth noting that manure application as a fertilizer is not a common practice in Argentina (Vázquez Amabile et al., 2015). In our approach, all of the manure is applied. This would potentially lead to a distinct seasonal cycle compared with a case in which all of the manure is stored for the whole year.

Figure 12 Monthly regional NH₃ emissions ($\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$). The CAMEO emissions accounting for natural and agricultural emissions aggregated with other sources are represented by the solid green line (CAMEO₊), and CAMEO emissions accounting for agricultural emissions are shown using the dotted green line (CAMEO_agri). The agricultural sector of CEDS alone and CEDS aggregated with other sources are represented by black dotted (CEDS_agri) and solid lines (CEDS₊), respectively. The IASI^inv product is shown in red (IASI^inv). The agricultural emissions from FANv2 aggregated with other sources are shown in blue (FAN₊). Other sources include biomass burning from van der Werf et al. (2010) and industrial and waste sectors from CEDS. The regions are defined in Fig. 5.

The temporal correlation scores between the reference (here IASI^inv) and CAMEO₊, CEDS₊, and FAN₊ calculated over the monthly time series for the 2011–2015 period are plotted in Fig. 13. The results highlight excellent month-to-month agreement with respect to variability between CAMEO₊ and IASI^inv in most regions of the globe. In the main hot spot regions, such as China, India, Europe, and the USA, correlations are between 0.7 and 0.9, while correlations between IASI^inv and CEDS₊ and between IASI^inv and FAN₊ hardly exceed 0.5. This means that our modeling approach enables a satisfying representation of the seasonal cycle in terms of agricultural and natural emissions in comparison with CEDS agricultural emissions, where a forced seasonal profile (two respective high volatilization peaks in May and September) is used, and with the FANv2 model, accounting for a more realistic representation. However, in the southeastern part of the USA and the Chaco region in Latin America, we observe a degradation of the seasonal pattern in our emissions compared with both CEDS₊ and FAN₊, which show high correlations. The Chaco region is one of the main hot spot in terms of natural soil emissions in our model, as highlighted in Fig. S6. It is characterized by savanna with grasslands, thorn forests, a mosaic of woods with savanna, and shrubs, and coarse grasses predominate (Berry et al., 1995). However, most of the natural emissions computed in ORCHIDEE originate from temperate broad-leaved summer green PFTs. Many studies have demonstrated the importance of the biomass burning events with respect to the emission quantities mainly occurring during the dry season in September–October (Castro Videla et al., 2013; Pereira et al., 2022). IASI^inv depicts an important peak in the NH₃ emissions (See Fig. S6) during this time of the year which can be attributed to fire events. We observe that natural and agricultural emissions have very similar patterns and are in the same range. However, the fire contribution in the CEDS₊ and CAMEO₊ datasets appears to be very low, supporting the limitation of using the GFEDs4 inventory in bottom-up NH₃ emission estimates for comparison with IASI^inv. The degradation of the correlation between IASI^inv and CAMEO₊ compared with the correlation with CEDS₊ in the Chaco region is explained by the fact that there is almost no temporal variability in the CEDS₊ at the annual scale.

In India, there is an interesting pattern with a clear longitudinal delimitation – western negative and eastern positive correlations – in CAMEO₊ and CEDS₊. In the northwestern part of India, CAMEO₊ performs better with respect to capturing the IASI^inv seasonality than CEDS₊ and FAN₊.

Figure 13 Temporal correlation scores between IASI^inv and CAMEO emissions aggregated with other sources (a) and between IASI^inv and total emissions from CEDS for agricultural, industrial, and waste sources aggregated with biomass burning (b) over the 2011–2015 period. Temporal correlation scores between IASI^inv and FANv2 aggregated with other sources are shown in panel (c). Other sources include biomass burning from van der Werf et al. (2010) and industrial and waste sectors from CEDS. Gray grid cells correspond to a standard deviation in the monthly IASI^inv time series lower than 0.2 $\mathrm{g}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$ in order to avoid biased correlation scores.

Based on the comparison with the seasonality of IASI^inv, there is a strong limitation in using CEDS as NH₃ emission information for a CTM in order to study its impact on atmospheric chemistry. More specifically, we demonstrate that using ORCHIDEE land-based emissions has the potential to improve the seasonal signal of the resulting ammonia concentration in the atmosphere.

However, using IASI^inv to evaluate our model results also has limitations due to the uncertainties associated with the satellite product and the derivation method. For example, many studies using IASI data do not consider winter observations in the USA nor Europe (Marais et al., 2021) due to a potential degradation of the data because of atmospheric conditions (e.g., cloud cover and low temperatures). The use of the NH₃ lifetime simulated by LMDZ-INCA in the inversion method is also associated with uncertainty. Moreover, the interpolation method used to regrid the IASI observations onto the LMDZ-INCA grid can be a source of uncertainties, as demonstrated in Evangeliou et al. (2020).

4 Conclusions

In this study, we implement a new module dedicated to global NH₃ emissions from agricultural practices including livestock waste management and mineral fertilizer application within the ORCHIDEE land surface model. Our development allowed us to consider dynamical variables (such as surface temperature and humidity) through different physical soil processes for the calculation of the NH₃ fluxes. This aspect, which is often neglected in bottom-up approaches, is key for a realistic seasonal representation of the emissions. In contrast to other emission models, our module interacts with the ORCHIDEE model for vegetation variables such as the biomass productivity. It allows the calculation of a grazing index, a global indicator of the pressure exerted on the vegetation. We estimate global agricultural emissions of about 44 Tg N yr⁻¹, with soil volatilization from fertilizer and manure application accounting for 78 % of emissions and indoor emissions (from livestock housing, yard and storage) accounting for the remainder. The spatial distribution of the calculated emissions is consistent with previous studies (the CEDS bottom-up inventory from McDuffie et al., 2020, and the process-based model by Vira et al., 2020) highlighting the most important emission regions, such as eastern China, northern India, the USA, and Europe, which are characterized by either a high N application rate or intensive livestock farming.

In order to evaluate the modeled emissions, different sensitivity simulations involving key parameter variations have been performed. The most important parameter driving the emissions is the pH of the N input, which induces an increase of about 74 % for NH₃ emissions when shifted from 7 to 7.5. Assuming a constant value for the pH simplifies our approach. Using a soil pH map would imply more complex processes, involving a change in the pH during N application. In addition to the ammonium content, the pH used is rarely available in the literature.

Manure management emissions are also associated with uncertainties due to the use of EFs. Even though EFs are calculated using an extensive definition of livestock and management systems as well as considering variations in climate and management practices, the use of European EFs for the whole globe is an important assumption. However, we demonstrate that the overall emissions are moderately sensitive to the EFs, with a global change of 6 % for emissions when the maximum range given by Sommer et al. (2019) is used. Similarly, regional parameters given by FAO (2018) were simplified to match the representation of the vegetation distribution in ORCHIDEE. Assuming livestock feed is composed of only grass and crop products neglects the use of the agro-industrial by-products, which is standard practice in Europe and the USA. Modeled emissions are also sensitive to the timing of the N application, especially in China, where a shortening of the fertilization period induces very high volatilization rates. Apart from this parameter, none of the other factors tested appear to be important drivers of the emission seasonality. Finally, the seasonality patterns have been further analyzed using satellite-derived emissions. The comparison suggests that ORCHIDEE simulates a very good representation of the seasonality of NH₃ emissions, with correlation scores larger than 0.7 in the most important emission regions.

In addition to the gain with respect to realistic seasonality, our approach fills the lack of estimates for emissions from natural soils, which are missing in almost every inventory. This is highly interesting for Africa and Latin America, where these sources are important (Hickman et al., 2018) and not well studied. These encouraging results prove the potential of coupling ORCHIDEE land-based emissions to CTMs, which are currently forced by bottom-up anthropogenic-centered inventories such as CEDS. This framework provides room to improve the representation of the emissions, as atmospheric variables are dynamically simulated by CTMs. For instance, the surface NH₃ concentrations used in the final calculation of the emissions could be updated at each time step instead of prescribing an external monthly file for a given climatology. In addition, a tight relationship exists between emission and deposition of NH₃, as NH₃ is particularly reactive, and deposition of ${\mathrm{NH}}_{\mathrm{4}}^{+}$ contributes to the reemission of NH₃ from natural and managed soils. By coupling emissions with the global LMDZ-INCA CTM through the dynamic calculation of wet and dry deposition, we plan to improve the representation of the emissions as well as the atmospheric concentrations. In addition, a further evaluation will consist of comparisons with atmospheric composition observations.