**Research article**| 17 Oct 2022

# High-resolution inverse modelling of European CH_{4} emissions using the novel FLEXPART-COSMO TM5 4DVAR inverse modelling system

_{4}emissions using the novel FLEXPART-COSMO TM5 4DVAR inverse modelling system

Peter Bergamaschi Arjo Segers Dominik Brunner Jean-Matthieu Haussaire Stephan Henne Michel Ramonet Tim Arnold Tobias Biermann Huilin Chen Sebastien Conil Marc Delmotte Grant Forster Arnoud Frumau Dagmar Kubistin Xin Lan Markus Leuenberger Matthias Lindauer Morgan Lopez Giovanni Manca Jennifer Müller-Williams Simon O'Doherty Bert Scheeren Martin Steinbacher Pamela Trisolino Gabriela Vítková and Camille Yver Kwok

^{1,},

^{2},

^{3},

^{3},

^{3},

^{4},

^{5,6},

^{7},

^{8},

^{9},

^{4},

^{10},

^{11},

^{12},

^{13,14},

^{15},

^{12},

^{4},

^{1},

^{12},

^{16},

^{8},

^{3},

^{17},

^{18},

^{4}

**Peter Bergamaschi et al.**Peter Bergamaschi Arjo Segers Dominik Brunner Jean-Matthieu Haussaire Stephan Henne Michel Ramonet Tim Arnold Tobias Biermann Huilin Chen Sebastien Conil Marc Delmotte Grant Forster Arnoud Frumau Dagmar Kubistin Xin Lan Markus Leuenberger Matthias Lindauer Morgan Lopez Giovanni Manca Jennifer Müller-Williams Simon O'Doherty Bert Scheeren Martin Steinbacher Pamela Trisolino Gabriela Vítková and Camille Yver Kwok

^{1,},

^{2},

^{3},

^{3},

^{3},

^{4},

^{5,6},

^{7},

^{8},

^{9},

^{4},

^{10},

^{11},

^{12},

^{13,14},

^{15},

^{12},

^{4},

^{1},

^{12},

^{16},

^{8},

^{3},

^{17},

^{18},

^{4}

^{1}European Commission Joint Research Centre (JRC), Ispra (Va), Italy^{2}Netherlands Organisation for Applied Scientific Research (TNO), Utrecht, the Netherlands^{3}Swiss Federal Laboratories for Materials Science and Technology (Empa), Dübendorf, Switzerland^{4}Laboratoire des Sciences du Climat et de l'Environnement (LSCE-IPSL), CEA-CNRS-UVSQ, Université Paris-Saclay, 91191 Gif-sur-Yvette, France^{5}National Physical Laboratory, Teddington, UK^{6}School of GeoSciences, University of Edinburgh, Edinburgh, UK^{7}Centre for Environmental and Climate Science (CEC), Lund University, Lund, Sweden^{8}Centre for Isotope Research (CIO), Energy and Sustainability Research Institute Groningen (ESRIG), University of Groningen, Groningen, the Netherlands^{9}Agence nationale pour la gestion des dechets radioactifs (Andra), DRD/OPE, Bure, France^{10}Centre for Ocean and Atmospheric Sciences, School of Environmental Sciences, University of East Anglia, Norwich, UK^{11}Netherlands Organisation for Applied Scientific Research (TNO), Petten, the Netherlands^{12}Deutscher Wetterdienst, Hohenpeissenberg Meteorological Observatory (MOHp), Hohenpeissenberg, Germany^{13}Cooperative Institute for Research in Environmental Sciences, University of Colorado Boulder, Boulder, CO, USA^{14}NOAA Earth System Research Laboratory, Global Monitoring Laboratory, Boulder, CO, USA^{15}University of Bern, Physics Institute, Climate and Environmental Division, and Oeschger Centre for Climate Change Research, Bern, Switzerland^{16}Atmospheric Chemistry Research Group, University of Bristol, Bristol, UK^{17}National Research Council of Italy, Institute of Atmospheric Sciences and Climate (CNR-ISAC), Bologna, Italy^{18}Global Change Research Institute of the Czech Academy of Sciences, Brno, Czech Republic^{}retired

^{1}European Commission Joint Research Centre (JRC), Ispra (Va), Italy^{2}Netherlands Organisation for Applied Scientific Research (TNO), Utrecht, the Netherlands^{3}Swiss Federal Laboratories for Materials Science and Technology (Empa), Dübendorf, Switzerland^{4}Laboratoire des Sciences du Climat et de l'Environnement (LSCE-IPSL), CEA-CNRS-UVSQ, Université Paris-Saclay, 91191 Gif-sur-Yvette, France^{5}National Physical Laboratory, Teddington, UK^{6}School of GeoSciences, University of Edinburgh, Edinburgh, UK^{7}Centre for Environmental and Climate Science (CEC), Lund University, Lund, Sweden^{8}Centre for Isotope Research (CIO), Energy and Sustainability Research Institute Groningen (ESRIG), University of Groningen, Groningen, the Netherlands^{9}Agence nationale pour la gestion des dechets radioactifs (Andra), DRD/OPE, Bure, France^{10}Centre for Ocean and Atmospheric Sciences, School of Environmental Sciences, University of East Anglia, Norwich, UK^{11}Netherlands Organisation for Applied Scientific Research (TNO), Petten, the Netherlands^{12}Deutscher Wetterdienst, Hohenpeissenberg Meteorological Observatory (MOHp), Hohenpeissenberg, Germany^{13}Cooperative Institute for Research in Environmental Sciences, University of Colorado Boulder, Boulder, CO, USA^{14}NOAA Earth System Research Laboratory, Global Monitoring Laboratory, Boulder, CO, USA^{15}University of Bern, Physics Institute, Climate and Environmental Division, and Oeschger Centre for Climate Change Research, Bern, Switzerland^{16}Atmospheric Chemistry Research Group, University of Bristol, Bristol, UK^{17}National Research Council of Italy, Institute of Atmospheric Sciences and Climate (CNR-ISAC), Bologna, Italy^{18}Global Change Research Institute of the Czech Academy of Sciences, Brno, Czech Republic^{}retired

**Correspondence**: Peter Bergamaschi (peter.bergamaschi@ext.ec.europa.eu) and Dominik Brunner
(dominik.brunner@empa.ch)

**Correspondence**: Peter Bergamaschi (peter.bergamaschi@ext.ec.europa.eu) and Dominik Brunner
(dominik.brunner@empa.ch)

Received: 18 Feb 2022 – Discussion started: 16 Mar 2022 – Revised: 27 Aug 2022 – Accepted: 04 Sep 2022 – Published: 17 Oct 2022

We present a novel high-resolution inverse modelling
system (“FLEXVAR”) based on FLEXPART-COSMO back trajectories driven by COSMO
meteorological fields at 7 km×7 km resolution over the European
COSMO-7 domain and the four-dimensional variational (4DVAR) data
assimilation technique. FLEXVAR is coupled offline with the global inverse
modelling system TM5-4DVAR to provide background mole fractions
(“baselines”) consistent with the global observations assimilated in
TM5-4DVAR. We have applied the FLEXVAR system for the inverse modelling of
European CH_{4} emissions in 2018 using 24 stations with in situ
measurements, complemented with data from five stations with discrete air
sampling (and additional stations outside the European COSMO-7 domain used
for the global TM5-4DVAR inversions). The sensitivity of the FLEXVAR
inversions to different approaches to calculate the baselines, different
parameterizations of the model representation error, different settings of
the prior error covariance parameters, different prior inventories, and
different observation data sets are investigated in detail. Furthermore, the
FLEXVAR inversions are compared to inversions with the FLEXPART extended
Kalman filter (“FLExKF”) system and with TM5-4DVAR inversions at $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$ resolution over Europe. The three inverse modelling
systems show overall good consistency of the major spatial patterns of the
derived inversion increments and in general only relatively small
differences in the derived annual total emissions of larger country regions.
At the same time, the FLEXVAR inversions at 7 km×7 km resolution
allow the observations to be better reproduced than the TM5-4DVAR simulations
at $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$. The three inverse models derive
higher annual total CH_{4} emissions in 2018 for Germany, France, and
BENELUX compared to the sum of anthropogenic emissions reported to UNFCCC
and natural emissions estimated from the Global Carbon Project CH_{4}
inventory, but the uncertainty ranges of top-down and bottom-up total
emission estimates overlap for all three country regions. In contrast, the
top-down estimates for the sum of emissions from the UK and
Ireland agree relatively well with the total of anthropogenic and natural
bottom-up inventories.

- Article
(9768 KB) -
Supplement
(12647 KB) - BibTeX
- EndNote

_{4}emissions using the novel FLEXPART-COSMO TM5 4DVAR inverse modelling system, Atmos. Chem. Phys., 22, 13243–13268, https://doi.org/10.5194/acp-22-13243-2022, 2022.

Atmospheric methane (CH_{4}) is the second most important anthropogenic
greenhouse gas (GHG) after carbon dioxide (CO_{2}), with an estimated
contribution of ∼16.3 *%* (0.520 W m^{−2}) to the direct
anthropogenic radiative forcing of all long-lived GHGs in 2020 (NOAA Annual
Greenhouse Gas Index (AGGI), evaluated relative to 1750; Butler
and Montzka, 2022). Including also additional indirect effects (e.g.
production of tropospheric ozone), however, the total radiative forcing of
CH_{4} is considerably higher, with current estimates of the
emission-based effective radiative forcing (ERF) of 1.21 (0.90 to 1.51)
W m^{−2} (Szopa et al.,
2021). The current global average CH_{4} mole fraction is 162 % higher
than preindustrial levels in 1750 (WMO, 2021) and continues to
increase, with recent growth rates (2014–2020: 10.1±3.2 ppb yr^{−1}) being again close to the high growth rates observed during the
1980s (1984–1989: 11.9±0.9 ppb yr^{−1}), while lower growth rates
were observed during the 1990s and almost zero growth rates during 2000–2006
(Dlugokencky, 2022).

Reducing CH_{4} emissions plays an essential role in mitigating climate
change, especially in the near term (Shindell et al., 2012, 2017; United Nations Environment Programme and Climate and Clean Air
Coalition, 2021), due to CH_{4}'s relatively short atmospheric lifetime of
around 10 years combined with its high radiative efficiency (resulting in a
global warming potential (GWP) around 80 times higher compared to CO_{2}
on a 20-year timescale; Forster
et al., 2021). The global emissions pathways to limit global warming to
1.5 ^{∘}C, compiled by IPCC (2018),
include significant reductions of CH_{4} emissions after 2020 (for
scenarios with no or limited overshoot of temperature above the
1.5 ^{∘}C target). The recognition of the importance of CH_{4}
emission reductions to mitigate climate change has also led to the recent
“Global Methane Pledge” (European Commission, 2021), with the
collective goal to reduce methane emissions by 2030 by at least 30 %
compared to 2020.

The development of emission reduction pathways as well as the control of
international climate agreements requires the accurate quantification of
current (and past) GHG emissions. For CH_{4}, however, the quantification
of emissions and sinks is particularly challenging, mainly owing to the
large spatial and temporal variability of emissions from many source sectors
and consequently large uncertainties in assumed mean emission factors (e.g.
for natural emissions from wetlands and anthropogenic emissions from
fugitive sources such as fossil fuels (coal, oil, gas) (e.g. Brandt et al.,
2014) or emissions from the waste sector).
Therefore, bottom-up inventories of CH_{4}, which are compiled by scaling
up emissions using activity data and emission factors, have significant
uncertainties. Complementary to bottom-up inventories, inverse modelling
provides top-down emission estimates using atmospheric measurements and
atmospheric transport models, by optimizing emissions from emission
inventories (used as prior estimates) to get an optimal agreement between
simulated and observed CH_{4} mole fractions, taking into account the
uncertainties of prior emission estimates, measurements, and model
simulations (e.g. Bergamaschi et al., 2018b; Houweling et al., 2017).
The Global Carbon Project CH_{4} (GCP-CH_{4}) provides synthesis
analyses of the global CH_{4} cycle based on comprehensive sets of
different bottom-up and top-down estimates from the international science
community working on this topic (Jackson et al., 2020; Saunois et al.,
2020; Stavert et al., 2021). The global inverse models used in these
analyses have generally relatively coarse horizontal resolution, in the
range of 2.5–6.0^{∘} (longitude) ×1.9–4.0^{∘} (latitude)
(Saunois et al., 2020). Therefore, such models are mainly suitable to
provide information on the global and larger regional scales.

In order to analyse regional emissions in more detail (and more accurately), specific regional inversions have been performed, employing regional atmospheric transport models at higher horizontal resolution (typically in the range of ∼20–100 km) and making use of the increasing number of regional in situ GHG measurements, which have become available in recent years in particular in Europe and North America (e.g. Bergamaschi et al., 2018a; Ganesan et al., 2015; Lunt et al., 2021; Manning et al., 2011; Miller et al., 2013). The regional models generally require global boundary conditions, which are usually provided from global inverse models. Alternatively, global models with zooming option for the specific region of interest have been employed (Bergamaschi et al., 2018a). A specific purpose of such regional inversions is to verify national bottom-up emission inventories reported to the United Nations Framework Convention on Climate Change (UNFCCC), finally aiming at an emission monitoring and verification system to support the international climate agreements using in situ and satellite observations (Bergamaschi et al., 2018b; Deng et al., 2022; National Research Council, 2010; Pinty et al., 2017, 2019). Within the European project VERIFY (https://verify.lsce.ipsl.fr/, last access: 4 October 2022) a pre-operational GHG verification system has been developed, employing various state-of-the-art global and regional atmospheric transport models (Petrescu et al., 2021a, b).

In order to further improve the atmospheric modelling, it is essential to
further increase the spatial resolution, aiming at further improving the
simulation of regional monitoring stations. A pioneering high-resolution
study has been reported by Henne et al. (2016), using the
FLEXPART-COSMO back trajectories driven by meteorological fields from the
Swiss national weather service (MeteoSwiss) at horizontal resolution of
approximately 7 km×7 km, analysing the CH_{4} emissions from
Switzerland using continuous measurements from six atmospheric monitoring
stations. The authors generated the FLEXPART-COSMO sensitivity fields (based
on the sampling of released particles) at even higher resolution (than the
resolution of the COSMO meteorological fields) of ${\mathrm{0.02}}^{{}^{\circ}}\times \mathrm{0.015}{}^{\circ}$ (≈1.7 km) over the Alpine domain (but coarser
horizontal resolution of $\mathrm{0.16}{}^{\circ}\times \mathrm{0.12}{}^{\circ}$ (≈13 km)
outside the Alpine domain). Since they solved the inverse problem
analytically (denoted by the authors as the “Bayesian method”), however, they
applied a reduced grid by merging model grid cells in areas with smaller
average source sensitivities in order to reduce the size of the inversion
problem to about 1000 unknowns. As an alternative to the Bayesian method,
Henne et al. (2016) also applied the extended Kalman filter method
described by Brunner et al. (2012), which – in
contrast to the Bayesian method – assimilates the observations
sequentially but for computational reasons also requires the application of
a reduced grid.

Aiming at high-resolution inversions of larger regions (such as the European domain), we have therefore developed a novel inversion framework (denoted “FLEXVAR”) based on the four-dimensional variational (4DVAR) data assimilation technique (Meirink et al., 2008; Talagrand and Courtier, 1987), which allows for the optimization of a much larger number of parameters and therefore also avoids the need to apply reduced grids. As in the Henne et al. (2016) study, the new system uses FLEXPART-COSMO back trajectories driven by COSMO meteorological fields at 7 km×7 km resolution but is optimizing emissions from individual grid cells over the whole European COSMO-7 domain with 393 (longitude) ×338 (latitude) grid cells. Furthermore, the new system uses background mole fractions (baselines) from global TM5-4DVAR inversions (using two different approaches); i.e. it is coupling the FLEXPART-COSMO inversions with the global/European TM5-4DVAR inversions offline.

The objective of this paper is to present the new FLEXVAR system and its
application to the inversion of European CH_{4} emissions for 2018 using a
comprehensive data set from 24 stations with in situ measurements,
complemented with data from 5 stations with discrete air sampling (and
additional stations outside the European COSMO-7 domain used for the global
TM5-4DVAR inversions). We analyse in detail the sensitivity of the FLEXVAR
inversions to internal parameterizations and model settings, as well as the
sensitivity to the main model input data, i.e. prior inventories and
observational data. Furthermore, we compare the FLEXVAR inversions with the
extended Kalman filter (FLExKF) method and with TM5-4DVAR inversions (at
$\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$ resolution over Europe). Finally, we present an
overall analysis of derived European CH_{4} emissions and comparison with
emissions reported to UNFCCC for some major countries (or group of
countries) which are best constrained by the available observations
(Germany, France, BENELUX, and UK and Ireland).

## 2.1 FLEXPART-COSMO back trajectories

FLEXPART is a Lagrangian atmospheric transport model that simulates the advective, turbulent, and convective transport by tracking the positions of a large number of infinitesimally small air parcels, so-called particles, either forward or backward in time (Pisso et al., 2019; Stohl et al., 2005). FLEXPART is an offline model that requires meteorological fields such as 3D wind fields from a numerical weather prediction (NWP) model as input. FLEXPART-COSMO is a version of the model that is driven by the output of the NWP model COSMO, which was jointly developed by a consortium of European weather services under the lead of the German meteorological service DWD (Baldauf et al., 2011). Different from all other FLEXPART versions, FLEXPART-COSMO operates on the native vertical grid of the driving model COSMO, which avoids potential loss of information and inaccuracies associated with the interpolation onto a different grid. More details on the model are provided in Henne et al. (2016) and Pisso et al. (2019).

In the backward mode, particles are released at the locations of individual observations and followed backwards in time over typically a few days. By sampling the near-surface residence times of the particles along their paths, a so-called source–receptor sensitivity matrix or “footprint” is computed, which describes the relationship between the change in mole fraction at the observation site and the fluxes discretized in space and time (Seibert and Frank, 2004). A time series of simulated mole fractions can be obtained by integrating the time series of source–receptor matrices with a discretized flux estimate. The simulations were driven by hourly output from the operational COSMO-7 analyses of the Swiss weather service MeteoSwiss at a horizontal spatial resolution of about 7 km×7 km, largely covering western and central Europe (Fig. 1). In the simulations used here, 50 000 virtual particles were released at all observation locations every 3 h (evenly distributed over each 3 h time interval) and traced backwards in time for 10 d (or until individual particles left the COSMO-7 domain). Despite the high spatial resolution, the orography of the COSMO-7 model is smoothed for complex terrain, leading to differences between the model surface altitude and the altitude of the observation site. When this difference is greater than 200 m, the release height of the particles has been chosen as the average between the measurement height above sea level and the model surface altitude. This approach has been used in previous studies (Henne et al., 2016) and was found to be the most representative release height (Brunner et al., 2013).

## 2.2 Coupled FLEXPART-COSMO TM5 4DVAR inverse modelling system FLEXVAR

### 2.2.1 Inversion framework

The new coupled FLEXPART-COSMO TM5 4DVAR inverse modelling system, denoted
FLEXVAR, allows for the optimization of emissions at grid scale using the
four-dimensional variational (4DVAR) data assimilation technique
(Meirink et al., 2008; Talagrand and Courtier,
1987), the FLEXPART-COSMO back trajectories described in Sect. 2.1, and
background mole fractions (“baselines”) from TM5-4DVAR (described in Sect. 2.2.2 and 2.4). The system follows the classical Bayesian approach
minimizing the cost function *J*(** x**):

where ** x** is the state vector,

*x*_{b}the prior estimate of the state vector (in data assimilation usually called the “background”),

**the set of observations (measurements) to be assimilated,**

*y*

*H***(**

*x***)**the observation operator (or model operator), representing the model simulation of the observations, and

**B**and

**R**the error covariance matrices of the prior estimate and the observations, respectively. For the regular inversions, a semi-lognormal probability density function (pdf) is applied for the emissions

**to be optimized (Bergamaschi et al., 2010), optimizing the emission deviation factors ${x}_{i,j,t}$,**

*e*
for each element ${e}_{i,j,t}$, representing the emissions of an individual
grid cell with longitude index *i*, latitude index *j*, and at emission time step
*t* (and ${e}_{\mathrm{b},\phantom{\rule{0.125em}{0ex}}i,j,t}$ the prior estimate of the emission
${e}_{i,j,t}$). In contrast, a linear expansion of the emission deviation
factors is used (resulting in a Gaussian pdf of the emissions) for
additional inversions for the evaluation of posterior uncertainties (as will
be described in more detail below):

For the prior estimate, the emission deviation factors ${x}_{i,j,t}$ are set to zero, i.e. ${e}_{i,j,t}\left({x}_{i,j,t}\right)={e}_{\mathrm{b},\phantom{\rule{0.125em}{0ex}}i,j,t}$ (both for the semi-lognormal and Gaussian pdf). In this study we present inversions optimizing the monthly total emissions of all FLEXPART-COSMO grid cells (at 7 km×7 km resolution) for 1 year (2018); hence the dimension of the state vector is

The background covariance matrix **B** is
parameterized as the product

where **S** is a diagonal matrix with the diagonal elements
containing the uncertainty of emissions (standard deviations), and
**C** is the correlation matrix, which is
parameterized as a Kronecker product of a horizontal correlation matrix
**C**_{hor} and a temporal correlation matrix
**C**_{t} (as in TM5-4DVAR, Meirink
et al., 2008):

The spatial covariance between two grid cells is parameterized using a Gaussian function:

where $d({i}_{\mathrm{1}},\phantom{\rule{0.125em}{0ex}}{i}_{\mathrm{2}},\phantom{\rule{0.125em}{0ex}}{j}_{\mathrm{1}},\phantom{\rule{0.125em}{0ex}}{j}_{\mathrm{2}})$ is the distance between two grid
cells (with longitude indices *i*_{1} and *i*_{2} and latitude indices
*j*_{1} and *j*_{2}) and *L*_{corr} a predefined correlation length
constant. The temporal correlation uses an exponential decay function:

where *d*_{t}(*t*_{1},*t*_{2}) is the temporal distance between two
emission time steps and *t*_{corr} a predefined temporal
correlation scale constant.

The observation error covariance matrix **R** considers only diagonal elements (i.e. ignores any error correlation
between different observations) and takes into account the uncertainties of
the measurements and the model representation error:

We use two different approaches to parameterize the model representation error, which are described in more detail in Sect. 2.2.3.

The observation operator *H***(***x***)** simulates the measurements
as a function of the state vector (i.e. as a function of emission deviation
factors). For a given measurement *m* in time interval *t*_{m} the simulation
is computed as

using the 3-dimensional FLEXPART-COSMO footprints ${\mathit{\varphi}}_{i,j,k}^{m}\phantom{\rule{0.125em}{0ex}}(t,t+\mathit{\delta}t$) (units: [1/((kg air) m^{−3} s^{−1})]), described in
Sect. 2.1. ${e}_{i,j,t}({x}_{i,j,t}$) denotes the CH_{4} emissions (computed according Eq. 2
or Eq. 3) in units of kilograms of methane per square metre per second (kg CH_{4} m^{−2} s^{−1}).
*M*_{air} and ${M}_{{\mathrm{CH}}_{\mathrm{4}}}$ are the molecular
masses of air and CH_{4}, respectively. *δ**h*_{k} is the layer
thickness of vertical layer *k* (here we use the two lowermost layers, each
with thickness of 50 m), and *w*_{k} is the weighting of layer *k* (here 0.5 for
the applied two layers). *T*_{m} represents the time interval of the applied
footprints (i.e. 10 d prior to the measurements) and *δ**t* the
averaging time (3 h) for the single footprints (computed for the time
interval between *t* and *t*+*δ**t*). *H*_{m}(** x**,

*t*

_{m}) represents the simulated enhancement of the CH

_{4}mole fraction above the baseline (which is evaluated using two different approaches as described in Sect. 2.2.2).

The minimization of the cost function Eq. (1) requires the evaluation of the gradient of the cost function with respect to the state vector:

where **H**^{T} is the adjoint model operator,
describing the sensitivity of the simulated observations with respect to
changes of the state vector. **H**^{T} can be
directly computed using the FLEXPART-COSMO footprints ${\mathit{\varphi}}_{i,j,k}^{m}\phantom{\rule{0.125em}{0ex}}(t,\phantom{\rule{0.125em}{0ex}}t+\mathit{\delta}t$).

In order to achieve better convergence of the minimization algorithm,
pre-conditioning is applied, transforming the state vector ** x** to
the control vector

**(similarly as, for example, in TM5-4DVAR, Meirink et al., 2008),**

*w*and in reverse direction,

The square root of the submatrix **C**_{t} is
calculated by eigenvalue decomposition. However, this is not possible for
the submatrix **C**_{hor} (in contrast to
TM5-4DVAR) due to the large size of this matrix of

Therefore, an Arnoldi factorization (Arnoldi, 1951) is used, computing
only the largest eigenvalues/eigenvectors. Consequently, the square root
matrix ${\mathbf{B}}^{\mathrm{1}/\mathrm{2}}$ and its inverse in Eqs. (13) and (12)
are replaced by corresponding factorizations. As default setting, we use a
fraction of 1 % of the eigenpairs, which is reducing the size of the
control vector to 1 % of the size of the state vector. Test inversions
with higher fractions of eigenpairs (up to 3 %) showed that a fraction of
1 % is generally sufficient for the range of correlation lengths
(*L*_{corr}=50 km 200 km) used in this study.

For the regular inversions, we use the limited memory quasi-Newton algorithm m1qn3 developed by Gilbert and Lemaréchal (1989), which also allows for the optimization of non-linear problems (as the semi-lognormal pdf (Eq. 2) introduces non-linearity to our optimization problem). In order to evaluate the posterior uncertainties, additional inversions are performed using a conjugate gradient algorithm (Fisher and Courtier, 1995; Lanczos, 1950; Meirink et al., 2008) for minimization and the linear expansion of the emission deviation factors (Eq. 3). The posterior covariance is then computed from the leading eigenvalues of the Hessian of the cost function:

with *υ*_{k} and *θ*_{k} the eigenvectors and eigenvalues of
the Hessian of the cost function. The inversions based on the conjugate
gradient algorithm yield in general rather similar posterior emissions as
the regular inversions using the m1qn3 algorithm. Derived aggregated annual
total CH_{4} emissions agree on average (average over all sensitivity
inversions) within −1.0 *%* to 2.1 % for the “country regions” discussed
in this paper (Germany, France, BENELUX, UK + Ireland; Fig. 2) and within
−6.4 *%* to 6.5 % for individual inversions.

### 2.2.2 Baselines

Since FLEXPART-COSMO is a limited domain model, providing the
source–receptor relationship due to emissions in the European COSMO-7
domain, it requires the provision of background mole fractions
(baselines), representing the CH_{4} mole fractions of the air masses
entering the COSMO-7 domain. Here, we simulate the baselines using the
global TM5-4DVAR inverse modelling system (described in Sect. 2.4), using
two different approaches to couple the regional FLEXPART-COSMO inversion
with TM5-4DVAR. The first approach applies the method described by
Rödenbeck et al. (2009) (denoted in the following as
“Rödenbeck baselines”), which is a flexible nesting scheme allowing for the
offline coupling of regional models with global inversions. For this
purpose, we use a manipulated version of TM5 which simulates only the “cis”
part of the CH_{4} fields (representing the CH_{4} enhancement due to
emissions in the area to the FLEXPART-COSMO model transported directly to
the European measurement stations, without leaving the COSMO-7 domain),
denoted Δ*c*_{mod 1, cis}. The baseline mole fractions,
*c*_{baseline}, are then computed as difference of the TM5-4DVAR
posterior simulation (*c*_{mod 1}) and the “cis” part
(Rödenbeck et al., 2009), both sampled at the location of
the corresponding station:

The second approach to couple FLEXPART-COSMO with TM5-4DVAR uses the
particle positions of the FLEXPART-COSMO back trajectories at the time of
termination, i.e. either at the applied maximum termination time (set to 10 d in our study) or when they leave the COSMO-7 domain (which can be well
before the 10 d). For each individual particle position, the CH_{4}
mole fraction is then extracted from the 3-dimensional TM5-4DVAR CH_{4}
fields. Since each 3 h average FLEXPART-COSMO footprint is based on the
release of 50 000 particles, the corresponding (3 h average) baseline
mole fraction is computed as the average of the CH_{4} mole fraction at the
termination points of all individual 50 000 particles. In the following, this
approach is denoted “particle position baselines”.

### 2.2.3 Model representation error

We applied two different approaches to estimate the model representation
error. The first approach, denoted “OBS”, is similar to the method described
by Henne et al. (2016), evaluating the residuals (difference between
observations and model simulations) as a function of CH_{4} enhancement.
However, some details of our method are different from Henne et al. (2016). Here, we use the following function to parameterize the model
representation error, Δ*y*_{MOD k,i}, as a function of
the absolute observed CH_{4} enhancement (i.e. observed CH_{4} mole
fraction minus CH_{4} background), |*y*_{k, i}|:

where *k* is the station index, *i* the index of the individual observational data
point of the time series of station *k*, and *f*_{0,k}, *ρ*_{k}, and
*a*_{k} the three fit parameters evaluated for each station. These fit
parameters are determined by calculating the fit curve of the absolute
residuals, $\left|{y}_{k,\phantom{\rule{0.125em}{0ex}}\phantom{\rule{0.125em}{0ex}}i}-{\left(\mathit{H}\mathbf{\left(}\mathit{x}\mathbf{\right)}\right)}_{k,i}\right|$ as a function of |*y*_{k, i}|. For
large *y*_{k, i}, the curve becomes linear with slope *ρ*_{k}; the
minimum value of the fit function is *f*_{0,k}, and *a*_{k} defines how fast
the function becomes linear. The evaluation of the fit function is generally
performed iteratively. In a first step, the fit parameters are computed
using the prior model simulations, while in subsequent iterations the
posterior simulations of the previous iteration are used. FLEXVAR includes
an outer loop system, which allows for an arbitrary number of iterations.
However, tests have shown that changes after the second iteration are
usually very small. Therefore, for the inversions presented in this paper,
we used generally only two iterations. Figure S1 in the Supplement (left column) illustrates
the fit functions for some selected stations.

As alternative to the model representation error OBS described above, a second approach, denoted “METEO”, has been developed, parameterizing the model representation error as a function of wind speed:

where *w*_{k, i} is the wind speed (in m s^{−1}) extracted from the COSMO
model corresponding to the data point *y*_{k, i}, and *f*_{0,k}, *ρ*_{k}, and *a*_{k} are the three fit parameters for station *k*. The rationale is
that if wind speed at the measurement location is low, the observed mole
fraction might be more influenced by local sources of emissions and
therefore less representative of the modelled mole fractions. Again, the fit
parameters are determined to get an optimal fit through the absolute
residuals, $\left|{y}_{k,\phantom{\rule{0.125em}{0ex}}\phantom{\rule{0.125em}{0ex}}i}-{\left(\mathit{H}\mathbf{\left(}\mathit{x}\mathbf{\right)}\right)}_{k,i}\right|$ and they are evaluated in two iterations. For large wind
speeds, Δ*y*_{MOD k,i} converges towards *f*_{0,k}, while
${f}_{\mathrm{0},k}+{\mathit{\rho}}_{k}$ represents the model representation error at wind speed
zero. The right column of Fig. S1 shows the fit functions for the METEO
model representation error for some selected stations.

## 2.3 FLEXPART extended Kalman filter inverse modelling system

FLExKF is an inverse modelling system based on an extended Kalman filter as described in detail in Brunner et al. (2012, 2017). Observations are assimilated sequentially (here day by day) to provide the best linear unbiased estimate of the emissions and their uncertainties based on measurements up to the present time of the assimilation. Rather than estimating monthly or annual mean emissions, the system adjusts the emissions continuously as the assimilation proceeds and in this way creates a smoothly varying emission field that later can be averaged to monthly or annual means. The filter includes a forecast step, which predicts the evolution of the emissions from one assimilation time step to the next. The simplest assumption is persistence (i.e. no change with time), but to incorporate seasonally varying prior emissions, a non-zero forecast update was implemented that follows the linear change in prior emissions from 1 month to the next. Since the forecast step is associated with an uncertainty, the posterior uncertainty can become larger than the prior uncertainty in regions that are poorly constrained by observations. In order to avoid an unrealistic growth of the uncertainties in these regions, the posterior uncertainties are reset to the prior uncertainty whenever they become larger.

The state vector consists of two components: (i) emissions on a reduced grid
with a total of 3608 elements for inversions with observation data set O1
and 6497 elements for inversions with observation data set O2 (Sect. 3.1 and
Table 1) and (ii) coefficients of an AR(1) autoregressive process describing
temporal correlations in the residuals at each individual site. As described
in Brunner et al. (2012), the reduced grid has high
spatial resolution near the measurement sites and lower resolution further
away, reflecting the reduced contribution of emissions at larger distance to
observed CH_{4} variations. It was constructed based on the combined total
annual footprint of all measurement sites. The state vector may contain the
emissions directly or the logarithm of the emissions. The latter option was
chosen here to enforce a positive solution, i.e. positive methane emissions
in each grid cell.

Another option in FLExKF is to optimize baseline mole fractions at each
observation site in addition to the gridded emission field. However, in the
configuration used here, the baselines (Rödenbeck baselines) described
in Sect. 2.2.2 were used directly without optimization. The results of
FLExKF should be readily comparable to those of FLEXVAR, since the same
FLEXPART-COSMO back trajectories, baselines, and observations were used.
Spatial correlations in the prior emission uncertainties were represented in
the prior error covariance matrix with a correlation length scale of 200 km
(exponential decay as in Eq. 8). The matrix was scaled such that the prior
uncertainty of the total domain emissions was 20 % (1*σ*).

## 2.4 TM5-4DVAR inverse modelling system

TM5-4DVAR is a global inverse modelling system based on the 4DVAR data
assimilation technique and has been described in detail by Meirink et al. (2008), while subsequent updates have been reported in
Bergamaschi et al. (2018a, 2010). TM5-4DVAR uses the Eulerian atmospheric
chemistry transport model TM5 (Krol et al., 2005), a
two-way nested zoom model, which allows the system to zoom in over specific
regions of interest. Here, we apply the $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$
zooming over the European domain −18^{∘} to 42^{∘} (longitude)
$\times \mathrm{32}{}^{\circ}$ to 64^{∘} (latitude) and an intermediate
3^{∘} (longitude) $\times \mathrm{2}{}^{\circ}$ (latitude) zooming over the
extended European domain −30^{∘} to 48^{∘} (longitude)
$\times \mathrm{26}{}^{\circ}$ to 70^{∘} (latitude), while the remaining
global domain is simulated at a horizontal resolution of 6^{∘}
(longitude) $\times \mathrm{4}{}^{\circ}$ (latitude). TM5 is an offline transport
model, driven by external meteorological data, preprocessed to provide
consistent meteorological fields for the different TM5 resolutions
(Krol et al., 2005). In this study, 3 h
interpolated meteorological fields from the European Centre for Medium-Range
Weather Forecasts (ECMWF) ERA-Interim reanalysis
(Dee et
al., 2011) have been applied, using 25 vertical layers (defined as a subset
of the 60 layers of the ERA-Interim reanalysis). About 5 layers represent
the boundary layer, 10 layers the free troposphere, and 10 layers the
stratosphere (Krol et al., 2005).

As for the regular FLEXVAR inversions, a semi-lognormal pdf has been used
(Eq. 2). Minimization of the cost function Eq. (1) is performed using the
m1qn3 algorithm (Gilbert and Lemaréchal, 1989) and the adjoint of
the tangent linear TM5 model (Krol et al., 2008; Meirink et al., 2008)
for evaluation of the gradient of the cost function Eq. (11). Four groups of
CH_{4} emissions are optimized independently: (1) wetlands, (2) rice, (3) biomass burning, and (4) all remaining sources (Bergamaschi et al.,
2018a). Uncertainties of 100 % per grid cell and month were applied for
each source group with a spatial correlation length scale of 200 km (Eq. 7). The temporal correlation timescales (Eq. 8) are set to 12 months for
the “remaining” CH_{4} sources (which are assumed to have no or only small
seasonal variations) and to zero for the other source groups (wetlands,
rice, and biomass burning), which have pronounced seasonal cycles. The model
representation error is parameterized as a function of local emissions (i.e.
emissions of the grid cell in which the corresponding monitoring station is
located) and 3-dimensional gradients of simulated mole fractions
(Bergamaschi et al., 2010). The photochemical
sinks of CH_{4} in the troposphere (OH) and stratosphere (OH, Cl, and
O(^{1}D)) are simulated as described in Bergamaschi et al. (2010). For the coupling of the FLEXPART-COSMO
inversions over the European COSMO domain with the global TM5-4DVAR
inversions, the baselines at the monitoring stations (within the COSMO-7
domain) have been computed using the two different approaches described in
Sect. 2.2.2.

## 3.1 Atmospheric observations

The atmospheric observations used in this study (within the COSMO-7 domain)
include ground-based CH_{4} data for 2018 from 24 stations with in situ
measurements, complemented with data from 5 stations with discrete air
sampling, as compiled in Table 1. Most of the in situ measurements (15
stations) are from the atmosphere network of the Integrated Carbon
Observation System (ICOS) (Heiskanen et al., 2022), a pan-European
infrastructure providing harmonized atmospheric measurements which are
rigorously standardized in terms of instrumentation, calibration, air
sampling, and quality control, including centralized data processing and data
evaluation at the ICOS Atmospheric Thematic Centre
(https://icos-atc.lsce.ipsl.fr/, last access: 5 October 2022) (Hazan et al., 2016; ICOS RI, 2020;
Yver-Kwok et al., 2021). Here, we use the ICOS Atmosphere Release of final,
quality-controlled data 2021-1 (ICOS RI, 2021). For station Lutjewad
(LUT), the ICOS data only start on 13 August 2018. Data before that date were
provided by the University of Groningen, with data processing very similar to
the ICOS data. Data from the ICOS station Ispra (IPR) have been further
processed using the robust extraction of baseline signal (REBS) spike
detection algorithm (El Yazidi et al., 2018; Ruckstuhl et al., 2012) in
order to filter out data affected by nearby farming activities. In addition
to the ICOS measurements, further in situ measurements have been used from
the UK Deriving Emissions linked to Climate Change (DECC) network (Bilsdale
(BSD), Tacolneston (TAC), Ridge Hill (RGL), Heathfield (HFD)), from the Advanced
Global Atmospheric Gases Experiment (AGAGE) (Mace Head (MHD)), from
the University of East Anglia (Weybourne (WAO)), from the Netherlands Organisation
for Applied Scientific Research (TNO) (Cabauw (CBW)), from Empa (Lägern
Hochwacht (LHW)), and from the University of Bern (Beromünster (BRM)). The
in situ measurements are complemented by discrete air samples (which are
usually collected weekly) from the NOAA Earth System Research Laboratory
(ESRL) global cooperative air sampling network (Dlugokencky et al., 2021), with five
stations within the COSMO-7 domain. Additional atmospheric data used for the
TM5-4DVAR inversions are compiled in Table S1 and include six further ICOS
stations with in situ measurements and 31 NOAA discrete air sampling sites
located outside the COSMO-7 domain. The selected stations outside the
European $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$ or $\mathrm{3}{}^{\circ}\times \mathrm{2}{}^{\circ}$ TM5 zoomed-in regions are mostly global background stations in
remote areas which can be reasonably well reproduced with the coarse global
TM5 resolution of $\mathrm{6}{}^{\circ}\times \mathrm{4}{}^{\circ}$.

The atmospheric CH_{4} data are reported on the WMO X2004A calibration
scale (Dlugokencky et al., 2005; NOAA,
2021), except the AGAGE MHD data which are reported on the Tohoku University
(TU) CH_{4} standard scale (Aoki et al., 1992; Prinn et al., 2000).
Comparison of parallel measurements by NOAA and AGAGE at five global sites over
more than 25 years showed that the two calibration scales are in close
agreement, with an average ratio of $\mathrm{1.0002}\pm \mathrm{0.0007}.$ Therefore, no
scale correction has been applied.

^{1} Data since 13 August 2018 from ICOS data release and before that date from
the University of Groningen.
^{2} Data filtered with REBS spike detection algorithm (see Sect. 3.1).
^{3} Observatoire Perenne de l'Environnement.
^{4} Centro de Investigación de la Baja Atmósfera.

For the in situ measurements (which are available quasi-continuously in
time), we assimilate only early afternoon data for stations in the boundary
layer and night-time data for mountain stations
(Bergamaschi et al., 2015),
selecting the 3 h time interval of the FLEXPART back trajectories (which
are provided for [00:00–03:00, 03:00–06:00, …] UTC), which is closest to
the time interval [12:00–15:00] LT for the stations in the boundary layer
and [00:00–03:00] LT for the mountain stations (indicated in Table 1 by
column “M”), respectively. This procedure ensures consistent averaging of
the FLEXPART back trajectories and the assimilated observations over the
same 3 h time intervals. Discrete air samples were taken as available,
i.e. without any temporal selection. The measurement uncertainty is set to
3 ppb for all observations (for observational part
**R**_{obs} of the observation error covariance
matrix Eq. 9).

In this study, we investigate two observation data sets (Table 1). The first data set, denoted O1, is considered the observational base data set and uses only the ICOS and NOAA data, while the second data set, O2, also includes all additional in situ measurements. The largest difference between the two data sets is the much better observational coverage of the British Isles in O2 with six in situ measurement stations located in that area, compared to only one station with discrete air sampling (MHD/NOAA) in O1.

## 3.2 Emission inventories

Three different emission inventories are used alternatively as prior
estimates of the major anthropogenic CH_{4} emissions (Table 2). The first
inventory is the Emissions Database for Global Atmospheric Research (EDGAR)
v6.0 (EDGAR v6.0, 2021), which provides monthly sector-specific
global grid maps of emissions at a horizontal resolution of $\mathrm{0.1}{}^{\circ}\times \mathrm{0.1}{}^{\circ}$ for 2000–2018. The second inventory,
TNO-VERIFYv3.0, is the third version of the TNO greenhouse gas and
co-emitted species (GHGco) emission database, developed by TNO within the
VERIFY project. TNO-VERIFYv3.0 provides annual European CH_{4} emissions
at a horizontal resolution of $\sim \mathrm{6}\phantom{\rule{0.125em}{0ex}}\mathrm{km}\times \mathrm{6}\phantom{\rule{0.125em}{0ex}}\mathrm{km}$ for the
years 2005–2018 but includes monthly emission profiles. The third emission
inventory has been provided by GCP-CH_{4} (Saunois et al., 2020)
globally at a horizontal resolution of $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$
for 2000–2017. In the absence of emission data for 2018 in the GCP-CH_{4}
inventory, we use here the 2017 data of this inventory. The resulting error
of this 1-year inconsistency, however, is considered to be much smaller
compared to the overall uncertainties of the emission inventories.

Natural CH_{4} emissions were used from the GCP-CH_{4} data set
(Saunois et al., 2020), providing estimates of the climatological mean
emissions of the major natural source categories. Furthermore, CH_{4}
emissions from biomass burning were taken from the Global Fire Emissions
Database (GFED) version 4.1 (Van Der
Werf et al., 2017). However, these were included only when using the EDGAR
v6.0 or TNO-VERIFYv3.0 inventories, while the GCP-CH_{4} (anthropogenic)
data set already includes emissions from biomass burning.

Using the above emission inventories, we have assembled the emission data sets E1, E2, and E3 as compiled in Table 2 and used as prior estimates for the different inversions described in Sect. 3.4. All emission data sets have been mapped on the COSMO-7 grid, using the Python package “emiproc” (Jähn et al., 2020), which has been integrated into the FLEXVAR inverse modelling system.

## 3.3 Post-processing of gridded emission data

In order to extract from gridded emission data (on COSMO-7 grid) total emissions of countries (or group of countries), country masks have been generated using the “Natural Earth dataset” (https://www.naturalearthdata.com/, last access: 30 July 2021), attributing each 7 km×7 km COSMO-7 grid cell to a certain country (or sea). Offshore emissions over the sea are not included in the country totals.

Since the COSMO-7 domain does not cover the upper northern part of the UK, a correction factor of 1.057 is applied to estimate the total emissions of the country region “UK + Ireland”; i.e. the UK + Ireland emissions extracted from the corresponding grid cells within the COSMO-7 domain are multiplied by this factor (for further details, see Sect. S1 in the Supplement). Furthermore, small correction factors are applied when extracting country total emissions from the gridded emissions of data set E3 (at horizontal resolution of $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$), since sampling of coastal $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$ grid cells with the corresponding 7 km×7 km COSMO-7 grid cells leads to a loss of emissions attributed to the countries, if the emissions of the coastal $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$ grid cell originate mainly from land (for further details, see Sect. S1).

## 3.4 Sensitivity inversions

Table 3 compiles the different FLEXVAR inversions presented in this paper.
INV-E1-O1 represents the base inversion, using the emission data set E1 as
prior estimates, the observation data set O1, the METEO model representation error
(Sect. 2.2.3), the Rödenbeck baselines (Sect. 2.2.2), and our default
settings for the prior error covariance. A first set of sensitivity
inversions investigates the impact of using alternatively the particle position baselines and the alternative parameterization OBS of the model
representation error (and the combination of both). In a further inversion
series, we analyse the sensitivity of the inversions to the main settings of
the prior error covariance matrix, i.e. for the spatial correlation length
constant, *L*_{corr}, the temporal correlation scale constant, *t*_{corr},
and the assumed uncertainty of emissions per grid cell and month.
Furthermore, we examine the sensitivity of the inversions to the use of the
alternative emission inventories E2 and E3 as prior estimates instead of E1 and the
use of the extended observational data set O2 instead of O1.

In addition to the FLEXVAR inversions compiled in Table 3, inversions with the FLExKF system (described in Sect. 2.3) and with TM5-4DVAR (described in Sect. 2.4) have been performed for comparison with FLEXVAR (and will be discussed in Sect. 4.3). These inversions have been made for both observational data sets, O1 and O2, using the emission inventory E3 as prior estimates. Furthermore, additional FLExKF inversions have been performed using alternatively E1 as prior estimates.

## 4.1 Sensitivity of FLEXVAR inversions to internal parameterizations and model settings

### 4.1.1 Sensitivity of FLEXVAR inversions to baselines

Figure 1 shows maps of European CH_{4} emissions derived for the base
inversion INV-E1-O1 and the sensitivity inversion INV-E1-O1-S1, in which the
particle position baselines were used instead of the Rödenbeck baselines. Both inversions display in general similar spatial patterns of
the inversion increments (i.e. difference between posterior and prior
emissions); however in most regions INV-E1-O1-S1 shows somewhat lower
CH_{4} emissions than INV-E1-O1, visible in the slightly larger areas with
negative inversion increments and slightly smaller areas with positive
inversion increments. Consequently, also the derived country total emissions
(shown in Fig. 2 and compiled in Table S4) are lower in INV-E1-O1-S1, e.g.
−6.6 *%* lower over Germany and −12.8 *%* lower over France compared to
INV-E1-O1.

Figure 3 illustrates the two different baselines at some example stations
during the 3-month period from 1 April until 1 July 2018. In general, both
baselines are rather similar, including their synoptic variability. However,
there are certain periods during which the particle position baselines
are somewhat higher than the Rödenbeck baselines, for example, at KIT, SAC,
and OPE during the period between day 140 and day 162. Consequently, the
observational forcing (i.e. the enhancement of the observations above
baseline) is lower during such periods for the particle position baselines, resulting in lower derived emissions. One major difference
between both approaches is that in the case of the Rödenbeck baselines, the
background mole fractions are transported to the stations by TM5, while in
the case of the particle position baselines, they are transported by FLEXPART.
In order to further investigate which baselines are more realistic, we have
compared model simulations and observations for “background conditions”,
defined as events when the contribution of European emissions (evaluated by
Eq. 10) is lower than a certain threshold (here set to 5 ppb). Figure S2
shows the comparison for eight stations, for which a sufficient number
(>20) of events with background conditions has been found. For
the Rödenbeck baselines, six of these eight stations show posterior biases
close to zero (<2 ppb), while PUY shows a small negative bias (−5.3 ppb) and CMN a small positive bias (3.7 ppb). In contrast, the particle position baselines result in a smaller negative bias at PUY (−3.3 ppb)
but larger positive biases at WAO (2.5 ppb), JFJ (5.5 ppb), and CMN (8.6 ppb). This analysis suggests that the performance of the Rödenbeck baselines is slightly better compared to the particle position baselines
under background conditions. However, we note that differences of the
baselines shown in Fig. 3 are mainly during periods of elevated CH_{4}
enhancements, for which it is more difficult to evaluate (based on the
observations) which baselines are more realistic.

### 4.1.2 Sensitivity of FLEXVAR inversions to parameterization of model representation error

Figure 1 illustrates the sensitivity of the derived emissions to the applied
parameterization of the model representation error. Inversion
INV-E1-O1-S2.1, for which the OBS model representation error has been
used, results in overall lower CH_{4} emissions compared to the base
inversion INV-E1-O1 with the METEO parameterization, again reflected in
the larger extension of the areas with negative inversion increments and
smaller extension (and magnitude) of the areas with positive inversion
increments. Accordingly, the annual total emissions derived in
INV-E1-O1-S2.1 are lower compared to INV-E1-O1 for all countries or
group of countries (denoted in the following as country regions) shown in
Fig. 2.

The OBS model representation error increases with increasing observed
CH_{4} enhancement (i.e. observed CH_{4} mole fraction minus CH_{4}
background) (Sect. 2.2.3 and Fig. S1) and shows a large dynamic range at
most stations, resulting in a generally relative low weighting in the
inversion of events with larger CH_{4} enhancements. In contrast, the
dynamic range of the METEO model representation error is smaller at most
stations, leading to a generally more equal weighting of all data points.
Using the METEO model representation error, the observations can be better
reproduced, achieving a higher average correlation coefficient (*r*=0.85)
and lower average root mean square difference (rms =30.0 ppb) compared to
the OBS model representation error (*r*=0.80; rms =35.4 ppb), as
shown in Fig. S3. Apart from the better statistical performance, the METEO
model representation is probably better at estimating the capability of the
model to reproduce the observations (which largely depends on the specific
meteorological situation), since wind speed might be a better indicator of
the representativeness of a certain data point than the observed CH_{4}
enhancement, as the latter not only depends on the meteorological situation,
but also on the regional CH_{4} emissions.

Given the relatively large impact of the parameterization of the model representation error and the baselines, we have also performed an inversion combining the OBS model representation error and the particle position baselines (inversion INV-E1-O1-S2.2), which yields further reduced country total emissions compared to INV-E1-O1-S2.1 and INV-E1-O1-S1 (Fig. 2).

### 4.1.3 Sensitivity of FLEXVAR inversions to model covariance settings

In the following, the sensitivity of the FLEXVAR inversions to the main
parameters of the prior covariance is investigated, i.e. horizontal
correlation length constant, temporal correlation scale constant, and
assumed uncertainties of emissions per grid cell and emission time step.
Figure S4 shows inversions for horizontal correlation length constants
*L*_{corr} (Eq. 7) of 50 km (INV-E1-O1-S3.1), 100 km
(default value; INV-E1-O1), and 200 km (INV-E1-O1-S3.2). As
expected, the spatial dimension of the inversion increments is increasing
with increasing *L*_{corr}. Despite these clearly visible
differences in the spatial patterns of the inversion increments, the impact
on the annual total emissions of the country regions shown in Fig. 2 is
relatively small, since apparently the differences in the smaller-scale
spatial patterns are largely averaged out over larger areas. Associated with
the increase of the horizontal correlation length constant is a significant
increase of the prior uncertainties of the annual total emissions per
country, since increasing horizontal correlation length constant implies
larger error correlations between neighbouring grid cells and hence
increasing aggregated uncertainties (as uncertainties per grid cell and
month were kept constant (at 100 %) in this sensitivity inversion series).
Analogously, the decrease in the temporal correlation scale constant,
*t*_{corr} (Eq. 8), results in a decrease of the aggregated
annual prior uncertainty, as illustrated by inversion INV-E1-O1-S5, in which
*t*_{corr} has been set to 1 month (instead of the default value of
12 months applied in all other inversions). Again, however, the effect on
the derived annual emissions of the country regions remains very small (Fig. 2).

Figure S5 shows the dependence of the inversions on the assumed uncertainties of prior emissions per grid cell and month for values of 50 % (INV-E1-O1-S4.1), 100 % (default value; INV-E1-O1), and 200 % (INV-E1-O1-S4.2). The increase of the assumed prior uncertainty leads to a significant increase of the derived regional inversion increments. This effect is most pronounced at larger distances from the monitoring stations, where observational constraints are relatively weak. Especially the large inversion increments visible in INV-E1-O1-S4.2 at the eastern domain boundary are probably an artefact, since the inversion may generate such patterns in regions far from the observations to compensate for systematic errors, for example, in model transport and with little penalty in the cost function in the case of prior uncertainties that are assumed very high.

Despite the dependence of the smaller-scale regional inversion increments on the assumed prior uncertainties, the impact on the derived annual total emissions remains again very small for the country regions shown in Fig. 2, since their emissions are relatively well constrained by the available observations and since differences in the smaller-scale inversion increments are averaged out over larger areas.

## 4.2 Sensitivity of FLEXVAR inversions to model input data

### 4.2.1 Sensitivity of FLEXVAR inversions to prior emission inventories

Figure 4 shows maps of the European CH_{4} emissions for INV-E1-O1,
INV-E2-O1, and INV-E3-O1, which use the three different emission data sets
E1, E2, and E3 (Sect. 3.2; Table 2) as prior emissions. While the major
patterns of the spatial prior emission distribution look relatively similar
for the three inventories (e.g. the high emissions over the BENELUX
countries and the Po valley), there are significant differences in the prior
country region total emissions (Fig. 2; Table S4). E2 has lower emissions
over Germany (16.1 %), France (15.1 %), and BENELUX (27.4 %) compared
to E1 (and 11.3 % lower over the whole COSMO-7 domain; Table 2), while E3
has higher total emissions over the COSMO-7 domain (15.1 % higher than
E1) and very high emissions in particular for UK + Ireland (42.1 % higher
than E1). Despite these considerable differences in the prior emissions, the
annual total posterior emissions of the country regions shown in Fig. 2 are
very similar for the three inversions. This indicates that the inversions
are largely driven by the observations. For UK + Ireland this is somewhat
surprising, since only one measurement station (MHD/NOAA) is located in
this country region in the applied observation data set O1, but apparently
the continental stations provide some constraints for the emissions from
UK + Ireland. We will see in the next section, however, that including
additional stations has a significant impact on the CH_{4} emissions
derived for UK + Ireland.

### 4.2.2 Sensitivity of FLEXVAR inversions to assimilated observations

While the base observation data set O1 uses only the ICOS in situ stations,
complemented by the NOAA discrete air sampling sites, nine further in situ
stations from other networks/institutions are added in observation data
set O2 (Table 1). Six of the additional stations are located on the British
Isles, two in Switzerland, and one in the Netherlands. Figure 5 displays the
inversions INV-E1-O1 and INV-E1-O2 using the two different observation data
sets. As expected, the largest differences are visible in the regions around
the additional stations. For UK + Ireland, the annual total emissions are
23.0 % higher in INV-E1-O2 (2.99 CH_{4} yr^{−1}) compared to INV-E1-O1
(2.43 CH_{4} yr^{−1}) (Fig. 2; Table S4). The significant additional
observational constraints for U K+ Ireland are also reflected in the
significantly lower posterior uncertainty for INV-E1-O2 (2*σ*
uncertainty: 0.6 Tg CH_{4} yr^{−1}) compared to INV-E1-O1 (2*σ*
uncertainty: 1.6 Tg CH_{4} yr^{−1}; Fig. 2; Table S4). For the BENELUX
country region, only a moderate change in the annual total emissions is
calculated (INV-E1-O1: 1.71 Tg CH_{4} yr^{−1}; INV-E1-O2: 1.82 Tg CH_{4} yr^{−1}; Fig. 2; Table S4), but the spatial distribution of
posterior emissions is somewhat different, with higher emissions around the
additional station CBW in INV-E1-O2 (Fig. 5). For Switzerland, a larger
(relative) difference of posterior emissions is calculated, with annual
total emission increasing from 0.15 Tg CH_{4} yr^{−1} (INV-E1-O1) to
0.22 Tg CH_{4} yr^{−1} (INV-E1-O2).

Using the extended observation data set O2, we have performed additional inversions, using alternatively the prior emission data sets E2 or E3 instead of E1. As for observation data set O1 (discussed in Sect. 4.2.1), the sensitivity of derived annual total emissions to the applied prior emission data set is relatively small (Fig. 2).

Furthermore, additional inversions (of observation data set O2) have been performed using alternatively the particle position baselines (INV-E1-O2-S1) or the alternative parameterization OBS of the model representation error (INV-E1-O2-S2.1). In a similar way, as shown with observation data set O1 (discussed in Sect. 4.1.1. and 4.1.2.), the use of these alternative parameterizations results in generally lower posterior emissions, with lowest posterior emission calculated in inversion INV-E1-O2-S2.2 (combining the OBS model representation error and the particle position baselines).

## 4.3 Model comparison and analysis of European CH_{4} emissions

In the following we compare the FLEXVAR inversions with inversions using the
extended Kalman filter (FLExKF) system (Sect. 2.3) and TM5-4DVAR (Sect. 2.4). Figure 6 shows the results of these three models using the emission
data set E3 as prior and the observation data set O2. Overall, all three
inverse models show relatively good consistency of the major spatial
patterns of the derived inversion increments, for example, the increase of
emissions over the BENELUX region and north-western France, the decrease of
emissions around Paris, and the decrease of offshore emissions over the
North Sea compared to the prior emissions. Since FLExKF uses the same
atmospheric transport as FLEXVAR, it is to be expected that the inversions
of these two models should give similar results. Nevertheless, there are
also some significant differences visible between the two models, especially
for the southern part of France, for which FLExKF yields overall lower
emissions than FLEXVAR. This difference is also clearly visible in the
derived country total emissions (Fig. 7; Table S4), with 10.3 % lower
annual total CH_{4} emission for France calculated by FLExKF (FLExKF
E3-O2: 3.82 CH_{4} yr^{−1}) compared to FLEXVAR (INV-E3-O2: 4.26 CH_{4} yr^{−1}). In contrast, FLExKF derives somewhat higher CH_{4}
emissions for BENELUX (6.3 %) and UK + Ireland (6.8 %) than FLEXVAR,
while emissions derived for Germany are very similar (within 1.4 %). One
major difference between FLExKF and FLEXVAR is the different
parameterization of the model representation error, leading to a different
weighting of the individual observational data points, which can cause
differences in the calculated regional inversion increments as shown for
FLEXVAR in Sect. 4.1.2. Another difference is the magnitude of the prior
uncertainties, though this was shown for FLEXVAR to have a rather small
impact on total emissions for the country regions presented in Fig. 2.
Furthermore, it is likely that the different inversion techniques have some
impact on the calculated solutions. For example, FLExKF yields generally
smoother seasonal variations of derived emissions, while FLEXVAR shows
larger month-to-month variability. The latter are, however, largely filtered
out by the use of 3-month running mean values for the seasonal variation of
the total emissions of country regions shown in Fig. 7 (left column).

The spatial emission patterns derived by TM5-4DVAR are in general similar to those calculated by FLEXVAR and FLExKF (Fig. 6) but show also some differences, for example, around the stations PUY and HPB, where TM5-4DVAR calculates higher emissions than FLEXVAR and FLExKF, probably related to the particular challenge to simulate mountain sites and sites in complex topography. Further differences between the models are the different derived seasonal variations of emissions, with larger variations calculated by TM5-4DVAR for Germany, France, and UK + Ireland compared to FLEXVAR and FLExKF (while the FLEXVAR inversions using the observation data set O2 show larger variations for BENELUX than the other models). In addition to the different model representation error in TM5-4DVAR, very likely the fundamentally different nature of the models (Eulerian vs. Lagrangian) and the related different simulation of transport play an essential role.

Nevertheless, the differences in the annual total emissions for the country regions are only moderate. For Germany, somewhat higher emissions are calculated by TM5-4DVAR compared to FLEXVAR and FLExKF, while the posterior emissions for France, BENELUX, and UK + Ireland derived by TM5-4DVAR are in the range of emissions calculated by FLEXVAR and FLExKF.

Figure 7 also includes inversions of the three models using the base observation data set O1. As discussed for FLEXVAR in Sect. 4.2.2., FLExKF and TM5-4DVAR also show higher emissions for UK + Ireland, when using O2 instead of O1 due to the six additional stations in data set O2 in that area. Furthermore, FLExKF inversions have also been performed using E1 instead of E3 as prior emissions (Fig. 7). As for FLEXVAR (Sect. 4.2.1), the impact on derived emissions is relatively small.

In order to evaluate the quality of the derived emissions, it is useful to
analyse how well the observations are reproduced by the models. Figure S6
compares the statistics (correlation coefficient and rms difference) for the
three models (using prior emission data set E3 and observation data set O2).
At most stations relatively high correlation coefficients and low rms
differences are obtained by all three models. However, stations with larger
regional emissions (e.g. LUT, CBW, BRM, IPR) or complex topography (e.g.
OXK, IPR) show generally poorer statistical performance. Figure S6 also
shows that the best statistical performance is achieved by FLEXVAR with a
mean correlation coefficient of *r*=0.86 (FLExKF: *r*=0.84, TM5-4DVAR:
*r*=0.81) and a mean rms difference of 28.21 ppb (FLExKF: 30.53 ppb,
TM5-4DVAR: 31.82 ppb). This finding demonstrates that the high spatial
resolution of FLEXVAR and FLExKF at 7 km×7 km allows the observations to be somewhat
better reproduced than the TM5-4DVAR simulations at
$\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$, although – beside the different
spatial resolution – other factors (such as fundamental differences in
the modelling of transport) are also likely to play a role. The slightly better
statistical performance of FLEXVAR compared to FLExKF could be due to the
higher degree of freedom to optimize the emissions in FLEXVAR but may also
be partly related to other factors, such as different weighting of
observations due to different parameterizations of the model representation
error and differences in the model covariance settings. For example, the FLEXVAR
inversion INV-E3-O2 used for the model comparison applies a smaller spatial
correlation length (*L*_{corr}=100 km) compared to FLExKF
(*L*_{corr}=200 km). Comparison of FLEXVAR inversions INV-E1-O1
and INV-E1-O1-S3.2 shows that increasing the correlation length from 100
to 200 km is indeed slightly deteriorating the statistical performance (mean
correlation coefficient and mean rms difference), but nevertheless FLEXVAR
(INV-E1-O1-S3.2) still performs slightly better compared to FLExKF
(inversion FLExKF E1-O1). On the other hand, FLExKF applies a higher prior
uncertainty than FLEXVAR (Table S4) in the model comparison discussed in the
paper. For FLEXVAR, increasing the prior uncertainty from 100 % to 200 %
(INV-E1-O1-S4.2 vs. INV-E1-O1), slightly improves the statistical
performance, i.e. partly compensating for the effect of a larger correlation
length (results not shown).

Figure S7 shows the time series of observed and simulated CH_{4} mole
fractions for all stations (inversion INV-E1-O2), illustrating that in
general the synoptic variability is well reproduced at most sites.
Furthermore, FLEXVAR also simulates the average diurnal cycle at most sites
realistically.

In the following, we compare the annual total CH_{4} emissions derived by
the inverse models with the anthropogenic CH_{4} emissions reported by the
countries to UNFCCC (UNFCCC, 2021). For a consistent comparison, it is
necessary also to take into account estimates of the natural CH_{4}
emissions, for which we use the bottom-up inventories of natural sources
from the GCP-CH_{4} data set (Saunois et al., 2020) (Table 2).
Furthermore, the comparison of top-down and bottom-up emission estimates
requires the inclusion of estimates of their uncertainties. For the
uncertainty estimate of the inverse models, we use the range of results from
the individual inversions (shown by the solid red rectangles in Fig. 7) and
the minimum–maximum values of the 2*σ* uncertainty ranges based on the
uncertainties computed for the individual inversions (shown by the error
bars). The total uncertainty ranges are evaluated separately (1) for the
whole set of FLEXVAR sensitivity inversions (as shown in Fig. 2) and (2) for
the whole set of all inversions, i.e. also including all FLExKF and
TM5-4DVAR inversions. The uncertainties of the UNFCCC emissions are based on
the uncertainties reported by the countries for the major CH_{4} source
categories, while estimates of the uncertainty of total CH_{4} emissions
are not provided by the countries. As in Bergamaschi et al. (2015), we estimate the total
uncertainties from the reported uncertainties per category, assuming – among
other things – uncorrelated uncertainties for the different major source
categories (for further details, see Sect. S2). The uncertainties of natural CH_{4} emissions from wetlands
were estimated from the ensemble of wetland models used for the GCP-CH_{4}
wetland emissions, taking the minimum–maximum range of the 11 individual
wetland models (Poulter et al., 2017). For other natural
CH_{4} emissions, we assume an uncertainty of 100 %.

Figure 7 shows that the CH_{4} emissions estimated by the inverse models
are higher than the sum of anthropogenic (UNFCCC) and natural bottom-up
inventories for Germany, France, and BENELUX, but the uncertainty ranges of
top-down and bottom-up estimates overlap for all three country regions. The
smallest overlap, however, is found for BENELUX. In contrast, the top-down
estimates for UK + Ireland agree relatively well with the total of
anthropogenic and natural bottom-up inventories. A tendency to higher
top-down emissions compared to the total (anthropogenic and natural)
bottom-up inventories for Germany, France, and BENELUX has also been found in
the analysis reported by Bergamaschi et al. (2018a) for the period
2006–2012, but also in that study uncertainty ranges of bottom-up and
top-down estimates were overlapping. Similar tendencies to higher top-down
emissions are apparent in the VERIFY analyses for the period 2005–2017
(VERIFY, 2021) using a larger ensemble of regional inversions, while
global inversions (with coarser resolution) showed in general lower
emissions, closer to the UNFCCC estimates for these country regions. Based
on the observation that several models showed clear seasonal cycles of the
derived emissions with maximum during summer, Bergamaschi et al. (2018a)
suggested that higher natural emissions could explain the difference between
top-down and bottom-up estimates. The FLEXVAR and FLExKF inversions analysed
in this study, however, show in general only relatively small seasonal
variations for Germany, France, and UK + Ireland compared to TM5-4DVAR. The
use of seasonal cycles to disentangle anthropogenic and natural sources is
further hampered by the fact that the seasonal cycles of major anthropogenic
sources are still not well characterized. Also, the anthropogenic emission
inventories used in this study show rather different seasonal cycles. Most
of the anthropogenic GCP-CH_{4} emission categories (which are largely
based on EDGARv4.3.2, except biomass burning, Saunois et al., 2020) have
no seasonality, except emissions from rice agriculture and biomass burning,
which however play only a minor role in Europe. EDGARv6.0 (used for E1) has
small seasonal variations of most energy-related source categories but
assumes constant emissions for the agricultural sources (except rice) and
for waste emission. In contrast, most sectors of the TNO-VERIFYv3.0
inventory (used for E2) show seasonal variations, including significant
seasonal variations of all agricultural sources, resulting in significant
seasonal variations of the total anthropogenic emissions with maximum
emissions in September (Fig. 2).

We have presented the novel inverse modelling system FLEXVAR based on the
4DVAR assimilation technique and FLEXPART-COSMO back trajectories driven by
COSMO meteorological fields at 7 km×7 km resolution over the
European COSMO-7 domain. A major advantage of the 4DVAR technique is that it
allows a much larger number of variables to be constrained (in our study about
1.6 million) compared to analytical inversion techniques. The offline
coupling with TM5-4DVAR ensures that the background mole fractions
(baselines) used in FLEXVAR are consistent with the global observations
assimilated in TM5-4DVAR. We have applied the FLEXVAR system for inversions
of European CH_{4} emissions in 2018 using 24 stations with in situ
measurements, complemented with data from 5 stations with discrete air
sampling (and additional stations outside the European COSMO-7 domain used
for the global TM5-4DVAR inversions).

We have investigated the sensitivity of the FLEXVAR inversions to internal
parameterizations, model settings, and main model input data. Using the
particle position baselines yields in general lower derived emissions
compared to inversions which apply the Rödenbeck baselines, resulting
in differences in the annual total emissions of 5 %–14 % for the analysed
country regions (Germany, France, BENELUX, and UK + Ireland). Furthermore we
found a significant impact of the applied parameterization of the model
representation error. Inversions using the OBS model representation error
derive, over large parts of the domain, somewhat lower emissions compared to
the METEO model representation error, with differences in the annual total
emissions of 0 %–15 % for the analysed country regions. Varying the main
parameters of the prior covariance (i.e. horizontal correlation length
constant, temporal correlation scale constant, and assumed uncertainties of
emissions per grid cell and month) has clearly visible effects on the
smaller-scale regional inversion increments, but the impact on the derived
annual total emissions remains very small for the analysed country regions,
since the differences in the smaller-scale spatial patterns are largely
averaged out over larger areas. Furthermore, the dependence of derived
emissions on the applied prior emission inventory has been found to be
relatively small for the country regions, which are well constrained by the
observations. Changing these observational constraints by including
additional sites, however, has a significant impact on the inversions,
especially in the vicinity of these sites. Using the extended observation
data set O2 (which includes six additional in situ stations located on the
British Isles) yields 23 %–28 % higher emissions for UK + Ireland
compared to inversions using only the base observation data set O1. At the
same time, the calculated uncertainty of the posterior emissions for
UK + Ireland is significantly reduced by these additional observational
constraints.
The comparison of the FLEXVAR inversions with inversions using the extended
Kalman filter (“FLExKF”) system (which both use the same atmospheric
transport model) shows overall good consistency of major spatial patterns of
the derived inversion increments but some difference (7 %–10 %) for the
derived total CH_{4} emission of France, probably mostly related to the
use of different parameterizations of the model representation error.
TM5-4DVAR also shows in general similar inversion increments and derives
posterior emissions for France, BENELUX, and UK + Ireland in the range of
emissions calculated by FLEXVAR and FLExKF. For Germany, however, TM5-4DVAR
estimates 5 %–11 % higher emissions than the other two models.

The FLEXVAR and FLExKF inversions at high spatial resolution of 7 km×7 km allow for a better reproduction of the observations compared to the TM5-4DVAR simulations at $\mathrm{1}{}^{\circ}\times \mathrm{1}{}^{\circ}$, reflected in the achieved higher correlation coefficient and lower rms difference between simulations and observations. Furthermore, the statistical performance of FLEXVAR is slightly better than that of FLExKF, which could be due to the higher degree of freedom to optimize the emissions in FLEXVAR but could be partly related also to other differences in the inversions, for example, the different parameterizations of the model representation error and differences in the model covariance settings.

The inverse models derive higher annual total CH_{4} emissions in 2018 for
Germany, France, and BENELUX compared to the sum of emissions reported to
UNFCCC and natural emissions (estimated from the GCP-CH_{4} inventory),
but the uncertainty ranges of top-down and bottom-up estimates overlap for
all three country regions. In contrast, the top-down estimates for
UK + Ireland agree relatively well with the total of anthropogenic and
natural bottom-up inventories.

Our study demonstrates that the new FLEXVAR system can be applied for
verification of reported emissions, as planned, for example, by Empa for its
quasi-operational system to estimate Switzerland's annual CH_{4} emissions
as a contribution to the Swiss National Inventory Reporting. FLEXVAR
inversions with the configuration presented in this paper could be performed
for the years 2002 to 2021, the period for which meteorological fields from
the COSMO-7 model at 7 km×7 km resolution are available. For
analysis periods after 2021, the use of different high-resolution
meteorological input fields could be considered, such as the
operational analysis data from the ECMWF IFS model at high resolution
($\mathrm{0.1}{}^{\circ}\times \mathrm{0.1}{}^{\circ}$) or the operational MeteoSwiss COSMO-1
analysis at horizontal resolution of 1 km×1 km. COSMO-1, however,
is limited to the larger Alpine area but can be nested into FLEXPART-IFS. A
FLEXPART-COSMO modelling system using COSMO-1 has already been developed by
Empa, including a modification of the turbulence parameterization
(Katharopoulos et al., 2022), which is required owing to the very high
resolution of 1 km×1 km.

While the relatively good agreement among the three models used in this study gives some confidence in the robustness of the inverse modelling results, further specific studies should be performed to assess the quality of the top-down estimates independently. Such assessments should include the comparison with further inverse models, comparison with independent regional emission estimates (e.g. based on aircraft or satellite measurements), and a more detailed validation of the applied atmospheric transport models (especially regarding the simulation of boundary layer height dynamics and vertical transport).

The code of the FLEXVAR inverse modelling system is available upon request. The atmospheric observations from ICOS are available at https://www.icos-cp.eu/data-products/atmosphere-release (ICOS RI, 2021). NOAA data are available at https://doi.org/10.15138/VNCZ-M766 (Dlugokencky et al., 2021), AGAGE data at https://doi.org/10.15485/1781803 (Prinn et al., 2021), and UK DECC data at https://archive.ceda.ac.uk/ (O'Doherty et al., 2022).

The supplement related to this article is available online at: https://doi.org/10.5194/acp-22-13243-2022-supplement.

PB led the FLEXVAR development. AS designed the FLEXVAR concept and performed the technical implementation of the FLEXVAR code. PB performed the FLEXVAR and TM5-4DVAR inversions. DB performed the FLExKF inversions and contributed to the integration of the FLEXPART-COSMO back trajectories into FLEXVAR. JMH generated the FLEXPART-COSMO back trajectories and contributed to the development of the FLEXVAR emission preprocessing. SH contributed to the development of the FLEXPART-COSMO modelling system. PB prepared the paper and figures with contributions from DB, AS, and JMH. SH, MR, TA, TB, HC, SC, MD, GF, AF, DK, XL, MLe, MLi, MLo, GM, JMW, SOD, BS, MS, PT, GV, and CYK provided atmospheric observation data.

The contact author has declared that none of the authors has any competing interests.

Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

The VERIFY project has received funding from the European Union's Horizon
2020 Research and Innovation programme under grant agreement no. 776810. We
are grateful to ECMWF for providing computing resources under the special
projects “Improve European and global CH_{4} and N_{2}O flux inversions
(2018–2020)” and “Extend and improve CH_{4} flux inversions at global
and European scale (2021)”. Furthermore, we thank Bradley Matthews for the
compilation of uncertainties of emissions reported to UNFCCC by EU member
states. ICOS Switzerland was funded by the Swiss National Science
Foundation, in-house contributions, and the State Secretariat for Education,
Research and Innovation.

This research has been supported by Horizon 2020 (grant no. VERIFY (776810)).

This paper was edited by Ilse Aben and reviewed by two anonymous referees.

Aoki, S., Nakazawa, T., Murayama, S., and Kawaguchi, S.: Measurements of atmospheric methane at the Japanese Antarctic Station, Syowa, Tellus B, 44, 273–281, 1992.

Arnoldi, W. E.: The principle of minimized iterations in the solution of the matrix eigenvalue problem, Q. Appl. Math., 9, 17–29, 1951.

Baldauf, M., Seifert, A., Förstner, J., Majewski, D., Raschendorfer, M., and Reinhardt, T.: Operational Convective-Scale Numerical Weather Prediction with the COSMO Model: Description and Sensitivities, Mon. Weather Rev., 139, 3887–3905, https://doi.org/10.1175/MWR-D-10-05013.1, 2011.

Bergamaschi, P., Krol, M., Meirink, J. F., Dentener, F., Segers, A., van
Aardenne, J., Monni, S., Vermeulen, A., Schmidt, M., Ramonet, M., Yver, C.,
Meinhardt, F., Nisbet, E. G., Fisher, R., O'Doherty, S., and Dlugokencky, E.
J.: Inverse modeling of European CH_{4} emissions 2001–2006, J. Geophys.
Res., 115, D22309, https://doi.org/10.1029/2010JD014180, 2010.

Bergamaschi, P., Corazza, M., Karstens, U., Athanassiadou, M., Thompson, R. L., Pison, I., Manning, A. J., Bousquet, P., Segers, A., Vermeulen, A. T., Janssens-Maenhout, G., Schmidt, M., Ramonet, M., Meinhardt, F., Aalto, T., Haszpra, L., Moncrieff, J., Popa, M. E., Lowry, D., Steinbacher, M., Jordan, A., O'Doherty, S., Piacentino, S., and Dlugokencky, E.: Top-down estimates of European CH_{4} and N_{2}O emissions based on four different inverse models, Atmos. Chem. Phys., 15, 715–736, https://doi.org/10.5194/acp-15-715-2015, 2015.

Bergamaschi, P., Karstens, U., Manning, A. J., Saunois, M., Tsuruta, A., Berchet, A., Vermeulen, A. T., Arnold, T., Janssens-Maenhout, G., Hammer, S., Levin, I., Schmidt, M., Ramonet, M., Lopez, M., Lavric, J., Aalto, T., Chen, H., Feist, D. G., Gerbig, C., Haszpra, L., Hermansen, O., Manca, G., Moncrieff, J., Meinhardt, F., Necki, J., Galkowski, M., O'Doherty, S., Paramonova, N., Scheeren, H. A., Steinbacher, M., and Dlugokencky, E.: Inverse modelling of European CH_{4} emissions during 2006–2012 using different inverse models and reassessed atmospheric observations, Atmos. Chem. Phys., 18, 901–920, https://doi.org/10.5194/acp-18-901-2018, 2018a.

Bergamaschi, P., Danila, A., Weiss, R. F., Ciais, P., Thompson, R. L., Brunner, D., Levin, I., Meijer, Y., Chevallier, F., Janssens-Maenhout, G., Bovensmann, H., Crisp, D., Basu, S., Dlugokencky, E., Engelen, R., Gerbig, C., Günther, D., Hammer, S., Henne, S., Houweling, S., Karstens, U., Kort, E., Maione, M., Manning, A. J., Miller, J., Montzka, S., Pandey, S., Peters, W., Peylin, P., Pinty, B., Ramonet, M., Reimann, S., Röckmann, T., Schmidt, M., Strogies, M., Sussams, J., Tarasova, O., Aardenne, J. v., Vermeulen, A. T., and Vogel, F.: Atmospheric monitoring and inverse modelling for verification of greenhouse gas inventories, EUR 29276 EN, Publications Office of the European Union, Luxembourg, 2018, ISBN 978-92-79-88938-7, https://doi.org/10.2760/759928, JRC111789, 2018b.

Brandt, A. R., Heath, G. A., Kort, E. A., O'Sullivan, F., Pétron, G., Jordaan, S. M., Tans, P., Wilcox, J., Gopstein, A. M., Arent, D., Wofsy, S., Brown, N. J., Bradley, R., Stucky, G. D., Eardley, D., and Harriss, R.: Methane Leaks from North American Natural Gas Systems, Science, 343, 733–735, https://doi.org/10.1126/science.1247045, 2014.

Brunner, D., Henne, S., Keller, C. A., Reimann, S., Vollmer, M. K., O'Doherty, S., and Maione, M.: An extended Kalman-filter for regional scale inverse emission estimation, Atmos. Chem. Phys., 12, 3455–3478, https://doi.org/10.5194/acp-12-3455-2012, 2012.

Brunner, D., Henne, S., Keller, C. A., Vollmer, M. K., and Reimann, S.: Estimating European Halocarbon Emissions Using Lagrangian Backward Transport Modeling and in Situ Measurements at the Jungfraujoch High-Alpine Site, in: Lagrangian Modeling of the Atmosphere, edited by: Lin, J. C., Gerbig, C., Brunner, D., Stohl, A., Luhar, A., and Webley, P., Vol. 200 of Geophysical Monographs V, AGU, Washington, DC, 207–221, https://doi.org/10.1029/2012GM001258, 2013.

Brunner, D., Arnold, T., Henne, S., Manning, A., Thompson, R. L., Maione, M., O'Doherty, S., and Reimann, S.: Comparison of four inverse modelling systems applied to the estimation of HFC-125, HFC-134a, and SF_{6} emissions over Europe, Atmos. Chem. Phys., 17, 10651–10674, https://doi.org/10.5194/acp-17-10651-2017, 2017.

Butler, J. H. and Montzka, S. A.: The NOAA annual greenhouse gas index (AGGI), https://www.esrl.noaa.gov/gmd/aggi/aggi.html, last access: 16 February 2022.

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A. C. M., Berg, L. v. d., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer, A. J., Haimberger, L., Healy, S. B., Hersbach, H., Hólm, E. V., Isaksen, L., Kållberg, P., Köhler, M., Matricardi, M., McNally, A. P., Monge-Sanz, B. M., Morcrette, J.-J., Park, B.-K., Peubey, C., Rosnay, P. d., Tavolato, C., Thépaut, J.-N., and Vitart, F.: The ERA-Interim reanalysis: configuration and performance of the data assimilation system, Q. J. Roy. Meteor. Soc., 137, 553–597, https://doi.org/10.1002/qj.828, 2011.

Deng, Z., Ciais, P., Tzompa-Sosa, Z. A., Saunois, M., Qiu, C., Tan, C., Sun, T., Ke, P., Cui, Y., Tanaka, K., Lin, X., Thompson, R. L., Tian, H., Yao, Y., Huang, Y., Lauerwald, R., Jain, A. K., Xu, X., Bastos, A., Sitch, S., Palmer, P. I., Lauvaux, T., d'Aspremont, A., Giron, C., Benoit, A., Poulter, B., Chang, J., Petrescu, A. M. R., Davis, S. J., Liu, Z., Grassi, G., Albergel, C., Tubiello, F. N., Perugini, L., Peters, W., and Chevallier, F.: Comparing national greenhouse gas budgets reported in UNFCCC inventories against atmospheric inversions, Earth Syst. Sci. Data, 14, 1639–1675, https://doi.org/10.5194/essd-14-1639-2022, 2022.

Dlugokencky, E.: Trends in Atmospheric Methane: https://www.esrl.noaa.gov/gmd/ccgg/trends_ch4/, last access: 16 February 2022.

Dlugokencky, E. J., Myers, R. C., Lang, P. M., Masarie, K. A., Crotwell, A.
M., Thoning, K. W., Hall, B. D., Elkins, J. W., and Steele, L. P.:
Conversion of NOAA atmospheric dry air CH_{4} mole fractions to a
gravimetrically prepared standard scale, J. Geophys. Res., 110,
D18306, https://doi.org/10.1029/2005JD006035, 2005.

Dlugokencky, E. J., Crotwell, A. M., Mund, J. W., Crotwell, M. J., and Thoning, K. W.: Atmospheric Methane Dry Air Mole Fractions from the NOAA GML Carbon Cycle Cooperative Global Air Sampling Network, 1983–2020, version: 2021-07-30, https://doi.org/10.15138/VNCZ-M766, 2021.

EDGAR v6.0: EDGAR – Emissions Database for Global Atmospheric Research, Global Greenhouse Gas Emissions, EDGAR v6.0, https://edgar.jrc.ec.europa.eu/dataset_ghg60 (last access: 8 August 2021), 2021.

El Yazidi, A., Ramonet, M., Ciais, P., Broquet, G., Pison, I., Abbaris, A., Brunner, D., Conil, S., Delmotte, M., Gheusi, F., Guerin, F., Hazan, L., Kachroudi, N., Kouvarakis, G., Mihalopoulos, N., Rivier, L., and Serça, D.: Identification of spikes associated with local sources in continuous time series of atmospheric CO, CO_{2} and CH_{4}, Atmos. Meas. Tech., 11, 1599–1614, https://doi.org/10.5194/amt-11-1599-2018, 2018.

European Commission: Launch by United States, the European Union, and Partners of the Global Methane Pledge to Keep 1.5C Within Reach, https://ec.europa.eu/commission/presscorner/detail/en/statement_21_5766, last access: 8 November 2021.

Fisher, M. and Courtier, P.: Estimating the covariance matrices of analysis and forecast error in variational data assimilation, ECMWF, Reading, UK, Technical Memorandum 220, 29 pp., 1995.

Forster, P., Storelvmo, T., Armour, K., Collins, W., Dufresne, J. L., Frame, D., Lunt, D. J., Mauritsen, T., Palmer, M. D., Watanabe, M., Wild, M., and Zhang, H.: The Earth's Energy Budget, Climate Feedbacks, and Climate Sensitivity, in: Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Masson-Delmotte, V., Zhai, P., Pirani, A., Connors, S. L., Péan, C., Berger, S., Caud, N., Chen, Y., Goldfarb, L., Gomis, M. I., Huang, M., Leitzell, K., Lonnoy, E., Matthews, J. B. R., Maycock, T. K., Waterfield, T., Yelekçi, O., Yu, R., and Zhou, B., Cambridge University Press, 923–1054, 2021.

Ganesan, A. L., Manning, A. J., Grant, A., Young, D., Oram, D. E., Sturges, W. T., Moncrieff, J. B., and O'Doherty, S.: Quantifying methane and nitrous oxide emissions from the UK and Ireland using a national-scale monitoring network, Atmos. Chem. Phys., 15, 6393–6406, https://doi.org/10.5194/acp-15-6393-2015, 2015.

Gilbert, J. C. and Lemaréchal, C.: Some numerical experiments with variable-storage quasi-Newton algorithms, Math. Programm., 45, 407–435, 1989.

Hazan, L., Tarniewicz, J., Ramonet, M., Laurent, O., and Abbaris, A.: Automatic processing of atmospheric CO_{2} and CH_{4} mole fractions at the ICOS Atmosphere Thematic Centre, Atmos. Meas. Tech., 9, 4719–4736, https://doi.org/10.5194/amt-9-4719-2016, 2016.

Heiskanen, J., Brümmer, C., Buchmann, N., Calfapietra, C., Chen, H., Gielen, B., Gkritzalis, T., Hammer, S., Hartman, S., Herbst, M., Janssens, I. A., Jordan, A., Juurola, E., Karstens, U., Kasurinen, V., Kruijt, B., Lankreijer, H., Levin, I., Linderson, M.-L., Loustau, D., Merbold, L., Myhre, C. L., Papale, D., Pavelka, M., Pilegaard, K., Ramonet, M., Rebmann, C., Rinne, J., Rivier, L., Saltikoff, E., Sanders, R., Steinbacher, M., Steinhoff, T., Watson, A., Vermeulen, A. T., Vesala, T., Vítková, G., and Kutsch, W.: The Integrated Carbon Observation System in Europe, B. Am. Meteor. Soc., 103, E855–E872 https://doi.org/10.1175/BAMS-D-19-0364.1, 2022.

Henne, S., Brunner, D., Oney, B., Leuenberger, M., Eugster, W., Bamberger, I., Meinhardt, F., Steinbacher, M., and Emmenegger, L.: Validation of the Swiss methane emission inventory by atmospheric observations and inverse modelling, Atmos. Chem. Phys., 16, 3683–3710, https://doi.org/10.5194/acp-16-3683-2016, 2016.

Houweling, S., Bergamaschi, P., Chevallier, F., Heimann, M., Kaminski, T., Krol, M., Michalak, A. M., and Patra, P.: Global inverse modeling of CH_{4} sources and sinks: an overview of methods, Atmos. Chem. Phys., 17, 235–256, https://doi.org/10.5194/acp-17-235-2017, 2017.

ICOS RI: ICOS Atmosphere Station Specifications V2.0, edited by: Laurent, O., ICOS ERIC, https://doi.org/10.18160/GK28-2188, 2020.

ICOS RI: ICOS Atmosphere Release 2021-1 of Level 2 Greenhouse Gas Mole
Fractions of CO_{2}, CH_{4}, N_{2}O, CO, meteorology and
^{14}CO_{2} (1.0), ICOS ERIC – Carbon Portal, https://doi.org/10.18160/wjy7-5d06, 2021.

IPCC: Summary for Policymakers, in: Global Warming of 1.5 ^{∘}C. An
IPCC Special Report on the impacts of global warming of 1.5 ^{∘}C
above pre-industrial levels and related global greenhouse gas emission
pathways, in the context of strengthening the global response to the threat
of climate change, sustainable development, and efforts to eradicate
poverty, edited by: Masson-Delmotte, V., Zhai, P., Pörtner, H.-O.,
Roberts, D., Skea, J., Shukla, P. R., Pirani, A., Moufouma-Okia, W.,
Péan, C., Pidcock, R., Connors, S., Matthews, J. B. R., Chen, Y., Zhou,
X., Gomis, M. I., Lonnoy, E., Maycock, T., Tignor, M., and Waterfield, T.,
Cambridge University Press, Cambridge, UK and New
York, NY, USA, 3–24, https://doi.org/10.1017/9781009157940.001, 2018.

Jackson, R. B., Saunois, M., Bousquet, P., Canadell, J. G., Poulter, B., Stavert, A. R., Bergamaschi, P., Niwa, Y., Segers, A., and Tsuruta, A.: Increasing anthropogenic methane emissions arise equally from agricultural and fossil fuel sources, Environ. Res. Lett., 15, 071002, https://doi.org/10.1088/1748-9326/ab9ed2, 2020.

Jähn, M., Kuhlmann, G., Mu, Q., Haussaire, J.-M., Ochsner, D., Osterried, K., Clément, V., and Brunner, D.: An online emission module for atmospheric chemistry transport models: implementation in COSMO-GHG v5.6a and COSMO-ART v5.1-3.1, Geosci. Model Dev., 13, 2379–2392, https://doi.org/10.5194/gmd-13-2379-2020, 2020.

Katharopoulos, I., Brunner, D., Emmenegger, L., Leuenberger, M., and Henne, S.: Lagrangian Particle Dispersion Models in the Grey Zone of Turbulence: Adaptations to FLEXPART-COSMO for Simulations at 1 km Grid Resolution, Bound.-Lay. Meteorol., 185, 129–160, https://doi.org/10.1007/s10546-022-00728-3, 2022.

Krol, M., Houweling, S., Bregman, B., van den Broek, M., Segers, A., van Velthoven, P., Peters, W., Dentener, F., and Bergamaschi, P.: The two-way nested global chemistry-transport zoom model TM5: algorithm and applications, Atmos. Chem. Phys., 5, 417–432, https://doi.org/10.5194/acp-5-417-2005, 2005.

Krol, M. C., Meirink, J. F., Bergamaschi, P., Mak, J. E., Lowe, D., Jöckel, P., Houweling, S., and Röckmann, T.: What can ^{14}CO measurements tell us about OH?, Atmos. Chem. Phys., 8, 5033–5044, https://doi.org/10.5194/acp-8-5033-2008, 2008.

Lanczos, C.: An iteration method for the solution of the eigenvalue problem, J. Res. Nat. Bur. Stand., 45, 255–282, 1950.

Lunt, M. F., Manning, A. J., Allen, G., Arnold, T., Bauguitte, S. J.-B., Boesch, H., Ganesan, A. L., Grant, A., Helfter, C., Nemitz, E., O'Doherty, S. J., Palmer, P. I., Pitt, J. R., Rennick, C., Say, D., Stanley, K. M., Stavert, A. R., Young, D., and Rigby, M.: Atmospheric observations consistent with reported decline in the UK's methane emissions (2013–2020), Atmos. Chem. Phys., 21, 16257–16276, https://doi.org/10.5194/acp-21-16257-2021, 2021.

Manning, A. J., O'Doherty, S., Jones, A. R., Simmonds, P. G., and Derwent, R. G.: Estimating UK methane and nitrous oxide emissions from 1990 to 2007 using an inversion modeling approach, J. Geophys. Res.-Atmos., 116, D02305, https://doi.org/10.1029/2010jd014763, 2011.

Meirink, J. F., Bergamaschi, P., and Krol, M. C.: Four-dimensional variational data assimilation for inverse modelling of atmospheric methane emissions: method and comparison with synthesis inversion, Atmos. Chem. Phys., 8, 6341–6353, https://doi.org/10.5194/acp-8-6341-2008, 2008.

Miller, S. M., Wofsy, S. C., Michalak, A. M., Kort, E. A., Andrews, A. E., Biraud, S. C., Dlugokencky, E. J., Fischer, M. L., Janssens-Maenhout, G., Miller, B. R., Montzka, S. A., and Sweeney, C.: Anthropogenic emissions of methane in the United States, P. Natl. Acad. Sci. USA, 110, 20018–20022, https://doi.org/10.1073/pnas.1314392110, 2013.

National Research Council: Verifying Greenhouse Gas Emissions: Methods to Support International Climate Agreements, The National Academies Press, Washington D.C., https://doi.org/10.17226/12883, 2010.

NOAA: Methane (CH_{4}) WMO Scale, Scale Update: WMO X2004A (Updated July
7, 2015), https://gml.noaa.gov/ccl/ch4_scale.html, last
access: 10 November 2021.

O'Doherty, S., Say, D., Stanley, K., Spain, G., Arnold, T., Rennick, C., Young, D., Stavert, A., Grant, A., Ganesan, A., Pitt, J., Wisher, A., Wenger, A., and Garrard, N.: UK DECC (Deriving Emissions linked to Climate Change) Network. Centre for Environmental Data Analysis, http://catalogue.ceda.ac.uk/uuid/f5b38d1654d84b03ba79060746541e4f (last access: 12 April 2021), 2022.

Petrescu, A. M. R., McGrath, M. J., Andrew, R. M., Peylin, P., Peters, G. P., Ciais, P., Broquet, G., Tubiello, F. N., Gerbig, C., Pongratz, J., Janssens-Maenhout, G., Grassi, G., Nabuurs, G.-J., Regnier, P., Lauerwald, R., Kuhnert, M., Balkovič, J., Schelhaas, M.-J., Denier van der Gon, H. A. C., Solazzo, E., Qiu, C., Pilli, R., Konovalov, I. B., Houghton, R. A., Günther, D., Perugini, L., Crippa, M., Ganzenmüller, R., Luijkx, I. T., Smith, P., Munassar, S., Thompson, R. L., Conchedda, G., Monteil, G., Scholze, M., Karstens, U., Brockmann, P., and Dolman, A. J.: The consolidated European synthesis of CO_{2} emissions and removals for the European Union and United Kingdom: 1990–2018, Earth Syst. Sci. Data, 13, 2363–2406, https://doi.org/10.5194/essd-13-2363-2021, 2021a.

Petrescu, A. M. R., Qiu, C., Ciais, P., Thompson, R. L., Peylin, P., McGrath, M. J., Solazzo, E., Janssens-Maenhout, G., Tubiello, F. N., Bergamaschi, P., Brunner, D., Peters, G. P., Höglund-Isaksson, L., Regnier, P., Lauerwald, R., Bastviken, D., Tsuruta, A., Winiwarter, W., Patra, P. K., Kuhnert, M., Oreggioni, G. D., Crippa, M., Saunois, M., Perugini, L., Markkanen, T., Aalto, T., Groot Zwaaftink, C. D., Tian, H., Yao, Y., Wilson, C., Conchedda, G., Günther, D., Leip, A., Smith, P., Haussaire, J.-M., Leppänen, A., Manning, A. J., McNorton, J., Brockmann, P., and Dolman, A. J.: The consolidated European synthesis of CH_{4} and N_{2}O emissions for the European Union and United Kingdom: 1990–2017, Earth Syst. Sci. Data, 13, 2307–2362, https://doi.org/10.5194/essd-13-2307-2021, 2021b.

Pinty, B., Janssens-Maenhout, G., Dowell, M., Zunker, H., Brunhes, T.,
Ciais, P., Dee, D., Denier van der Gon, H., Dolman, H., Drinkwater, M.,
Engelen, R., Heimann, M., Holmlund, K., Husband, R., Kentarchos, A., Meijer,
Y., Palmer, P., and Scholze, M.: An Operational Anthropogenic CO_{2}
Emissions Monitoring & Verification Support capacity – Baseline
Requirements, Model Components and Functional Architecture, European
Commission Joint Research Centre, EUR 28736 EN, https://doi.org/10.2760/08644, 2017.

Pinty, B., Ciais, P., Dee, D., Dolman, H., Dowell, M., Engelen, R.,
Holmlund, K., Janssens-Maenhout, G., Meijer, Y., Palmer, P., M. Scholze,
Gon, H. D. v. d., Heimann, M., Juvyns, O., Kentarchos, A., and Zunker, H.:
An Operational Anthropogenic CO_{2} Emissions Monitoring & Verification
Support Capacity – Needs and high level requirements for in situ
measurements, European Commission Joint Research Centre, EUR 29817 EN, https://doi.org/10.2760/182790, 2019.

Pisso, I., Sollum, E., Grythe, H., Kristiansen, N. I., Cassiani, M., Eckhardt, S., Arnold, D., Morton, D., Thompson, R. L., Groot Zwaaftink, C. D., Evangeliou, N., Sodemann, H., Haimberger, L., Henne, S., Brunner, D., Burkhart, J. F., Fouilloux, A., Brioude, J., Philipp, A., Seibert, P., and Stohl, A.: The Lagrangian particle dispersion model FLEXPART version 10.4, Geosci. Model Dev., 12, 4955–4997, https://doi.org/10.5194/gmd-12-4955-2019, 2019.

Poulter, B., Bousquet, P., Canadell, J. G., Ciais, P., Peregon, A., Saunois, M., Arora, V. K., Beerling, D. J., Brovkin, V., Jones, C. D., Joos, F., Gedney, N., Ito, A., Kleinen, T., Koven, C. D., McDonald, K., Melton, J. R., Peng, C., Peng, S., Prigent, C., Schroeder, R., Riley, W. J., Saito, M., Spahni, R., Tian, H., Taylor, L., Viovy, N., Wilton, D., Wiltshire, A., Xu, X., Zhang, B., Zhang, Z., and Zhu, Q.: Global wetland contribution to 2000–2012 atmospheric methane growth rate dynamics, Environ. Res. Lett., 12, 094013, https://doi.org/10.1088/1748-9326/aa8391, 2017.

Prinn, R. G., Weiss, R. F., Fraser, P. J., Simmonds, P. G., Cunnold, D. M., Alyea, F. N., O'Doherty, S., Salameh, P., Miller, B. R., Huang, J., Wang, R. H. J., Hartley, D. E., Harth, C., Steele, L. P., Sturrock, G., Midgley, P. M., and McCulloch, A.: A history of chemically and radiatively important gases in air deduced from ALE/GAGE/AGAGE, J. Geophys. Res.-Atmos., 105, 17751–17792, https://doi.org/10.1029/2000jd900141, 2000.

Prinn, R., Weiss, R., Arduini, J., Arnold, T., DeWitt, H. L., Fraser, P., Ganesan, A., Gasore, J., Harth, C., Hermansen, O., Kim, J., Krummel, P., Li, S., Loh, Z., Lunder, C., Maione, M., Manning, A., Miller, B., Mitrevski, B., Mühle, J., O'Doherty, S., Park, S., Reimann, S., Rigby, M., Saito, T., Salameh, P., Schmidt, R., Simmonds, P., Steele, P., Vollmer, M., Wang, H. R., Yao, B., Young, D., and Zhou, L: In-situ measurements of chemically and radiatively important atmospheric gases from the Advanced Global Atmospheric Gas Experiment (AGAGE) and affiliated stations, The Advanced Global Atmospheric Gases Experiment (AGAGE), ESS-DIVE repository, https://doi.org/10.15485/1781803, 2021.

Rödenbeck, C., Gerbig, C., Trusilova, K., and Heimann, M.: A two-step scheme for high-resolution regional atmospheric trace gas inversions based on independent models, Atmos. Chem. Phys., 9, 5331–5342, https://doi.org/10.5194/acp-9-5331-2009, 2009.

Ruckstuhl, A. F., Henne, S., Reimann, S., Steinbacher, M., Vollmer, M. K., O'Doherty, S., Buchmann, B., and Hueglin, C.: Robust extraction of baseline signal of atmospheric trace species using local regression, Atmos. Meas. Tech., 5, 2613–2624, https://doi.org/10.5194/amt-5-2613-2012, 2012.

Saunois, M., Stavert, A. R., Poulter, B., Bousquet, P., Canadell, J. G., Jackson, R. B., Raymond, P. A., Dlugokencky, E. J., Houweling, S., Patra, P. K., Ciais, P., Arora, V. K., Bastviken, D., Bergamaschi, P., Blake, D. R., Brailsford, G., Bruhwiler, L., Carlson, K. M., Carrol, M., Castaldi, S., Chandra, N., Crevoisier, C., Crill, P. M., Covey, K., Curry, C. L., Etiope, G., Frankenberg, C., Gedney, N., Hegglin, M. I., Höglund-Isaksson, L., Hugelius, G., Ishizawa, M., Ito, A., Janssens-Maenhout, G., Jensen, K. M., Joos, F., Kleinen, T., Krummel, P. B., Langenfelds, R. L., Laruelle, G. G., Liu, L., Machida, T., Maksyutov, S., McDonald, K. C., McNorton, J., Miller, P. A., Melton, J. R., Morino, I., Müller, J., Murguia-Flores, F., Naik, V., Niwa, Y., Noce, S., O'Doherty, S., Parker, R. J., Peng, C., Peng, S., Peters, G. P., Prigent, C., Prinn, R., Ramonet, M., Regnier, P., Riley, W. J., Rosentreter, J. A., Segers, A., Simpson, I. J., Shi, H., Smith, S. J., Steele, L. P., Thornton, B. F., Tian, H., Tohjima, Y., Tubiello, F. N., Tsuruta, A., Viovy, N., Voulgarakis, A., Weber, T. S., van Weele, M., van der Werf, G. R., Weiss, R. F., Worthy, D., Wunch, D., Yin, Y., Yoshida, Y., Zhang, W., Zhang, Z., Zhao, Y., Zheng, B., Zhu, Q., Zhu, Q., and Zhuang, Q.: The Global Methane Budget 2000–2017, Earth Syst. Sci. Data, 12, 1561–1623, https://doi.org/10.5194/essd-12-1561-2020, 2020.

Seibert, P. and Frank, A.: Source-receptor matrix calculation with a Lagrangian particle dispersion model in backward mode, Atmos. Chem. Phys., 4, 51–63, https://doi.org/10.5194/acp-4-51-2004, 2004.

Shindell, D., Kuylenstierna, J. C. I., Vignati, E., van Dingenen, R., Amann, M., Klimont, Z., Anenberg, S. C., Muller, N., Janssens-Maenhout, G., Raes, F., Schwartz, J., Faluvegi, G., Pozzoli, L., Kupiainen, K., Höglund-Isaksson, L., Emberson, L., Streets, D., Ramanathan, V., Hicks, K., Oanh, N. T. K., Milly, G., Williams, M., Demkine, V., and Fowler, D.: Simultaneously Mitigating Near-Term Climate Change and Improving Human Health and Food Security, Science, 335, 183–189, https://doi.org/10.1126/science.1210026, 2012.

Shindell, D., Borgford-Parnell, N., Brauer, M., Haines, A., Kuylenstierna, J. C. I., Leonard, S. A., Ramanathan, V., Ravishankara, A., Amann, M., and Srivastava, L.: A climate policy pathway for near- and long-term benefits, 356, 493–494, https://doi.org/10.1126/science.aak9521, 2017.

Stavert, A. R., Saunois, M., Canadell, J. G., Poulter, B., Jackson, R. B., Regnier, P., Lauerwald, R., Raymond, P. A., Allen, G. H., Patra, P. K., Bergamaschi, P., Bousquet, P., Chandra, N., Ciais, P., Gustafson, A., Ishizawa, M., Ito, A., Kleinen, T., Maksyutov, S., McNorton, J., Melton, J. R., Müller, J., Niwa, Y., Peng, S., Riley, W. J., Segers, A., Tian, H., Tsuruta, A., Yin, Y., Zhang, Z., Zheng, B., and Zhuang, Q.: Regional trends and drivers of the global methane budget, Glob. Change Biol., 28, 182–200, https://doi.org/10.1111/gcb.15901, 2021.

Stohl, A., Forster, C., Frank, A., Seibert, P., and Wotawa, G.: Technical note: The Lagrangian particle dispersion model FLEXPART version 6.2, Atmos. Chem. Phys., 5, 2461–2474, https://doi.org/10.5194/acp-5-2461-2005, 2005.

Szopa, S., Naik, V., Adhikary, B., Artaxo, P., Berntsen, T., Collins, W. D., Fuzzi, S., Gallardo, L., Kiendler-Scharr, A., Klimont, Z., Liao, H., Unger, N., and Zanis, P.: Short-Lived Climate Forcers. In Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Masson-Delmotte, V., Zhai, P., Pirani, A., Connors, S. L., Péan, C., Berger, S., Caud, N., Chen, Y., Goldfarb, L., Gomis,M. I., Huang, M., Leitzell, K., Lonnoy, E., Matthews, J. B. R., Maycock, T. K., Waterfield, T., Yelekçi, O., Yu, R., and Zhou, B., Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA, 817–922, 2021.

Talagrand, O. and Courtier, P.: Variational assimilation of meteorological observations with the adjoint vorticity equation. I: Theory, Q. J. Roy. Meteor. Soc., 113, 1311–1328, 1987.

UNFCCC: National Inventory Submissions 2021, https://unfccc.int/ghg-inventories-annex-i-parties/2021, last access: 23 November 2021.

United Nations Environment Programme and Climate and Clean Air Coalition: Global Methane Assessment: Benefits and Costs of Mitigating Methane Emissions, Nairobi, United Nations Environment Programme, ISBN: 978-92-807-3854-4, 2021.

van der Werf, G. R., Randerson, J. T., Giglio, L., van Leeuwen, T. T., Chen, Y., Rogers, B. M., Mu, M., van Marle, M. J. E., Morton, D. C., Collatz, G. J., Yokelson, R. J., and Kasibhatla, P. S.: Global fire emissions estimates during 1997–2016, Earth Syst. Sci. Data, 9, 697–720, https://doi.org/10.5194/essd-9-697-2017, 2017.

VERIFY: VERIFY – FactSheets, http://webportals.ipsl.jussieu.fr/VERIFY/FactSheets/, last access: 23 November 2021.

WMO: WMO Greenhouse Gas Bulletin (GHG Bulletin) – No.17: The State of Greenhouse Gases in the Atmosphere Based on Global Observations through 2020, World Meteorological Organization, Geneva, 10, 2021.

Yver-Kwok, C., Philippon, C., Bergamaschi, P., Biermann, T., Calzolari, F., Chen, H., Conil, S., Cristofanelli, P., Delmotte, M., Hatakka, J., Heliasz, M., Hermansen, O., Komínková, K., Kubistin, D., Kumps, N., Laurent, O., Laurila, T., Lehner, I., Levula, J., Lindauer, M., Lopez, M., Mammarella, I., Manca, G., Marklund, P., Metzger, J.-M., Mölder, M., Platt, S. M., Ramonet, M., Rivier, L., Scheeren, B., Sha, M. K., Smith, P., Steinbacher, M., Vítková, G., and Wyss, S.: Evaluation and optimization of ICOS atmosphere station data as part of the labeling process, Atmos. Meas. Tech., 14, 89–116, https://doi.org/10.5194/amt-14-89-2021, 2021.

_{4}emissions in 2018. The new system combines a high spatial resolution of 7 km x 7 km with a variational data assimilation technique, which allows CH

_{4}emissions to be optimized from individual model grid cells. The high resolution allows the observations to be better reproduced, while the derived emissions show overall good consistency with two existing models.