High-resolution hybrid inversion of IASI ammonia columns to constrain US ammonia emissions using the CMAQ adjoint model

Ammonia (NH₃) emissions have large impacts on air quality and nitrogen deposition, influencing human health and the well-being of sensitive ecosystems. Large uncertainties exist in the “bottom-up” NH₃ emission inventories due to limited source information and a historical lack of measurements, hindering the assessment of NH₃-related environmental impacts. The increasing capability of satellites to measure NH₃ abundance and the development of modeling tools enable us to better constrain NH₃ emission estimates at high spatial resolution. In this study, we constrain the NH₃ emission estimates from the widely used 2011 National Emissions Inventory (2011 NEI) in the US using Infrared Atmospheric Sounding Interferometer NH₃ column density measurements (IASI-NH₃) gridded at a 36 km by 36 km horizontal resolution. With a hybrid inverse modeling approach, we use the Community Multiscale Air Quality Modeling System (CMAQ) and its multiphase adjoint model to optimize NH₃ emission estimates in April, July, and October. Our optimized emission estimates suggest that the total NH₃ emissions are biased low by 26 % in 2011 NEI in April with overestimation in the Midwest and underestimation in the Southern States. In July and October, the estimates from NEI agree well with the optimized emission estimates, despite a low bias in hotspot regions. Evaluation of the inversion performance using independent observations shows reduced underestimation in simulated ambient NH₃ concentration in all 3 months and reduced underestimation in NH${}_{\mathrm{4}}^{+}$ wet deposition in April. Implementing the optimized NH₃ emission estimates improves the model performance in simulating PM_2.5 concentration in the Midwest in April. The model results suggest that the estimated contribution of ammonium nitrate would be biased high in a priori NEI-based assessments. The higher emission estimates in this study also imply a higher ecological impact of nitrogen deposition originating from NH₃ emissions.

1 Introduction

Ammonia (NH₃) emissions play a major role in ambient aerosol formation and reactive nitrogen deposition (Stevens, 2019: Houlton et al., 2013). However, our understanding of NH₃ sources and sinks is limited by the large uncertainties present in the NH₃ emissions inventories (Xu et al., 2019; McQuilling and Adams, 2015). In chemical transport models, uncertainties in NH₃ emissions propagate into the dynamic modeling of the atmospheric transport, chemistry, and deposition of NH₃, other reactive nitrogen species, and other key atmospheric constituents associated with NH₃ (Heald et al., 2012; Paulot et al., 2013; Kelly et al., 2014; Zhang et al., 2018b), hindering an accurate assessment of the various NH₃-related environmental impacts and the associated sources. The large uncertainties in the NH₃ emission inventories are partially due to a lack of sufficient in situ NH₃ measurements that could be used to constrain emission estimates (Zhu et al., 2015).

Emerging satellite observations of gaseous NH₃ provide a unique opportunity to better constrain the bottom-up NH₃ emission estimates for both their spatial distribution and seasonality. Bottom-up inventories calculate the NH₃ emissions based on estimated activity levels and corresponding emission factors, both of which are subject to high uncertainties, particularly for agricultural sources, the major contributor (Cooter et al., 2012; McQuilling and Adams, 2015). Several studies have utilized NH₃ column density retrieved from the Infrared Atmospheric Sounding Interferometer (IASI) (Clarisse et al., 2009; Van Damme et al., 2015b) or the Atmospheric Infrared Sounder (AIRS; Warner et al., 2016) as well as the inferred surface mixing ratio of NH₃ from the Cross-track Infrared Sounder (CrIS; Shephard and Cady-Pereira, 2015; Shephard et al., 2020) to characterize the spatiotemporal distribution of NH₃. These satellite measurements are useful for supplementing emission inventories to identify and quantify underestimated or missing emission hotspots, especially in intensive agricultural zones (Van Damme et al., 2018; Dammers et al., 2019; Clarisse et al., 2019). These studies find that the satellite-derived emission estimates are often twice as much as the bottom-up estimates on a regional scale and can be over 10 times higher over hotspots. However, the NH₃ retrievals from satellites are also subject to large uncertainties when the signal-to-noise ratio is low, which limits their ability to accurately measure NH₃ columns in low-emission areas (Clarisse et al., 2010; Van Damme et al., 2015a).

Inverse-modeling-based optimization combines the information from a priori emission inventories and observations and allows us to use the information from both. As one of the inverse modeling methods, the four-dimensional variational assimilation (4D-Var) method seeks the best emission estimate by minimizing a cost function that measures the differences between observations and model predictions, as well as the differences between a prior and adjusted emission estimates. 4D-Var can be computationally expensive at fine model resolutions or with a large set of observations to be assimilated (Brasseur and Jacob, 2017). Recent studies have taken advantage of the implementation of the adjoint technique in the chemical transport models to conduct 4D-Var for optimizing emissions estimation (Zhu et al., 2013; Paulot et al., 2014; Zhang et al., 2018c). The adjoint-based inversion method calculates the gradients of a cost function analytically and searches for the solution using a steepest-descent optimization algorithm through iterating (Brasseur and Jacob, 2017). By testing the performance of the inverse modeling method using artificial observational data, Li et al. (2019) proposed that a two-step optimization process, which combines the iterative mass balance (IMB) method and the 4D-Var method, can further reduce the computational cost. The IMB method assumes a linear relationship between the NH₃ column density and local NH₃ emission and searches for the emission scaling factors iteratively until the simulated NH₃ column density converges to the observations. At a coarse ($\mathrm{2}{}^{\circ }\times \mathrm{2.5}{}^{\circ }$) resolution, the IMB method is as effective as the 4D-Var method and requires two-thirds less computational time. In the second step, emission scaling factors obtained from the IMB method with a coarser resolution are used as an initial starting point for the 4D-Var optimization process to reduce the overall computational time (Li et al., 2019).

This work utilizes satellite observations from the IASI-NH₃ column density measurements (IASI-NH₃) (Clarisse et al., 2009;Van Damme et al., 2017), to provide a high-resolution, optimized NH₃ emission inventory for the US developed using an adjoint inverse modeling technique (Li et al., 2019), the robustness of which is demonstrated by evaluation against multiple independent in situ measurements. The IASI-NH₃ dataset was applied to optimize NH₃ emission estimates from the 2011 National Emissions Inventory (2011 NEI) using the Community Multiscale Air Quality Modeling System (CMAQ) and its adjoint model at a 36 km × 36 km resolution. The multiphase adjoint model for CMAQ v5.0 was developed recently, including full adjoints for gas-phase chemistry, aerosols, cloud process, diffusion, and advection (Zhao et al., 2020). Both process-by-process and full adjoint model evaluations show reasonable accuracy based on agreements between the adjoint sensitivities and forward sensitivities (Zhao et al., 2020). Previous inversion-based NH₃ emission constraints using in situ measures are limited by the spatial coverage and representativeness of the measurements (Gilliland et al., 2006; Henze et al., 2009; Paulot et al., 2014;). Zhu et al. (2013) first attempted to optimize the NH₃ emission inventory using NH₃ derived from the Tropospheric Emission Spectrometer satellite at $\mathrm{2}{}^{\circ }\times \mathrm{2.5}{}^{\circ }$ resolution (Zhu et al., 2013). Inverse modeling at such a coarse resolution is limited to refining regional emissions. Similar to the inversion using CrIS NH₃ measurements (Cao et al., 2020), inversion with the IASI-NH₃ dataset allows us to perform the optimization at a finer resolution with its daily global spatial coverage. Furthermore, the hybrid inversion approach adopted in this study allows us to calculate full adjoint sensitivities online instead of using approximated sensitivities from the offline simulations (Zhu et al., 2013; Cao et al., 2020). The performance of our optimized estimates and the 2011 NEI are evaluated and compared based on in situ observed ambient NH₃ concentrations and NH${}_{\mathrm{4}}^{+}$ wet deposition. Finally, by substituting the a priori NH₃ emissions with the optimized emissions, we assess the subsequent changes in simulated ambient PM_2.5 concentrations and nitrogen deposition exceedances.

2 Materials and methods

2.1 IASI-NH₃ observations

NH₃ column densities retrieved from IASI on board the Metop-A satellite are assimilated to constrain spatially resolved NH₃ emissions using the 2011 NEI as the a priori inventory (Clarisse et al., 2009; Van Damme et al., 2014; USEPA, 2014). The polar sun-synchronous satellite has a 12 km diameter footprint at nadir and a bidaily global coverage. Only observations from the morning pass around 09:30 local standard time (LST) are used due to more favorable thermal contrast and smaller errors as compared to the night pass around 21:30 LST. A comparison between the IASI-NH₃ data and ground-based Fourier transform infrared (FTIR) observations shows a correlation between the two with r=0.8 and the slope = 0.73, indicating a tendency of IASI-NH₃ to underestimate the FTIR observations (Dammers et al., 2016). A comparison between IASI-NH₃ and airborne measurements also indicated an underestimation in California, while the comparison between IASI-NH₃ and ground observation from the Ammonia Monitoring Network (AMoN) indicated an overestimation (Van Damme et al., 2015a; NADP, 2014). Overall, the evaluations show broad consistency between IASI-NH₃ and other independent measurements, with no consistent biases identified. These evaluations were based on previous datasets. Here we use a new version that relies on another retrieval algorithm, which among other things has better performance for measurements under unfavorable conditions (Whitburn et al., 2016; Van Damme et al., 2017).

Specifically, the NH₃ products for 2011 from ANNI-NH₃-v2.2R-I datasets were used (Van Damme et al., 2017). The algorithm relies on the conversion of hyperspectral range indices to NH₃ column density using a neural network that takes into account 20 input parameters, characterizing temperature, pressure, humidity, and NH₃ vertical profiles. A relative uncertainty estimate is provided along with each of the NH₃ vertical column densities in the dataset. Small negative columns are possible – and these are valid observations, needed to reduce overall biases in the dataset. As the retrieval is unconstrained, no averaging kernels are calculated. We therefore directly compare the IASI-NH₃ column density with the simulated column density in CMAQ. Such comparison may be biased because the sensitivity of retrieved NH₃ column densities to NH₃ concentrations is height-dependent (typically peaks around 700–850 hPa) (Dammers et al., 2017; Shephard et al., 2015). Although the CMAQ-simulated NH₃ columns are also most sensitive to NH₃ concentration changes between 700 to 900 hPa (Fig. S1), we cannot quantify the relating uncertainties without knowing the averaging kernels. Without information on averaging kernels, differences between NH₃ vertical profiles in CMAQ and the ones used for retrieval may also contribute to the bias between retrieved and modeled column densities, depending on the magnitude of differences (Whitburn et al., 2016).

The retrieved NH₃ columns densities are regridded to the 36 km by 36 km CMAQ grid for 4D-Var data assimilation and 216 km by 216 km resolution (a six-grid-by-six-grid CMAQ simulation grid matrix) for iterative mass balance (Fig. 1). The mean column density (Ω_o) is calculated as the arithmetic mean of all retrievals with their centroids falling in the same grid cell, following the recommendation that the unweighted mean is preferred for the updated version of IASI-NH₃ as error-weighting can lead to biases (Van Damme et al., 2017). The error (molec cm⁻²) corresponding to the mean column density in each grid is calculated as

$$\begin{array}{}\text{(1)}& \stackrel{\mathrm{‾}}{\mathit{\sigma }}=\sqrt{{\displaystyle \frac{\sum {\left({\mathit{\sigma }}_i\times {\mathrm{\Omega }}_i\right)}^{\mathrm{2}}}{n-\mathrm{1}}}},\end{array}$$

where $\stackrel{\mathrm{‾}}{\mathit{\sigma }}$ is the mean error (molec cm⁻²), Ω_i is the ith retrieval of NH₃ column density from IASI-NH₃ Level 2 data, σ_i is the relative error associated with each Ω_i as reported, and n is the number of retrievals within each grid cell during the defined time period. For 4D-Var inversion and IMB inversion, daily and monthly means and errors are calculated, respectively.

Figure 1 IASI monthly average NH₃ column density in April, July, and October 2011 at 36 km by 36 km (a, b, c) and 216 km by 216 km (d, e, f) resolutions within the model simulation domain of this study. The average relative error associated with the column density is shown in the corner of each plot.

The observations from April, July, and October are used to constrain the monthly NH₃ emission estimates in corresponding months from 2011 NEI. Limited by the high computational cost of adjoint-model-based inversion, the optimization is only performed for the 3 months selected instead of a full year. Observations from winter months are not used because they are too noisy when the thermal contrast is low (Dammers et al., 2016).

2.2 NH₃ emission from 2011 NEI

The EPA 2011 NEI is used for a priori emission estimates. Major NH₃ sources include livestock waste management, fertilizer application, mobile sources, fire, and fuel combustion, with the majority being emitted by the first two sources. Specifically, the emissions from livestock waste management are estimated based on county-level animal population data and process-based daily emission factors. Emissions from fertilizer applications are estimated based on county-level fertilizer quantities and fixed emission factors, following the CMU Ammonia Model (USEPA, 2015). The NH₃ emissions over Mexico and Canada are derived from the simulation results of a fully coupled bi-directional agroecosystem and chemical transport model (FEST_C_EPIC_CMAQ_BIDI) (Shen et al., 2020). Emissions for other species also come from the 2011 NEI.

2.3 CMAQ and its adjoint

We use CMAQ v5.0 (Byun and Schere, 2006; USEPA, 2012) and its adjoint (Zhao et al., 2020), driven by meteorological fields produced from the Weather Research and Forecasting (WRF) Model v3.8.1 with grid nudging using the North American Regional Reanalysis dataset (NOAA, 2019). The simulated meteorological fields show good agreement with surface observations (Fig. S2) (NOAA, 2020). The CB05 chemical mechanism was adopted for gas-phase chemistry (Yarwood et al., 2005). The model implements ISORROPIA-II in the aerosol module (AERO06) to calculate the gas–particle partitioning of NH₃ and NH${}_{\mathrm{4}}^{+}$ (Fountoukis and Nenes, 2007). The simulation domain covers the contiguous US (CONUS) and part of Canada and Mexico with a 36 km by 36 km horizontal resolution and 13 vertical layers extending up to 14.5 kPa (∼16 km) (Fig. 1). To evaluation CMAQ model performance, the simulated gas–particle partitioning ratio of NH₃–NH${}_{\mathrm{4}}^{+}$ and NH${}_{\mathrm{4}}^{+}$ deposition is compared with observations from AMoN, the Clean Air Status and Trends Network, and the National Atmospheric Deposition Program (NADP) (Figs. S3 and S4). CMAQ captures the overall spatial pattern of these governing processes for atmospheric NH₃ abundance, considering the uncertainties in emissions, model parameters, and meteorological fields. Expanded evaluation of CMAQ model performance in simulating gas–particle partitioning and nitrogen deposition has been conducted in previous studies (Chen et al., 2019, 2020). Monthly simulations are conducted for April, July, and October in 2011 with a 10 d spin-up for each month.

2.4 Hybrid inversion approach

We chose the hybrid inversion approach to combine the advantage of the faster computational speed of the mass balance method and the better optimization performance of the 4D-Var method. The first step is to apply the IMB approach to adjust the a priori (2011 NEI) NH₃ emission at 216 km by 216 km resolution (referred to as the coarse grid hereafter) based on the ratio between the monthly averaged observed (Ω_o) and simulated (Ω_a) NH₃ column density at the satellite overpassing time, iteratively. At each iteration, the emission in each grid cell is scaled by the ratio following the equation below:

$$\begin{array}{}\text{(2)}& E_{\mathrm{t}}={\displaystyle \frac{{\mathrm{\Omega }}_{\mathit{o}}}{{\mathrm{\Omega }}_a}}\times E_{\mathrm{a}},\end{array}$$

where E_t and E_a are the new and a priori emission estimates, respectively. The method has been described in detail in previous studies (Li et al., 2019; Cooper et al., 2017; Martin et al., 2003). The IMB is applied at the coarse grid so that the NH₃ column will be dominated by the local emissions instead of transport from neighboring grids (Li et al., 2019). The coarse resolution also reduces the uncertainty associated with IASI-NH₃ as the number of retrievals increases in each grid cell. For grid cells with mean relative error larger than 100 %, the satellite observations are considered to be too noisy to provide useful constraints and the a priori emission estimates are retained. The iteration stops when the normalized mean square error either decreases by less than 10 % or begins to increase. The final scaling factor (ε₀) for each grid cell is the multiplication of the scaling factors derived at each iteration and downscaled to 36 km by 36 km resolution by assigning the same value to the six-by-six grid matrix. This scaling factor is applied to the 2011 NEI emissions to create the revised a priori estimate for the 4D-Var inversion.

Next, the 4D-Var inversion is performed. The solution of the optimization problem is sought iteratively by minimizing the cost function (J) defined as the combination of error-weighted, squared difference between the emission scaling factor and unity and the error-weighted, squared difference between IASI-NH₃ and the simulated column density, as below:

$$\begin{array}{}\text{(3)}& \begin{array}{rl}& J=\mathit{\gamma }{\left(\mathit{\epsilon }i-{\mathit{\epsilon }}_{\mathrm{0}}\right)}^T{S}_{\mathbf{a}}^{-\mathrm{1}}\left(\mathit{\epsilon }-{\mathit{\epsilon }}_{\mathrm{0}}\right)+{\left({\mathrm{\Omega }}_{\mathit{o}}-F\left(\mathit{\epsilon }\right)\right)}^T\\ & S_{\mathbf{o}}^{-\mathrm{1}}({\mathrm{\Omega }}_{\mathit{o}}-F(\mathit{\epsilon }\left)\right).\end{array}\end{array}$$

ε is the monthly emission scaling factor to be optimized at each iteration where $\mathit{\epsilon }=\log \left(E_{\mathrm{t}}/E_{\mathrm{a}}\right)$ on the 36 km by 36 km CMAQ grid, consisting of 6104 overland grid cells in the CONUS. S_a and S_o are error covariance matrices for the a priori emission estimates and IASI-NH₃ retrievals, respectively. With limited information on the spatial correlation of the error covariance, the two matrices are assumed to be diagonal (Paulot et al., 2014; Zhu et al., 2013). For S_o, the grid average absolute error is used to represent the observational error. Our test shows that negative Ω_o will lead to a continuous decrease in the adjusted emission for the grid cell because modeled column density cannot become negative. To limit the influence of these negative Ω_o, their original weights are multiplied by 0.01. For S_a, the uncertainty in each grid cell is assumed to be 100 % of the a priori emissions. F(ε) is CMAQ-simulated NH₃ column density sampled at the satellite passing time if there is at least one IASI-NH₃ retrieval in that grid cell; γ is the regularization factor balancing the relative contribution of the a priori emission inventory and IASI-NH₃ retrievals to the J value. γ is chosen to be 800 for April and 500 for July and October based on the L-curve criteria (Hansen, 1999) (Fig. S5).

The gradients of the cost function to NH₃ emissions are calculated by the CMAQ adjoint model. In each iteration, the emission-weighted monthly averaged sensitivities in each grid cell are supplied to the L-BFGS-B optimization routine contained in the “optimr” package in R to find the scaling factors that will achieve the minimum of the cost function (Zhu et al., 1997; Byrd et al., 1995). NH₃ column density is re-simulated using adjusted emissions by the new set of scaling factors. The iteration process is terminated when the decrease in J is less than 2 % or the local minimum is reached (Li et al., 2019; Zhu et al., 2013).

2.5 Posterior evaluation

The posterior emissions are evaluated by comparing the model simulation from optimized emissions with observations. Simulated results are compared with ambient NH₃ concentrations from AMoN (NADP, 2014) and the NH${}_{\mathrm{4}}^{+}$ wet deposition from NADP (NADP, 2019). The simulated NH₃ concentration in ppmv is converted to micrograms per cubic meter (µg m⁻³) using local temperature and pressure from the model meteorological inputs. For evaluation against the NH${}_{\mathrm{4}}^{+}$ wet deposition, the simulated deposition is scaled by the ratio between measured and simulated precipitation to eliminate the bias introduced by precipitation fields (Appel et al., 2011).

3 Results

3.1 Optimization performance evaluation

The optimized NH₃ emissions reduce the bias in the NH₃ columns between the satellite observation and the model prediction as shown by the decrease in the values of normalized root mean square error (NRMSE) and normalized mean biases (NMBs) in Fig. 2. There are negative biases using 2011 NEI in all 3 months, especially in areas with high emission rates. Although the IMB inversion can lower the NRMSE, it tends to over-adjust and introduce a positive bias, likely because of the coarse resolution and neglect of the impact of transport. The 4D-Var inversion effectively decreases the positive bias and further reduces the NRMSE. The cost function value reduces by 85, 46, and 38 % with the 4D-Var inversion in April, July, and October, respectively. We find that it is more challenging to adjust the emissions in April than in the other 2 months because of the greater differences in the magnitude and the spatial distribution of the emissions. The optimized NH₃ emission successfully captures the high NH₃ column density in the Southern States (Texas and Oklahoma), reducing the NRMSE by half in that region. Despite the general improvement in the model performance, negative biases in July increase in California's San Joaquin Valley. Scaling up the emission in the San Joaquin Valley will result in high NH₃ concentrations downwind even when the local NH₃ emissions downwind are zeroed, whereas the IASI-NH₃ observed concentrations downwind are low. The transported hotspot downwind of the San Joaquin Valley in CMAQ only occurs in July, suggesting near-field removal may not be captured at the current resolution, and warrants further investigation. Grid-by-grid comparison between model-simulated NH₃ column density using the a priori and optimized estimates with IASI-NH₃ shows improved agreement in both high- and low-emission grid cells after optimization (Fig. S6). It shows that the hybrid inversion approach can alleviate the weakness of direct 4D-Var inversion, which tends to over-adjust high-emission regions and under-adjust low-emission regions, mainly because the IMB inversion provides a better initial state.

Figure 2 CMAQ-simulated monthly average NH₃ column density for April, July, and October 2011 using the a priori emissions (a, b, c), the emissions adjusted by IMB (d, e, f), and the final optimized emissions using the hybrid approach (g, h, i). For comparison with the IASI-NH₃ retrievals, simulated NH₃ columns at the passing time were derived when there were observations in that grid cell. Normalized root mean square error (NRMSE) and normalized mean bias (NMB) between the simulated values and IASI-NH₃ are provided. Residue map (IASI-NH₃ − simulated NH₃ column densities) is shown in the corner of each plot.

The IMB inversion took three iterations to achieve the convergence condition for each month, and subsequently the 4D-Var inversion took 10, 4, and 6 iterations for April, July, and October, respectively. Fewer iterations are needed with the hybrid approach than the direct 4D-Var inversion, which typically takes up to 15 to 20 iterations of adjoint simulation (Paulot et al., 2014; Zhang et al., 2018a). The CPU time of a forward simulation is only one-fifth of an adjoint simulation. In total, the CPU time required by the hybrid approach is expected to be one-third to two-thirds lower than the direct 4D-Var inversion approach.

3.2 Optimized estimate of NH₃ emissions

The monthly total NH₃ emission in the CONUS increases by 35 % in April, 18 % in July, and 10 % in October for the optimized estimates. Spatially, the distribution for high-emission regions shifts from the Midwest in the 2011 NEI to the Southern States in the optimized estimates in April, whereas the hotspot regions remain consistent in July and October (Fig. 3). Regional total emissions are summarized according to the USDA farm production regions, which define the areas with similar crop production activities (Cooter et al., 2012). In general, the regional variation of NH₃ emissions in April is dominated by fertilizer application. The optimized estimates in the Corn Belt and Lake States regions are lower than the 2011 NEI, where high contributions from fertilizer applications were estimated. In contrast, the optimized estimates are 2–3 times higher than the 2011 NEI estimates in the Delta States and Southern States, where the a priori estimates for NH₃ emission from fertilizer application are low. The higher NH₃ emission estimates in the Southern States are driven by the enhanced NH₃ column densities from IASI over that region. IASI-NH₃ column densities are higher in 2011 than those in adjacent years (Fig. S7), which coincides with the higher surface temperature observed in 2011 (NOAA 2019) (Fig. S8). NH₃ emission will increase due to enhanced NH₃ volatilization from agricultural lands under warmer conditions (Bash et al., 2013; Shen et al., 2020). In fact, the optimized NH₃ emission pattern in April is more consistent with the spatial pattern of inorganic nitrogen fertilizer estimated based on plant demand (Cooter et al., 2012). NH₃ emission in 2011 estimated by CMAQ with a NH₃ bidirectional exchange model also predicted higher NH₃ emission in the Southern States (Shen et al., 2020). The ratio between NH₃ emission estimates in Southern States and those within the CONUS is 26 and 18 % in the optimized estimates and estimates including NH₃ bidirectional exchange, respectively. In comparison, the ratio is only 10 % in the a priori NEI estimates, suggesting a potential low bias in 2011 NEI. In July, regional differences are smaller except for the Northern Plain and Mountain region. In the Northern Plain, the NH₃ emission is 66 % higher in the optimized estimates, driven by the emission increase in hotspot areas with concentrated animal feeding operations (CAFOs) (USDA, 2012; Van Damme et al., 2017, Clarisse et al., 2019). The potential bias in different sectors suggests the need for sectoral inversion when a larger observational dataset becomes available in the future. In October, the relative difference is less than 10 % in most of the regions, indicating that the 2011 NEI appropriately reflects the NH₃ emission pattern. There is a significant increase in the NH₃ emissions in Mexico during all 3 months. Such an emission increment is crucial to improving the model performance in both Mexico and the southwestern US. However, it was not a goal of this study to determine emissions biases in Mexico given the limited information on NH₃ emissions.

Figure 3 The spatial distribution of monthly total NH₃ emission from the a priori (a, b, c) and optimized (d, e, f) estimates in April, July, and October. The total emission based on the a priori and optimized estimates is summarized for each USDA farm production region (g, h, i). The source contributions to total emission are shown for the a priori estimates.

The total NH₃ emissions in the optimized estimates are 623, 564, and 335 Gg per month in April, July, and October, respectively. In comparison, the emission estimates in the 2011 NEI are 462, 475, and 304 Gg per month for the 3 months. Similar to a bottom-up agricultural NH₃ emission inventory (MASAGE_NH₃) and two inverse-model-optimized estimates based on NH${}_{\mathrm{4}}^{+}$ wet deposition, we find a higher emission in the spring season (Paulot et al., 2014; Gilliland et al., 2006), while others, including the NEI, estimate a summertime peak (Zhu et al., 2013; USEPA, 2015; Cooter et al., 2012; Cao et al., 2020). The large variation between different inventories warrants both improved information on bottom-up inventories and more observations to support inverse model optimization in the spring season. Better knowledge about agricultural activities and more independent ground and space observations are needed. Besides the a priori emission inventory and observational constraints, the inversion performance will also be affected by other processes (e.g., gas–particle partition, transport, cloud and precipitation, and dry and wet deposition) governing the atmospheric abundance of NH₃. Future works refining the pertinent processes will also help improve the optimized NH₃ emission estimates. It should also be noted that there are interannual variations in emission inventories developed for different years. The good spatial agreement with IASI-NH₃ indicates that the 2011 NEI captures the NH₃ emission pattern in general in these 2 months. Although the inversion is only applied for the selected 3 months, the simulated NH₃ column densities using the a priori inventory are consistently lower than the IASI-NH₃ observations in 2011 (Fig. S9), suggesting that the NH₃ emission estimates in 2011 NEI may be biased low in other months, too.

3.3 Evaluation of the optimized emission estimates against independent datasets

The robustness of the NH₃ emission optimization is evaluated by comparing the model outputs using both the a priori and optimized emission estimates with independent observations. The bias and uncertainties inherited in the CMAQ forward model and its adjoint, as well as the assumptions made about the uncertainties of the a priori emission inventory and IASI-NH₃ observations, will all influence the robustness. Here, we choose to evaluate the outputs against (1) biweekly average ambient NH₃ concentrations measured by AMoN and (2) weekly average NH${}_{\mathrm{4}}^{+}$ wet deposition measured by NADP (Fig. 4).

Figure 4 Evaluation of the simulated NH₃ surface concentration (a, b, c) and NH${}_{\mathrm{4}}^{+}$ wet deposition (d, e, f) against biweekly NH₃ concentration observations from AMoN and weekly NH${}_{\mathrm{4}}^{+}$ wet deposition observations from NADP, respectively. The orange circles and blue dots represent comparison using the a priori and optimized NH₃ emission estimates, respectively. Summary statistics including sample size (N), normalized mean bias (NMB), normalized root mean square error (NRMSE), least square error regression slope and intercept, and R square (R²) for all comparisons are listed below the plots.

In general, the optimized NH₃ emission reduces the negative NMB when comparing the CMAQ outputs with AMoN NH₃ concentration for all 3 months. There is a greater improvement at the high-concentration end than the low-concentration end because both the IASI satellite and the passive samplers at the AMoN sites have higher uncertainties in areas with low NH₃ abundance (Van Damme et al., 2015a; Puchalski et al., 2011). Yet the NRMSE gets higher and R² gets lower in April, indicating a higher spatial variation in the residuals. There is an over-adjustment for sites in Pennsylvania in April, where there is a hotspot observed by IASI on 14 and 15 April. The hotspot possibly came from a large transported plume at a higher altitude from the central US to Pennsylvania (Figs. S10 and S11), which is not measured by ground observations at AMoN sites at biweekly resolution. If that is the case, the hybrid inverse modeling framework would have difficulties in reproducing the long-range transport contribution for two reasons. First, local emissions in Pennsylvania would be enhanced in the IMB inversion, and inter-grid transport were neglected at 216 km by 216 km resolution. Second, the following 4D-Var inversion very likely reached a local optimal by adjusting emissions from local and surrounding grid cells near the observed hotspot rather than grid cells at distance. Furthermore, the IASI-NH₃ column densities may be overestimated because vertical profiles with the highest concentrations near the surface were assumed in the retrieval process (Whitburn et al., 2016).

For evaluation against NADP observations, there is a noticeably improved agreement in April, with reduced negative NMB and reduced discrepancies for most of the data pairs. For July, the emission optimization only slightly improved the model performance. For October, the optimization increased the NMB from −1.8 to 4.8 %. This indicates that NH₃ emission is not the dominant explanatory factor for bias in simulated NH${}_{\mathrm{4}}^{+}$ wet deposition that is commonly observed in chemical transport models (Appel et al., 2011; Paulot et al., 2014). A better representation of the cloud, precipitation, and deposition processes in the WRF Model and CMAQ is needed to close the gap between simulated and observed NH${}_{\mathrm{4}}^{+}$ deposition amount. Overall, the improved model operational performance for ambient NH₃ suggests that the inverse model optimization applied in this study provides improvements in the NH₃ emission estimates during all 3 months in most of the CONUS, except in Pennsylvania and surrounding regions in April. The hybrid inverse modeling technique may over-adjust local emissions in hotspots dominated by long-range transport.

4 Implications

4.1 Ambient aerosol concentration

As a major precursor of ambient aerosol formation, the NH₃ emission inventory is believed to be a major source of uncertainty in PM_2.5 assessment in several parts of the CONUS (Henze et al., 2009; Schiferl et al., 2014; Heald et al., 2012), which can further bias the source contribution assessments on PM_2.5-related health impacts (Lee et al., 2015; Zhao et al., 2020). Comparison of the simulated PM_2.5 mass concentration using the a priori and optimized NH₃ emission estimates shows that the NH₃ emission bias in April is a major factor for bias in the modeled PM_2.5 concentration leading to high or low bias in ammonium nitrate (NH₄NO₃) formation (Fig. 5). The relative change of the monthly average PM_2.5 concentration is over 5 % in one-fifth of the CONUS, including an increase in the Northeast, the Pacific West, the Rocky Mountains, part of Texas, and the Gulf Coast region, and a decrease in the Midwest. For most of these regions, over 90 % of the change is driven by the change in concentration of NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$.

Figure 5 The changes in monthly average PM_2.5, NH${}_{\mathrm{4}}^{+}$, and NO${}_{\mathrm{3}}^{-}$ mass concentration in April due to the NH₃ emission adjustment in the optimized estimates. The change is defined as conc_optimized−conc_a priori, where conc_optimized and conc_a priori represent the simulated monthly average mass concentration using the optimized and a priori NH₃ emission estimates, respectively. The difference between the observed NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ mass concentration and simulated concentrations using the a priori NH₃ emission (conc_obs−conc_a priori, where conc_obs represents the observed monthly average mass concentration) is overlaid using colored dots with the same color scheme.

Comparison of the simulated monthly average NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ concentration using the a priori estimates against ambient monitoring network data (USEPA, 2018) shows that there is a high bias in the Midwest region and Pennsylvania state, and a low bias for the rest of the sites (Table 1). Simulations using the optimized NH₃ emission estimates reduce the high bias in the Midwest region but exacerbate the high bias in Pennsylvania state and the surrounding areas. For the other sites, the impact of optimization is mixed but minor in general.

Table 1 Statistical summary of the correlation between simulated monthly average NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ concentrations and observations in April^*.

^* The correlation between observed concentrations and simulated ones based on a priori and optimized NH₃ emission estimates is compared. The sites are grouped as the Midwest region, Pennsylvania state and surrounding areas, and other areas.

For the Midwest, our optimized NH₃ emission is 12 % lower than the 2011 NEI, leading to a 5–30 % decrease in NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ concentration. Overestimation of NO${}_{\mathrm{3}}^{-}$ in the Midwest has been recognized in previous model evaluations. Previous studies have attempted to moderate the high bias by lowering the nitric acid (HNO₃) concentration through either lowering both the daytime and nighttime HNO₃ formation rate or raising the deposition removal rate (Heald et al., 2012; Zhang et al., 2012; Walker et al., 2012). It was found that such modification in the model parameterization cannot fully account for the overestimation (Heald et al., 2012; Zhang et al., 2012; Walker et al., 2012). Our study implies that the springtime overestimation can partly be explained by the overestimation in NH₃ emissions which drives the high bias in NH₄NO₃ formation.

The large increase of the NH₄NO₃ concentration in Pennsylvania state and the surrounding areas is due to the over-amplified local NH₃ emissions in the optimized estimates to match the high NH₃ column density in IASI-NH₃ 2011, as discussed earlier. It leads to a higher magnitude of biases in NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ concentration as compared to ground measurements. The fact that the simulated ambient NH₃ concentration, NH${}_{\mathrm{4}}^{+}$ concentration, and NH${}_{\mathrm{4}}^{+}$ wet deposition using the optimized NH₃ estimates are biased high in comparison with independent ground measurements suggests the enhanced NH₃ abundance observed from IASI is possibly driven by long-range transport at higher altitudes instead of local surface emissions.

For the rest of the CONUS, there is only a slight impact of the optimization on simulated NH₄NO₃ formation. For example, although the NH₃ emission is doubled in the San Joaquin Valley in California, the modeled NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ concentrations are still biased low using the optimized estimates. A sensitivity test using GEOS-Chem shows that the San Joaquin Valley region is nitric-acid-limited instead of ammonia-limited (Walker et al., 2012), suggesting that there is an underestimation in HNO₃ formation. A comparison of the simulated and measured speciated PM_2.5 shows that there is a low bias in non-volatile cation concentrations at the sites in the San Joaquin Valley, limiting the formation of NH₄NO₃ through gas–particle partitioning (Chen et al., 2019). Thus, attempts to close the gap between the simulated and monitored NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$ concentrations by scaling NH₃ emission alone are ineffective and might lead to an overestimation in local NH₃ emissions.

For July and October, there is a very limited difference between the simulated PM_2.5 concentration using the optimized and a priori NH₃ emission estimates, as expected, because the change in NH₃ emission is minor. There are only 1 and 4 % of the CONUS with a relative change in PM_2.5 concentration over 5 %, respectively. This result shows that the uncertainty in NH₃ emission estimates is moderate and is not a major contributor to biases in modeled PM_2.5 in July and October.

4.2 Reactive nitrogen deposition

The uncertainties in NH₃ emission inventory also impact the reactive nitrogen (Nr) deposition assessment, which informs the ecosystem impacts evaluation and effective mitigation actions (Ellis et al., 2013). To evaluate the impact of the NH₃ emission optimization on simulated Nr deposition, the Nr deposition amount simulated using optimized and a priori emission estimates is analyzed in all biodiversity-protected areas designated by the USGS (Fig. S12) within the CONUS (USGS, 2018). In total, the Nr deposition increased by 27, 9, and 5 % on average in these protected areas in April, July, and October, respectively. A regional comparison based on the Level I ecoregions (Pardo et al., 2015) shows that the deposition increment is the highest in the Tropical Wet Forests (+64 %), followed by the Great Plain region (+46 %), in April (Fig. 6). Although the overall increase is small in July and October, the increment can be high in individual ecoregions, including Southern Semiarid Highlands (+95 % in July) and Temperate Sierras (+62 % in July). In addition to the increment in deposition amount, higher NH₃ emission, especially in intensive agriculture regions, may indicate higher source contribution from agricultural NH₃ than previous estimates (Lee et al., 2016).

Figure 6 The changes in the simulated monthly reactive nitrogen (Nr) deposition amount in protected areas for biodiversity conservation caused by the emission adjustment in April, July, and October. For each month, the left bar is for the a priori deposition amounts and the right bar is for the optimized deposition amounts. The deposition is grouped for 10 Level I ecoregions defined by the Commission for Environmental Cooperation, including Northern Forests (NF), Great Plains (GP), Northwestern Forested Mountains (NFM), Marine West Coast Forest (MWCF), North American Deserts (NAD), Mediterranean California (MC), Southern Semiarid Highlands (SSH), Temperate Sierras (TS), and Tropical Wet Forests (TWF).

Driven by the increase in the reduced form of Nr (NH₃ and NH${}_{\mathrm{4}}^{+})$ deposition, a higher share of the reduced form of Nr to the total Nr deposition is found in most of the ecoregions for all 3 months than in the NEI-based estimates. More detrimental impacts on sensitive species and biodiversity are expected when this change in dominant Nr form is considered in addition to the increase in magnitude because the growth of many sensitive plant species will be inhibited by a high NH${}_{\mathrm{4}}^{+}$-to-NO${}_{\mathrm{3}}^{-}$ ratio in soil and water (Bobbink and Hicks, 2014).

5 Conclusions

We apply the newly developed multiphase adjoint of the CMAQ v5.0 chemical transport model and NH₃ column observations from the satellite-borne IASI to optimize NH₃ emissions estimates in the CONUS using a hybrid inversion modeling approach. The approach consists of a coarse-resolution iterative mass balance inversion (216 km by 216 km) and a fine-resolution 4D-VAR inversion (36 km by 36 km) and is performed using IASI-NH₃ observations in April, July, and October. The hybrid approach overcomes the over-adjusting problem for high-emission areas in the direct 4D-Var method and reduces the computational cost, but it may introduce over-adjustment in special cases where the NH₃ abundance is dominated by transport instead of local emissions.

We use the NH₃emission from 2011 NEI, commonly used in regional and national simulations and assessments as the a priori emission. We find that the optimized NH₃ emission inventory differs greatly with the 2011 NEI in April. The emission in the Midwest is overestimated and the emission in the Southern States is underestimated in the 2011 NEI. Overall, the optimized emission is 35 % higher in April. The optimized emission estimates in July and October are also higher (18 and 10 %) than the 2011 NEI estimates, but the spatial distribution agrees well. The IASI-NH₃ observations indicate a consistent underestimation of NH₃ emissions in California's San Joaquin Valley in all 3 months; however, the inverse modeling fails to properly scale up the emissions in July. The evaluation of simulation outputs against ground measurements including ambient NH₃ concentrations from AMoN and NH${}_{\mathrm{4}}^{+}$ wet deposition from NADP shows that the optimized NH₃ emission estimates reduce the NMB between model outputs and independent observations, especially in April. The NRMSE remains high, indicating (1) the potential to further optimize NH₃ emission estimates when more representative observations of ambient NH₃ abundance become available and (2) the need to address the uncertainties in other processes affecting the NH₃ abundance, such as gas–particle partitioning, dry and wet deposition, and in-cloud processes.

Application of the optimized NH₃ emission estimates also yields a better agreement between the simulated and observed PM_2.5 concentration in April in the Midwest region by improving the model performance on simulated NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{3}}^{-}$. This is consistent with previous findings that the uncertainty in NH₃ emission is a key factor limiting the model performance of PM_2.5. The optimized NH₃ emission estimates in general increase the Nr deposition amount and the relative importance of reduced-form Nr, highlighting the importance of constraining NH₃ emission estimates for accurately assessing nitrogen deposition and ecosystem health over sensitive regions.