Inverse modelling of NO x emissions over eastern China : uncertainties due to chemical non-linearity

Satellite observations of nitrogen dioxide (NO2) have often been used to derive nitrogen oxides (NOx =NO+NO2) emissions. A widely used inversion method was developed by Martin et al. (2003). Refinements of this method were subsequently developed. In the context of this inversion method, we show that the local derivative (of a first-order Taylor expansion) is more appropriate than the “bulk ratio” (ratio of emission to column) used in the original formulation for polluted regions. Using the bulk ratio can lead to biases in regions of high NOx emissions such as eastern China due to chemical non-linearity. Inverse modelling using the local derivative method is applied to both GOME-2 and OMI satellite measurements to estimate anthropogenic NOx emissions over eastern China. Compared with the traditional method using bulk ratio, the local derivative method produces more consistent NOx emission estimates between the inversion results using GOME-2 and OMI measurements. The results also show significant changes in the spatial distribution of NOx emissions, especially over high emission regions of eastern China. We further discuss a potential pitfall of using the difference of two satellite measurements to derive NOx emissions. Our analysis suggests that chemical non-linearity needs to be accounted for and that a careful bias analysis is required in order to use the satellite differential method in inverse modelling of NOx emissions.


Introduction
Nitrogen oxides (NO x = NO + NO 2 ) play an important role in tropospheric photochemistry of ozone and secondary aerosol formation.They are normally emitted from both anthropogenic (e.g.fossil fuel combustion) and natural sources (e.g.lightning, soil and wild fire).The rapid economic growth in East Asia has led to a significant increase of energy consumption, thereby increasing anthropogenic NO x emissions, during last two decades (Boersma et al., 2008;Ghude et al., 2009;Gu et al., 2013;Ma et al., 2006;Mijling et al., 2013;Richter et al., 2005;Stavrakou et al., 2008;van der A et al., 2006;Zhang et al., 2007).
The traditional bottom-up emission inventory relies on detailed information about sources and emission factors and, therefore, can have large uncertainties, especially in countries such as China where the emission information is incomplete (Streets et al., 2003).In recent years, satellite measurements of nitrogen dioxide (NO 2 ) from multiple instruments, including Global Ozone Monitoring Experiment (GOME), Scanning Imaging Absorption Spectrometer for Atmospheric Chartography (SCIAMACHY), Ozone Monitoring Instrument (OMI) and GOME-2, have been widely used to derive top-down estimation of NO x emissions as an alternative to the bottom-up emission inventory (Gu et al., 2014;Jaegle et al., 2005;Lamsal et al., 2011;Lin et al., 2010;Stavrakou et al., 2008;Zhang et al., 2012).
Inverse modelling of NO x emissions, incorporating model simulations and satellite observations, can help reduce uncertainties in bottom-up inventories, especially in China (Lin et al., 2012;Mijling and van der A, 2012;Zhao and Wang, 2009).Martin et al. (2003) developed the first monthly inversion method by scaling the bottom-up emission inventory with top-down constraints from GOME NO 2 column measurements.Zhao and Wang (2009) improved the method by carrying out the emission inversion iteratively on a daily ba-sis.Gu et al. (2014) further refined this method by including column NO 2 retrieval in inverse modelling, which removes the biases introduced by inconsistency between NO 2 profiles in the retrieval and inverse modelling.Lin et al. (2010) developed a method of inverse modelling using column NO 2 difference between GOME-2 and OMI.
In the emission inversion problem, the top-down emission is estimated by modifying the a priori emission with a term proportional to the difference between the simulated and observed column density, which can be expressed as follows: where E is the top-down emission, E a is the a priori emission used in model simulation, s is the satellite-observed column, m is the model-simulated column and α is the correction rate.In the original formulation by Martin et al. (2003) and all subsequent studies using this or some variants of this method (e.g.Jaegle et al., 2005), α is calculated as the ratio between the a priori emission and simulated column (referred to as the bulk ratio hereafter).
Note that this formulation has an implicit assumption that the emission is linearly proportional to the column density.However, previous studies of the emission and column trends suggested that the non-linearity between NO x emission and tropospheric NO 2 column is non-trivial (Gu et al., 2013;Lamsal et al., 2011;Lu and Streets, 2012;Stavrakou et al., 2008).This non-linear relationship is mainly due to the nonlinear photochemical feedbacks between NO x and OH (the increase in NO x promotes OH production and reduces its lifetime in the low-emission condition, but suppresses OH production and increases its lifetime in the high-emission condition).Figure 1 shows tropospheric NO 2 column densities as a function of surface emissions at GOME-2 and OMI overpass time over eastern China simulated in the Regional chEmical and trAnsport Model (REAM) (model details are given in Sect.2.2).Clearly, tropospheric NO 2 column is not linearly proportional to NO x emission, which was assumed in Eq. ( 2).Another illustration of this problem for inverse modelling over China can be found in Gu et al. (2013).
The difference in the non-linearity in the ratio of emission to column NO 2 between GOME-2 and OMI overpass time reflects in part the difference in the relative importance between transport and chemistry at a different time of day.
In this study, we examine whether considering the nonlinear ratio between NO x emissions and NO 2 columns can improve the inverse modelling results and reduce the discrepancies in emission estimates using different satellite observations.Applying the Taylor expansion to Eq. ( 1), we define a local derivative in place of the bulk ratio (Eq.2).We implement the local derivative of emission to column NO 2 ratio in the REAM model and examine its effects on the inverse modelling estimates of surface anthropogenic NO x emissions with GOME-2 and OMI measurements over China in August 2007.The new inversion results are compared to those using the bulk ratio method and the satellite differential method (Lin et al., 2010) and the implications for NO x inverse modelling are discussed.
2 Satellite data and inverse modelling method

Satellite data
We use the measurements from GOME-2 and OMI instruments in this study.Both instruments are nadir-viewing spectrometers (Boersma et al., 2004(Boersma et al., , 2007(Boersma et al., , 2011)).The OMI instrument was launched on board the Aura satellite in July 2004, and it has a spatial resolution of 24 × 13 km 2 at nadir (Levelt et al., 2006).The GOME-2 instrument used in this study was launched onboard the MetOp-A satellite in June 2006 with ground pixel size of 80 × 40 km 2 (Irie et al., 2012).The local overpass time across the equator is around 09:30 LT for GOME-2 and around 13:30 LT for OMI.
We use the Royal Netherlands Meteorological Institute (KNMI) OMI (DOMINO2 v2.0) and GOME-2 (TM4NO2A v2.3) tropospheric NO 2 vertical column density (VCD) products for this study.We excluded the data flagged with row anomalies from the OMI measurements (http://www.knmi.nl/omi/research/product/rowanomaly-background.php).To reduce the cloud interference, we only use measured NO 2 column data when cloud fraction is < 20 % for both measurements.The NO 2 VCD products were gridded to the REAM model domain with a 36 × 36 km 2 horizontal resolution.
The error in the retrieval of NO 2 tropospheric VCD is determined by those in total slant column density (SCD), stratospheric SCD and tropospheric air mass factor (AMF) estimation.In this study, we use total and stratospheric SCD errors from KNMI DOMINO2 and TM4NO2A products and compute the uncertainty of tropospheric AMF estimation following the KNMI algorithm.The details of error analysis were described by Boersma et al. (2004Boersma et al. ( , 2007Boersma et al. ( , 2011) ) and Hains et al. (2010).In general, the uncertainties of total and stratospheric SCD estimations are small (< 0.7 × 10 15 mol cm −2 ) relative to high tropospheric VCDs over eastern China (Zhao and Wang, 2009).The uncertainty in tropospheric AMF comes from surface albedo, cloud fraction, cloud pressure and profile shape.The uncertainty from an a priori profile can lead to ∼ 10 % error in tropospheric VCD retrievals (Martin et al., 2003;Zhao and Wang, 2009).The total uncertainty of an individual retrieval is up to 50 % over highly polluted eastern China for both OMI and GOME-2.More detailed discussion of the retrieval uncertainties and their effects on NO x emission inversion can be found in the studies by Boersma et al. (2004, 2007, 2011), Martin et al. (2003) and Zhao and Wang (2009).

REAM model
The 3-D REAM has been applied in a number of tropospheric chemistry and transport studies over East Asia, North America and polar regions (Choi et al., 2005(Choi et al., , 2008a, b;, b;Gu et al., 2013Gu et al., , 2014;;Jing et al., 2006;Liu et al., 2010Liu et al., , 2012a, b;, b;Wang et al., 2006Wang et al., , 2007;;Yang et al., 2011;Zeng et al., 2003Zeng et al., , 2006;;Zhao et al., 2009aZhao et al., , b, 2010)).The model has a horizontal resolution of 36 × 36 km 2 with 30 vertical layers in the troposphere.Transport is driven by WRF assimilated meteorological fields constrained by the NCEP reanalysis products (http://www.esrl.noaa.gov/psd/).The chemistry mechanism in the REAM is adopted from the GEOS-Chem model (Bey et al., 2001) with updates of kinetics data (http://jpldataeval.jpl.nasa.gov/).The anthropogenic NO x and VOCs emissions are from Zhang et al. (2009).The biomass burning emissions are taken from the Global Fire Emissions Database, Version 2 (GFEDv2.1;available at http://daac.ornl.gov/).The algorithm of soil NO x emissions follows Yienger and Levy (1995) as adopted by Wang et al. (1998).The lightning NO x emission is parameterised as in Zhao et al. (2010).The simulation was conducted in August 2007 in this study.For consistence with both OMI and GOME-2 retrievals, the same model-simulated NO 2 VCDs on each day (at corresponding satellite overpass time) are collected with an averaging kernel for inverse modelling with both satellite products.

Inverse modelling method
In the bulk ratio inverse modelling, we use the same framework by Martin et al. (2003) outlined above.As discussed in the introduction, the bulk ratio (α) used in the traditional method is based on the assumption of a linear relationship between NO x emission and tropospheric column NO 2 , which is only accurate under a low emission condition.As shown in Fig. 1, using the bulk ratio for inverse modelling overlooks the non-linear chemical feedback of NO x lifetime at different satellite overpass times, which not only affects the accuracy of inverse modelling results but also leads to inconsistencies between emission inversions using different satellite measurements.To account for the non-linearity, we apply the Taylor expansion to Eq. ( 1) and obtain a local derivative ratio (α * ) to replace Eq. ( 2): where E a is the change of the a priori emission and m is the change of model-simulated column.
The Taylor expansion formulation of Eqs. ( 1) and ( 3) is equivalent to the previous works by Vinken et al. (2014a, b) and Castellanos et al. (2014), who accounted for the nonlinear chemistry effect by adding the ratio of relative emission to relative column NO 2 changes, first introduced by Lamsal et al. (2011), to the formulation of Eqs. ( 1) and (2).For inversion of ship emissions, Vinken et al. (2014b) computed averaged non-linear factors for selected regions with perturbations proportional to the model-observation column difference.For polluted regions, the effect of profile change is small (Zhang et al., 2016) and is not included in this study.In the previous studies, model sensitivities are computed with domain-wide emission perturbations.This approach introduces uncertainties since the column change in a given grid cell is affected by emission changes in both that grid cell and upwind grid cells; the effects of upwind grid cells are a complex function of emissions, chemistry and transport.As in Gu et al. (2013), we carried out single-cell-based emission perturbation in the 3-D REAM model to compute the value of α * .We first archived the 3-D influxes of all model tracers at each time step in the standard simulation.In the perturbation simulation (15 % of anthropogenic NO x emissions), the influxes of all model tracers for all grid cells were replaced with the archived values at each time step.Consequently, the emission perturbation only affects NO x chemistry, out flux and concentration in the same grid cell.Using the perturbation and standard simulation results, we computed the value of α * and applied it to Eq. ( 1) to estimate inverted NO x emission over eastern China in August 2007 with GOME-2 and OMI observations.The results from the local derivative method are compared with those using the bulk ratio method in later sections.
The uncertainties of the a posteriori emissions come from those in a priori and top-down emission estimates (Martin et al., 2003).Uncertainties in top-down emission estimates are derived from those in tropospheric NO 2 VCD retrievals and model simulations.The retrieval uncertainty is discussed in Sect.2.1.The uncertainty of model simulation is estimated at 30 % and that of the bottom-up inventory is ∼ 60 % over China (Zhao and Wang, 2009)  the a posteriori emission is typically in the range of 20-40 % over polluted eastern China (Gu et al., 2014;Zhao and Wang, 2009).
3 Results and discussion

Comparisons between α and α *
Figure 2 compares the relative difference between local derivative ratio (α * ) and bulk ratio (α) values over eastern China for August 2007, at GOME-2 and OMI overpass times respectively.For both satellites, α * is higher (> 20 %) than α over most low-emission rural regions but lower (−20 ∼ −60 %) than α over most high-emission regions including eastern coastal areas and Sichuan province.For emission estimates, the inversion biases over high emission regions tend to be more important than rural regions.The high bias of α relative to α * implies that top-down NO x emission estimates tend to over-correct for a given difference between observed and simulated column NO 2 .We note here that the inversion bias can be either positive or negative depending on the column difference.

NO x emission inversion consistency between OMI and GOME-2
Another way to look at the effects of correction biases using the bulk ratio relative to local derivative method is to compare the inversion estimated NO x emissions using OMI and GOME-2 measurements.Gu et al. (2014) showed that inversion results using standard DOMINO products of GOME-2 tend to be higher than OMI possibly due to a bias in the TM4 NO 2 profiles used in GOME-2 retrievals.Coupling this tendency of a high retrieval bias of GOME-2 with the different sensitivities discussed in the previous paragraph implies that the bulk ratio formulation would tends to overcorrect and estimate higher NO x emissions using GOME-2 than OMI measurements in polluted regions.The bias expectation is confirmed in Fig. 3.While the inversion results using DOMINO GOME-2 and OMI products generally show higher NO x emissions in the former, the difference between the inversions using two satellites is smaller using the local derivative than bulk ratio method, particularly over high emission regions.For example, we compare the emission estimates between the bulk ratio and local derivative methods in three largest megacities in China (e.g.Beijing, Shanghai and Guangzhou).Using the bulk ratio method, the relative difference of inversion emission estimates of using two satellites products relative to OMI results are 54.2, 70.5 and 55.6 % for the three megacities.The local derivative method reduces the corresponding relative differences down to 9.0, 5.3 and 13.5 % (Fig. 3 and Table 1).These results demonstrate that a significant fraction of the discrepancy in inverted NO x emissions between different instruments can be attributed to neglecting the chemical non-linearity in the traditional bulk ratio method and that improvement can be achieved with the local derivative method that we proposed in this study.The remaining discrepancy between inversion emission estimates using GOME-2 and OMI observations is likely due to NO 2 profiles used in the retrieval and possible systematic biases between the two satellite instruments (Gu et al., 2014).

Comparison of inversion results between using the bulk ratio method and the satellite differential approach
Inversion can be applied to the column NO 2 measurements from two satellites (Lin and McElroy, 2010).Lin et al. (2010) Figure 3. Relative difference of inversion emission estimates for NO x with GOME-2 relative to OMI products using the local derivative and the bulk ratio methods: regional distributions over eastern China for August 2007 (left) and the differences in three megacities (right).developed a method to utilise information from two different satellites (referred to as the satellite differential method in this study) to improve the emission estimates.The formulation of the method can be simplified as follows, where E j is the a posteriori emission, E a is the a priori emission, OMI is OMI-observed NO 2 column, GOME-2 is GOME-2-observed NO 2 column, m1330 is model-simulated NO 2 column at OMI overpass time, m0930 is modelsimulated NO 2 column at GOME-2 overpass time, t is the time gap between two satellite overpass time and τ is the lifetime of NO x .The derivation of this method is analogous to the bulk ratio method (Lin et al., 2010) and is not directly comparable to the local derivative method.We, therefore, compare the inversion results using the satellite differential method with those using the bulk ratio method.
In Table 1, the inversion emission estimates using the satellite differential method are compared with those using the bulk ratio method over eastern China in August 2007.
As discussed in the previous section, the bulk ratio method leads to consistently high emission estimates using GOME-2 products compared to using OMI over high emission ratios, e.g.50-70 % in the three megacities.The estimates using the satellite differential method are even lower than the bulk ratio inversion estimates using OMI products.The reason is that in our analysis GOME-2 / m0930 is consistently larger than OMI / m1330 over eastern China (Fig. 4).The mathematical analysis in Appendix A shows that, consequently, a posteriori emission estimates using the satellite differential method are consistently lower than those using the bulk ratio method with either GOME-2 or OMI products.Lin and McElroy (2010) also noticed that the original Lin et al. (2010) method can lead to systematic inversion biases as we demonstrated here.They corrected the biases by not using the inversion results if the retrieved OMI and GOME-2 are greater than 120 % of the simulated corresponding values.When the differences between the two satellite products are not well characterised, we find that the local derivative method is more robust relative to the bulk ratio and satellite differential methods (Table 1).

Conclusions
We show in this study that the non-linearity of NO x chemistry implies that the local derivative (of a first-order Taylor expansion) is better suited in the inversion formulation developed by Martin et al. (2003) than bulk ratio, in agreement with Vinken et al. (2014a, b) and Castellanos et al. (2014).In this study, single-grid-cell-based perturbation sensitivity calculation was used instead of previous domain-wide perturbations such that upwind emission changes do not affect local derivative estimates.The latter effects can be more appropriately accounted for using the iterative method by Zhao and Wang (2009) and Gu et al. (2014).In the context of the inversion formulation by Martin et al. (2003), we compared the bulk ratio and local derivative methods in inverse modelling of anthropogenic NO x emissions over eastern China for August 2007.At the observation time of OMI and GOME-2, the local derivative ratios (α * ) are smaller (−20 ∼ −60 %) than the bulk ratios (α) over most high-emission regions, but are higher (> 20 %) than bulk ratios (α) over most lowemission rural regions.Over high emission regions, the inversion emission estimates using the local derivative method produces more consistent results between OMI and GOME-2 products than the bulk ratio method.In our work, the observed to simulated tropospheric column NO 2 are consistently higher for GOME-2 than OMI over high emission regions of eastern China, leading to a consistent low bias in a posteriori emission estimates by the satellite differential method relative to those by the other methods.Computationally, the local derivative ratio is more complex to compute than bulk ratio.The iterative method in the previous work by Zhao and Wang (2009) and Gu et al. (2014) can largely mitigate the biases introduced by the bulk ratio method, although the inversion convergence time will likely be reduced when using local derivative than bulk ratios.For certain applications such as deriving the timeline of emission reduction during the Beijing Olympics (Yang et al., 2011), reducing the inversion convergence time is critically important.Further studies are needed to quantify the improvements.
Table 1 shows that the inversion results using the satellite differential method are consistently lower than that of the bulk ratio method using OMI measurements, which are lower than the bulk ratio method using GOME-2 measurements over high emission regions.Lin and McElroy (2010) suggested it may be introduced by systematic errors in satellite retrievals.Here we show the mathematical analysis demonstrating the reasons for the consistent difference.In our analysis, the following inequality is found over high emission regions of eastern China (Fig. 4): Given a priori emission estimation of Eqs. ( 2) and (4), Eq. (A2) implies that E j < E b,OMI < E b,GOME-2 , (A3) where E j is the inversion emission estimate by the satellite differential method, E b,OMI is the inversion emission estimate by the bulk ratio method using OMI measurements, and E b,GOME-2 is the inversion emission estimate by the bulk ratio method using GOME-2 measurements.

Figure 1 .
Figure 1.Simulated NO 2 column density as a function of surface NO x emission at GOME-2 (black) and OMI (yellow) overpass time over eastern China for August 2007.The dots are grid average data from REAM simulations binned by a column NO 2 interval of 6 × 10 14 mol cm −2 and the solid lines are from least-square polynomial regression results.

Figure 2 .
Figure 2. Relative difference between monthly mean local derivative ratio α * and bulk ratio α, defined as 1 − α/α * , for GOME-2 (left) and OMI (middle) in REAM simulations.The right panel shows NO x emission distribution estimated from GOME-2 observations by using the local derivative method over eastern China for August 2007.

Table 1 .
NO x inversion emission estimates for August 2007 using the local derivative, bulk ratio and satellite differential inversion methods.