The manuscript by Chen et al. presents the estimate of methane emissions for 2019 in China by using the TROPOMI data, an inversion method, and GEOS-Chem model simulations. Overall, paper is interesting and could be published after the following concerns are addressed. 

1. Uncertainty in the abstract and the text. Should the uncertainty range be centered at the best estimate?  In some cases, it is; CH4 from coal is 16.6 with uncertainty range of 15.6-17.6. In some case, it is not; total anthropogenic emission is 65 with uncertainty range of 57.7 - 68.4. Please explain why it is so, and how the uncertainty is defined? Is it because the errors from the prior emission are assumed to be log-normal? but this still doesn't explain why in some cases, the best estimate is indeed centered at the uncertainty range. 
2. The motivation of this paper should be strengthened in the introduction section. It appears that Qu et al (2021) already done global estimate of CH4 emission. How does the results from that paper compare with bottom-up estimate in China and motivate this work?  Is it because Qu et al didn't do sector attribution at fine spatial resolution? More articulation above L90 is needed here.
3. Bias in boundary conditions (text around L200). While the importance of the boundary conditions is recognized here, does the optimization here attribute the innovation to the emissions within China only or also emissions in elsewhere as well? How are the boundaries defined and will the results be sensitive to the location of boundaries? Methane is long lived, and model-observation difference in CH4 can be attributed to the upper wind regions that are far from the domain of interest, even at the seasonal scale.
4. Section 2.4. Does the state vector include sector resolved emissions or simply spatially resolved total emissions? 'total' could be added here. Why number of Gaussian functions is 600 ?  Later, 'the inversion optimizes the emissions for each Gaussian function along with 16 boundary conditions ... for a total of .. 616 elements".  Should the inversion optimize the state vector x? Consider to add 'in the state vector' at the end of this sentence. Also, should all three Gaussian parameters be optimized, and so there would be 3*600 +16 = 1816 elements? Otherwise, will it be possible that the new source of CH4 in the regions where the prior emissions are low can not be accurately located in the domain with GMM? The fine-resolution TROPOMI data has the potential to accurately identify the local new sources. GMM appears to have the limitation to locate a new source spatially in areas that the prior emission informs the possibility for large-spatial aggregation. True?
5. Aggregation vs. smoothness. From Figure 4c, it appears there are still distinct boundaries in post/prior ratios, which could suggest an artifact as a result of aggregation. I would recommend to plot it with a different color bar that has more smooth transition of the colors and discuss further if there are still distinct boundaries. 
6. Gaussian representation. 'a mean location, spatial standard deviation, and emission magnitude'. How are those parameters defined by using the a prori emission? What does the mean location mean? Is it emission-weighted geographic center? 
7. L240. 'the relationship between methane emissions and concentrations in the nested GEOS-Chem simulation is linear'. Is this an assumption? In what time scale? I would think even there is no methane emission in China, global transport would still result in the methane in the air in China.
8. Section 2.5.  The sectoral attribution. How is the summation matrix W defined?  In addition, from what can tell here, it seems that the attribution yields the single scaling factor for the national emission from each sector. In other words, it didn’t spatially resolve the scaling factor for each sector in every aggregated area. will it make sense that the sectoral attribution is made in proportion of the relative weights of each sector from the prior emission in every aggregated area?  Description of the assumption or reasons for the feasibility of the attribution method is needed here. 
9. Section 4.1. I could miss something here. but, is there any map showing the change of emission by sector? so it can support that 'our downward correction of .... is driven by the Shanxi province and southwestern China'? This is also relevant to my comments in #4 and #8 as well as the abstract - 'our higher livestock emissions are attributed ... to northern China'. A map showing the posterior/prior ratio for each sector can be useful in the main text. 

Minor.

In abstract. mention the year for the prior emission. 

The manuscript by Chen et al. presents valuable and timely estimates of sector-specific CH4 emissions in China inferred from TROPOMI data using an inversion with the GEOS-Chem model.  Overall, the topic, scope, and presentation are well suited for ACP, and the results are of value to the community.  The paper emphasizes the importance of using log-normal distributions to describe CH4 emissions, and so I focused on this a bit in my review.  In this regard, though, I’m concerned about one issue described in detail below.  Given the importance of this part of their approach, as it related to identification of sector-specific emissions, I think it warrants thoughtful revisions. I’ve provided comments on some other points as well, again focusing mostly on aspects of the inversion method. 

Major comments: 

I do have a question about the overall framework as presented here, and in prior work from this group applying the same approach.  To start with, I would appreciate seeing what equation is being minimized by Equation 5.  In other words, can the cost function for the non-Gaussian case be written out, explicitly?  There is a bit of literature on this topic, and results one could adopt readily, starting from this point, but without seeing what cost function they are actually minimizing, I’m inferring a bit.   

It does though appear that equation (5) is derived from minimization of a cost function based on ln(x) but not the log of y.  This is potentially concerning for the following reason.  The authors claim that they can substitute Sa’ and Ka’ for Sa and K, and then directly use the equations derived in Gaussian space for the posterior error (Eq 3) and averaging kernel (Eq 4).  Now, I’ve seen and used this sort of Gaussian anamorphosis before, but only in the case wherein both the variable and observations are transformed from log-normal to Gaussian space. In this situation, indeed Eq 3 and 4 are applicable, as the solution is the BLUE solution (see equations 5.58 - 5.60 of Cohn 1997).  Also see for example the application studies of Brioude 2011, Saide 2012, Cui 2019, or the theoretical works of Fletcher 2010 and Bocquet 2010, none of which, I note, are cited in the present manuscript. 

However, with only x transformed but not y, I’m not sure they have arrived at a BLUE solution. What is the basis for the assertion that Eq. (3) gives the posterior error of the solution to (5)? I believe it is up to the authors to show that is indeed the case, or find a paper that shows this, or if not then to demonstrate that the approximation is in some way tolerable.  Citing Lu 2022 is not sufficient, as that paper makes the same assumption without discussion.  Alternatively, the authors might choose instead an ensemble approach to calculating their posterior error, as we did in Cui 2019 (using ~1000 samples). 

That being said, I appreciate the authors’ ensemble / sensitivity approach to error estimation, and it could be the distinction I’m making has little impact on total error given that their error ranges from the ensemble of sensitivity runs was often larger than their analytic posterior error.  But the interpretation of the error correlations  and  averaging kernel is rather critical to the results in the manuscript regrading sector-specific emissions, and it seems these are based entirely on the analytical calculations.  Though that brings me to my last point on this topic..

Lastly, and separately, that the ensemble standard deviation is typically larger than the analytic standard deviation (line 331 / Fig S1) also makes me wonder about the adequacy of the analytic posterior and Ak, even if calculated correctly.  Did the authors consider updating \hat{S} using information from the ensemble, rather than just the diagonal? If they in fact did do this (i.e., if Fig 6 incorporated the spread from the ensemble), it might be worth stating more clearly, as I didn’t catch that. 

Minor comments: 

97: I know what you mean, but technically it isn’t “no additional” cost, though indeed it is quite minimal relative to additional CTM runs. Still, the language could be more precise.

201: The previous global inversion of Qu 2021 used for the boundary conditions was based off an earlier version of TROPOMI data. Are the biases in that version (1.02) different from the one currently used (2.01)? I know the paper discusses how biases with respect to GOSAT have changed substantially between these two, by about 20 ppb, and that is in part impacted by the regional scope and higher spatial resolution of this analysis. Still, I wonder if the use of boundary values optimized by v1.02 data that is generally biased higher than the v2.02 data contributes to in part to the emissions increases found within this region inversion. Or is this compensated for by the boundary values mostly being scaled down (Table S1)? 

215: What is the reason for choosing 600 basis functions? This comes across as a bit arbitrary.  Later we learn the inversion DOFs is 167, so 600 does seem like a safe dimension for the inversion. Still, there’s no guarantee with this framework that all modes of variability in the emissions at 0.25 degree resolution captured by the observations are spanned by the 600 element GMM basis.  

230: The use of the word “optimal” here is ambiguous.  More precisely, the equation presented below provides the maximum likelihood estimate of the posterior (which for normally distributed errors is also the mean), and is better described as such. 

238: More precisely, additional regularization is used because both the observation and prior error covariance may not be well known.  It’s not as if gamma only corrects for unknowns in the former.  In fact, your results somewhat support this.  If gamma were really only compensating for misspecification of So, then we would wonder why your estimate of So is off be nearly two orders of magnitude.  In contrast, if we assume some of this is owing to underestimation of Sa, that helps explain why the increments found for individual posterior x values relative to xa can exceed the prior uncertainty of 50%. 

241: Actually, the relationship between x and concentrations is affine, not linear, since you don’t include the initial concentrations in x. But I’m guessing that the timescale of the inversion is long enough that the role of the initial conditions is minimal. 

269:  Since your emissions don’t actually follow a normal distribution and thus (x-xa)^2 isn’t chi-squared, does this still hold? Shouldn’t it be based on ln(x) and Sa’ to be chi-squared? I doubt it impacts your results though, and as you’ve already checked your sensitivity to gamma the overall conclusions likely remain unchanged. Still, in this case one might consider other methods for estimating gamma. By construction of So, I think your observation term is order m to begin with, so that method wouldn’t work.   It looks to me like the value identified in Lu 2021 Fig 4b might be at the corner of an L-curve, as the red line started to flatten again where the blue line rises dramatically. 

291: It’s not clear to me how the authors arrived at k = 10. The cited work of Lu 2022 does not explain either, and only cites Rodgers 2000.  Levenberg-Marquardt schemes usually adjust k dynamically, letting it be smaller when the minimization is progressing well (i.e., approaching taking Gauss-Newton steps) and making it be larger (i.e., approaching steepest-descent steps) if the minimization starts to become unstable.  How many iterations were typically required? 

281 - 305: I note that discussion of Equation 5 (lines 281 - 289, then again 295-305) conceptually follows Lu 2022 sentence by sentence, which while paraphrased, is not original content. Would it be better to just cite Lu 2022 here rather than repeat?  Otherwise we have about a dozen instances in sequence like the following:

Lu 2022: The boundary conditions are still optimized with normal error distributions, assuming an error standard deviation of 10 ppb.

Chen 2022: The normal error assumption is retained for the boundary conditions elements of the state vector, with a prior error standard deviation of 10 ppb. 

Perhaps this is a moot point though, if the authors do revise this section given some of the questions raised earlier. 

Throughout: could use some more proofreading, for example

14: as a prior 

23: the uncertainty range

34: unaccounted for (?)

84: extra space

85: missing spaces

References:

Bocquet, M., Pires, C. A., & Wu, L. (2010). Beyond Gaussian statistical modeling in geophysical data assimilation. Monthly Weather Review, 138(8), 2997–3023. https://doi.org/10.1175/2010MWR3164.1.

Brioude, J., et al. (2011), Top down estimate of anthropogenic emission inventories and their interannual variability in Houston using a mesoscale inverse modeling technique, J. Geophys. Res., 116, D20305, doi:10.1029/2011JD01621.

Cohn, S. E., 1997: An introduction to estimation theory. J. Meteor. Soc. Japan, 75 , 257–288.

Cui, Y., D. K. Henze, J. Brioude, W. M. Angevine, Z. Liu, N. Bousserez, J. Guerrette, S. A. McKeen, J. Peischl, B. Yuan, T. Ryerson, G. Frost, M. Trainer (2019), Inversion estimates of lognormally distributed methane emission fluxes from the Haynesville–Bossier oil and gas production region using airborne measurements, J. Geophys. Res., 124, 3520–3531, https://doi.org/10.1029/2018JD029489.

Fletcher, S. J. (2010). Mixed Gaussian lognormal four dimensional data assimilation. Tellus, 62 (3), 266–287

Saide, P., Peterson, D. A., da Silva, A., Anderson, B., Ziemba, L. D., Diskin, G., Sachse, G., Hair, J., Butler, C., Fenn, M., Jimenez, J. L., Campuzano-Jost, P., Perring, A. E., Schwarz, J. P., Markovic, M. Z., Russell, P., Redemann, J., Shinozuka, Y., Streets, D. G., Yan, F., Dibb, J., Yokelson, R., Toon, O. B., Hyer, E., & Carmichael, G. R. (2015). Revealing important nocturnal and day-to-day variations in fire smoke emissions through a multiplatform inversion. Geophysical Research Letters, 42, 3609–3618. https://doi.org/10.1002/2015GL063737