Direct estimates of biomass burning NO x emissions and lifetimes using daily observations from TROPOMI

. Biomass burning emits an estimated 25% of global annual nitrogen oxides (NO x ), an important constituent that participates in the oxidative chemistry of the atmosphere. Estimates of NO x emission factors, representing the amount of NO x per mass burned, are primarily based on field or laboratory case studies, but the sporadic and transient nature of wildfires 10 makes it challenging to verify whether these case studies represent the behavior of the global fires that occur on earth. Satellite remote sensing provides a unique view of the earth, allowing the study of emission and downwind evolution of NO x from a large number of fires. We describe direct estimates of NO x emissions and lifetimes for fires using an exponentially modified Gaussian analysis of daily TROPOspheric Monitoring Instrument (TROPOMI) retrievals of NO 2 tropospheric columns. We update the a priori profile of NO 2 with a fine-resolution (0.25˚) global model simulation from NASA’s GEOS Composition 15 Forecasting System (GEOS-CF), which largely enhances NO 2 columns over fire plumes. We derive representative NO x emission factors for six fuel types globally by linking TROPOMI derived NO x emissions with observations of fire radiative power from Moderate Resolution Imaging Spectroradiometer (MODIS). Satellite-derived NO x emission factors are largely consistent with those derived from in-situ measurements. We observe decreasing NO x lifetime with fire intensity, which we infer is due to the increase in both NO x abundance and hydroxyl radical production. Our ﬁndings suggest promise for applying 20 space-based observations to track the emissions and chemical evolution of reactive nitrogen from wildfires.

Satellite instruments observe fire NOx plumes as a mixture of fresh and aged smoke. NOx is a short-lived species, and its 65 concentration will decay in the plume due to the formation of nitric acid (HNO3) and other reservoir species. The relationship between satellite observed NOx concentration and emissions depends on the loss rate of NOx. The chemical processes governing the lifetime of NOx in the fire plumes are poorly understood (Alvarado et al., 2010). Previous studies assume a constant chemical NOx lifetime of 2 hours (Mebust et al., 2011;Mebust and Cohen, 2014). Laughner and Cohen (2019) provide space-based evidence of changing NOx lifetime over U.S. cities as NOx emissions decline. As fire intensity varies by several 70 orders of magnitude, assuming constant NOx lifetime for all fires will likely introduce errors in the derived NOx emissions (De Foy et al., 2014). The improved spatial resolution of TROPOMI allows direct measurements of the length scale of NO2 decay.
By analyzing the plume evolution downwind, we derive an effective NOx lifetime. Beirle et al. (2011) first proposed an exponentially modified Gaussian (EMG) approach to directly estimate NOx emissions and lifetime from satellite observations, which has been widely used to derive NOx emissions from anthropogenic sources (Beirle et al., 2011;Lu et al., 2015;Goldberg 75 et al., 2019;Laughner and Cohen, 2019). Our study is the first to apply the EMG approach to simultaneously estimate NOx emissions and lifetime from biomass burning plumes. The resulting emission estimates provide a large ensemble with which to evaluate current emission models and also provide detailed constraints on the chemical evolution of NOx. The resulting lifetimes provide insights into hydroxyl radical abundances in the plume and thus constraints on the lifetime of other chemicals emitted from fires. 80

Datasets
TROPOMI is a nadir-viewing hyperspectral spectrometer launched on October 13, 2017 by the European Space Agency (ESA) for the European Union's Copernicus Sentinel 5 Precursor (S5p) satellite mission. TROPOMI provides afternoon global observations in the UV−visible−near infrared−shortwave spectra with a fine spatial resolution of 7 × 3.5 km 2 at nadir (increased to 5.5 × 3.5 km 2 since August 2019). We obtain the daily Level-2 TROPOMI retrievals of NO2 tropospheric column density 85 data from April 2018 to June 2020 from NASA Goddard Earth Sciences (GES) Data and Information Services Center (DISC).
The retrieval of the NO2 tropospheric vertical column includes three steps: (1) retrieval of the total slant column density along the optical path using differential optical absorption spectroscopy (Boersma et al., 2011).; (2) subtraction of the total slant column density from stratospheric NO2 slant column based on information from a data assimilation system (Boersma et al., 2018); (3) conversion of the tropospheric slant column density to vertical column density using air mass factors (AMFs), which 90 are obtained from radiative transfer calculations that account for the viewing geometry, cloud fraction, surface properties, and the a priori vertical profile of NO2 (Boersma et al., 2018). We include TROPOMI observations with the quality assurance value higher than 0.5, which filter out problematic retrievals but still keep good quality retrievals over cloud (or aerosols). In addition to NO2, we obtain TROPOMI aerosol layer height (ALH) or pressure (ALP) data, which provides height information of aerosol layer in the troposphere. Retrieval of ALH or ALP is based on the O2 absorption band at near-infrared wavelengths 95 between 759 and 770 nm (Graaf et al., 2019). Details of the aerosol layer retrieval algorithms can be found at Graaf et al., (2019) and Nanda et al. (2019).
We use the Moderate Resolution Imaging Spectroradiometer (MODIS) Active Fire products (Collection 6) to provide information on the intensity and location of fires (Giglio et al., 2016), which are available from NASA's Fire Information for Resource Management System (FIRMS). We include daytime MODIS measurements from the Aqua satellite to match with 100 the overpass of TROPOMI. Fire detection from MODIS is performed using a contextual algorithm that measures the infrared radiation from fires (Giglio et al., 2016). Each hotspot is recorded as the center of a ~1 × 1 km 2 pixel that contains one or more fires, and the FRP is estimated via an empirical relationship using the 4µm band brightness temperatures (Kaufman et al., 1998). We group fire pixels whose distances are within 20 km as a single fire event, and the center of the fire is calculated as the mean of fire pixel locations weighted by pixel FRP. To assess the sensitivity to the choice of FRP product, we also process 105 the daytime Suomi NPP Visible Infrared Imaging Radiometer Suite (VIIRS) observations of active fire accessed from NASA's FIRMS. VIIRS fire product uses a similar algorithm as MODIS for fire detection (Schroeder et al., 2014). To assess potential effects of aerosol from plumes on satellite retrieval of NO2, we acquire the Multi-angle Implementation of Atmospheric Correction (MAIAC) Aerosol Optical Depth (AOD) Level-2 1-km daily gridded product (MCD19A2) from NASA's Earth Observing System Data and Information System (EOSDIS). Details of the retrieval of AOD can be found at Lyapustin et al. 110 (2012Lyapustin et al. 110 ( , 2018. We classify the fire episodes based on the fuel classification in the Global Fire Emission Database (GFED), which is estimated using the MODIS land cover type product and University of Maryland classification scheme (Friedl et al., 2010;van der Werf et al., 2017). We assign the fuel type to grid cells with mixed fuel types based on the dominant fuel type. We follow the definition of GFED, grouping savanna, grassland and shrubland fires as a single herbaceous fuel type. To assess if the NOx EF 115 varies among these three herbaceous types, we use the 500-m yearly MODIS land cover product (v5) to classify the herbaceous fires based on the dominant land cover type (Friedl et al., 2010). We use wind fields from the hourly ERA-5 reanalysis data developed by the European Center for Medium-range Weather Forecast (ECMWF), which provides meteorological variables at 0.25˚ × 0.25˚resolution with 37 pressure levels from 1000 hPa to 1 hPa from 1979 to present (Hersbach et al., 2020). We sample ERA-5 wind data closest to the center of each fire episode at the TROPOMI overpass time (~1 PM local time). 120

An improved fire a priori for TROPOMI NO2
The a priori vertical profiles of NO2 used in the standard TROPOMI products are obtained from global daily model simulations (TM5) with coarse resolution (1˚) and monthly average biomass burning emissions (Williams et al., 2017). Fires are intrinsically episodic and occur over land areas that are often as small as a few kilometers. Here we re-compute the tropospheric 125 AMFs using the vertical NO2 profiles provided by the NASA GEOS-CF simulations with 0.25˚ resolution. GEOS-CF system combines GEOS weather analysis and forecasting system with GEOS-Chem chemistry scheme version 12.0.1 (Bey et al., https://doi.org/10.5194/acp-2021-381 Preprint. Keller et al., 2014Keller et al., , 2021Long et al., 2015). GEOS-CF includes detailed gas-phase and aerosol chemistry (Knowland et al., 2020;Keller et al., 2021). The near-real-time satellite-based Quick Fire Emission Database (QFED v2.5) is used to provide daily biomass burning emissions (Darmenov and Silva, 2015). In the GEOS-CF system, 65% of biomass burning emissions 130 are distributed within the boundary layer, and the other 35% is distributed evenly between 3.5 and 5.5 km (Fischer et al., 2014).
GEOS-CF provides hourly global vertical profiles of NO2 at 23 pressure levels from 1000 hPa to 10 hPa since 2018.
We sample GEOS-CF products at the time and location of all fire episodes identified. For each episode, we spatially interpolate the GEOS-CF simulated NO2 profiles to match the resolution of TROPOMI products. The AMF (AMFGC, clear) for clear sky conditions can be calculated following Eq. (1): 135 where ml is scattering weight, which is a function of satellite viewing geometry, surface pressure and reflectivity etc.; xGC,l is the GEOS-CF sub-column for layer l. We acquire averaging kernels (AK) from TROPOMI Level-2 products, and we interpolate GEOS-CF vertical profiles to the 34 vertical layers that provide information on AKs. AK is equal to the ratio of the scattering weight to the tropospheric AMFs computed from the a priori profile (AMFa) (Eskes and Boersma, 2003): 140 Given that AMFa is the ratio of slant columns to vertical columns, we can relate the vertical columns with GEOS-CF simulated profile (ΩGC, clear) to the originally retrieved vertical columns (Ωa) as Eq. (4): 145 For partly cloudy scenes, the air mass factor can be written as a linear combination of a clear and cloudy air mass factor (Boersma et al., 2004): where fcloud is radiance weighed cloud fraction, and AMFGC, cloud is essentially the above-cloud component of Eq.

Estimation of emissions and lifetimes of wildfires
We apply an EMG approach to estimate NOx emissions and lifetimes from each fire episode. For each fire episode, we first rotate TROPOMI NO2 swath data along the wind direction in the range of 200 km around the fire centre, and map the rotated 155 (2020). Next, we integrate the TROPOMI NO2 columns in the across-wind direction within ±100 km, which gives reduced one-dimensional line densities. The NO2 line densities (L) are then fitted with an EMG model that is a convolution of a gaussian shaped emission and an exponential decay function (Beirle et al., 2011;Lu et al., 2015;Laughner and Cohen, 2019) following Eq. (7): 160 where x0 is the e-folding distance that represents the length scale of the NO2 decay; μx is the location of the apparent source relative to the fire centre; σx represents the Gaussian smoothing length scale; is a cumulative distribution function; a is a scale factor that represents the observed total number of NO2 molecules in the fire plumes, and B represents the background NO2. We use the best guesses for initial values following (Laughner and Cohen, 2019). The effective NO2 lifetime ( EMG) and 165 the estimated NOx emissions (EEMG) can be calculated from the fitted x0 and a following Eq. (8) and Eq. (9): where w is the wind speed, and is the ratio of NOx to NO2, which is assumed to be 1.32. Previous studies either use the averaged wind speeds of the first several layers (Beirle et al., 2011;Lu et al., 2015) or choose a constant layer such as 900 hPa 170 (Mebust et al., 2011), but injection height of wildfires varies significantly, especially for large fires which inject emissions into high altitudes (Val Martin et al., 2010). To account for varying injection height, we use TROPOMI ALH as an approximation of the fire injection height instead of assuming a constant layer. We vertically interpolate ERA-5 wind data to the pressure level of aerosol layer. For the fires without valid ALH (~36% of the selected fires), we use 900hPa, as the ALP level for the majority of selected fires is near 900 hPa (see Sect. 4.1). We assume a constant because O3 varies little in the plume and the 175 time scale for NO and NO2 to reach steady state is of order 100s (Alvarado and Prinn, 2009). Mebust et al. (2011) suggest the uncertainty of this value is ~20%.

Idealized plume model
To understand the factors that control the NOx lifetime, we employ a one-dimensional (1-D) multi-box plume model based on the Python Editable Chemical Atmospheric Numerical Solver (PECANS). PECANS is a flexible idealized atmospheric 180 chemistry modelling framework that allows for one box to three-dimensional multi-box simulations of atmospheric chemistry with idealized transport. In this study, we set the model to be 1-D with 600 km domain size and 2.5 km resolution, which is analogous to the integrated 1-D NO2 line density along wind direction. The wind speed is fixed at 5 m/s, and the diffusion coefficients are also fixed at 100 m 2 /s. We assume a simplified set of reactions to represent a chemical condition within the two groups; the first group (hereafter RVOC) do not contribute PAN formation; the second group is modelled as an immediate PAN precursor (hereafter OVOC), specifically acetaldehyde. We include a Gaussian shaped NOx emission source (expressed in NO) at x=200 km with 6 km in half width. The concentrations for O3, hydroxyl radical production rate P(HOx), VOC 190 reactivities, and alkyl nitrate branching ratio are given as model input. We run PECANS repeatedly with varying NOx emissions, P(HOx), RVOC, and OVOC. The O3 concentration is fixed at 65 ppbv, and a fixed branching ratio to form RONO2 in RO2 + NO reaction of 0.05 is used. Each model run outputs the concentration of NOx and its major sinks along the wind direction, which are then used to calculate both EMG fitted and chemical NOx lifetime.

Selection of wildfires 195
Since not all fire plumes are detectable from space, and the EMG approach works best for single sources with clear plume patterns, we apply the following four criteria to select candidate fires from fire plumes identified from MODIS FRP observations: First, fires should be large enough so that an apparent enhancement of NO2 is observed from TROPOMI. Considering the detection limit of TROPOMI NO2 and the greater uncertainty of MODIS FRP for small fires (Kaufman et al., 1998), fire 200 episodes with MODIS FRP higher than 200 MW are selected. We then select fires where TROPOMI NO2 tropospheric columns on the fire day are at least one standard deviation higher than the mean TROPOMI NO2 columns 3 to 30 days before and after the fire day (56 days in total, defined as ΩNO2_B).
Second, nearby fire plumes should not contribute to the NO2 line density of the selected fires. A major challenge of applying the EMG approach is that it only applies to single source, but isolated wildfires are rare in nature. To reduce the influence of 205 nearby plumes, we develop an algorithm that automatically detects and filter out the surrounding plumes. We first identify plume-affected grid cells defined as ΩNO2 higher than background NO2 (ΩNO2_B). Next, pixels are grouped to separate plumes based on their connections with surrounding pixels, assuming that plume pixels belonging to the same fire event should be connected. We then exclude the plumes that do not belong to the centre fire plumes. The filtered areas are filled with background ΩNO2_B. The ability of clustering depends on the choice of ΩNO2_B: high ΩNO2_B may truncate plumes as edges are 210 counted as background, while low ΩNO2_B may lead to background being counted as fire plumes, so that nearby plumes are connected with centre plumes. To optimize the performance, we repeat the clustering and filtering steps with different ΩNO2_B, and select the ΩNO2_B that maximizes the filtered size of nearby plumes while retains the centre plume. This filtering algorithm, however, does not apply to the case where fire plumes are overlapped. Therefore, we further exclude the fire plumes where comparable or larger fires are detected (i.e., total FRP of the nearby plumes is greater than one-third of the selected fire) over 215 the region after applying the filtering.
Third, the fire plumes should align with the wind direction. We define a rotation bias as the angle between the wind direction and the observed apparent plume direction. From previous step, we obtain an approximate region of fire plumes, whose coordinates in x and y directions can be fitted with a line using linear regression, where the slope of the line can be converted to degrees (Fig. S1). We only select fire episodes with rotation biases between -30˚ and 30˚.
Fourth, the fire plumes should give good fitting statistics that satisfy the following criteria: 1). R 2 > 0.5; 2) σx < x0, meaning that emission width is smaller than the e-folding distance, which could prevent the case in which emission shape confounds with lifetime; 3) |μx| < 50 km, meaning that the apparent source centre is not too far from the fire centre. Fires with more than 50% missing TROPOMI NO2 values are excluded. We also exclude fires in which TROPOMI NO2 line densities are monotonically increasing or decreasing within the region. The outcome of EMG function is sensitive to the initial condition 225 of each fitting parameter. To test the sensitivity of the fitting results to initial conditions, we repeat the fitting with varying initial conditions 50 times, and we exclude fires where the standard deviation of resulting emissions is more than 50% of the emissions. After excluding the fires sensitive to initial conditions, the uncertainty of the emission due to initial conditions is ~5%.

Characteristics of selected fires 235
Applying above selection criteria, we identified 3248 fire episodes globally between April 2018 and June 2020 suitable for the EMG approach (Fig. 1a). The majority of the fires (77%) occur in the savanna, grassland and shrubland ecosystems. We identified 573 (18%) forest fires, including 227 over the boreal forest, 153 over temperate forest and 193 over the tropical forest fires. Twenty fire episodes are classified as peatland, which occurred in equatorial Asia. We also identified 158 fires (5%) due to agriculture waste burning distributed across different regions. Figures 1b to 1d show the distributions of MODIS 240 FRP, TROPOMI ALP and MAIAC AOD of the selected fire episodes. The MODIS FRP is below 10000 MW for 95% of the selected fires, and 34 fires have FRP larger than 20,000 MW. The mean TROPOMI ALP is 828 hPa (or 1836 m for ALH), with 1σ standard deviation being 118 hPa (or 1751 m). Assuming TROPOMI ALH is indicative of fire injection height, more than 80% of selected fires inject fire emissions to a pressure level between 700 and 950 hPa. About 83% of the fire episodes show MAIAC AOD less than 0.3 near the fire centre, and only 64 (3%) fire episodes have AOD higher than 1.0. In summary, 245 most of our selected fires can be characterized as median to large fires with relatively low injection height and small AOD.  (Fig. 2a). Second, we update the a priori profile of NO2 to improve the estimate of NO2 column (Sect. 3.1), which leads to an enhancement of NO2 250 gradient near the fire centre (Fig. 2b). Third, we filter two nearby fire plumes, and the nearby plumes are filled with background NO2 (Fig. 2c). Fourth, we integrate across the wind direction to obtain a 1-D line density and fit with the EMG function ( Fig.   2d and Eqs. (7) to (9)). The EMG model captures the variation of the line density (R 2 = 0.98). The lifetime is estimated to be 1.6 hours, and the total NOx emissions are estimated as 7870 g/s. shows the fitted line density using the EMG function (Eq. (7)). The lifetime is estimated from Eq. (8), and the total NOx emissions 260 are estimated using Eq. (9). FRP is calculated as the sum of the FPR of all fire pixels detected by MODIS shown on (c).

Satellite-derived fire NOx emissions
Deriving NOx emissions for the large ensemble of fires, we investigate what drives the variation of NOx emissions and lifetimes among these fires. As MODIS FRP is considered as a good indicator of fire intensity (Ichoku and Kaufman, 2005), the emission coefficient (EC), defined as the mass of pollutant emitted per unit of radiative energy (i.e., Emissions/FRP), has been used to 265 derive the emissions of chemical species from fires (Ichoku and Kaufman, 2005;Mebust et al., 2011;Mebust and Cohen, 2014). Figure 3 shows the relationship between TROPOMI derived NOx emissions with MODIS FRP for six different fuel types. Overall, we find TROPOMI derived fire NOx emissions are positively correlated with MODIS FRP. We assess an overall EC by fitting a line with intercept fixed at zero. FRP explains 39% to 78% variance in emissions with the highest R 2 for tropical forest fires and lowest for agricultural fires. The variability not accounted for may be related to the uncertainty of 270 satellite retrieval of NO2, errors with the EMG approach, uncertainties with FRP (see Sect. 5), and/or true differences in NOx  (Fig. S2). Only 20 fires are classified as peatland fires, and we find a relatively good correlation between NOx emissions and FRP with R 2 =0.62 and EC = 0.75 [0.47, 1.03] g/MJ. For agricultural fires, there is a large scatter between NOx emissions and FRP with R 2 = 0.39, and the estimated EC is 1.10 [0.88, 1.31] g/MJ. Updating the a priori profile 280 of NO2 enhances the spatial gradient of NO2 (Fig. 2b), allowing for better estimates of fire NOx emissions. Indeed, we find that using TROPOMI standard products gives a weaker correlation between FRP and NOx emissions, and the ECs decrease by 46% on average (Fig. S3).

Comparison with previous studies
To compare with previous studies, we convert the ECs to emission factors (EFs) assuming a constant ratio of fuel consumption to FRP, Kr = 0.41 kg MJ -1 , (Vermote et al., 2009;Mebust et al., 2011;Mebust and Cohen, 2014). Figure 4 shows the TROPOMI derived NOx EFs and the associated uncertainties (95% CI) compared to previous studies. We find that satellite-derived NOx 295 EFs (hereafter EFssat) are largely consistent with the mean reported in the Andreae (2019, hereafter EFsandreae), which represent an up-to-date compilation of field and laboratory measurements over the last two decades. In most fuel types except for temperate forest, EFssat are largely consistent with EFsandreae to within 30% difference. Our derived NOx EFsat for tropical forest (3.17 g/kg) is nearly twice as large as that in the boreal forest (1.70 g/kg), consistent with Andreae (2019), which also shows larger NOx EFs over the tropical forest (2.81 g/kg) than boreal forest (1.18 g/kg). However, for the temperate forest, the derived 300 NOx EFsat (1.36 g/kg) is less than half of EFandreae (3.02 g/kg). There is a large spread of NOx EFs for the temperate forests in literature, ranging from 0.49 g/kg (Liu et al., 2017) to 7.44 g/kg (Yokelson et al., 2007). Our derived NOx EFsat, however, is close to the in-situ estimates of NOx EFs (1.56 g/kg) from the recent aircraft campaign (i.e., WE-CAN) over the western US during summer 2018 (Lindaas et al., 2021). In non-forest fuel types, we find the smallest NOx EFsat (1.83 g/kg) over peatland, followed by grassland (2.49 g/kg), and agriculture (2.68 g/kg), which are close to EFandreae. Our derived NOx EFs are nearly 3 305 times larger than a previous study based on OMI observations, which suggest NOx EFs are lower than 1g/kg in all fuel types (Mebust and Cohen, 2014). The discrepancy appears to be primarily due to less accurate representation of biomass burning emissions in the a priori profile of NO2. Using the standard TROPOMI NO2 products without updating the a priori profile, the derived NOx EFs are similar to those developed by Mebust and Cohen (2014), 44 to 66% of EFsat, and 26 to 68% of EFsandreae, suggesting that updating the a priori profile is the dominant source of bias in the previous research using satellite 310 observations to derive EFs.

Figure 4 Comparison of the TROPOMI-derived NOx emission factors with previous studies and those used in global biomass burning emission inventories. We include two estimates of NOx emission factors: one using the original TROPOMI NO2 (purple); the other with updated the a priori profile for AMF calculation (black). The error bars of TROPOMI NOx EFs represent 95% CI calculated based on Student's t-distribution test. The red dots show the mean NOx EFs reported in previous studies compiled by Andreae (2019), and the error bars represent the standard deviation. The error bars of Mebust and Cohen (2014) are calculated using nonparametric bootstrap resampling. The error bars of Lindaas et al. (2021) indicate the overall uncertainty of measurements.
As the development of biomass burning emission inventories is done by separate groups with different approaches, NOx EFs used in these inventories also differ. We compare our derived NOx EFs with those used in four commonly used global biomass 320 burning emission inventories (Fig. 4) and (4) QFED (Darmenov and Silva, 2015). We find that satellite-derived EFssat best agree with those used in FINN for forest and grassland in terms of absolute magnitude and variations among fuel types. GFED and FINN use smaller EFs over boreal (0.9 g/kg and 1.8 g/kg) and temperate (1.9 g/kg and 1.3 g/kg) forest than tropical forest (2.6 k/kg), which is contrary to GFAS 325 and QFED that use higher NOx EFs of 3.4 g/kg and 3.0 g/kg for the boreal and temperate forest than that for tropical forest (2.3 and 1.6 g/kg respectively). Our derived NOx EFsat for tropical forest fires is larger than those used in emission inventories.
For grassland, our derived NOx EFsat of 2.49 g/kg is closest to that used in FINN (2.8 g/kg) and GFAS (2.1 g/kg), and smaller than that used in GFED and QFED (3.9 g/kg). For peatland fires, both GFED and GFAS use EF of 1.0 g/kg, which is smaller than our estimated EF of 1.83 g/kg, but we acknowledge a large uncertainty in our derived EF for peatland given the small 330 number of samples. For agricultural fires, our derived NOx EFsat (2.68 g/kg) is slightly higher than that used in GFAS (2.3 g/kg), but smaller than that used in GFED (3.11 g/kg) and FINN (3.5 g/kg).

Satellite-derived NOx lifetime and its driven factors
We find a large variation of NOx lifetime in fire plumes. Figure 5 shows the variation of mean NOx lifetime as a function of NOx emissions at different wind speeds. We find an overall negative relationship between NOx lifetime and emissions: NOx 335 lifetime decreases from over 5 hours for fires with emissions less than 500 g/s to less than 2 hours for fires higher than 5000 g/s (Fig. 5). The varying lifetime with emissions suggests that the assumption of constant lifetime used in previous studies leads to an overestimate in emissions for small fires, while an underestimate of emissions for big fires. At low emission levels (< 2000 g/s), NOx lifetime tends to decrease with increasing wind speed, which is due to the vertical and horizontal diffusion effects that dilute the concentration of NOx and thus alters the rate of NOx removal due to the feedback on OH and the rate of 340 NOx removal (Valin et al., 2013). However, as NOx emissions further increase (> 2000 g/s), the chemical removal becomes fast compared to dilution, and NOx lifetime no longer depends on wind speed (Fig. 5). We note that the EMG derived NOx emissions depend on both NOx abundance and lifetime (Eq. (9)), and thus the observed negative emission-lifetime relationship may partially reflect that emissions are estimated from derived NOx lifetime. However, we find a similar negative relationship between NOx lifetime and TROPOMI NO2 column density near the fire centre (Fig. S4), indicating that chemical feedback of 345 NOx is the primary driver of the derived NOx lifetime. The NOx chemical lifetime, in theory, is determined by its loss to HNO3 and RONO2. We use the 1-D PECANS model to 350 simulate NOx evolution downwind fire plumes and calculate a lifetime by fitting model simulated NOx concentration with the EMG function (Eq. (7)). Figure 6 shows the dependence of the EMG fitted NOx lifetime on NOx emissions rate and P(HOx).
At low NOx emissions, the NOx lifetime decreases rapidly with increasing NOx emissions, while almost independent of P(HOx), indicating a NOx-limited regime. At the NOx-limited regime, increasing NOx facilitates the conversion from HO2 to OH, and thus faster loss of NOx through formation of HNO3 (e.g., Valin et al., 2014;Romer et al., 2020). The loss of NOx through 355 formation of RONO2 also increases with NOx emissions (e.g., Romer et al., 2020). As NOx emissions further increase, the NOx lifetime shows a strong dependence on P(HOx), and the NOx lifetime increases slightly with NO emissions, indicating a NOxsaturated regime. In the NOx-saturated regime, as the loss of NOx through the formation of HNO3 consumes both NOx and OH, increasing NOx leads to decreasing oxidative capacity and thus a longer NOx lifetime. If the NOx lifetime is driven entirely by changes in NOx concentration, the derived NOx lifetime should first decrease and then increase with NOx emissions, which 360 is not found from the observed lifetime-emission relationship. Therefore, we infer that it is likely that P(HOx) increases with fire intensity in fire plumes, which combined with increasing NOx abundances, leads to an overall decrease of NOx lifetime with NOx emissions. If we assume VOC reactivities and branching ratio are fixed, we could use TROPOMI retrieved NOx emissions and lifetime to infer an approximate level of P(HOx). Figure 6 labels the satellite retrieved mean NOx lifetime at a given emission level, and the corresponding P(HOx). To match the observed negative lifetime-emission relationship, P(HOx) 365 should also increase by near a factor of 4 from 15 × 10 6 molec/cm 3 for fires with lifetime longer than 4 hours to 60 × 10 6 molec/cm 3 for fires with lifetime smaller than 1 hour. The increase in P(HOx) may be related to increasing emissions of HONO that generate OH. A recent study suggests that previously underestimated HONO emissions from fires are responsible for twothirds of the HOx production from fresh fire plumes . The changes of P(HOx), however, should be of secondary importance compared to NOx emissions in driving the observed variations in NOx lifetime, which can be evidenced 370 from the slower rate of the increase in inferred P(HOx) and small changes of the NOx lifetime at high NOx emissions (Fig. 5). In addition to NOx emissions and P(HOx), VOC reactivity is the third factor that affects the fitted NOx lifetime. At low NOx emissions, increasing RVOC reactivity leads to a shorter NOx lifetime, but its impacts become smaller and reverse with 380 increasing NOx emissions (Fig. 7a). At low NOx emissions, increasing RVOC reactivity facilitates the loss of NOx through the formation of RONO2. At high NOx emissions, as both RVOC and NOx consumes OH, increasing VOC leads to a longer NOx lifetime. The formation of PAN acts as a temporal reservoir of NOx, which also affects the evolution of NOx. The formation of PAN depends on the concentration of OVOCs, OH level and temperature. Figure 7b shows the dependence of EMG fitted NOx lifetime as a function of NOx emission rate and the reactivity of PAN's immediate precursor (OVOC). We find that 385 increasing PAN formation through increasing OVOC reactivity will lead to an overall increase in EMG fitted NOx lifetime.
The impact of OVOC is especially evident at low levels of NOx emissions (Fig. 7b). Without PAN formation, the fitted NOx lifetime will be shorter than that derived from TROPOMI observations, suggesting that PAN formation plays a non-negligible role in determining the evolution of NOx and also the effective lifetime of NOx.

Uncertainties of satellite retrieved NO2 columns 395
Satellite retrievals of NO2 columns are subject to uncertainties with spectral fitting, and uncertainties from the a priori NO2 profile shape and scattering weights needed for calculating AMFs and stratospheric NO2 columns. Boersma et al. (2018) suggest an overall uncertainty of 35% to 45% for single-pixel OMI NO2 retrieval, which should be smaller for TROPOMI NO2 given its improved signal-to-noise ratio (Veefkind et al., 2012). Due to the sporadic nature of fires, there is no validation of TROPOMI NO2 columns with fire NOx with ground-based measurements. Validation with ground-based Differential Optical 400 Absorption Spectroscopy (DOAS) measurements in urban stations show an overall good agreement between TROPOMI retrieved and ground-based NO2 columns, and an overall negative bias of -22% is reported, but the biases vary with stations (Lambert et al., 2020).
Over polluted regions, the uncertainty and bias of single-pixel satellite retrieval of NO2 columns is dominated by uncertainties of the AMF (Boersma et al., 2018). We replace the a priori profile shape of NO2 using NASA GEOS-CF simulations that 405 include daily biomass burning emissions, which leads to higher NOx emissions factors that are more consistent with in situ measurements. Bousserez (2014) similarly suggests that using a fire-specific NO2 profile shape can lead to a near 60% reduction in AMF. In the NASA GEOS-CF, 65% biomass burning emissions are distributed within the boundary layer, which will lead to negative biases in AMF for large fires with plumes that rise well above the boundary layer. In our study, since the ALH is below 2000m for the majority of fires, model uncertainties of the injection height should have small impacts on the 410 retrieval of NO2 columns. However, we notice that fire events with low ALP (< 700 hPa) tend to show higher NOx emissions per FRP, and most outliers in Fig. 3 are associated with high aerosol layers. Since satellite instruments are more sensitive to NO2 at higher altitudes, the inaccurate assumption of NO2 concentrated within the boundary layer will lead to an underestimate of the AMF and thus an overestimate of the retrieved vertical NO2 columns if the majority of NO2 is injected to the free troposphere. 415 In addition, the large amounts of aerosols emitted from wildfires may also impact the retrieval of NO2 columns from space.
The impacts of aerosols depend on the magnitude, location, and optical properties (absorbing vs. scattering) of aerosols (Bousserez, 2014;Lin et al., 2015). In the TROPOMI NO2 retrieval algorithm, the effects of aerosols are implicitly accounted for through modifying the cloud properties (Boersma et al., 2011). Bousserez (2014) suggests an explicit correction is needed in the presence of clouds and scattering aerosols, and the effects of aerosol correction can be as large as 100% when cloud 420 fraction is 30% and AOD is higher than 1. In our study, the mean cloud fraction of the selected fire plumes is 9%, and the mean AOD is 0.22, corresponding to a small uncertainty of less than 20% based on Bousserez (2014).

Uncertainties of MODIS FRP
Since we derive NOx EFs from the linear regression between TROPOMI NOx emissions and MODIS FRP, uncertainties in MODIS FRP should also contribute to uncertainties of the satellite derived NOx EFs. Detection of MODIS FRP may be 425 obscured by cloud, aerosols or canopy cover. The uncertainty of FRP, however, is lower than 5% for fires that aggregate 30 or more active fire pixels together (Freeborn et al., 2014). Our selected fire events are aggregated by 37 active fire pixels on average (SD = 31). To evaluate if our results are robust with the choice of FRP data, we conduct similar analysis with FRP measurements from Suomi NPP VIIRS sensors. In general, MODIS and VIIRS FRP are in good agreement (R 2 = 0.72, Fig.   S5). VIIRS FRP is lower than MODIS FRP by 5.8% on average. Deriving NOx emission factors using VIIRS FRP, we find a 430 slight increase of NOx EFs for forest and agricultural fires, but a decrease for peatland (Fig. S6). For herbaceous fires, where a large number of fires are sampled, the derived NOx EF remains almost unchanged, suggesting that the difference in MODIS and VIIRS FRP should diminish as we increase the sample size. The overall relationship between NOx emissions and FRP is similar for both VIIRS and MODIS, though stronger correlation is found for MODIS FRP. Overall, we estimate the difference in NOx EFs using MODIS versus VIIRS FRP is ~20%. The choice of the mass-to-energy conversion factor (Kr) is another 435 source of uncertainty. While previous studies suggest the relationship between FRP and fuel consumption is highly linear and significant (Freeborn et al., 2008;Wooster et al., 2005), Kr value is 0.368 g/MJ in Wooster et al. (2005), but 0.453g/MJ in Freeborn et al. (2008), suggesting an uncertainty of order 10% for the value of Kr.

Uncertainties in the EMG approach
The first step of the EMG approach is to rotate TROPOMI observations along the wind direction. The derived NOx lifetime 440 and emissions are therefore subject to uncertainties of the wind direction and speed due to uncertain plume heights that cross wind shear, or the plume thermodynamics that are not captured by ERA5 wind data. In the case the fire plume does not align with wind direction, calculating line density along the wind direction should lead to an underestimate of the e-folding distance.
In this study, we only select fires with less than 30˚ rotation bias, and the mean rotation bias is 10˚, which corresponds to less than 2% underestimate of the e-folding distance. Uncertainty and variance of the wind speed, however, should lead to errors 445 in the derived NOx lifetime. Here we determine the wind speed by interpolating the wind profile to TROPOMI derived aerosol layer. Comparison of TROPOMI ALH with plume height from the Cloud-Aerosol LIdar with Orthogonal Polarization (CALIOP) and the Multi-angle Imaging SpectroRadiometer (MISR) measurements suggest that TROPOMI ALH is overall 500 m lower (Griffin et al., 2020;Nanda et al., 2020). We estimate that an increase of 500 m ALH corresponds to ~22% increase of the wind speed on average, meaning that NOx lifetime will decrease by ~18%. 450 Here we use the EMG approach to derive an effective NOx lifetime, which is a combination of both the chemical lifetime of NOx and the lifetime due to mixing such as those plume movement in different directions that reduces the line density. In addition, chemical nonlinearities can result in an effective chemical lifetime that is averaged over the plume where at each point in the plume evolution a different chemical lifetime occurs. To assess if EMG fitted NOx lifetime is indicative of the chemical lifetime, we use the PECANS model to calculate an EMG fitted lifetime and a chemical lifetime of NOx from two 455 permanent losses of NOx through the formation of HNO3 and RONO2. The chemical lifetime of NOx varies with location ( Fig.   8a), which reaches a minimum near the source centre at low NOx emissions (NOx-limited regime), but shows maximum at high NOx emissions (NOx-saturated regime). The effect of varying lifetime on the emission estimates is not considered with the EMG approach, which instead gives an overall effective NOx-lifetime of the plume. At low NOx, the EMG fitted lifetime is higher than the chemical lifetime at the fire centre, but the reverse occurs at high NOx (Fig. 8a). We notice that the EMG 460 fitted lifetime is overall more consistent with the chemical lifetime around 40 to 50 km downwind from the fire centre, which should vary depending on the model input. Figure 8b compares the EMG fitted lifetime with overall chemical lifetime over downwind regions under different NOx emissions, P(HOx) and VOC reactivities. We find a moderate correlation between between the EMG fitted and chemical lifetime (R 2 = 0.93). PAN is generally considered as a sink of NOx near the fire, but a 465 source of NOx over downwind region, which deepens the gradient of NO2 near the source but flattens the gradient downwind (Valin et al., 2013). The EMG approach is unable to capture the effects of PAN formation on the evolution of NOx as it assumes NOx decays exponentially. In the presence of PAN formation, we find the EMG approach tends to overestimate NOx lifetime at low NOx emission (< 5000 g/s), while underestimating NOx lifetime at high NOx emissions and low P(HOx). Overall, we estimate the overestimate at low NOx emissions (< 5000 g/s) cause around 33% negative biases to the derived emissions, while 470 18% positive biases at high NOx emissions (> 5000 g/s). and VOC reactivities, where the colours represent different emission levels, and the symbols represent different levels of P(HOx). The chemical lifetime is calculated as the mean NOx concentration divided by the mean chemical loss of NOx through the formation of HNO3 and RONO2 over downwind region. We define the regional mean as the region between fire centre and the downwind area where NOx concentration is higher than background, where background value is estimated from EMG function (B in Eq. (7)).

6 Conclusions
We estimate NOx emissions and lifetime from over 3000 fires globally using daily TROPOMI retrievals of NO2 tropospheric columns, and derive NOx emission factors by linking TROPOMI derived NOx emissions with MODIS FRP. The overall derived NOx emissions factors are 1.70 g/kg, 1.36 g/kg and 3.17 g/kg for the boreal, temperate and tropical forest, 2.49 g/kg for herbaceous (grassland, savannas and shrubland combined) fires, 1.83 g/kg for peatland, and 2.67 g/kg for agricultural 485 burning. Satellite-derived NOx emission factors are largely consistent with the mean NOx emission factors reported by previous field studies (Andreae, 2019). By studying a large number of fires globally, we provide more representative NOx emission factors. These top-down emission estimates of NOx could be used to assess biomass burning emission inventories in terms of both emission factors and fuel consumption, which could help diagnose the causes of discrepancies among different emission inventories. 490 Our study features three improvements over previous studies that use satellite measurements to derive NOx emissions (Mebust et al., 2011;Mebust and Cohen, 2014;Schreier et al., 2015). First, we use observations from the newly launched TROPOMI with finer spatial resolution and improved signal-to-noise ratio. Second, we replace the a priori profile of TROPOMI NO2 retrieval with a high-resolution global model simulation from NASA GEOS-CF simulations to calculate AMF. This update result in steeper gradients between the plumes and the background, and more accurate description of NO2 vertical shape, 495 reducing the discrepancy between satellite and in-situ derived estimates of NOx emission factors. Third, we relax the assumptions of constant NOx lifetime by directly estimating lifetime through fitting the evolution of NOx downwind with the EMG approach.
We observe decreasing NOx lifetime with increasing fire NOx emissions, which is indicative of NOx-limited chemistry, where increasing NOx emissions makes the chemical loss of NOx more efficient. However, for the largest fires with high NOx, a 500 regime transition from a NOx-limited to NOx-saturated regime is expected, where increasing NOx emissions leads to a longer NOx lifetime. Using a 1-D idealized plume model to interpret the factors affecting the NOx lifetime, we infer that P(HOx) must also increase with fire intensity, consistent with observations that indicate a large source of HONO in fires Peng et al., 2020). The formation of PAN also impacts the NOx lifetime, but the evolution of NOx due to the formation of PAN and its thermal decomposition over downwind areas is not well captured by the EMG approach that assumes exponential decay 505 with NOx downwind. Future studies will benefit from the integrative analysis of satellite retrievals of NO2, HONO and PAN to more completely describe the chemical evolution of reactive nitrogen from wildfires, thus allowing for better prediction of the air quality impacts of fires. The newly launched or upcoming geostationary satellite instruments such as GEMS and TEMPO will offer an unprecedented opportunity to continuously observe the emissions and chemical evolution of NOx from fires that will no longer be limited to a single snapshot (Chance et al., 2013;Kim et al., 2019). 510