Estimating power plant CO₂ emission using OCO-2 XCO₂ and high resolution WRF-Chem simulations

Our ability to predict future climate change requires an improved understanding of carbon dioxide (CO₂) emissions and uptake. Top-down inversion approaches that couple atmospheric observations of CO₂ with chemistry transport models have been widely used to estimate biospheric and oceanic CO₂ flux, with the fossil fuel emission often prescribed as a known quantity [2]. In such approaches, errors in the prescribed fossil fuel emission will lead to biases in the inversely estimated biospheric and oceanic fluxes [3–5], and large uncertainties exist in the fossil fuel emission inventory data [6, 7].

For instance, inverse modeling of anthropogenic CO₂ emissions has been applied to large cities [8–10], and power plants [1, 11]. These studies can be broadly grouped into two categories. The first category solves for long term average fossil fuel CO₂ emission by assimilating in situ CO₂ mixing ratio measurements. For instance, Nickless et al [10] solved for weekly average emissions on a 101 × 101 grid over Cape Town, South Africa, and Lauvaux et al [8] solved for 5-day average emissions on a 87 × 87 grid over Indianapolis, IN, USA.

The second category solves for short term emission rates through assimilating remotely sensed column average CO₂ (XCO₂) measurements. This category includes the Observation System Simulation Experiments (OSSE) studies by Pillai et al [12] and Broquet et al [13] using the simulated CarbonSat XCO₂ data [14], by Nassar et al [1] using OCO-2 [15] XCO₂, and by Krings et al [11] using airborne XCO₂ measurements. These studies utilize the spatial coverage of XCO₂ to distinguish the CO₂ plume caused by localized anthropogenic sources from the background. This is done either by identifying the background XCO₂ before extracting the plume enhancement [1, 11], or through simultaneous optimization of anthropogenic emissions and the background XCO₂ [12, 13].

Because XCO₂ enhancement from localized anthropogenic sources (large cities and power plants) decreases rapidly with increasing pixel size [14], high spatial resolution imaging is critical for estimating the emission rates of these sources. Although OCO-2 is not an optimal platform for estimating localized anthropogenic emissions due to large gaps between its narrow swaths, a few studies have shown it is possible in certain situations [1, 16].

The objective of this manuscript is to estimate the power plant emission rates through OCO-2 XCO₂ data and high resolution WRF-Chem simulations [17, 18]. Our method is similar to that of Pillai et al [12] and Broquet et al [13], but differs in a key aspect: instead of Bayesian inversion, we use a maximum likelihood estimation, which obtains power plant CO₂ emissions from OCO-2 XCO₂ data alone without prior emission information.

In this section, we first describe the WRF-Chem modeling system, followed by the XCO₂ observation operator, Jacobian matrix calculation, OCO-2 XCO₂ data, and finally the maximum likelihood estimation approach.

2.1. WRF-Chem modeling system

We use WRF-Chem (version 3.6.1) [18] as the atmospheric transport model. For each inversion case, WRF-Chem is set up with one-way four-level nested domains with horizontal grid spacing of 36, 12, 4, and 1.33 km from the outermost to the innermost domains. Each inner domain is located approximately at the center of its parent domain. The innermost domain (d04) has 97 grid points in the north-south direction and 61 grid points in the east-west direction, and is centered approximately at the power plant location for each inversion case.

For all WRF-Chem simulations, all four of the model domains have 52 vertical levels between the surface and the model top of 100 hPa. For a column where the surface is at the sea level, 19 of the 52 vertical levels are below 2 km. WRF physics packages used for the simulations are: WSM3 microphysics [19], Dudhia shortwave [20] and RRTM longwave radiation [21], the Unified Noah Land Surface Model, the Kain–Fritsch cumulus scheme [22] (only applied to the 36 and 12 km domains), and ACM2 boundary layer scheme [23]. Positive-definite transport is applied to the advection of scalars and CO₂. For each inversion case, the WRF-Chem simulation runs for 36 h, including a 12 h spin-up. Simulation results are saved at 15-min intervals. The simulation result that is closest to the OCO-2 overpass time for a given case is used for its inversion. The WRF-Chem simulation start time, time used for inversion, and OCO-2 overpass time for the seven inversion cases are given in table 1.

Table 1.  OCO-2 data, power plants, and WRF-Chem simulation information. Number of footprints is of OCO-2 v9 data. The numbers from v7 and v8 are slightly different.

Site	Date	OCO-2 overpass time	Simulation used for inversion	State vector	Number of OCO-2 footprints	Power plant
Kansas	2015-12-04	19:33	19:30	${[{e}_{1}b]}^{T}$	444	Westar
Ohio	2015-07-30	18:37	18:30	${[{e}_{1}b]}^{T}$	322	Gavin/Kyger
Kentucky	2015-08-13	18:50	18:45	${[{e}_{1}b]}^{T}$	377	Ghent
India	2014-10-23	07:38	07:45	${[{e}_{1}{e}_{2}b]}^{T}$	439	Sasan
						Vindhyachal
India	2014-11-10	07:50	07:45	${[{e}_{1}{e}_{2}b]}^{T}$	354	Sasan
						Vindhyachal
South Africa	2014-11-07	11:50	11:45	${[{e}_{1}b]}^{T}$	458	Matimba
South Africa	2016-10-11	11:50	11:45	${[{e}_{1}b]}^{T}$	424	Matimba

Meteorological initial and lateral boundary conditions are from the National Centers for Environmental Prediction (NCEP) Climate Forecast System Version 2 (CFSv2) 1° × 1° six hourly products [24]. CO₂ initial and lateral boundary conditions are from the CarbonTracker version 2017 (CT2017) [25] posterior mole fraction dataset. The CT2017 1° × 1° three hourly posterior fluxes for biosphere, fossil fuel, biomass burning, and ocean are horizontally interpolated to the WRF-Chem domains using bicubic interpolation. They are processed into three hourly emissions files which are ingested by the WRF-Chem model during its integration. In our simulations, CO₂ is treated as a chemically inert tracer and its uptakes by the biosphere is through using the CarbonTracker biospheric fluxes.

2.2. OCO-2 XCO₂ observation operator

WRF-Chem simulated three-dimensional CO₂ mixing ratios are used by the observation operator to calculate the modeled XCO₂. We implement the XCO₂ observation operator following Connor et al [26] and O'Dell et al [27] as:

$$\begin{eqnarray}&&{X}_{{\mathrm{CO}}_{2}}^{{\rm{m}}}={X}_{{\mathrm{CO}}_{2}}^{{\rm{a}}}+\displaystyle \sum _{j}{h}_{j}\,{a}_{{\mathrm{CO}}_{2},j}{\left({x}_{{\rm{m}}}-{x}_{{\rm{a}}}\right)}_{j},\end{eqnarray} \tag{ 1 }$$

where ${X}_{{{\rm{CO}}}_{{\rm{2}}}}^{m}$ is the modeled XCO₂, ${X}_{{{\rm{CO}}}_{2}}^{a}$ is the prior XCO₂, ${a}_{{{\rm{CO}}}_{2},j}$ are elements of the column averaging kernel, h_j is the pressure weighting function, x_a is the prior CO₂ vertical profile, and x_m is the WRF-Chem simulated CO₂ profile interpolated to the OCO-2 XCO₂ pressure levels assuming CO₂ varies linearly with atmospheric pressure [27]. The prior CO₂ profile, prior XCO₂, column averaging kernel, and pressure weighting function are taken directly from the OCO-2 XCO₂ data. An analysis of the column averaging kernel on emission estimations is given in section 1 of the supplementary material available online at stacks.iop.org/ERL/14/085001/mmedia.

2.3. Jacobian matrix $K$ of the forward model

The forward model f(x), consists of the WRF-Chem and the observation operator, and relates power plant emission rates to the modeled XCO₂. The forward model's Jacobian matrix K quantifies the sensitivity of the modeled XCO₂ to the variation in power plant emission rates. We calculate ${\bf{K}}$ using a one-sided finite difference:

$$\begin{eqnarray}&&\displaystyle \frac{\partial f}{\partial {\bf{e}}}=\displaystyle \frac{f({\bf{e}}+{\rm{\Delta }}{\bf{e}})-f({\bf{e}})}{{\rm{\Delta }}{\bf{e}}}.\end{eqnarray} \tag{ 2 }$$

Here, f(e) is a base run and f(e + Δe) is a perturbed run of WRF-Chem. A base run is a 24 h simulation with all its CO₂ flux prescribed using the three hourly CT2017 CO₂ fluxes as described in section 2.2. A perturbed run is identical to its base run counterpart, except that a small extra CO₂ emission ${\rm{\Delta }}e$ is added to the WRF-Chem grid at the location of a given power plant. Note that while all the CT2017 based CO₂ fluxes are placed at the lowest level of the WRF-Chem grid, the perturbation Δe for power plant emissions is added to the vertical level which corresponds to the smoke stack height. Multiple stacks of a given power plant are treated as a single point source (see Nassar et al [1] for the details of the smoke stack heights). We modified the WRF-Chem code to add the power plant emissions to the appropriate vertical levels based on the model's terrain-following hydrostatic-pressure vertical coordinate. The details of this modification are given in section 6 of the supplementary material.

2.4. Maximum likelihood estimation

To obtain power plant emission values from the OCO-2 XCO₂ data alone without prior information, we use the maximum likelihood approach. A maximum likelihood estimation is obtained by minimizing the cost function:

$$\begin{eqnarray}&&J({\bf{x}})={\left({\bf{Hx}}-{\bf{y}}\right)}^{{\rm{T}}}{{\bf{R}}}^{-1}({\bf{Hx}}-{\bf{y}}),\end{eqnarray} \tag{ 3 }$$

where the state vector ${\bf{x}}={[{\bf{e}}b]}^{{\rm{T}}}$ consists of two components: the power plant emission rates vector e (section 2.3), and the background XCO₂ b. The background XCO₂ is similar to that used in Pillai et al [12] and Broquet et al [13]. y is the observation vector whose elements are XCO₂ values from OCO-2 footprints. In this work we use three versions of XCO₂ data: v7 [28], v8 [29], and v9 [30]. H relates the state vector x to the observation vector y, and is formed by the forward model Jacobian matrix K (section 2.3) and the unity:

$$\begin{eqnarray}{\bf{H}}=\left[\begin{array}{cc}{\bf{K}} & 0\\ 0 & 1\end{array}\right]\end{eqnarray} \tag{ 4 }$$

where the unity value represents the background XCO₂ influence on the modeled XCO₂. R is the error covariance matrix (section 2.5). When the model-observation residual error is assumed to follow a Gaussian distribution, the maximum likelihood estimation $\hat{{\bf{x}}}$ is obtained as

$$\begin{eqnarray}&&\hat{{\bf{x}}}={\left({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{H}}\right)}^{-1}({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{y}}),\end{eqnarray} \tag{ 5 }$$

and its covariance matrix A is

$$\begin{eqnarray}&&{\bf{A}}={\left({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{H}}\right)}^{-1}.\end{eqnarray} \tag{ 6 }$$

Section 2 of the supplementary material includes an analysis of the residual error of the maximum likelihood estimation.

2.5. Error covariance matrix $R$

The error covariance matrix R should include observational uncertainties and transport model uncertainties:

$$\begin{eqnarray}&&{\bf{R}}={{\bf{R}}}_{{\rm{obs}}}+{{\bf{R}}}_{{\rm{model}}}.\end{eqnarray} \tag{ 7 }$$

We set the observation uncertainty covariance matrix R_obs as

$$\begin{eqnarray}{R}_{\mathrm{obs}}[i,j]=\left\{\begin{array}{lr}\mathrm{corr}(i,j){\sigma }_{{\mathrm{xco}}_{2},i}\,{\sigma }_{{\mathrm{xco}}_{2},j} & i\ne j\\ {\sigma }_{{\mathrm{xco}}_{2},i}\,{\sigma }_{{\mathrm{xco}}_{2},j} & i=j\end{array}\right.,\end{eqnarray} \tag{ 8 }$$

where ${\sigma }_{{{\rm{XCO}}}_{2}}$ is the XCO₂ posterior error (uncertainty) provided by the OCO-2 data. The correlation between the uncertainty of two OCO-2 footprints is calculated using

$$\begin{eqnarray}&&\mathrm{corr}(i,j)=\displaystyle \frac{{t}_{\max }-{t}_{i,j}}{{t}_{\max }}\times {\mathrm{corr}}_{\max },\end{eqnarray} \tag{ 9 }$$

where t_i,j is the time difference between two individual OCO-2 footprints, t_max is the maximal time difference, and ${\mathrm{corr}}_{\max }$ is the maximal error correlation. This configuration sets the XCO₂ error correlation to corr_max for two footprints adjacent in time, and linearly scales to zero when the time difference increases to t_max. We set t_max = 10.0 s and corr_max = 0.3 based on Worden et al [31].

Transport model uncertainty R is much more difficult to quantify [32, 33], so we adopt a simple approach used in Deng et al [34] of inflating the observation error.

$$\begin{eqnarray}&&{\bf{R}}=c\times {{\bf{R}}}_{\mathrm{obs}},\end{eqnarray} \tag{ 10 }$$

where c is a scaling factor. We set c = 1.5 and examined its impacts on the estimation uncertainty in section 5.1 of the supplementary material.

Nassar et al [1] identified seven OCO-2 overpasses that captured emission plumes from middle-to-large-sized coal power plants. We apply the method described in section 2 to these cases to estimate the power plants' daily mean CO₂ emission rate.

When equation (5) is used for the estimation, there are different choices for H, y, and R. First, the K component of H represents the forward model Jacobian matrix. Because the WRF-Chem CO₂ transport process is not linear, the value of K varies with the choice of emission vector e when equation (2) is used. Second, there are different choices for the observation vector y from different OCO-2 data versions. Third, R varies with the scaling factor c of equation (10).

Because all of these changes will influence the maximum likelihood estimation results, we assess the emission estimation sensitivity to the above variation in H, y, and R. For the sensitivity to transport model non-linearity, we calculate the K (and thus H) for a set of four emission vectors e: for cases with one power plant (table 1), e is set to [1], [25], [50], and [100]; for the cases with two power plants, e is set to ${[11]}^{{\rm{T}}}$, ${[2525]}^{{\rm{T}}},{[5050]}^{{\rm{T}}}$, and ${[100100]}^{{\rm{T}}}$. We note that calculating K at different emission vectors e in our maximum likelihood estimation is similar to choosing different emission priors for a Bayesian inversion [12, 13]. For a Bayesian inversion that seeks to minimize the cost function $J{({\bf{x}})={({\bf{Hx}}-{\bf{y}})}^{{\rm{T}}}{{\bf{R}}}^{-1}({\bf{Hx}}-{\bf{y}})+({\bf{x}}-{{\bf{x}}}^{b})}^{{\rm{T}}}{{\bf{B}}}^{-1}({\bf{x}}-{{\bf{x}}}^{b})$, its posterior is ${{\bf{x}}}^{a}={({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{H}}+{{\bf{B}}}^{-1})}^{-1}({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{y}}+{{\bf{B}}}^{-1}{{\bf{x}}}^{b})$ [35] . When the prior uncertainty covariance matrix B is set to infinity (or very large values), the prior x^b becomes non-informative, and the posterior x^a is numerically equivalent to the maximum likelihood estimate $\hat{{\bf{x}}}={({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{H}})}^{-1}({{\bf{H}}}^{{\rm{T}}}{{\bf{R}}}^{-1}{\bf{y}})$ [36]. In this regard, using K calculated at different e is similar to using different prior emissions in a Bayesian estimation. This results in a set of four H for each of the seven cases. For the sensitivity to the variation in observation data, we use the three latest versions (v7, v8, and v9) of OCO-2 XCO₂ to form the observation vector y. Between the v7 and v8 data, there are a large number of changes in the OCO-2 XCO₂ retrieval algorithm related to spectroscopy, aerosol treatment, the surface model, and prior meteorology [29]. The differences between the v8 and v9 data are mainly from an update in surface pressure estimation [30]. The three versions of OCO-2 XCO₂ used for the seven cases are shown in figures 1–7 of the supplementary material. This results in a set of three y for each of the seven cases. With the scaling factor set as c = 1.5, the above combinations of the different H and y results in 12 estimation results for each case. We list these estimations for the seven cases in tables 2–4.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The 2015-12-04 OCO-2 observation over the Westar Jeffrey Energy Center in Kansas has ideal conditions for emission estimation: the wind direction is aligned with the OCO-2 orbit track and the OCO-2 data were acquired in the nadir mode. Figure 1(a) shows the outermost (d01) WRF-Chem simulation domain for the Kansas case, where the three red color boxes represent the three nested inner domains (d02, d03, and d04), respectively. Figure 1(b) is a zoomed-in view of the 1.33 × 1.33 km innermost domain. There are 444 OCO-2 footprints within this innermost domain and they are used for the estimation (the observation vector y has 444 elements). To show more details, figures 1(c), (d), and (e) are plotted over the the area of the innermost domain around the power plant and its CO₂ plume. The spatial range of these three figures are marked by the red colored box in figure 1(b). Figure 1(c) shows the modeled XCO₂ calculated by the forward model (WRF-Chem and observation operator), with horizontal wind fields at the power plant smoke stack height level superimposed. Figure 1(e) shows the sensitivity of the modeled XCO₂ to the variation of power plant emission rate, which is represented by the forward model Jacobian matrix K. As explained in the last paragraph, we calculate K at a set of four emission rates, and figures 1(c) and (e) are the results from using e = [25] ktCO₂/day. We note that the modeled XCO₂ (as depicted in figure 1(c)) is not used by the maximum likelihood estimation (equation (5)), but K is used (as part of H matrix).

A visual comparison of figures 1(c)–(e) shows that the WRF-Chem wind field (figure 1(c)) and the resulting advection generated a power plant CO₂ plume that matches well with the OCO-2 observed plume (figure 1(d)). Our estimation results (table 2) put the Westar daily mean emission rates at 32.48–33.87 ktCO₂/day using v7 data, 26.23–27.27 ktCO₂/day using v8 data, and 22.42–23.39 ktCO₂/day using v9 data. In comparison, the reported emission value at Westar is 26.67 ktCO₂/day. Using the v7 data, Nassar et al [1] estimated the Westar emission rate at 31.21 ± 3.71 ktCO₂/day, which is close to our estimation using v7 data. The small difference (<10°) between the simulated wind direction and the OCO-2 plume could contribute to the emission estimation error.

Among the estimations using the three OCO-2 data versions, there exists a fair amount of differences. Because all estimations are based on the same set of WRF-Chem simulations, the differences in the estimations are attributed to the differences among different OCO-2 data versions, including the differences in XCO₂ values, uncertainty, retrieval pressure levels, pressure weighting function, and averaging kernel. For the results using each of the OCO-2 data versions, there are only small differences (∼1 ktCO₂/day) in the estimation, indicating our estimation approach is only weakly sensitive to the transport model non-linearity. In the Kentucky 2015-08-13 case, the OCO-2 track is about 8 km west of the Ghent power plant (figure 2). Our estimated emission values (9.84–20.40 ktCO₂/day) are substantially lower than the reported value of 29.17 ktCO₂/day. These underestimations are primarily caused by the complex and inaccurate horizontal wind field simulation by WRF-Chem. While the plume was located largely north of the river (figure 2(d)), the WRF-Chem simulations put the plume south of the river (figures 2(c) and (e)). In comparison, Nassar et al [1] achieved an estimation of 29.46 ± 15.58 ktCO₂. Their estimation applied a 22.5° wind direction adjustment to optimize the match between the Gaussian plume model and the plume detected by the OCO-2 data. We note that the simulated low wind speed (∼2 m s⁻¹) in this case is more likely to be associated with variability in wind direction, as smaller scale features that are difficult to resolve exert greater influence in the absence of stronger winds.

In the Ohio 2015-07-30 case, two power plants, Gavin and Kyger Creek, are located within 2.5 km of each other (figure 3). We initially attempted to estimate the emission of the two power plants separately, but the resulting emissions show large oscillations with K calculated using different emission vectors e (results not shown). This is most likely because the K vector of the two power plants largely overlap, causing the estimation to be under-constrained. We redesigned the system to estimate the combined emission rates of two plants by using a single control variable in the state vector for both the plants. However, our estimations are still substantially lower than the reported value (50.54 ktCO₂/day): 13.49–14.17 ktCO₂/day using the v7 data, 21.18–22.29 ktCO₂/day using the v8 data, and 30.95–32.52 ktCO₂/day using the v9 data. That the v9 data results in the closest estimation to the reported value may indicate the improved data quality in this latest OCO-2 data version. Similar to the Kansas case, the spread in the emission estimation among those using the same OCO-2 data version are small (∼1 ktCO₂/day). The estimation of Nassar et al [1] using the v7 data is 48.66 ktCO₂/day, substantially higher than our results using v7 (13.49–14.17 ktCO₂/day). There are a large number of OCO-2 footprints with elevated XCO₂ values south of the power plant plume in the v7 data, which are not in the v8 or v9 data (figure 3 of the supplementary material). These elevated XCO₂ values contribute to our higher background estimation (398.03 ppm) than that of Nassar et al [1] (396.64 ppm), thus leading to our lower emission estimation.

There are two OCO-2 overpasses (2014-10-23 and 2014-11-10, figures 4 and 5) over the Indian site, which has five power plants in its vicinity (figure 4(b)). Like Nassar et al [1], we found that only emissions from the Sasan and Vindhyachal power plants were captured by the OCO-2 observations. Accordingly, the state vector is [e₁, e₂, b]^T with e₁ and e₂ for the CO₂ emission rates of Sasan and Vindyachal, respectively (table 1). From the 2014-10-23 data (figure 4), our estimations for the Sasan plant are 43.56–46.27 ktCO₂/day using the v7 data, 42.13–44.59 ktCO₂/day using the v8 data, and 46.78–49.37 ktCO₂/day using the v9 data. Our estimations for the Vindhyacal plant are 85.90–90.53, 66.75–70.07, and 68.54–72.35 ktCO₂/day using the three data versions, respectively. Overall, these estimations underestimate for Sasan (reported value 66.23) but fare better for Vindhyachal (reported value 88.76 ktCO₂/day). In comparison, Nassar et al [1] use the reported emission value of Vindhyachal to remove its plume from the OCO-2 XCO₂ before estimating the emission rate of the Sasan plant. Their study estimated the Sasan plant emission to be 63.93 ± 9.98 ktCO₂/day, and attributed the large estimation uncertainty to the CO₂ from the Vindhyachal plant. The ability of our method to estimate emissions of multiple power plants is an advantage that can be expanded to estimate urban emissions in the future. The India 2014-11-10 observation was acquired in glint mode and with wind largely perpendicular to the OCO-2 ground track (figure 5(c)). Both factors limit the observation data constraints for emission estimation. Our estimations from the 2014-11-10 data are lower than from the 2014-10-23 nadir mode data: 19.62–21.10, 34.15–37.59, and 30.68–34.19 ktCO₂/day for the Sasan plant and 50.45–53.15, 51.86–53.27, and 51.05–52.78 ktCO₂/day for the Vindhyachal plant using the three data versions.

There are two glint overpasses (2014-11-07 and 2016-10-11) over the Matimba power plant in South Africa. Our estimation for the 2014-11-07 case (figure 6) are quite similar among the results using the three data versions: 39.98–41.58 (v7), 39.69–42.4 (v8), and 40.85–42.35 (v9) ktCO₂/day. Our estimation for the 2016-10-11 case are close to that for the 2014-11-07 case: 39.98–40.73 (v7), 46.37–47.04 (v8), and 40.47–41.55 (v9) ktCO₂/day. While these estimations are lower than the reported value of 66.25 ktCO₂/day, they are consistent with each other. In addition, the estimations from Nassar et al [1] are also substantially lower than the reported value in both cases: 33.05 ± 10.57 and 33.66 ± 3.42, respectively.

We developed and tested a method for estimating power plant emissions using OCO-2 XCO₂ data and WRF-Chem simulations. This method is similar to the Bayesian inversion approach used by Pillai et al [12] and Broquet et al [13], in that both use a high resolution transport model to link emissions to XCO₂ observations. However, we use the maximum likelihood approach to estimate power plant emissions directly from the XCO₂ without prior power plant emission information. In this regard, our method is similar to that of Nassar et al [1], but we use WRF-Chem instead of the Gaussian plume model for CO₂ transport simulation. Compared with Nassar et al [1], our method provides a formal framework to address dynamic meteorological fields and possible expansion to spatially more complex anthropogenic emission sources like urban areas. One disadvantage of our method is that it is less tolerant of observation data errors as its emission results using OCO-2 v7 data are less accurate compared with those of the plume model approach, which can avoid the visually detectable observation data error through separating background and plume enhancement footprints. It should also be pointed out that the narrow swath (≤10.3 km) and column averaged measurements (XCO₂) of OCO-2 pose challenges for accurate emission estimations.

To effectively utilize OCO-2's fine footprint and minimize representation errors, we use 1.33 ×1.33 km horizontal grid spacing in the innermost domain of the WRF-Chem simulations to calculate the Jacobian matrix and modeled XCO₂. However, high resolution modeling does not necessarily increase simulation accuracy, as demonstrated in Feng et al [37]. Following Deng et al [34], we inflate the OCO-2 XCO₂ observation errors to account for the transport model errors. A comprehensive treatment would require much more rigorous analyses, such as those conducted in [32, 33, 38, 39]. Although a rigorous quantification of transport model error will give a more realistic representation of the true uncertainty of the posterior estimates, a complementary approach is used to improve the transport model accuracy. Such approaches include establishing a suitable WRF-Chem model configuration for each particular site through validations against meteorological measurements [40], assimilating meteorological observations during simulations [41, 42], and using ensemble transport model simulations. In addition, variational assimilation approaches can be used to address the sensitivity of the estimation to the transport model non-linearity [43, 44].

We have shown that the maximum likelihood estimation using high resolution WRF-Chem simulation is a viable method for quantifying power plant CO₂ emissions from OCO-2 under certain conditions. The ability to quantify CO₂ emissions using existing and future spaceborne XCO₂ observations at the scale of an individual facility, an urban area or at the national level is currently an area of significant scientific interest, but also has important implications for future climate change mitigation policies including the transparency framework of the United Nations Paris Agreement. Efforts are currently underway by Europe [45, 46] and multiple individual nations to understand the requirements and thus design spaceborne systems for anthropogenic CO₂ monitoring. The Committee on Earth Observation Satellites (CEOS) is reviewing the state of science in the field, assessing the technical requirements for monitoring, and advising space agencies on the way forward to develop an internationally coordinated constellation for this objective [47]. While at present, the observational coverage of space-based XCO₂ data is the greatest limitation prohibiting routine CO₂ emission monitoring, with recent new missions and more to come in the future as space agencies work toward establishing the system envisioned in the CEOS report [47], it is important for modeling methods to advance in parallel to demonstrate the value of the improved observations and make optimal use of them when they become available.

We thank NASA and the OCO-2 project for making OCO-2 data available for this work. We thank NOAA for making the NCEP/CFSR and CarbonTracker data available for the public. We thank the WRF-Chem development team for making the code available for the public. We thank Dr Chin-I Cheng of the Central Michigan University Statistical Consulting Center for her helpful discussion. This work was partially supported by a Central Michigan University Research Incentive Fund. This is contribution 124 of the Central Michigan University Institute for Great Lakes Research.