Anonymous Referee #1

Comments

This manuscript presents a Bayesian inversion of the global methane budget based on a harmonized 35-year dual-isotope dataset (δ13C–CH4 and δD–CH4). By combining a two-box model with a discrete parameter tuning (DPT) approach, the study demonstrates the added value of δD–CH4 in constraining methane sources and sinks. The work is methodologically innovative and addresses the ongoing debate on the drivers of global methane growth since 2006. The manuscript is clearly structured and merits publication after consideration of the points raised below.

Response: We thank the reviewer for their comments and suggestions. They have been addressed below and incorporated into the revised manuscript.

Major Comments

1.1 Tropical and low-latitude regions dominate global methane emissions and oxidation. The division of the globe into two hemispheric boxes at the equator (Fig. 1), together with the markedly different hemispheric source isotopic signatures listed in Table 1, implies that the same source type in the tropical region may be assigned substantially different isotopic values north and south of the equator. Given the large contribution of tropical methane to the global budget (some studies suggest on the order of 60–65%), the authors should discuss whether this structural simplification could introduce biases or otherwise influence the inferred global emissions and their temporal trends.

Response: We agree that two hemispheric boxes represent a structural simplification, and we clarify here why this approach remains appropriate for the objectives of this study. Our analysis is explicitly designed to constrain global- and to a lesser extent hemispheric-scale changes in the methane budget, rather than regional source attribution. The isotopic records used are derived from high-latitude background stations, which are intentionally selected because they are relatively well mixed, minimally influenced by local or regional emissions, and representative of hemispheric background air. These data are homogenised following Dasgupta et al. (2025), ensuring that the retained signal reflects hemispheric-scale variability rather than local source effects. Latitudinal gradients are not ignored; instead, they are aggregated into a hemispheric representation, ensuring that each hemisphere is correctly constrained. Specifically, the hemispheric gradients applied in the model are derived from long-term NOAA observations spanning a wide range of tropical, mid-latitude, and polar stations. The methodology for constructing and applying these gradients is described in detail in the Supporting Information (section 5, Fig. S4). We emphasise that this framework is sufficient for distinguishing, for example, between global biogenic and fossil contributions, and for identifying hemispheric asymmetries in growth rates, which are central to this paper. However, we explicitly acknowledge that more detailed spatial attribution (e.g., separating tropical from high-latitude wetlands) requires three-dimensional transport inversions and denser isotope observations. These limitations are discussed in the manuscript, and our conclusions are framed accordingly.

1.2 The posterior emissions shown in Figures 3 and 4 exhibit substantial interannual variability, particularly for wetlands, where year-to-year changes can exceed 50 Tg. It would be helpful to clarify whether these fluctuations are interpreted as reflecting true physical variability (e.g., climate-driven wetland dynamics), compensatory adjustments among sources within an underdetermined inversion, or artefacts arising from the two-box framework or the adopted prior uncertainty settings.

Response: The reviewer notes a pronounced interannual variability in the posterior wetland emissions. We clarify that these oscillations primarily arise from compensatory adjustments among sources (and lifetime) within an underdetermined inversion, rather than from artefacts arising from the two-box framework or the adopted prior uncertainty settings. Importantly, wetlands dominate the global methane budget (~250 Tg yr⁻¹), and under the assumed 30% prior uncertainty, this corresponds to an absolute uncertainty of ~75 Tg yr⁻¹. In contrast, pyrogenic emissions are much smaller (~28 Tg yr⁻¹), so the same fractional uncertainty corresponds to only ~8 Tg yr⁻¹. As a result, compensatory adjustments required to balance atmospheric growth and lifetime variability are most strongly expressed in the wetland category, simply because it is the largest and most flexible source. We also show that tightening prior uncertainties artificially amplifies oscillatory behaviour, while the error-scaled and DPT ensemble approaches reduce extreme variability without altering the inferred long-term growth. This demonstrates that the observed interannual variability is a predictable outcome of uncertainty scaling and source dominance, not evidence of instability, and does not undermine the main conclusions. These oscillations therefore reflect the inversion’s attempt to reconcile year-to-year changes in atmospheric methane growth with modest variations in sink strength, particularly OH-driven changes in lifetime. Only a minor fraction of this oscillation, if any, may be caused by teleconnections such as ENSO as seen in pyrogenic and wetland posterior emissions corresponding to El Niño and La Niña years, respectively. Smoothening of the posterior emissions by weighting the preceding years’ budgets in the inversion's cost function is expected to improve with newer versions of the model, which the authors are currently developing.

1.3 The DPT methodology is a key strength of the study. However, the physical basis for the chosen acceptance threshold (mean normalized RMSE < 0.1) requires further clarification. Please indicate what typical absolute residuals (in ‰) for δ13C–CH4 and δD–CH4 this threshold corresponds to, and discuss whether modest adjustments to this cut-off value would materially affect the main conclusions regarding source partitioning and long-term trends.

Response: The DPT ensemble uses a normalised RMSE threshold of 0.1 to identify model configurations that reproduce the atmospheric isotope observations within observational uncertainty. We clarify this point with an additional figure (Fig. S7b) showing that the threshold is well justified by the data distribution. As shown in this normalised RMSE plot, 92.8% of all isotope observations are ≲ 0.1, indicating that this cut-off retains the vast majority of observationally consistent solutions while excluding clearly unrealistic combinations of source signatures and lifetimes. The normalisation accounts for differences in absolute uncertainty between δ¹³C and δD, ensuring that neither tracer dominates the selection criterion. In absolute terms, retained solutions correspond to mean residuals on the order of ~0.05‰ for δ¹³C and ~0.2‰ for δD, well within observational and inter-laboratory uncertainty. Sensitivity tests using slightly stricter or looser thresholds produce nearly identical ensemble medians for global source trends, confirming that our conclusions are not sensitive to the precise numerical value of the threshold.

1.4 The adopted fossil fuel δD signatures (–192‰ for the Southern Hemisphere and –191‰ for the Northern Hemisphere) are consistent with a global average, yet regional variations are known to exist (e.g., with latitude, and among conventional gas, coal, and shale gas). Given that a key conclusion of the study is that including δD–CH4 reduces the inferred growth in fossil fuel emissions, a brief discussion of the sensitivity of this conclusion to potential systematic shifts in the assumed fossil δD signature would be valuable. For example, could a systematically heavier fossil δD signature permit a larger fossil fuel contribution while still reproducing the observed atmospheric δD trend?

Response: We acknowledge the reviewer’s concern regarding uncertainty and potential variability in fossil methane δD source signatures. In the DPT framework, fossil δD values are explored within ±5‰ of literature-reported values, allowing for realistic variability while remaining consistent with available measurements. The manuscript explicitly discusses scenarios in which shale gas methane, which can exhibit isotopically heavier δD signatures than conventional fossil methane due to its formation processes and production characteristics, contributes more strongly to the fossil fuel signal Uveges et al, 2025). We therefore do not assume a single immutable fossil δD endmember; instead, the DPT ensemble allows fossil δD to vary within bounds that encompass both conventional and shale-related values reported in the literature. A similar discussion follows Riddell-Young et al. (2025). Despite this flexibility, the inversion computes limited fossil growth across the ensemble. Even when fossil δD values are shifted toward the heavier end of the explored range, consistent with shale gas signatures, the atmospheric δD constraint still requires substantial increases in biogenic emissions to explain, for example, the observed post-2006 methane growth. This reflects the large and persistent isotopic separation between fossil and biogenic methane in hydrogen isotopes, as well as the sensitivity of our conclusions to potential systematic shifts in the assumed fossil δD signature.

1.5 δD–CH4 is highly sensitive to oxidation kinetics, yet the two-box model aggregates tropospheric oxidants (e.g., OH and Cl) and does not resolve their latitudinal gradients. While these limitations are acknowledged in Section 4.6, the manuscript would be strengthened by more clearly distinguishing which major conclusions (such as the dominance of biogenic methane growth after 2006) are likely robust at the global scale, and which may be more sensitive to the simplified representation of hemispheric transport and oxidation processes.

Response: As stated in the manuscript, methane lifetime interpretations should primarily be considered at the global scale. Given the 2-box setup, global-scale source totals and growth trends are expected to be robust, whereas hemispheric attribution is more sensitive to structural assumptions, including interhemispheric exchange time and prior uncertainty weighting (see Section 4.4). At the global scale, the dominance of biogenic growth after 2006 is robust: it holds across all five inversion scenarios, across the full DPT ensemble, and is insensitive to the precise lifetime partitioning between hemispheres (Table S4, Fig. 5). The global wetland growth estimate (~33 Tg) and the reduced fossil growth (~18 Tg vs. +44 Tg prior) are consistent regardless of whether OH is spatially uniform or hemispherically weighted. At the hemispheric scale, conclusions are more sensitive to structural assumptions. The NH/SH wetland split, the hemispheric lifetime divergence, and the hemispheric fossil attribution are all sensitive to interhemispheric lifetime (Fig. S5) and prior uncertainty weighting (Section 4.4). These should be treated as indicative rather than definitive. The aggregation of OH and Cl into a single tropospheric sink term adds further uncertainty to hemispheric lifetime estimates, specifically, though the global mean lifetime (~9.1 yr) is stable across all scenarios.

Minor Comments

1. Please clarify the usage and grouping of the terms thermogenic, pyrogenic, and fossil

throughout the manuscript to ensure consistency (e.g., specifying whether “fossil” is treated as a subset of “thermogenic”).

Response: Fossil fuel emissions are a subset of Thermogenic emissions. In our study, the total Thermogenic category as derived from the emission inventories (from EDGAR) includes Fuel exploration, Oil refineries, Chemical processes, and Power industry (Table S1). This has been clarified in the revised manuscript.

1. Line 252: “GIANS” should be corrected to “GAINS”.

Response: Rectified

1. Line 366: “DTP” should be corrected to “DPT”.

Response: Rectified

1. Upon first use, please spell out “RMSE” before using the abbreviation.

Response: Revised

1. Given that the isotope dataset is dominated by high-latitude monitoring sites, a brief

discussion of how representative these constraints are for tropical emissions would further

strengthen the interpretation.

Response: The 2-box approach, by design, collapses latitudinal structure into two hemispheric reservoirs, which removes regional signals, does not resolve mid-latitude or tropical gradients, and importantly, does not include intrahemispheric transport.

Anonymous Referee #2

General Comments

This study presents a joint inversion framework that assimilates atmospheric CH4 mole fraction, carbon isotopic composition (δ13C-CH4), and hydrogen isotopic composition (δD-CH4) to better constrain methane source partitioning and sink processes. The main objective was to assess whether combining δ13C-CH4 and δD-CH4 provides stronger constraints on fossil versus biogenic CH4 sources compared to single-isotope or mole-fraction-only inversions. However, several aspects of the methodology, interpretation, and robustness of the conclusions require further clarification and stronger justification before the main claims can be fully supported.

Response: We thank the reviewer for their comments and suggestions. They have been addressed below and incorporated into the revised manuscript.

2.1 The description of the box-model framework, particularly with respect to sink processes are briefly presented. The governing continuity equations for both CH4 and its isotopes should be explicitly presented. For isotope-based budget analyses, the representation of sinks is critically important, and the current description lacks insufficient details. The assumed interhemispheric exchange time also appears unrealistically short and may strongly influence the inferred hemispheric source partitioning. Previous three-dimensional modeling studies (e.g., TransCom) suggest exchange times on the order of ~1 year. Given the sensitivity of inversion results to this parameter, the use of more realistic values are strongly recommended.

Response:

The governing continuity equation and its numerical implementation are given in SI Section 7. This has been added to the code/data availability section. Sink-specific KIEs for all three isotopologues and all three sink processes are listed in Table 2. The pseudo-KIE derivation for the stratospheric sink is detailed in SI Section 9. On the exchange time: tau = 0.75 yr is shorter than typical 3-D model values (~1 yr), but this is expected in a two-box framework, which lacks vertical stratification and explicit interhemispheric transport. A shorter tau compensates for these structural simplifications and produces realistic hemispheric CH₄ gradients that are consistent with Naus et al. (2019). The sensitivity of results to tau (0.65 to 0.95 yr) is shown in Fig. S5 and has been independently verified with SF₆ inversion. Global totals are robust; however, hemispheric partitioning is more sensitive to tau. This is why we are confident in the global-scale inclusions and present the hemispheric-scale conclusions as indicative rather than definitive (Section 4.4).

2.2 In the comparison of inversion results, the species-specific inversions appear to reproduce their respective observations more closely (e.g., the δ13C-based inversion better matches δ¹³C observations, while the δD-based inversion better reproduces δD-CH4; Fig. 2). However, after approximately 2005, the differences between the dual-isotope inversion and the single-isotope inversions become difficult to distinguish. This raises an important question regarding the added value of the dual-isotope framework and whether it provides a significantly stronger constraint than single-isotope approaches.

Response:

As documented in SI Section 10 (Fig. S10), the dual-isotope inversion uniquely achieves a good simultaneous fit for both isotopic tracers. For δ¹³C-CH₄, the dual-isotope inversion achieves a mean residual of −0.05 ‰ and a median of −0.06 ‰, closely matching the carbon-only result (mean −0.07 ‰). For δD-CH₄, the dual-isotope inversion achieves a mean residual of −0.18 ‰ (median −0.09 ‰), compared to −0.55 ‰ (median −0.47 ‰) for the hydrogen-only inversion and −1.74 ‰ (median −1.18 ‰) for the carbon-only inversion. By definition, a single-isotope inversion constrained only by δ¹³C cannot reproduce the deuterium record (and vice versa). The dual-isotope setup is therefore the only configuration that simultaneously satisfies both isotopic constraints, necessary for a physically consistent global budget.

The more consequential effect of δD-CH₄ is on source attribution rather than on tracer fit in the later period. As shown in Fig. 3d and discussed in Section 4.1, the temporal evolution of fossil fuel emissions diverges substantially between carbon-only (+28.7 Tg growth) and dual-isotope/hydrogen-only (+19–21 Tg) inversions. This divergence is most apparent before ~2005 and reflects δD-CH₄'s strong sensitivity to the enriched deuterium signature of fossil CH₄ (δD ≈ −191 ‰), which is approximately 110–150 ‰ heavier than biogenic sources (Table 1). The tropospheric deuterium KIE of 1.313 produces a fractionation roughly 46 times larger than the carbon KIE of 1.0068, making δD-CH₄ particularly diagnostic of fossil/biogenic partitioning.

This temporal redistribution is not visible in the tracer fit after 2005 in Fig. 2. By this period, all inversions have converged toward a similar emission trajectory, consistent with both isotopic records. The dual-isotope framework primarily constrains the early-period fossil baseline and the temporal shape of the fossil emission trajectory, reflected in the growth estimates over the full analysis period.

2.3 Another concern is the large interannual variability (IAV) inferred for agriculture and waste emissions. These sectors are generally expected to exhibit relatively smooth, monotonic trends rather than strong short-term fluctuations. The pronounced variability appears somewhat artificial and may reflect compensatory adjustments among source categories within an underdetermined inversion system. This raises concerns about the robustness of the inferred emissions and the extent to which artifacts from the simplified two-box framework may influence the results.

Response: The large interannual variability (IAV) in posterior agriculture and waste emissions primarily arises from compensatory adjustments among sources (and lifetime) within an underdetermined inversion, rather than from artefacts arising from the two-box framework or the adopted prior uncertainty settings. Importantly, wetlands dominate the global methane budget (~250 Tg yr⁻¹), and under the assumed 30% prior uncertainty, this corresponds to an absolute uncertainty of ~75 Tg yr⁻¹. In contrast, pyrogenic emissions are much smaller (~28 Tg yr⁻¹), so the same fractional uncertainty corresponds to only ~8 Tg yr⁻¹. As a result, compensatory adjustments required to balance atmospheric growth and lifetime variability are most strongly expressed in the wetland category, simply because it is the largest and most flexible source. We also show that tightening prior uncertainties artificially amplifies oscillatory behaviour, while the error-scaled and DPT ensemble approaches reduce extreme variability without altering the inferred long-term growth. This demonstrates that the observed interannual variability is a predictable outcome of uncertainty scaling and source dominance, not evidence of instability, and does not undermine the main conclusions. These oscillations therefore reflect the inversion’s attempt to reconcile year-to-year changes in atmospheric methane growth with modest variations in sink strength, particularly OH-driven changes in lifetime. Only a minor fraction of this oscillation, if any, may be caused by teleconnections such as ENSO as seen in pyrogenic and wetland posterior emissions corresponding to El Niño and La Niña years, respectively. Smoothening of the posterior emissions by weighting the preceding years’ budgets in the inversion's cost function is expected to improve with newer versions of the model, which the authors are currently developing. This is discussed in Section 4.4 of the manuscript.

Section 4.2 of the manuscript discusses how agriculture and waste occupy intermediate isotopic space and are susceptible to compensatory adjustments. For this reason, Fig. 4 presents 5-year moving averages alongside raw posteriors, and growth statistics in Table S4 use 3-year endpoint averages. All five inversion scenarios agree on global totals to within a few Tg (Table S4), confirming robustness of the growth estimates despite year-to-year noise.

2.4 Furthermore, the inversion results are not fully consistent with the observed isotopic trends in both hemispheres. For example, the dual-isotope simulation suggests a more rapid decline in δ13C over the Southern Hemisphere than indicated by observations (Fig. 3). Similarly, for δD-CH4, all inversion cases show broadly similar temporal behavior, yet the model underestimates the observed trend in the Northern Hemisphere (Fig. 3). These discrepancies suggest that the proposed source attribution may not fully capture the underlying processes.

Response: The reviewer cites Fig. 3, which shows posterior emissions. The tracer fits are shown in Fig. 2. In Fig. 2b, the SH posterior δ¹³C does sit systematically below the observations across most of the analysis period. This small negative bias is quantified in Fig. S10 (mean δ¹³C residual ~−0.05 to −0.07 ‰) and falls within observational uncertainty. It most likely reflects residual uncertainty in the SH biogenic source signatures and the latitudinal gradient correction applied to convert high-latitude station records to hemispheric averages (SI Section 5). In addition, the DPT optimisation already adjusts SH wetland δ¹³C toward slightly heavier values (Table 1), reducing but not fully eliminating this bias. For NH δD-CH₄ (Fig. 2c), a modest underestimation in the early record is similarly small and quantified in Fig. S10. Both residuals are consistent with known uncertainties in tropical source signatures and the structural limitations of the two-box framework, as discussed in Section 4.6.

2.5 In particular, the interpretation of increasing microbial emissions—especially from Southern Hemisphere wetlands—appears somewhat overstated given the available evidence. The observed δ13C trend could also be consistent with alternative explanations, such as a slower increase in thermogenic emissions or different sink dynamics. Therefore, the strength of the conclusions regarding source attribution should be moderated, and uncertainties more clearly acknowledged.

Response: The observed atmospheric δ¹³C trend is, in principle, also consistent with alternative explanations, such as a slower rise in thermogenic emissions, changes in source signatures, or variations in sink strength. This is precisely why we present five distinct inversion approaches (Prior, Carbon-only, Hydrogen-only, Dual-isotope, DPT-ensemble, Error-scaled) rather than a single preferred solution. All five scenarios are cost-minimised solutions to the same observational constraints under different prior assumptions, and all converge on the conclusion that biogenic growth dominates (Table S4, Fig. 5). This cross-approach consistency is the basis for our confidence in the qualitative conclusion, not any single inversion run. Our inversion framework explicitly accounts for uncertainty through: (a) ±30% prior emission uncertainties propagated into posterior standard deviations (reported in Table S4 as ± values); (b) the DPT ensemble, which aggregates over >13 million parameter combinations and retains all runs with normalised RMSE < 0.1; this ensemble explicitly quantifies parametric uncertainty in source signatures, KIEs, and error scaling; and (c) the error-scaled scenario, which tests sensitivity to prior weighting. The uncertainty ranges are presented throughout the results, and we do not claim that the inversion yields a unique solution; rather, we note that the range of acceptable solutions consistently favours biogenic-dominant growth.

We agree that the language can be moderated and will revise to clarify that the inversion yields a cost-optimised solution given the observations and priors, not a unique attribution. Alternative interpretations involving slower thermogenic growth or sink variability are discussed in Sections 4.1 and 4.3, and the two-box limitations are acknowledged in Section 4.6.

Specific Comment on Figure 1

2.6 In Figure 1, the second-last pictorial element appears to represent fossil-fuel combustion. However, this may not be the most appropriate visualization for methane emissions, as fossil-fuel CH4 emissions predominantly arise from leakage during extraction, processing, and transport rather than combustion itself. Additionally, the clarity of the figure would be improved by explicitly labeling all source categories.

Response: Following the reviewer’s suggestion, we have revised Figure 1. The Fossil Fuel pictorial representation (second from bottom) has been updated. Source labels have been added.