Global distribution of methane emissions, emission trends, and OH concentrations and trends inferred from an inversion of GOSAT satellite data for 2010–2015

2.1 GOSAT observations

The TANSO-FTS instrument on board the Greenhouse Gases Observing Satellite (GOSAT) observes column-averaged dry-air methane mixing ratios by solar backscatter in the SWIR with near-unit sensitivity down to the surface (Butz et al., 2011). The satellite is in polar sun-synchronous orbit. Observations are made at around 13:00 local time for circular pixels of 10 km diameter. In the default observation mode, the pixels are separated by ∼ 250 km along track and cross track, with repeated observation of the same pixels every 3 days. Denser observations are also made in target mode over features of interest. GOSAT spectra have shown no significant drift or degradation of data quality since the beginning of the record (Kuze et al., 2016). We use the University of Leicester version 7 CO₂ proxy retrieval over land (Parker et al., 2011, 2015) from January 2010 to December 2015 in order to have even observations of all seasons. The single-observation precision is 13 ppb, and the relative (regional) bias is 2 ppb compared to ground-based column-averaged dry-air mole factions from the Total Carbon Column Observing Network (TCCON; Buchwitz et al., 2015). Other retrievals of GOSAT data are consistent with the University of Leicester product (Buchwitz et al., 2015). Figure 1 illustrates the GOSAT data ingested in our inversion, representing a total of 1 211 468 retrievals. Glint data over the oceans and data poleward of 60^∘ are not included because of seasonal sampling biases (Turner et al., 2015).

Figure 12010–2015 average of the GOSAT methane dry column mixing ratios used in our inversion. Data are from the University of Leicester version 7 CO₂ proxy retrieval (Parker et al., 2015), excluding glint observations over the oceans and observations poleward of 60^∘. GOSAT pixels are of 10 km circular diameter and are inflated here to 0.5^∘ for visibility. The red stripes are an averaging artifact as these retrievals are from towards the end of the 2010–2015 time period when methane was higher.

Download

2.2 Prior estimates

The inversion requires prior estimates and error statistics for all components of the state vector including methane emissions on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ GEOS-Chem grid (1009 ice-free land-containing grid cells with prior emissions larger than $\mathrm{8}\times {\mathrm{10}}^{-\mathrm{3}}$ Mg km⁻² a⁻¹, covering 99 % of global emissions), 2010–2015 linear emission trends on the same grid, and global mean OH concentrations for individual years 2009–2015 (2009 is only used for initialization), for a total of 2025 state vector elements.

Figure 2Prior estimates of methane emissions from wetlands, livestock, coal mining, oil and gas, wastewater and landfills, and other sources. Values are 2010–2015 averages and are shown on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ GEOS-Chem grid used for the inversion. Global totals for each source type are given in Table 1.

Download

Table 1 gives our global prior inventory with the contributions from different source types, and Fig. 2 shows the spatial distributions. Monthly wetland emissions for individual years are from the WetCHARTS v1.0 extended ensemble mean (Bloom et al., 2017). For anthropogenic emissions we use the EDGAR v4.3.2 global emission inventory for 2012 (https://edgar.jrc.ec.europa.eu/, last access: 1 December 2017; Janssens-Maenhout et al., 2019) as worldwide default, including additional information from EDGAR to subset the “fuel exploitation” emissions category into oil and gas and coal mining. Over the continental United States, we replace EDGAR v4.3.2 with a gridded version of the US EPA greenhouse gas inventory (Maasakkers et al., 2016). In Canada and Mexico, we use the oil and gas emissions from Sheng et al. (2017). Anthropogenic emissions are assumed as aseasonal for lack of better prior information except for manure management and rice cultivation. Seasonal scaling of manure management emissions is done using the temperature dependence of Maasakkers et al. (2016). Seasonal scaling of rice cultivation emissions is based on Zhang et al. (2016). Daily global open fire emissions are from QFED (Darmenov and da Silva, 2013). Termite emissions are from Fung et al. (1991). Emissions from geological macroseeps (oil and gas seeps and mud volcanoes) are based on Etiope (2015) and Kvenvolden and Rogers (2005). For areal seepage, we use the sedimentary basins (microseepage) and potential geothermal seepage maps from Kvenvolden and Rogers (2005) with the emission factor previously used by Lyon et al. (2015). Over the United States, we use the sedimentary basin map from the Energy Information Administration (EIA, 2016) and basin-specific emission factors from Etiope and Klusman (2010). While global geological emissions have previously been estimated to be over 50 Tg a⁻¹ (Kirschke et al., 2013), Petrenko et al. (2017) showed that based on ice core measurements they should be no higher than 15 Tg a⁻¹.

Table 1Prior global estimates of methane sources and sinks (mean 2010–2015 values).

^* Including fossil fuel combustion, industrial processes, and agricultural field burning.

Download Print Version | Download XLSX

Construction of the prior error covariance matrix S_a requires estimates of error variances for the prior emissions on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid. For wetland emissions, we use the standard deviation of $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ annual emissions from the WetCHARTs ensemble members (Bloom et al., 2017). The error variance averages 58 % on the grid level. For US anthropogenic emissions and oil and gas emissions in Canada and Mexico, we use the scale-dependent error variances from Maasakkers et al. (2016). For lack of better information, we assume 50 % error standard deviation for EDGAR v4.3.2 emissions (Turner et al., 2015) and 100 % for non-wetland natural emissions. The diagonal terms of S_a are then constructed by adding the error variances of individual source types for $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells in quadrature, capping total errors at 50 %. We assume no error spatial covariance on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid so that S_a is diagonal. This is a reasonable assumption for anthropogenic emissions (Maasakkers et al., 2016), though errors on wetland emissions may still be correlated on that scale (Bloom et al., 2017).

Our state vector in the inversion includes linear emission trends for $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells over the 2010–2015 period, superimposed on interannual variability in the case of wetlands and fires. Our global prior estimate of mean methane emissions for the 2010–2015 period exceeds the sinks by 13 Tg a⁻¹ (Table 1), which drives a 5 ppb a⁻¹ increase in methane concentrations over that period, even in the absence of an emission trend. Therefore our prior estimate of linear emission trends for individual $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells is zero, with an absolute error standard deviation of 10 % of the local prior emissions over the 2010–2015 time period (1.7 % a⁻¹). This error standard deviation is based on trend estimates for North America inferred from GOSAT data (Turner et al., 2016; Sheng et al., 2018a).

The prior estimate of the global tropospheric OH concentration is based on a GEOS-Chem full-chemistry simulation (Wecht et al., 2014) that yields a methane lifetime ${\mathit{\tau }}_{{\mathrm{CH}}_{\mathrm{4}}}^{\mathrm{OH}}$ of 10.6 years, consistent with the best estimate inferred from the methyl chloroform proxy (Prather et al., 2012) and the 9.7±1.5 years estimate from the ACCMIP model ensemble (Naik et al., 2013). Here and elsewhere, ${\mathit{\tau }}_{{\mathrm{CH}}_{\mathrm{4}}}^{\mathrm{OH}}$ is defined as the ratio between the total mass of atmospheric methane (including the stratosphere) and the annual loss rate from oxidation by OH below the tropopause. The uncertainty in the methane lifetime is about 10 % (Prather et al., 2012), but the uncertainty on OH interannual variability is less, about 3 % (Holmes et al., 2013). We assume a 3 % error standard deviation in the global annual mean OH concentration for our standard inversion but also conduct a sensitivity study with 10 % error standard deviation. We further conduct an inversion taking the OH trend over the 2010–2015 period as linear and assuming in that case error standard deviations of 10 % for the mean global OH concentration and 5 % a⁻¹ (absolute) for the linear trend. Scaling of global OH concentrations in the inversion is done without modifying the spatial or seasonal OH distribution. Zhang et al. (2018) found that inversions of atmospheric methane data using the 3-D GEOS-Chem OH fields give consistent results with inversions using other global OH distributions from the ACCMIP model ensemble (Naik et al., 2013).

2.3 Forward model

We use the GEOS-Chem CTM v11-01 at $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid resolution (Wecht et al., 2014; Turner et al., 2015) as forward model for the inversion. The model is driven with 2009–2015 MERRA-2 meteorological fields (Bosilovich et al., 2016) from the NASA Global Modeling and Assimilation Office (GMAO). Atmospheric methane concentrations are initialized on January 2009 using the previous GOSAT inversion results of Turner et al. (2015), shown in that work to be unbiased compared to surface and aircraft background data including for the tropospheric meridional gradient.

The loss from oxidation by tropospheric OH is computed with archived 3-D monthly fields of OH concentrations from a GEOS-Chem full-chemistry simulation as described by Wecht et al. (2014). Local tropopause information is from the MERRA-2 data. The global loss rate for individual years is optimized in the inversion by uniform scaling of the OH concentrations. Other minor loss terms include stratospheric oxidation computed with archived monthly loss frequencies from the NASA Global Modeling Initiative model (Murray et al., 2012), tropospheric oxidation by Cl atoms computed using archived Cl concentration fields from Sherwen et al. (2016) and the reaction rate constant from Allan et al. (2007), and soil uptake as described by Fung et al. (1991) with temperature-based seasonality based on Ridgwell et al. (1999). The loss from oxidation by Cl totals 9 Tg a⁻¹, intermediate between the 12–13 Tg a⁻¹ estimated by Hossaini et al. (2016) using the TOMCAT chemical transport model and 5.3 Tg a⁻¹ estimated by Wang et al. (2019) in a GEOS-Chem simulation with full accounting of tropospheric chlorine. These minor sinks are not optimized in the inversion.

The GEOS-Chem simulation of GOSAT methane columns features a latitude-dependent background bias that needs to be corrected (Turner et al., 2015). This bias likely reflects a model overestimate of methane in the extratropical stratosphere (Saad et al., 2016), which is common across global models due to excessive meridional transport in the stratosphere (Patra et al., 2011) and was first seen in a SCIAMACHY inversion using the TM5 chemical transport model (Bergamaschi et al., 2007). Stanevich (2018) found a significant difference in methane columns simulated by GEOS-Chem at $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ compared to $\mathrm{2}{}^{\circ }\times \mathrm{2.5}{}^{\circ }$ resolution, but we find that this difference is mainly in the stratosphere (Appendix A). We remove the background bias by applying the latitudinal correction based on background grid cells from Turner et al. (2015), recomputed with the University of Leicester v7 GOSAT proxy retrieval (Parker et al., 2015) and the MERRA-2 meteorological fields. The mean model–GOSAT difference in column mean mixing ratio for background $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells is fitted to a second-order polynomial of latitude:

where θ is the latitude in degrees, and ξ is the model correction in ppb. This correction is similar to Turner et al. (2015), who used $\mathit{\xi }=\left(\mathrm{5}{\mathit{\theta }}^{\mathrm{2}}-\mathrm{5}\mathit{\theta }\right)\times {\mathrm{10}}^{-\mathrm{3}}-\mathrm{0.5}$. A seasonal bias remains after application of this correction, and we fix it by removing the zonal monthly mean concentration differences averaged over rolling 12^∘ latitudinal bands. This seasonal bias may be due to errors in the seasonality of emissions or atmospheric transport (Saad et al., 2016; Bader et al., 2017; Stanevich, 2018). We find that the seasonal correction does not affect the inversion results significantly, as shown in Appendix B, where we optimize emissions for individual seasons separately without applying a seasonal correction.

2.4 Observational error covariance matrix

The observational error covariance matrix S_O includes contributions from random instrument and forward model errors. We construct it applying the residual error method of Heald et al. (2004) using the 2010–2015 time series of local methane column differences $\mathrm{\Delta }=y_{\text{GEOS-CHEM, \hspace{0.17em}prior}}-y_{\text{GOSAT}}$ for individual $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells between the GEOS-Chem model with prior estimates (emissions and OH concentrations) and the GOSAT observations after background bias correction. The mean difference $\stackrel{\mathrm{‾}}{\mathrm{\Delta }}=\stackrel{\mathrm{‾}}{y_{\text{GEOS-CHEM, \hspace{0.17em} prior}}-y_{\text{GOSAT}}}$ is to be corrected in the inversion, while the residual error ${\mathrm{\Delta }}^{\prime }=\mathrm{\Delta }-\stackrel{\mathrm{‾}}{\mathrm{\Delta }}$ is taken as the observational error. Statistics of Δ^′ define the observational error variance (diagonal of the observational error covariance matrix). The same method was previously used in the satellite-based methane inversions by Wecht et al. (2014) and Turner et al. (2015). The resulting observational error standard deviation averages 13 ppb. The mean instrument error standard deviation is 11 ppb (Parker et al., 2015), implying that most of the observational error is generally from the instrument rather than from the forward model. This would indeed be expected for the random error of individual measurements. For a given measurement, if the local error standard deviation computed by the residual error method is smaller than the reported measurement precision, then we use the latter instead; this is the case for 10 % of retrievals. All observational error standard deviations are set to be at least 10 ppb (this threshold affects 8 % of retrievals). S_O is taken to be diagonal for lack of better information, but the general effect of error correlation in the observations is accounted for in the inversion by a regularization factor (Sect. 2.5).

2.5 Inversion procedure

We perform inversions with two different specifications of prior error variance statistics: normal and lognormal. Assumption of normally distributed errors enables a linear optimization problem with an analytical solution including closed-form error characterization (Rodgers, 2000). Assumption of lognormal errors may be more appropriate for modeling the high tail of the probability density function (Zavala-Araiza et al., 2015) and also has the advantage of enforcing positive solutions (Miller et al., 2014), but the optimization problem is then nonlinear. By comparing the two approaches we can evaluate consistency in results.

Both inversions minimize the Bayesian cost function J(x) (Rodgers, 2000):

where x is the state vector, x_a is the prior estimate, S_a is the prior error covariance matrix, F(x) is the simulation of observations y by the GEOS-Chem model, S_O is the observational error covariance matrix, and γ is a regularization factor (Brasseur and Jacob, 2017). The variances in S_O are underestimated because of correlation in the observational error that is missing in the diagonal formulation of S_O and is difficult to quantify. We use γ to scale the original diagonal S_O to get an optimal covariance matrix to be used in the inversion. Zhang et al. (2018) showed in an observing system simulation experiment (OSSE) for inversion of methane satellite data that a regularization factor γ=0.05 adjusts the variances optimally and prevents overfitting. This was done by calculating the likelihood at $\widehat{\mathit{x}}$ for a range of values of γ. Diagnosis of overfit and optimization of γ is readily done in an OSSE such as in Zhang et al. (2018) where the “true” solution is known. Here we find that using γ=1 (as in the pure Bayesian statement of the optimization problem) produces checkerboard patterns in the solution that are likely spurious. We choose γ=0.05 consistent with Zhang et al. (2018) for our base inversion as providing the best balance between prior and observational terms in the posterior value of the cost function. We examine the sensitivity to the choice of γ by conducting a sensitivity inversion with γ=0.1.

Further balancing of the cost function is needed because the global OH concentration and its interannual variability are represented by only seven state vector elements, while the emissions on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid are represented by 1009 elements. To provide equal weight to OH and emissions for explaining global methane trends, we increase the weight of the OH terms in the cost function (through the OH components of S_a) by the ratio of the number of state vector elements 1009∕7 so that from a cost-function perspective, a change in OH and global methane emissions is equally expensive. The sensitivity inversion assuming 10 % prior error standard deviation on OH instead of 3 % is equivalent to decreasing this weighting by a factor of 11.

The GEOS-Chem forward model y=F(x) relating methane column concentrations y to the state vector x is essentially linear. There is a small nonlinearity from the optimization of OH concentrations because changes in the methane concentrations affect the loss rate (Houweling et al., 2017), which we neglect because changes in methane concentrations are small, and methane is well mixed globally. We therefore express the forward model as $F\left(\mathit{x}\right)=\mathbf{K}\mathit{x}+\mathit{c}$, where $\mathbf{K}=\partial y/\partial x$ is the Jacobian matrix of the model, and c is an initialization constant (January 2009 concentrations taken from Turner et al., 2015). Replacing F(x)=Kx in Eq. (2) and subtracting the initialization constant c from the observations, the minimization problem $\text{d}J\left(\mathit{x}\right)/\text{d}\mathit{x}=\mathbf{0}$ has an analytical solution for the optimal posterior solution $\widehat{\mathit{x}}$ (Rodgers, 2000):

The posterior error covariance matrix $\widehat{\mathbf{S}}$ describing the error statistics of $\widehat{\mathit{x}}$ is given by

and the averaging kernel matrix ($\mathbf{A}=\partial \widehat{\mathit{x}}/\partial \mathit{x}$) defining the sensitivity of the solution to the true state is given by

The trace of the averaging kernel matrix defines the degrees of freedom for signal (DOFS) of the inversion, that is the number of pieces of information on the state vector that can be gained from the observing system.

The analytical solution as described by Eqs. (3)–(5) requires the explicit construction of the Jacobian matrix K characterizing the GEOS-Chem model. We do this column by column, with GEOS-Chem simulations perturbing each element of the state vector independently. This is readily achievable, even for 2025 state vector elements as a massively parallel computation. Sparse matrix algebra is used where possible in solving Eqs. (3)–(5), taking advantage of the diagonal structure of the error covariance matrices.

The analytical solution to the Bayesian optimization problem requires assumption of Gaussian errors, but this allows for the possibility of negative values of state vector elements. Negative emissions could conceivably be attributed to locally strong soil uptake or oxidation by Cl atoms but may also be unphysical (Miller et al., 2014). We can enforce positivity in the Bayesian solution by optimizing for ln (x) instead of x, with normal Gaussian errors specified for ln (x) (corresponding to lognormal errors for x). The model is then nonlinear, so that the solution and the corresponding error statistics must be found iteratively with an updated Jacobian matrix ${\mathbf{K}}_N^{\prime }=\partial \mathit{y}/\partial \ln \mathit{x}$ at each iteration N. This recomputation is immediate using the previously derived Jacobian matrix K for the linear problem, since the individual scalar elements $\partial y_i/\partial \ln \left(x_i\right)$ of K^′ are related to those of K by $\partial y_i/\partial \ln \left(x_j\right)=x_j\partial y_i/\partial x_j$. Thus only a simple scaling of the linear Jacobian matrix is required at each iteration. This conversion to log space is done only for the emissions component of x. Emission trends and global OH concentrations are still optimized with normal error distributions, and no scaling is applied to those rows of the Jacobian.

Optimizing emissions in log space means that the best posterior estimate is for the median of emissions instead of the mean. The mean and the median of the lognormal distribution are not equal, so results cannot be summed over grid squares to provide a best estimate of the mean. For this reason, analysis of aggregate and global emissions and sinks will be done with the inversion using normal errors.

The iterative solution for the inverse problem with lognormal errors is obtained with the Levenberg–Marquardt method (Rodgers, 2000) for each iteration N:

where ${\mathit{x}}^{\mathbf{\prime }}=\ln \mathit{x}$, the initial guess ${{\mathit{x}}^{\mathbf{\prime }}}_{\mathrm{0}}$ is the prior estimate, and κ is a coefficient for the iterative approach to the solution that is set to 100 to start and is gradually decreased as the solution is approached. The prior error covariance matrix ${\mathbf{S}}_{\text{a}}^{\prime }$ (diagonal elements $s_{\text{A}}^{\prime }$) defining error variances for ln x_a is derived from the previously described prior error covariance matrix S_a (diagonal elements s_A) by scaling the error variances for the individual elements:

2.6 Error correlations between global estimates of sources and sinks

Inversion results for the spatial distributions of emissions and trends on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid are mainly informed by local and regional patterns of methane concentration. However, implied inversion results for the global methane emission and its trend may be significantly correlated with those for the global tropospheric OH concentration and its trend. Some separation is expected because sources of methane have a different spatial and seasonal imprint on the global methane distribution than the OH sink (Zhang et al., 2018), but it is important to quantify the error correlation, i.e., the extent to which adjustments to the global methane emission and its trend may be aliased by adjustments to the global OH concentration and its trend.

To do this we reduce the dimensionality of the inverse analysis by collapsing global emissions and trends into one state vector element each. Following Calisesi et al. (2005), if the state vector can be transformed using a summation matrix W as

then the averaging kernel matrix of the reduced system (A_red) is given by

where ${\mathbf{W}}^{*}=({\mathbf{W}}^T\mathbf{W}{)}^{-\mathrm{1}}{\mathbf{W}}^T$ is the generalized pseudo-inverse of W. Our original state vector x in this case includes mean 2010–2015 emissions and their linear trends on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid and the global mean tropospheric OH concentration for 2010–2015 and its linear trend. Again, the minor sinks in Table 1 are not optimized and are maintained instead at their prior values. We apply the summation matrix W to the emission terms and thus reduce the state vector to four elements defining the global methane budget (global mean emission, global mean OH concentration, global emission trend, global OH trend). The off-diagonal terms of the reduced averaging kernel matrix A_red then measure the extent to which differences relative to the true state are aliased between sources and sinks in the optimization of this global budget. The advantage of this summation approach, as compared to a global inversion including just four elements, is that the distribution of methane emissions and its trends is still optimized.

Figure 3Comparisons of observed methane concentrations to the GEOS-Chem forward model using either prior or posterior (optimized) estimates of 2010–2015 emissions and OH concentrations. Panels (a) and (b) show differences between the model and GOSAT observations for 2010–2015 means on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid. Panels (c) and (d) show the monthly time series of the differences averaged over latitude bands. Panels (e–g) show independent 2010–2015 comparisons to global observations from NOAA surface stations, HIPPO aircraft meridional cross sections over the Pacific (2010 and 2011, with the model sampled along the flight tracks), and TCCON. Reduced major axis (RMA) regressions are as shown along with the 1 : 1 line (in grey). HIPPO observations are averaged over GEOS-Chem grid cells. The NOAA surface stations and HIPPO aircraft measure local methane dry-air mole fractions, while TCCON measures column-averaged dry-air mole factions. We apply the same latitudinal and seasonal corrections to TCCON that we applied to GOSAT.

Download

3.1 Spatial distribution and source attribution of methane emissions

Figure 4 shows the global distribution of mean 2010–2015 posterior emissions from the base inversion and from the sensitivity inversion assuming lognormal errors in the prior emission estimates. Correction patterns are very similar between the two inversions. Small negative emissions are found in the base inversion for 6 of the 1009 optimized grid cells. The inversion assuming lognormal errors does not allow these negative emissions. Downward corrections tend to be smaller in the inversion assuming lognormal errors, while positive corrections are larger and more concentrated in a few grid cells, as would be expected from the different shapes of the error standard deviation distributions.

Figure 4b shows the diagonal terms of the averaging kernel matrix for the base inversion (averaging kernel sensitivities), measuring the ability of the observations to constrain the inversion. The trace of the averaging kernel matrix (DOFS = 128) measures the number of independent pieces of information constrained by the inversion. A Bayesian inversion without correcting for overfit (γ=1 in Eq. 3) would erroneously produce much higher DOFS. We find that the inversion provides strong constraints on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid for source regions in East Asia, central Africa, and South America. Averaging kernel sensitivities are generally weaker over North America and in Europe, indicating that the inversion provides more diffuse spatial information in these regions.

We find that the EDGAR v4.3.2 inventory prominently overestimates anthropogenic emissions over eastern China, likely from coal production, and around the Persian Gulf, likely from oil and gas production. The finding of an EDGAR overestimate in China is consistent with previous global inversions of GOSAT data using EDGAR v4.1, v4.2, and v4.2FT2010 as prior estimates (Monteil et al., 2013; Thompson et al., 2015; Alexe et al., 2015; Turner et al., 2015; Pandey et al., 2016) and a regional inversion using EDGAR v4.3.2 (Miller et al., 2019). The overestimate of coal mining emissions may be because standard IPCC emission factors used by EDGAR v4.2 were too high for Chinese coal mines, and recovery of coal mine methane is not sufficiently taken into account (Peng et al., 2016). Emission factors were decreased in EDGAR v4.3.2 (Janssens-Maenhout et al., 2019), but we still find an overestimate. We find that EDGAR underestimates emissions over Japan and Southeast Asia, where rice cultivation is the largest anthropogenic source, but there are also large wetland emissions. There are also large corrections in wetland areas of central Africa, South America, and North America.

We do not find large correction factors over the United States, except for the southeastern coast which is likely due to an overestimate of methane emissions from coastal wetlands in the prior WetCHARTs inventory. This overestimate of US coastal wetland emissions in WetCHARTs is consistent with a previous inversion of aircraft observations over the Southeast United States by Sheng et al. (2018b) and may be explained by low soil organic carbon in these ecosystems (Holmquist et al., 2018) and/or the overestimated impacts of partial wetland land-cover classes predominant in the southeastern United States (Lehner and Döll, 2004; Bloom et al., 2017). Previous inversions found factor of 2 underestimates of EDGAR v4.2 emissions of the South Central United States (Miller et al., 2013; Turner et al., 2015), but we do not find such an underestimate here and attribute this to our use of the gridded version of the US EPA inventory as prior estimate (Maasakkers et al., 2016). EDGAR v4.2 allocated oil and gas emissions mainly according to population, which greatly underestimates emissions in oil and gas production regions in the South Central United States (Maasakkers et al., 2016).

Figure 5Global methane emissions by source type in the prior estimate for the inversion (Table 1, “other” includes fossil fuel combustion, industrial processes, and agricultural field burning) and in the posterior estimate. Values are 2010–2015 means. The attribution to source types in the posterior estimate is done by assuming that the relative contributions of different source types in individual $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells are correct in the prior estimate. Posterior estimates are from the base inversion, and error bars show the ranges of results from the inversion ensemble.

Download

Improved estimates of global methane emissions for the individual source types of Table 1 can be inferred from our results by assuming that the relative contributions from different source types in a given $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cell are correct in the prior inventory. The global posterior estimate for a given source type is then obtained by applying the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ posterior / prior ratios from Fig. 4 to the distribution of source types in Fig. 2. Results in Fig. 5 indicate little change to 2010–2015 average emissions compared to the global prior inventory by source type, even though there are large regional reallocations. Coal mining emissions decrease by 29 %, mainly due to China, and rice cultivation and livestock increase by 15 % and 8 % respectively, mainly driven by the tropics.

There has been particular interest in quantifying emissions from oil and gas exploitation because of the potential for large reductions of these emissions through simple control measures (Zavala-Araiza et al., 2015; Alvarez et al., 2018). The EDGAR v4.3.2 national oil and gas emission totals can differ greatly from the national (spatially unresolved) totals reported by individual countries to the United Nations Framework Convention on Climate Change (UNFCCC, 2017). This is shown in Fig. 6 with national oil and gas emissions from the top 10 countries in either the EDGAR v4.3.2 or UNFCCC inventories. We can estimate national oil and gas emission totals from our inversion by again assuming that the relative contributions of oil and gas to total emissions in individual $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells are correct and by further mapping the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ correction factors to the $\mathrm{0.1}{}^{\circ }\times \mathrm{0.1}{}^{\circ }$ EDGAR emission grid. The emission-weighted scaling factor is then used with the national oil and gas totals reported by EDGAR. Russia is the largest national source, but the inversion is limited in its ability to constrain oil and gas emissions there because a third of these emissions are north of 60^∘ N in EDGAR v4.3.2 (Fig. 2).

Figure 6National estimates of methane emissions from the oil and gas industry for countries in the top 10 of either the EDGAR v4.3.2 or UNFCCC inventories. Values reported by individual countries to the UNFCCC for 2012 (Annex I countries) or the closest year (non-Annex I countries: Nigeria (1994), Venezuela (1999), Algeria (2000), Iran (2000), India (2010), Saudi Arabia (2010), and China (2012)) are compared to 2012 emissions from EDGAR v4.3.2 national oil and gas totals and to the posterior values from our base inversion as described in the text. Black lines are ranges for the ensemble of inversions. A large part of Russian emissions are too far north to be effectively constrained by the inversion.

Download

Figure 72010–2015 methane emission trends and global tropospheric OH trends as optimized by the inversion of GOSAT data and corresponding averaging kernel sensitivities (diagonal terms of the average kernel matrix). The degrees of information for signal or DOFS (trace of the averaging kernel matrix) is shown inset. Panel (c) gives the global attribution of the emission trends to individual source types, with ranges from the inversion ensemble. Shaded sections of the bars indicate the contribution from the tropics (24^∘ S–24^∘ N). The vertical bars in the OH trend panel are the posterior error standard deviations from the base inversion. The 2010–2015 decreasing trend in OH concentrations is not statistically significant (95 % confidence level).

Download

Results in Fig. 6 show that the inversion generally pushes the prior EDGAR v4.3.2 estimates of oil and gas emissions toward the UNFCCC values. One would expect the UNFCCC national reports to provide better estimates than EDGAR v4.3.2 because of their use of local information (Scarpelli et al., 2018) as compared to the more generic estimates used by EDGAR on the basis of IPCC Tier 1 methodology (IPCC, 2006). Thus we find that EDGARv4.3.2 greatly underestimates emissions in Uzbekistan, which are high because of leaky infrastructure (Scarpelli et al., 2018). For Iran, Algeria, Nigeria, Saudi Arabia, and Qatar we find much lower emissions than EDGAR v4.3.2 that are more consistent with the UNFCCC data. For China we are in better agreement with EDGAR v4.3.2 than with the UNFCCC estimate, which relies on anomalously low emission factors (Larsen et al., 2015). In Venezuela we find higher emissions than both EDGAR v4.3.2 and UNFCCC. The latest available report from Venezuela to the UNFCCC dates back to 1999.

3.2 Spatial distribution and source attribution of methane emission trends

Figure 7 shows base inversion results for the linear emission trends on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid for 2010–2015 and the associated averaging kernel sensitivities. Also shown in panel d is the 2010–2015 time series of posterior OH concentrations with error standard deviations from the posterior error covariance matrix. We find no significant OH trend over the period, although uncertainties are large. The information on the spatial distribution of emission trends originates from local and regional gradients of atmospheric methane observed by GOSAT, and we find from the posterior error covariance matrix of the inversion that it is not correlated with information on OH concentrations. Thus the large posterior uncertainty in global OH concentrations does not induce any significant correlated error in the spatial distribution of emission trends. This may be expected in view of the long lifetime of methane relative to the relevant timescales for atmospheric transport.

The GOSAT data provide seven independent pieces of information (DOFS) on the spatial distribution of the emission trend. Again, a Bayesian inversion without correcting for overfit (γ=1) would erroneously indicate much higher DOFS. We find increasing emissions in the tropics and little change at higher latitudes. There are well-defined anthropogenic positive trends over China, India, and the Persian Gulf. Trends in China are in areas with dominant emissions from coal mining but also significant contributions from livestock and waste. In an inversion of surface observations, Thompson et al. (2015) previously found an increasing trend over China for 2000–2011, which they attributed to coal mining. Miller et al. (2019) found that this trend continued up to 2015 using GOSAT in a regional inversion. Trends over India are in areas of rice production but may also reflect waste management and livestock. The trend over India is 0.4 (0.3–0.5) Tg a⁻¹ (range of the inversion ensemble), consistent with the 2010–2015 trend of 0.7±0.5 Tg a⁻¹ from a regional GOSAT inversion by Miller et al. (2019). Ganesan et al. (2017) found a nonsignificant trend (0.2±0.7 Tg a⁻¹) over India for 2010–2015 using an ensemble of GOSAT, commercial aircraft (CARIBIC), and surface station methane data, but our estimate is not incompatible with their range. EDGAR v4.3.2 predicts a 0.4 Tg a⁻¹ increase in anthropogenic emissions from India between 2010 and 2012, mainly from livestock, coal, and waste based on increasing activity data (this trend is not included in our prior estimate). The trend over the United States is less well defined and not well constrained but suggests an increase over the eastern part of the country where multiple source types could contribute (Sheng et al., 2018a, b).

Figure 7c shows the attribution of the global increasing trend in emissions to individual source types, following the same assumption that was used in Fig. 5 to attribute emissions to source types. We further separate tropical and extratropical contributions. Boreal wetland trends cannot be constrained by our inversion effectively (no observations north of 60^∘ N). 43 % of the 5 Tg a⁻¹ global emission trend found in the inversion for 2010–2015 is driven by wetlands (mainly tropical), 16 % by livestock, and 11 % by oil and gas. No source type shows a global decrease. Our source attribution of the methane trend is consistent with isotopic evidence, suggesting that the increase in methane over the past decade has been driven by biogenic sources outside the Arctic (Nisbet et al., 2016; Schwietzke et al., 2016; Schaefer et al., 2016), including tropical wetlands (McNorton et al., 2016). Worden et al. (2017) previously found a decrease in biomass burning from 2001–2007 to 2008–2014 but no significant change for the 2010–2015 period investigated here. Their argument that a decrease in the biomass burning emissions would have masked the effect of an increase in fossil fuel emissions on the isotope signature of methane would not apply for our time period.

Table 2Global 2010–2015 methane budget^a.

^a From the inversion optimizing (1) mean 2010–2015 methane emissions on the $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid, (2) linear methane emission trends on that same grid, (3) global mean 2010–2015 tropospheric OH concentration, and (4) linear trend in global OH concentrations. Expected values and error standard deviations are shown. The prior estimates are described in Sect. 2.2. The posterior global emission and its trend are obtained by summing the contributions from all $\mathrm{4}{}^{\circ }\times \mathrm{5}{}^{\circ }$ grid cells, and the error standard deviations are computed accounting for posterior error correlation. Minor methane sinks totaling 61 Tg a⁻¹ are not optimized in the inversion. ^b Methane lifetime against oxidation by tropospheric OH, computed as the ratio between the total atmospheric mass of methane (including the stratosphere) and the annual loss rate from oxidation by OH in the troposphere.

Download Print Version | Download XLSX

3.3 Global methane budget and trends

The previous sections showed that our inversion of the GOSAT data is able to provide relatively fine information on the spatial distribution of methane emissions (DOFS = 128) as well as some information on the spatial distribution of 2010–2015 emission trends (DOFS = 7). This information on the spatial distribution originates from local and regional gradients of atmospheric methane observed by GOSAT. We now examine to what extent error correlations may limit our ability to independently quantify the global emission of methane, the global tropospheric OH concentrations, and their respective trends.

Figure 8Constraints on the global 2010–2015 methane budget from our inversion of GOSAT data. The lines show the rows of the averaging kernel matrix A_red (Eq. 9) for the reduced four-element state vector consisting of the 2010–2015 mean emission, the linear emission trend, the 2010–2015 mean tropospheric OH concentration, and the linear OH trend.

Download

Figure 9Joint probability density functions (pdfs) for the global methane budget as constrained by the 2010–2015 GOSAT data. Panel (a) shows the joint pdfs of the 2010–2015 global mean methane emission and methane lifetime against oxidation by tropospheric OH. Panel (b) shows the joint pdfs of the 2010–2015 global emission trend and the OH trend. Contours show confidence ranges from 0.1 to 0.9. The error correlation coefficients are shown inset. The tilt of the ellipse indicates the extent of error correlation.

Download

To analyze the constraints from the inversion on the global budget of methane, we collapse the inversion to the reduced four-member global state vector of 2010–2015 mean values described in Sect. 2.6 (global methane emission, global emission trend, global tropospheric OH concentration, global OH trend). We use normal errors for all state vector elements (using lognormal errors could bias the mean). Table 2 compares the prior and posterior values for this global budget. The uncertainty in global emissions and trends is likely underestimated because of the lack of prior error covariance assumed between the 1009 grid cells. The global mean tropospheric OH concentration is expressed in terms of the corresponding methane lifetime ${\mathit{\tau }}_{{\mathrm{CH}}_{\mathrm{4}}}^{\mathrm{OH}}$. Figure 8 shows the averaging kernel rows for this reduced global state vector (A_red in Sect. 2.6), measuring the sensitivity of the inversion results to the true values (diagonal terms) and the aliasing due to error correlations (off-diagonal terms). We find that the mean 2010–2015 global methane emission and OH concentration are strongly and independently constrained, with averaging kernel sensitivities near unity and little error correlation. On the other hand, there is strong negative error correlation between emission trends and OH trends, and the OH trend can only be weakly constrained. This is illustrated further in Fig. 9 with the joint probability density function (pdf) plots of the posterior estimates, where the confidence levels measure the probability of a given value, and the tilts of the ellipses measure the error correlations.

A major implication of being able to constrain the global methane emission and the global OH concentration independently is that satellite observations of atmospheric methane can provide an independent proxy for quantifying the global mean tropospheric OH concentration. Our posterior estimate of the methane lifetime ${\mathit{\tau }}_{{\mathrm{CH}}_{\mathrm{4}}}^{\mathrm{OH}}$ is 10.8±0.4 years. It is strongly constrained by the inversion, as shown by the averaging kernel sensitivity near unity, and is thus largely independent of the prior estimate of 10.6±1.1 years. So far the main method for estimating global OH has been through the methyl chloroform budget (Prather et al., 2012), but this is becoming problematic as methyl chloroform concentrations decrease, and previously minor potential sources like ocean outgassing may become significant (Wennberg et al., 2004; Liang et al., 2017). Satellite observations of methane could provide an alternative. Our inversion confirms the best estimate of global OH from the methyl chloroform budget (Prather et al., 2012) but reduces its uncertainty from 10 % to 4 %. The magnitude of reduction may be overoptimistic because of the idealized treatment of error statistics, the assumption that the global 3-D OH distribution in the forward model is correct, and the assumption that the minor sinks (Table 1) are correct. Zhang et al. (2018) present a more thorough error analysis of this potential of methane satellite observations as proxy for global OH concentrations.

We find on the other hand that there is large error correlation between our estimates of global 2010–2015 emission trends and OH trends and limited ability to constrain the OH trend. We find that most of the increase in methane is explained by increasing emissions. Our posterior estimates for the 2010–2015 trends are $+\mathrm{0.84}\pm \mathrm{0.04}$ % a⁻¹ (4.6±0.2 Tg a⁻¹) for emissions and $-\mathrm{0.2}\pm \mathrm{0.8}$ % a⁻¹ ($-\mathrm{1.0}\pm \mathrm{3.8}$ Tg a⁻¹) for OH. The joint pdf in Fig. 9 illustrates the error correlation between the two. Other factors driving the 2010–2015 atmospheric methane trend are the initial imbalance in the 2010 budget, which we can derive from the posterior estimates of the mean 2010–2015 budget imbalance and trends, and the interannual variability of wetlands emissions as represented by WetCHARTS. Figure 10 shows the contributions of these different terms to the observed 2010–2015 methane growth. 2010 was a relatively high year for tropical wetlands emissions according to WetCHARTS, which acts to dampen the overall trend. We can state with some confidence that increasing tropical emissions (Fig. 7) made an important contribution to the 2010–2015 methane trend, but any conclusion about the effect of an OH trend is highly uncertain, including in its sign. Our 2010–2015 growth rate averages 6.8 ppb a⁻¹, compared to 7.3 ppb a⁻¹ in the NOAA record (https://www.esrl.noaa.gov/gmd/ccgg/trends_ch4/, last access: 27 April 2019). The increase in the NOAA record is higher because of especially strong growth in 2014 (12.8 ppb), which is not fully captured by the linearized optimization used here. In our base inversion, this anomaly is explained by a reduced sink from OH.

Figure 10Attribution of the 2010–2015 increase in the atmospheric burden of methane. The grey bars show the trend imposed by the 2010 imbalance between sources and sinks combined with the interannual variability (IAV) of the prior estimate (mainly from wetlands). This trend decreases over the 2010–2015 period because the methane sink rises in response to the increasing methane concentration and also because wetland emissions in 2010 are higher than in other years. Purple and orange show the contributions of the 2010–2015 methane emission trends and OH trends. The apportionment of the emission trend by source region and source type is shown in Fig. 7. The OH trend has high uncertainty as discussed in the text.

Download