Cities are greenhouse gas emission hot spots, making them targets for emission reduction policies. Effective emission reduction policies must be supported by accurate and transparent emissions accounting. Top-down approaches to emissions estimation, based on atmospheric greenhouse gas measurements, are an important and complementary tool to assess, improve, and update the emission inventories on which policy decisions are based and assessed. In this study, we present results from 9 research flights measuring CO2 and CH4 around New York City during the nongrowing seasons of 2018–2020. We used an ensemble of dispersion model runs in a Bayesian inverse modeling framework to derive campaign-average posterior emission estimates for the New York–Newark, NJ, urban area of (125 ± 39) kmol CO2 s–1 and (0.62 ± 0.19) kmol CH4 s–1 (reported as mean ± 1σ variability across the nine flights). We also derived emission estimates of (45 ± 18) kmol CO2 s–1 and (0.20 ± 0.07) kmol CH4 s–1 for the 5 boroughs of New York City. These emission rates, among the first top-down estimates for New York City, are consistent with inventory estimates for CO2 but are 2.4 times larger than the gridded EPA CH4 inventory, consistent with previous work suggesting CH4 emissions from cities throughout the northeast United States are currently underestimated.
Introduction
Cities are large sources of greenhouse gas emissions, with an estimated 36.8% of total global emissions coming from within urban extents in the year 2000 (Marcotullio et al., 2013). In North America, this urban share of greenhouse gas emissions was estimated to be even larger (49.2%; Marcotullio et al., 2013).
To reduce emissions and mitigate climate change, many cities are implementing new policy practices and regulations (Trencher et al., 2016; Sethi et al., 2020). New York City has legally binding emission reduction targets of 40% by 2030 and 80% by 2050, relative to 2005 levels (New York City Mayor’s Office of Sustainability, 2016). As part of the legislative effort to meet these targets, the city council passed the Climate Mobilization Act, which introduces emission limits on large- and medium-sized buildings under Local Law 97 (Climate Mobilization Act, 2019).
Progress toward emission targets and the efficacy of emission reduction policies can be assessed through the compilation of city-specific inventories, often referred to as self-reported inventories. These inventories estimate greenhouse gas emissions using a bottom-up approach, where the total emission from each source category is calculated by multiplying activity data (e.g., fuel sales) by an emission factor (e.g., greenhouse gas emissions per unit of fuel sold). Self-reported inventories typically provide emission estimates at the whole-city level, broken down by emission source. These estimates can include direct emissions within the city limits (scope 1 emissions), emissions associated with electricity used within the city but generated elsewhere (scope 2 emissions) and other indirect emissions that occur elsewhere (e.g., landfill CH4) as a result of activity (e.g., waste production) within the city (scope 3 emissions; Nangini et al., 2019). The New York City emissions inventory (MacWhinney and Barnett, 2019) is compiled using a mixed-scope approach known as the city-induced framework (Fong et al., 2014), which is designed to account for emissions attributable to activities within the city limits, regardless of emission location. This mixed-scope approach includes some emissions from each of scopes 1, 2, and 3 but does not include all the emissions from any single scope. Quantifying indirect emissions is important from a policy perspective, but it is challenging to assess mixed-scope inventories with atmospheric measurements, which are only sensitive to scope 1 emissions.
Spatially resolved inventories of scope 1 emissions have been compiled at the national level, either using bottom-up approaches (Gately and Hutyra, 2017; Gurney et al., 2020) or using proxy data to spatially disaggregate national emission totals (Maasakkers et al., 2016; Oda et al., 2018; Janssens-Maenhout et al., 2019). These can be used to assess and improve the self-reported inventories or adopted as the primary method by which cities tailor and assess their emission reduction policies (Hutyra et al., 2014; Gurney and Shepson, 2021). A recent study by Gurney et al. (2021) compared scope 1 emissions from the national Vulcan inventory with self-reported inventories from 48 U.S. cities and found that most cities underreported emissions. The spatial information in gridded national inventories allows them to be combined with atmospheric transport models to estimate atmospheric mole fractions, allowing direct comparison to measurements and thus facilitating top-down approaches to estimating emissions.
Many top-down studies have used tower-based measurements of greenhouse gas mole fractions in an inverse modeling approach to estimate urban emissions (Bréon et al., 2015; Lamb et al., 2016; Lauvaux et al., 2016; 2020; Staufer et al., 2016; Deng et al., 2017; Huang et al., 2019; Nickless et al., 2019; Yadav et al., 2019). Tower-based sites allow for continuous year-round observations, enabling top-down fluxes to be estimated at annual timescales. On the other hand, an individual tower has limited spatial sensitivity, so multiple towers may be required to provide the necessary spatial coverage to accurately estimate whole-city emissions. Aircraft-based measurements cannot provide the continuous temporal coverage of tower observations. However, by sampling a wide area in a short period of time, well-tailored flight tracks enable aircraft to provide useful snapshot estimates of whole-city emissions, with the ready ability to define and change the area of study.
Aircraft measurements can be used in an inverse modeling framework to estimate urban emissions. This approach incorporates information regarding the prior fluxes, transport model, and measurements, including the error covariance structure for each of these elements. Several previous studies have applied such a framework to estimate urban fluxes using aircraft measurements (Brioude et al., 2011; Brioude et al., 2013; Pisso et al., 2019; Lopez-Coto et al., 2020). However, this approach has yet to be applied to estimate New York City greenhouse gas emissions. In this study, we use measurements from 9 research flights (during the nongrowing season) and an ensemble of dispersion model runs in an inverse modeling approach to estimate carbon dioxide and methane emission rates for the New York–Newark, NJ, urban area and the five boroughs of New York City.
Methods
Aircraft measurements
Measurements were made on board the Purdue University Airborne Laboratory for Atmospheric Research (ALAR), a modified Beechcraft Duchess. Full details of the aircraft payload and sampling configuration are provided by Cambaliza et al. (2014). Carbon dioxide (CO2) and methane (CH4) mole fractions were measured using a Picarro Cavity Ringdown Spectrometer (Crosson, 2008), calibrated in flight using 3 calibration cylinders traceable to the National Oceanic and Atmospheric Administration (NOAA) reference scales for CO2 (X2007; Tans et al., 2017) and CH4 (X2004A; Dlugokencky et al., 2005). The cylinders were prepared, filled, and certified by NOAA; details are provided in SI Table S1.1.
Three Picarro analyzers were used over the course of the campaign: 2 different G2301-f analyzers and 1 G2401-m analyzer. The 2 G2301-f analyzers were operated in low-flow mode, with each species measured at approximately 1.2 s intervals on 1 analyzer and 1.4 s intervals on the other. The G2401-m analyzer took a measurement of each species at approximately 2.3 s intervals. These analyzers exhibited typical precisions of 0.2 ppm (µmol mol–1) for CO2 and 3 ppb (nmol mol–1) for CH4, based on analysis of data gathered during in-flight and ground-based calibrations. A linear calibration curve was derived for each flight day using in-flight measurements of all 3 calibrations cylinders on board the aircraft. In all cases, we obtained r2 values greater than 0.999, demonstrating good instrument linearity.
Flights were performed on 9 separate days between November 2018 and March 2020, with all flights taking place in either November, February, or March. Flight tracks are shown in Figure S4.1, which also gives the dates and times of the flights. The total flight duration for each day (as measured from the first usable measurement to the final usable measurement) varied between 2.5 and 5.5 h.
Transport modeling
The surface influence, or footprint, of the sampled air was determined using HYSPLIT v5.0.0 (Hybrid Single Particle Lagrangian Integrated Trajectory Model; Draxler and Hess, 1998; Stein et al., 2015; Loughner et al., 2021), a Lagrangian particle dispersion model developed by the NOAA Air Resources Laboratory. This version of HYSPLIT incorporates elements of STILT (Stochastic Time-Inverted Lagrangian Transport Model; Gerbig et al., 2003; Lin et al., 2003; Fasoli et al., 2018), a related model that branched off an earlier version of HYSPLIT. This combines features from STILT, including options for a different vertical turbulence scheme that prevents well-mixed particles from accumulating in low-turbulence regions (Thomson et al., 1997), a varying vertical Lagrangian timescale (Hanna, 1982) and an additional boundary layer turbulence parameterization (Hanna, 1982), with recent model development and bug fixes within HYSPLIT (Loughner et al., 2021).
Each of our HYSPLIT runs was repeated using 4 different choices of input meteorology and 2 different turbulence parameterizations (Hanna, 1982; Kantha and Clayson, 2000) giving us an 8-member ensemble of dispersion model runs (hereafter referred to as the transport model ensemble). The input meteorology was taken from 4 publicly available archived products: the European Centre for Medium Range Weather Forecasts Fifth Reanalysis (ERA5), the NOAA Global Forecast System (GFS) model, the NOAA North American Mesoscale Forecast System (NAM) model, and the NOAA High-Resolution Rapid Refresh (HRRR) model (see SI Section S2 for more details on each product).
A separate HYSPLIT model run was conducted for every minute of each flight that included aircraft sampling within the boundary layer. The mean measured mole fraction for each minute was calculated as an average over all measurement points where the aircraft was within the boundary layer, with the corresponding model particle release consisting of approximately 1,000 particles distributed equally across these measurement locations (covering a horizontal distance of approximately 3–4 km at typical aircraft speeds). In this way, we ensured that the model run was representative of the measurements we included in our average.
Tiered multi-resolution Bayesian inversion framework
Fluxes were optimized separately for each flight using a Bayesian inversion framework, after first subtracting the background mole fraction from the aircraft measurements (where the background is defined as the mole fraction of CO2 or CH4 within the measured air prior to it entering the domain; see SI Section S9 for details). This process consisted of multiple stages:
Optimize fluxes on a large domain at coarse spatial resolution (0.08° × 0.08°).
Use the posterior flux map from step 1 to estimate spatial variability in the mole fractions of CO2 and CH4 flowing into a smaller, high resolution (0.02° × 0.02°) domain, nested within the large domain and centered on New York City.
Optimize fluxes on the smaller, high-resolution domain, taking into account the background variability calculated in step 2.
The large domain (henceforth referred to as the d01 domain) was bounded by 34.4°N, 44.4°N, 83.7°W, and 69.7°W, while the smaller nested domain (henceforth the d03 domain) was bounded by 39.2°N, 42.0°N, 75.7°W, and 72.1°W (see SI Figure S3.1 for map). This approach takes account of upwind sources that influenced our measurements and provides an optimized representation of their influence. This allows us to better isolate the specific contribution from the area of interest (e.g., the New York–Newark urban area) by providing an optimized background for each minute of the flight.
A gridded footprint (representing the influence of surface fluxes on the sampled air) was calculated on both the d01 and d03 domains for each minute of the flight. Modeled mole fraction enhancements (in units of µmol mol–1; commonly referred to as ppm) were derived by multiplying these gridded footprints (in units of ppm [µmol m–2 s–1]–1) by surface fluxes (in units of µmol m–2 s–1) of CO2 and CH4. The inventories used to provide prior fluxes are described in the following section. For both the initial d01 inversion and the nested d03 inversion, optimized posterior fluxes were obtained by minimizing the cost function J(x):
This cost function assumes normally distributed uncertainties and has the following analytical solution (Enting, 2002; Tarantola, 2005):
Here x is a vector containing the fluxes that are to be optimized, xa contains the posterior (optimized) fluxes, xb contains prior fluxes derived from an emissions inventory and H is a matrix containing the modeled footprint for each minute of the flight. The vector y contains the measured mole fraction enhancements for each minute of the flight, relative to a background value representing the mole fraction present in the sampled air when it entered the domain. For the d01 inversion, this background value was taken to be constant, while for the d03 inversion, a spatially varying background was calculated based on the posterior fluxes within the d01 domain (see SI Section S9 for more details, including details of a sensitivity analysis to determine the impact of the nested inversion approach). The error covariance matrices Pb and R represent uncertainty in the prior fluxes and model–measurement mismatch, respectively.
Error covariances
Following previous studies (Lauvaux et al., 2016; Lauvaux et al., 2020; Lopez-Coto et al., 2017; Lopez-Coto et al., 2020), the off-diagonal components of the error covariance matrix Pb were populated using an exponential covariance model (in space). The diagonal elements of Pb were populated on the assumption that the 1σ flux uncertainty in each grid cell equaled 100% of the prior flux within that grid cell, following the values used by Lopez-Coto et al. (2017; Lopez-Coto et al., 2020). For comparison, Andres et al. (2016) reported 2σ grid-cell-level uncertainties for a global CO2 inventory ranging from 4% to more than 190%, averaging about 120%. Gately and Hutyra (2017) reported 50%–250% differences between ACES (Anthropogenic Carbon Emissions System; version 1) and 3 global disaggregated CO2 inventories at city scales, with grid-cell-level differences of over 100% for half of the grid cells (urban and rural) in the northeast United States. For CH4, there are no good estimates of grid-cell-level uncertainty, but considering recent work indicating inventories significantly underestimate urban CH4 emissions, and the lack of knowledge about CH4 sources at urban scales, we can only assume that the relative uncertainties are at least as large as for CO2. A lower limit uncertainty in each grid cell was set to 1 µmol m–2 s–1 and 1 nmol m–2 s–1 for CO2 and CH4, respectively, allowing some correction in grid cells with very low prior fluxes.
The off-diagonal elements of the prior flux error covariance matrix, representing the correlation between prior flux errors in different grid cells, were calculated according to an exponential decay model based on the distance between grid cells (Lauvaux et al., 2012; Lauvaux et al., 2016; Lauvaux et al., 2020; Lopez-Coto et al., 2017; Lopez-Coto et al., 2020). The correlation length was set to 10 km following Lopez-Coto et al. (2017) who found this value to be appropriate for studies at urban scales. A variogram analysis of the spread of the prior flux ensemble indicated that the exponential covariance model and the correlation length used here are appropriate for both species in both domains (see SI Section S11). The same approach to constructing Pb was used for both the d01 and d03 inversions.
A double exponential covariance model (in space and time) was selected for the model–measurement mismatch error covariance matrix R, with a correlation length of 1 km and a correlation time of 1 h based on the short correlation length and time scales reported for atmospheric trace gases in urban environments (Shusterman et al., 2018; Turner et al., 2020). The diagonal elements of R represent the combined uncertainty due to random error in the measurements, model and background for each minute of the flight. This was calculated by combining background uncertainty (see below) with the variance across the eight transport model ensemble members for a given minute and the variance in measured mole fraction during that minute. In all cases, low limits for the uncertainty of 0.2 ppm (µmol mol–1) and 3 ppb (nmol mol–1) were used for CO2 and CH4 respectively, based on an analysis of instrument precision from calibration data.
The background uncertainty for the d01 inversion was derived by calculating the respective variances in measured and modeled mole fractions across the background data points, then summing these variances (see SI Section S9 for the definition of the background points). For the d03 inversion, these terms were combined with the variance of the spatially varying component of the background (i.e., the outside contribution derived from the posterior d01 fluxes), taken across the 8 transport model ensemble members.
Prior fluxes
Separate inversions were performed using 3 different CO2 flux priors and 4 different CH4 flux priors. For CO2, we used fluxes from Vulcan v3.0 (Gurney et al., 2020), ACES v2.0 (Gately and Hutyra, 2017), and EDGAR v5.0 (Emission Database for Global Atmospheric Research; European Commission Joint Research Centre (JRC)/Netherlands Environmental Assessment Agency (PBL), 2019; Crippa et al., 2020). For CH4, we used fluxes from EDGAR v5.0, EDGAR v4.2 (European Commission Joint Research Centre (JRC)/Netherlands Environmental Assessment Agency (PBL), 2011) and the gridded Environmental Protection Agency (GEPA) inventory (Maasakkers et al., 2016). Because EDGAR v4.2 uses a population proxy to distribute fluxes, it contains much larger emissions from urban areas compared to EDGAR v5.0 and the GEPA inventory (conversely it is known to underestimate emissions in oil and natural gas production regions). To bridge this large emissions gap between inventories, we also used a composite CH4 prior containing the average of the GEPA fluxes and the EDGAR v4.2 fluxes in each grid cell.
In all cases, we used annual average emissions for the most recent year available: 2015 for Vulcan, 2017 for ACES, 2018 for EDGAR v5.0 CO2, 2015 for EDGAR v5.0 CH4, 2012 for the GEPA, and 2008 for EDGAR v4.2. The use of annual average fluxes (even for those inventories for which hourly fluxes are available) was motivated by the fact that we did not have inventory fluxes for the actual flight days (which were more recent than the latest available inventory data). Furthermore, holding the prior emission rates constant for each flight ensures that flight-to-flight differences in posterior emission rates are not a consequence of prior temporal variability. While prior constraints are important for solving ill-posed problems, it is known that the choice of prior flux map can influence the posterior fluxes (Lauvaux et al., 2020). By using multiple priors for each gas, we reduce the dependency of our mean posterior fluxes on one specific prior.
A conservative method was used to regrid each prior to the extent and resolution of the d01 and d03 domains (Zhuang, 2020). Canadian emissions for the U.S.-specific priors (Vulcan, ACES, and GEPA) were taken from EDGAR v5.0. A more detailed discussion of prior fluxes is given in SI Section S3.
Sensitivity analysis
In addition to our base case ensemble of inversions, we performed a sensitivity analysis to better understand the impact of different background choices and prior error covariance parameters in our posterior results, as well as to quantify uncertainties in the methodology. A full description of the sensitivity test cases is provided in SI Sections S10 and S11.
Results and discussion
Total posterior emission rates were calculated for both the 5 boroughs of New York City (NYC) and the wider New York–Newark urban area (see SI Section S3 for area definitions). The wider urban area better represents the area of sensitivity for our flights (see SI Figure S4.1 for footprint maps); therefore, these wider urban area emission totals are the focus of analysis in this section. Emission rate estimates for the 5 boroughs of NYC are also presented due to the high policy relevance of this area.
Urban area posterior emission rates
Figure 1 shows total emission rates within the New York–Newark urban area boundary, broken down by flight, prior, and transport model. The mean posterior emission rate for CO2 was (125 ± 39) kmol s–1, where the reported uncertainty represents the 1σ variability across the 9 flights (using ensemble-average totals for each flight). Also shown (in purple) are inventory emissions from the hourly ACES inventory and hourly Vulcan inventory, calculated using only dates and hours that were representative of our flights (see SI Section S5 for details). The mean value of these representative inventory emissions was (145 ± 21) kmol s–1 for ACES (15.5% larger than our mean posterior estimate) and (124 ± 20) kmol s–1 for Vulcan (0.9% lower than our mean posterior estimate). The inventory variability reported here represents 1σ across 45 representative days in each inventory. Both our mean posterior emission rate and these representative inventory emission rates are larger than the corresponding annual average emission rates (116 kmol s–1 for ACES and 105 kmol s–1 for Vulcan) that we used as priors. While we do not expect these representative inventory values (from previous years) to agree perfectly with our mean posterior emissions rates, the agreement shown here is well within the observed flight-to-flight variability of our estimates as well as the expected daily variability of the inventories.
The posterior CO2 emission estimates represent total emissions and are influenced by any contribution from biospheric respiration and photosynthesis, while the inventories only include anthropogenic fluxes. However, the flights in this study were conducted during the nongrowing season and the estimated biospheric contribution was very small, approximately 2% on average (see SI Section S6 for details).
The posterior CH4 emission rates calculated for the New York–Newark urban area using all 4 priors are much larger than emission rates from the GEPA and EDGAR v5.0 inventories but much smaller than the emission rate from EDGAR v4.2. The mean posterior CH4 emission rate was (0.62 ± 0.19) kmol s–1 (1σ variability across flights). This mean posterior total is 2.4 times larger than the annual value from the GEPA inventory.
The temporal variability of urban CH4 emissions on diurnal, weekly, and seasonal timescales is poorly understood. None of the CH4 inventories used in this study provide emission estimates on hourly timescales, so it is not possible to repeat the calculation of representative inventory totals presented above for CO2. Recent studies in Los Angeles (Yadav et al., 2019) and Washington, DC–Baltimore (Huang et al., 2019) found large seasonal cycles in urban CH4 emissions. Floerchinger et al. (2021) observed significant seasonal variability in the biogenic fraction of CH4 emissions for several cities, including New York. More flights are required to provide the annual coverage necessary to confirm if strong seasonality also exists in total CH4 emissions for the New York–Newark urban area.
It is important to emphasize that our mean posterior emission rates (for both CO2 and CH4) do not represent estimates of annual average emissions but can be regarded as the current best top-down estimate for the nongrowing-season, daytime value. However, the finding that emissions are underestimated in the GEPA inventory is likely to hold at the annual timescale, as emissions during much of the rest of the year would have to be almost zero in order to offset the 2.4 times larger CH4 emissions observed during November, February, and March in this study.
Sources of variability and uncertainty analysis
The flight-to-flight variability in urban area posterior emission rate was 31% for both CO2 and CH4 (1σ, relative to the mean posterior emission rate). These values for relative flight-to-flight variability are similar to those reported by Lopez-Coto et al. (2020) for the Washington, DC–Baltimore urban area. Table 1 gives the 1σ spread in posterior emission rate across base case ensemble members (transport models and priors) and sensitivity test ensemble members (prior error covariance parameters and background definitions—see SI for more details), including the average spread for a single flight and the spread in campaign-average emissions. For both species, flight-to-flight variability in posterior emission rate was larger than the combined spread across the ensembles (calculated by adding individual ensemble spreads in quadrature). It is noteworthy that the combined ensemble spread for CH4 is much smaller than the difference between the mean posterior emission rate and the inventory values.
. | Single-Flight Ensemble Spread . | Campaign-Average Ensemble Spread . | ||
---|---|---|---|---|
Ensemble . | CO2 (%) . | CH4 (%) . | CO2 (%) . | CH4 (%) . |
Transport model | 13 | 12 | 3.8 | 3.6 |
Prior | 8.2 | 18 | 7.9 | 18 |
Prior flux covariance parameters | 4.6 | 4.2 | 3.7 | 2.8 |
Background | 4.6 | 8.6 | 4.3 | 7.7 |
Combined | 16 | 24 | 10 | 20 |
. | Single-Flight Ensemble Spread . | Campaign-Average Ensemble Spread . | ||
---|---|---|---|---|
Ensemble . | CO2 (%) . | CH4 (%) . | CO2 (%) . | CH4 (%) . |
Transport model | 13 | 12 | 3.8 | 3.6 |
Prior | 8.2 | 18 | 7.9 | 18 |
Prior flux covariance parameters | 4.6 | 4.2 | 3.7 | 2.8 |
Background | 4.6 | 8.6 | 4.3 | 7.7 |
Combined | 16 | 24 | 10 | 20 |
To provide a rough estimate of expected variability in true CO2 emissions from flight to flight, we calculated the variance in inventory emissions for the urban area using dates and hours corresponding to each flight (see SI Section S5 for details). The standard deviation in inventory emissions across 45 representative days in ACES was equal to 15% of the average emissions across these days, while for Vulcan, this value was 16% (boxplots showing these representative emissions for both inventories are shown in purple in Figure 1).
There are multiple potential reasons why the flight-to-flight variability (31%) of our posterior emissions is larger than that calculated using representative inventory emissions. Differences in sampling pattern for each flight (shown in SI Figure S4.1) mean that different parts of the urban area were sampled on different flights. Because the measurements are not sensitive to the whole urban area for each flight, emissions in unsampled areas default to the emission rate of the prior. In addition, different areas were sampled at different times within each flight, which for sources with large hourly variability in emission rate can result in apparent flight-to-flight variability in posterior emissions. Thus, this irregular sampling combined with large spatiotemporal variability in fluxes can become aliased as apparent temporal variability, as demonstrated by Lopez-Coto et al. (2020). In addition, for a given ensemble (e.g., transport or background), errors in the ensemble mean can differ in magnitude and direction for each flight, thus increasing the observed flight-to-flight variability.
Finally, the variability in inventory emissions may not be truly representative of real emission variability because the temporal emission profiles in the inventories are partly based on interpolation and/or downscaling using proxy data (resulting in smoothed temporal emission profiles). It is also worth considering that inventory emissions are for previous years, which may have had lower variability than the flight days, and that the flights were conducted in 2 different winters (2018–2019 and 2019–2020) while the inventories are only for a single year each, so potential interannual variability is not accounted for.
Quantifying the accuracy of the emission rate estimates is challenging since the true value is not known. On one hand, the combined campaign-average ensemble spreads (10% for CO2, 20% for CH4) could be larger than the uncertainty of the mean posterior emission rate because this mean value represents an average over the base case ensemble members. Conversely, the potential methodological sources of flight-to-flight variability discussed above (e.g., emission aliasing) are not accounted for in the campaign-average ensemble spread but nonetheless contribute toward the overall uncertainty. In general, we focus on flight-to-flight variability when reporting results, rather than the campaign-average ensemble spread, but it is important to note that the 1σ flight-to-flight variability includes both real variability in emissions and methodological uncertainty.
NYC posterior emission rates
We also calculated mean posterior emission rates for the 5 boroughs of NYC: (45 ± 18) kmol CO2 s–1 and (0.20 ± 0.07) kmol CH4 s–1 (1σ variability across flights). As was the case for the urban area totals, inventory CO2 emission rates calculated for NYC using dates and hours representative of our flights are consistent with our posterior values. The representative NYC emission rate for ACES was (45 ± 9) kmol CO2 s–1, while for Vulcan it was (52 ± 10) kmol CO2 s–1. The fact that the mean NYC posterior emission rate for CO2 is in line with expectations based on representative inventory emissions (within 15% for both inventories) is an encouraging sign that the inversion is able to estimate emissions at spatial scales smaller than the wider urban area. Further details of the spatial structure of prior and posterior fluxes in central NYC are given in SI Section S8.
It is difficult to compare our posterior emission rates directly to the NYC self-reported inventory (SRI) because that inventory includes a mixture of scope 1, 2, and 3 emissions and does not provide scope 1 totals for each emitted species. The 2019 SRI total for CH4 is 0.21 kmol s–1 (New York City Mayor’s Office of Sustainability, 2020), but this is dominated by waste emissions (0.17 kmol s–1), which are mainly scope 3 (as landfill waste is exported from the city; MacWhinney and Klagsbald, 2017). The SRI CH4 emission estimate for natural gas distribution leakage is a scope 1 value; for 2019, this was calculated to be 0.041 kmol s–1 (21% of our posterior emission total). Considering that previous studies (Plant et al., 2019; Floerchinger et al., 2021) have shown that the majority of New York CH4 emissions are from thermogenic (i.e., natural gas) sources, it is likely that the SRI underestimates scope 1 thermogenic emissions. More measurements of NYC C2H6 emissions could enable a direct quantitative assessment of thermogenic SRI CH4 emissions. CO2 emissions in the NYC SRI are also comprised of multiple scopes, but a recent study (Gurney et al., 2021) has shown that total emissions from scope-1-only source categories in the NYC SRI are lower than the corresponding total Vulcan emissions for these categories by 19%.
Posterior–prior spatial differences
The spatial distributions of the difference between posterior and prior fluxes are shown in Figure 2. Prior emissions were adjusted upward throughout the urban area for Edgar v5.0 CO2, Edgar v5.0 CH4, and GEPA CH4. Conversely, the large prior CH4 emissions in Edgar v4.2 (distributed according to population) were adjusted downward throughout the urban area. Posterior fluxes estimated using the Vulcan and ACES CO2 priors showed an upward adjustment throughout most of the urban area, but a downward adjustment in central New York City (specifically in the southern part of Manhattan—see SI Figure S8.1 for a zoomed in map). The spread across posterior fluxes using 3 different priors for this area is smaller than the spread across the prior fluxes themselves. However, further research is required, including higher resolution transport models, to understand whether small-scale spatial features observed in these posterior emission maps reflect real spatial patterns in emissions or artifacts of the relatively coarse resolution of the transport models used (see SI Section S2 for more details regarding the transport model configuration). New York City is a very complex land/urban/water interface, and increased resolution would allow for better representation of smaller spatial features that could alter the local circulations.
Emission ratio
The average posterior CH4:CO2 emission ratio for the New York–Newark urban area is (4.9 ± 0.7) nmol µmol–1 (1σ variability across flights), and for New York City it is (4.4 ± 1.1) nmol µmol–1. As for the posterior CO2 and CH4 emission rate estimates, this ratio should be interpreted as an average over the dates and times (daytime, nongrowing season) of our flights, not an annually representative value. Nonetheless, it is interesting to compare our emission ratio against enhancement ratios for the New York–Newark urban area reported by previous studies. Plant et al. (2019) reported an average CH4:CO2 enhancement ratio of 7.2 (+1.4/–1.0) nmol µmol–1 (95% confidence interval), taken over 10 flight days during April and May 2018. A lower enhancement ratio of 3.7 (+0.3/–0.2) nmol µmol–1 (95% confidence interval) was reported by Floerchinger et al. (2021) for a single (“winter”) flight in March 2018 (see SI Figure S7.1 for a comparison plot).
Plant et al. (2019) also estimated a CH4 emission rate for the urban area, which they calculated as the product of their CH4:CO2 enhancement ratio and an inventory estimate of CO2 emissions. They derived a total CH4 emission rate that was 2.7 times larger than the GEPA inventory value, in relatively close agreement to the results of our study (2.4 times larger than the GEPA inventory). At first sight, this appears at odds with the fact that the Plant et al. (2019) enhancement ratio is significantly larger than our estimated emission ratio. This apparent discrepancy is reconciled by considering that the CO2 inventory emission rate used by Plant et al. (2019) is lower than our mean posterior CO2 emissions. This is not unexpected, as CO2 emission inventories suggest emissions during April and May are typically lower than emissions during November, February, and March (as shown in SI Figure S5.1).
Conclusions
Our posterior emission rate estimates show good agreement with representative values from the ACES and Vulcan CO2 inventories. In contrast, none of the available CH4 inventories fell within a factor of 2 of our posterior estimates. It is of particular policy relevance that our posterior CH4 emission rate estimate for the New York–Newark urban area was 2.4 times larger than the GEPA inventory. Taken together with previous studies in other cities, there is strong evidence that the GEPA inventory widely underestimates urban CH4 emissions.
Underestimation of urban CH4 emissions could impede urban emission reduction policies if cities neglect to address important (but difficult to quantify) emission sources such as natural gas leaks. Accurate, spatially resolved bottom-up estimates of U.S. CH4 emissions would allow cities to develop better informed emission reduction policies, but uncertainty in the magnitude and location of key sources (e.g., emissions from natural gas distribution) makes it difficult to compile such an inventory with the required accuracy.
Aircraft measurements enable top-down estimates of urban greenhouse gas emissions at spatial scales representative of the whole urban area. To track changes in annual emissions, it is necessary to design aircraft sampling campaigns to include repeat flights covering weekly and seasonal timescales of emission variability. In addition, as suggested by Lopez-Coto et al. (2020), increasing the number of simultaneous measurements might at least partially alleviate the flight-to-flight variability caused by irregular sampling. In the case of New York City, increasing the seasonal coverage of flights to better quantify the poorly understood seasonal cycle of CH4 emissions is a particular priority. Regular aircraft sampling could be especially valuable in cases such as New York that lack an extensive network of urban tower measurements.
Inverse modeling can play a key role going forward by providing independent quantification of emissions in near-real time. It can be used to evaluate emission estimates in the past (if measurements exist), and it can be used to update inventories and tailor policy in the present. There is an inherent positive feedback loop between top-down inverse modeling emission estimates and inventory development. Top-down emission estimates can be used to improve the models and data on which the inventories rely, resulting in better inventory products. These improved inventories then provide more accurate priors for subsequent inverse modeling studies, resulting in improved emission estimates and reduced posterior uncertainties. Together, top-down and bottom-up emission estimation approaches can effectively guide emission reduction policy in urban areas.
Data accessibility statement
Aircraft measurement data used in this study are available from the Harvard Dataverse at https://doi.org/10.7910/DVN/UECQXX.
Supplemental files
The supplemental files for this article can be found as follows:
Appendix A. Further details regarding dispersion model configuration, domains, footprints, inventory emissions, biospheric influence, emission ratios, posterior spatial distribution, inverse modeling configuration and sensitivity tests. (PDF)
Appendix B. Prior and posterior emission rates from individual ensemble members, and representative inventory emission rates. (XLSX)
Acknowledgments
The authors thank FAA air traffic control in and around NYC for their support with our flight patterns and the Jonathan Amy Facility for Chemical Instrumentation at Purdue University for continued assistance in maintaining equipment in ALAR. They also thank Róisín Commane for the loan of her Picarro analyzer during the November 2019 flights.
Funding
We thank the National Institute of Standards and Technology (NIST) for financial support of this work, via Award Nos. 70NANB16H262, 70NANB19H167, and 70NANB21H021. K. R. Gurney acknowledges the support from the National Aeronautics and Space Administration (NASA) Grant NNX14AJ20G and the NASA Carbon Monitoring System program: “Understanding User Needs for Carbon Monitoring Information” (subcontract 1491755).
Competing interests
The authors declare that no competing interests exist. Paul B. Shepson is an associate editor at Elementa. He was not involved in the review process of this article. Certain commercial equipment, instruments, or materials are identified in this paper in order to specify the experimental procedure adequately. Such identification is not intended to imply recommendation or endorsement by NIST nor is it intended to imply that the materials or equipment identified are necessarily the best available for the purpose.
Author contributions
Contributed to conception and design: JRP, ILC, KDH, CRF, AK, JRW, PBS.
Contributed to acquisition of data: JRP, ILC, KDH, JT, RK, TJ, BHS, CPL, CKG, LRH, KRG, GSR, JL, SG, PBS.
Contributed to analysis and interpretation of data: JRP, ILC, KDH, CRF, CKG, GSR, SG, AK, PBS.
Drafted and/or revised the article: JRP, ILC, KDH, CPL, LRH, GSR, SG, AK, PBS.
Approved the submitted version for publication: JRP, ILC, KDH, JT, RK, TJ, BHS, CRF, CPL, CKG, LRH, KRG, GSR, JL, SG, AK, JRW, PBS.
References
How to cite this article: Pitt, JR, Lopez-Coto, I, Hajny, KD, Tomlin, J, Kaeser, R, Jayarathne, T, Stirm, BH, Floerchinger, CR, Loughner, CP, Gately, CK, Hutyra, LR, Gurney, KR, Roest, GS, Liang, J, Gourdji, S, Karion, A, Whetstone, JR, Shepson, PB. 2022. New York City greenhouse gas emissions estimated with inverse modeling of aircraft measurements. Elementa: Science of the Anthropocene 10(1). DOI: https://doi.org/10.1525/elementa.2021.00082
Domain Editor-in-Chief: Detlev Helmig, Boulder AIR LLC, Boulder, CO, USA
Associate Editor: Isobel J. Simpson, Department of Chemistry, University of California, Irvine, CA, USA
Knowledge Domain: Atmospheric Science