# Modeled production, oxidation, and transport processes of wetland methane emissions in temperate, boreal, and Arctic regions

Masahito Ueyama, Sara H Knox, Kyle B Delwiche, Sheel Bansal, William J Riley, Dennis Baldocchi, Takashi Hirano, Gavin Mcnicol, Karina Schafer, Lisamarie Windham-myers, et al.

► To cite this version:

Masahito Ueyama, Sara H Knox, Kyle B Delwiche, Sheel Bansal, William J Riley, et al.. Modeled production, oxidation, and transport processes of wetland methane emissions in temperate, boreal, and Arctic regions. Global Change Biology, 2023, 29 (8), pp.2313-2334. 10.1111/gcb.16594. hal-03955531

HAL Id: hal-03955531

https://hal.science/hal-03955531

Submitted on 25 Jan 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
Modeled production, oxidation and transport processes of wetland methane emissions in temperate, boreal, and Arctic regions

Masahito Ueyama¹, Sara H. Knox², Kyle B. Delwiche³, Sheel Bansal⁴, William J. Riley⁵, Dennis
Baldocchi³, Takashi Hirano⁶, Gavin McNicol⁷, Karina Schafer⁸, Lisamarie Windham-Myers⁹,
Benjamin Poulter¹⁰, Robert B. Jackson¹¹, Kuang-Yu Chang⁵, Jiquen Chen¹², Housen Chu¹³, Ankur
R. Desai¹⁴, Sébastien Gogo¹⁵, Hiroki Iwata¹⁶, Minseok Kang¹⁷, Ivan Mammarella¹⁸, Matthias
Peichl¹⁹, Oliver Sonnentag²⁰, Eeva-Stiina Tuittila²¹, Youngryel Ryu²², Eugénie S. Euskirchen²³,
Mathias Göckede²⁴, Adrien Jacotot²⁵, Mats Nilsson²⁶, Torsten Sachs²⁷

¹ Graduate School of Agriculture, Osaka Metropolitan University, Sakai, Japan

² Department of Geography, The University of British Columbia, Vancouver, Canada

³ Department of Environmental Science, Policy & Management, UC Berkeley, Berkeley, CA, USA

⁴ U.S. Geological Survey, Northern Prairie Wildlife Research Center, Jamestown, ND 58401 USA

⁵ Climate and Ecosystem Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA, USA

⁶ Research Faculty of Agriculture, Hokkaido University, Sapporo, Japan

⁷ Department of Earth and Environmental Sciences, University of Illinois Chicago, Chicago, IL, USA

⁸ Department of Earth and Env Science, Rutgers University Newark, NJ, USA

⁹ Water Mission Area, U.S. Geological Survey, Menlo Park, CA, USA

¹⁰ NASA Goddard Space Flight Center, Biospheric Sciences Laboratory, Greenbelt, MD, USA

¹¹ Department of Earth System Science, Stanford University, Stanford, CA, USA

¹² Department of Geography, Environment, and Spatial Sciences, Michigan State University, East Lansing, MI, USA

¹³ Climate and Ecosystem Sciences Division, Lawrence Berkeley National Lab, Berkeley, CA, USA

¹⁴ Dept of Atmospheric and Oceanic Sciences, University of Wisconsin-Madison, Madison, WI, USA

¹⁵ ECOBIO (Écosystèmes, Biodiversité, Évolution), Université Rennes 1, CNRS, UMR 6553, 35000, Rennes, France
32 ¹⁶ Department of Environmental Science, Faculty of Science, Shinshu University, Matsumoto, 33 Japan

34 ¹⁷ National Center for Agro Meteorology, Seoul, South Korea

35 ¹⁸ Institute for Atmospheric and Earth System Research/Physics, Faculty of Science, University of 36 Helsinki, Helsinki, Finland

37 ¹⁹ Department of Forest Ecology and Management, Swedish University of Agricultural Sciences, 38 901 83 Umeå, Sweden

39 ²⁰ Université de Montréal, Département de géographie, Université de Montréal, Montréal, QC 40 H2V 0B3, Canada

41 ²¹ School of Forest Sciences, University of Eastern Finland, Joesnuu, Finland

42 ²² Department of Landscape Architecture and Rural Systems Engineering, Seoul National 43 University, Seoul, South Korea

44 ²³ University of Alaska Fairbanks, Institute of Arctic Biology, Fairbanks, AK, USA

45 ²⁴ Max Planck Institute for Biogeochemistry, Department of Biogeochemical Signals, Jena, 46 Germany

47 ²⁵ INRAE, UMR 1069 SAS, 65 rue de Saint-Brieuc, 35042 Rennes, France

48 ²⁶ Department of Forest Ecology and Management, Swedish University of Agricultural Sciences, 49 901 83 Umeå, Sweden

50 ²⁷ GFZ German Research Centre for Geosciences, Telegrafenberg, 14473 Potsdam, Germany

51 Corresponding author: M. Ueyama, Graduate School of Agriculture, Osaka Metropolitan 52 University, Sakai, Japan (Tel: +81-72-254-9432) (mueyama@omu.ac.jp)

54 ORCiD

56 M. Ueyama : https://orcid.org/0000-0002-4000-4888

57 S. Knox : https://orcid.org/0000-0003-2255-5835

58 K. B. Delwiche : https://orcid.org/0000-0002-5981-2500

59 S. Bansal : https://orcid.org/0000-0003-1233-1707

60 W. J. Riley : https://orcid.org/0000-0002-4615-2304

61 D. Baldocchi : http://orcid.org/0000-0003-3496-4919

62 T. Hirano : https://orcid.org/0000-0002-0325-3922
63 G. McNicol : https://orcid.org/0000-0002-6655-8045

64 K. Schafer : No ID

65 L. Windham-Myers : https://orcid.org/ 0000-0003-0281-9581

66 B. Poulter : https://orcid.org/ 0000-0002-9493-8600

67 R. B. Jackson : https://orcid.org/0000-0001-8846-7147

68 K.-Y. Chang : https://orcid.org/ 0000-0002-7859-5871

69 J. Chen : https://orcid.org/0000-0003-0761-9458

70 H. Chu : https://orcid.org/0000-0002-8131-4938

71 A. R. Desai : https://orcid.org/0000-0002-5226-6041

72 S. Gogo : https://orcid.org/0000-0002-0867-497X

73 H. Iwata : https://orcid.org/0000-0002-8962-8982

74 M. Kang : https://orcid.org/0000-0003-4901-4465

75 I. Mammarella : https://orcid.org/0000-0002-8516-3356

76 M. Peichl : https://orcid.org/0000-0002-9940-5846

77 O. Sonnentag : https://orcid.org/0000-0001-9333-9721

78 E.-S. Tuittila : https://orcid.org/0000-0001-8861-3167

79 Y. Ryu : https://orcid.org/0000-0001-6238-2479

80 E. S. Euskirchen : https://orcid.org/0000-0002-0848-4295

81 M. Göckede : https://orcid.org/0000-0003-2833-8401

82 A. Jacotot : https://orcid.org/0000-0002-0126-7597

83 M. Nilsson : https://orcid.org/0000-0003-3765-6399

84 T. Sachs : https://orcid.org/0000-0002-9959-4771
86 **Abstract**

87 Wetlands are the largest natural source of methane (CH₄) to the atmosphere. The eddy covariance 88 method provides robust measurements of net ecosystem exchange of CH₄, but interpreting its 89 spatio-temporal variations is challenging due to the co-occurrence of CH₄ production, oxidation, 90 and transport dynamics. Here we estimate these three processes using a data-model fusion 91 approach across 25 wetlands in temperate, boreal, and Arctic regions. Our data-constrained model 92 – iPEACE – reasonably reproduced CH₄ emissions at 19 of the 25 sites with normalized root mean 93 square error of 0.59, correlation coefficient of 0.82, and normalized standard deviation of 0.87.

94 Among the three processes, CH₄ production appeared to be the most important process, followed 95 by oxidation in explaining inter-site variations in CH₄ emissions. Based on a sensitivity analysis, 96 CH₄ emissions were generally more sensitive to decreased water table than to increased gross 97 primary productivity or soil temperature. For periods with leaf area index (LAI) of ≥20% of its 98 annual peak, plant-mediated transport appeared to be the major pathway for CH₄ transport. 99 Contributions from ebullition and diffusion were relatively high during low LAI (<20 %) periods.

100 The lag time between CH₄ production and CH₄ emissions tended to be short in fen sites (3 ± 2 101 days) and long in bog sites (13 ± 10 days). Based on a principal component analysis, we found that 102 parameters for CH₄ production, plant-mediated transport, and diffusion through water explained 103 77% of the variance in the parameters across the 19 sites, highlighting the importance of these 104 parameters for predicting wetland CH₄ emissions across biomes. These processes and associated 105 parameters for CH₄ emissions among and within the wetlands provide useful insights for 106 interpreting observed net CH₄ fluxes, estimating sensitivities to biophysical variables, and 107 modeling global CH₄ fluxes.

108

109 Running title: Estimating processes for wetland CH₄ emissions

110 **Keywords:** Methane emissions, Eddy covariance, Bayesian optimization, Multi-site synthesis,

111 Methane model, Data-model fusion
## 1. Introduction

Wetlands are the largest natural source of methane ($CH_4$) – a potent greenhouse gas contributing to climate warming. Methane emissions from wetlands contribute approximately 20% of total annual $CH_4$ emissions (Saunois et al., 2020). Despite their importance, estimates of wetland $CH_4$ emissions are highly uncertain (Bohn et al., 2015; Melton et al., 2013) because direct measurements of $CH_4$ emissions (Delwiche et al., 2021) are far fewer than those of carbon dioxide ($CO_2$) fluxes (Pastorello et al., 2020). In particular, the variability in $CH_4$ emissions appears high across spatial and temporal scales (Knox et al., 2019; Delwiche et al., 2021). As a result of the associated uncertainties, current estimates of the global $CH_4$ budget contain large discrepancies between top-down and bottom-up approaches (Saunois et al., 2020; Jackson et al. 2020).

Methane emissions from wetlands also exhibit a wide range of magnitudes and responses to biophysical variables. Because $CH_4$ is primarily produced by anaerobic methanogens and oxidized by aerobic bacteria (Bridgham et al., 2013; Conrad 2009), water table depth (WTD) has been identified as an important thermodynamic boundary and thus potential predictor of wetland $CH_4$ emissions (Brown et al., 2014; Moore and Roulet, 1993; Rinne et al., 2018). Methanogens produce $CH_4$ using substrates both from carbon recently fixed through photosynthesis (Whiting and Chanton, 1993) and previously fixed carbon (Glaser et al., 2004; Karofeld and Tõnisson, 2014). Thus, $CH_4$ emissions are often correlated with plant primary production and/or soil respiration (Turetsky et al., 2014; Villa et al., 2020; Whiting and Chanton, 1993). Because temperature affects $CH_4$ production kinetics, soil temperature is typically correlated with $CH_4$ emissions (Knox et al., 2019; Yvon-Durocher et al., 2014), albeit substantial seasonal hysteresis was reported to occur in many sites, likely due to substrate-temperature driver interactions (Chang et al. 2020, 2021). In addition to production and oxidation, transport pathways are also crucial in modeling $CH_4$ emissions. Because $CH_4$ in soils is transported through plant aerenchyma, ebullition bubbles through standing water, and/or diffusion, $CH_4$ emissions were shown to be often correlated with leaf area index (LAI), latent heat flux, and/or barometric pressure (PA) (Kwon et al., 2017; Sturtevant et al., 2016; Tokida et al., 2005; Ueyama et al. 2020).

To better understand wetland $CH_4$ emissions, the eddy covariance (EC) method has been widely used at various wetlands along with measurements of other ancillary covariates such as WTD and soil temperature (Morin, 2018; Knox et al., 2019; Delwiche et al., 2021). The EC method provides quasi-continuous measurements of $CH_4$, $CO_2$, and energy exchanges between the land
surface and the atmosphere (Baldocchi, 2014). The direct measurements have been used to evaluate magnitudes of CH₄ emissions, their interannual variations, and their responses to various biophysical variables (Chang et al., 2021; Chu et al., 2014; Knox et al., 2019; Rinne et al., 2018; Yuan et al., 2022). Previous studies have identified biophysical variables such as soil and air temperature and WTD as the primary drivers for the temporal and spatial variations in CH₄ emissions (Knox et al., 2019; Turetsky et al., 2014; Yuan et al., 2022), but their importance varies substantially among wetlands and across time scales (Koebsch et al., 2015; Knox et al., 2021). Furthermore, complex interactions hinder the use of simple correlation analyses for disentangling responses of CH₄ emissions to biophysical variables, leading to large uncertainties when interpreting observations (Chang et al. 2020; Sturtevant et al., 2016; Knox et al., 2021). Recently, the FLUXNET-CH₄ database was curated for supporting synthesis of wetland CH₄ emissions using the EC methods (Knox et al., 2019; Delwiche et al., 2021) and, for example, was used to evaluate inter-site variations in CH₄ emissions (Chang et al. 2021; Knox et al., 2021, 2021; Yuan et al., 2022).

To improve the mechanistic understanding and accurate modeling of CH₄ emissions, the relative contributions of CH₄ emission pathways have been measured or estimated with various field measurements (Table 1). These measurements include chamber techniques (Korrensalo et al., 2022; Tokida et al., 2007a, b), bubble traps (Stanley et al., 2019), isotope techniques (Dorodnikov et al., 2011), and dissolved CH₄ concentrations in pore water (McNicol et al., 2017). Recently, a wavelet analysis of EC measurements examined the contribution of ebullition to total CH₄ emissions (Göckede et al., 2019; Iwata et al., 2018; Hwang et al., 2020; Richardson et al., 2022; Schaller et al., 2018). These analyses revealed that plant-mediated transport was the most important pathway for wetland CH₄ emissions (up to 98% of the total emissions), but the other two pathways were also important under environmental conditions such as flooded wetlands without emergent vegetation and shallow ponds. Many process-based models (Table 1) have also shown that CH₄ emissions occur mostly through plant-mediated transport (Castro-Morales et al., 2018; Ma et al., 2017; Peltola et al., 2018; Susiluoto et al., 2018; Wania et al., 2010), although one model found ebullition was the dominant pathway (Ito, 2019). Although previous studies conducted across relatively few wetland sites are useful for understanding CH₄ transport pathways, comparisons of transport mechanisms across multiple wetlands remain challenging. The challenge lies in uncertainties in measurement techniques, spatial representation of measured processes in
175 the field, and different model structures in process-based models.

176 Data-model fusion approaches have recently been used for evaluating wetland CH₄ emissions
177 (Ma et al., 2017; Müller et al., 2015; Salmon et al., 2022; Susiluoto et al., 2018; Ueyama et al.,
178 2022). These methods use observed data for constraining process-based models that are often
179 difficult to calibrate, and can be used to evaluate processes of CH₄ emissions and their sensitivity
180 to biophysical drivers. To reduce the uncertainties in a process-based model, Müller et al. (2015)
181 used observed data for constraining a model for CH₄ dynamics and found that detailed process-
182 based models were not well constrained owing to the complexity of the model. Susiluto et al.
183 (2018) calibrated a detailed model using nine years of EC-based CH₄ flux data in a northern fen.
184 Their results suggested that CH₄ production was the most important factor responsible for the
185 interannual variations in CH₄ emissions, whereas plant-mediated transport was the most important
186 CH₄ transport pathway. Data-model fusion approaches to study CH₄ emissions have been applied
187 only for a limited number of individual sites; thus, their applicability should be evaluated across
188 wide arrays of wetland sites and biomes.

189 Recently, Ueyama et al. (2022) developed a process-based model (i.e., inferring Processes for
190 Ecosystem-Atmosphere CH4 Exchange - iPEACE) for partitioning CH4 emissions using a data-
191 model fusion approach for a cool temperate bog in Japan. Their approach constrained the model
192 using CH4 emissions and associated biophysical variables from the EC tower with the goal to
193 determine a parameter set for reproducing daily CH4 emissions under various environmental
194 conditions. These conditions included growing and dormant seasons, wet and dry conditions, high
195 and low LAI, and various ranges of gross primary production (GPP), soil temperature, and PA.
196 The model reasonably identified processes that were qualitatively consistent with previous field
197 experiments to shed light on processes in the bog. Findings include: 1) ebullition and plant-
198 mediated transport as the important CH4 transport pathways, 2) high contributions of the deep
199 organic layer (i.e., <30 cm) to total CH4 emissions due to very low CH4 concentrations in the
200 surface organic layer (Tokida et al., 2007a), and 3) gaseous-bubble accumulation in deep organic
201 layer (Tokida et al., 2005, 2007a, b). A chamber-based study further suggested that contributions
202 of bubble transport to total CH4 emissions ranged from 67-95% during the snow-free season in the
203 bog (Tokida et al., 2005, 2007a, b), which was close to the iPEACE model estimates (64%).

204 Here, we modified iPEACE to simulate CH₄ fluxes and infer processes related to CH₄
205 emissions (i.e., production, oxidation, and transport pathways) from 25 wetlands across mid- to
high-latitudes included in the FLUXNET-CH₄ database. Applying the data-model fusion method (Ueyama et al., 2022) across these wetland sites spanning temperate, boreal, and Arctic regions, our objectives were to: 1) evaluate the model's suitability for simulating CH₄ emissions across wetland types, 2) quantify inter-site variations in estimated processes related to CH₄ emissions, 3) evaluate the sensitivities of CH₄ emissions to GPP, soil temperature, and WTD, and 4) examine inter-site variations in parameters for improved predictions of wetland CH₄ emissions.

## 2. Materials and methods

### 2.1 Dataset and model inputs

We used daily EC CH₄ flux data archived in the FLUXNET-CH₄ database (Delwiche et al., 2021). We selected all mid- to high-latitude freshwater wetland sites from the database (Table 2) that contained all relevant forcing variables (i.e., soil and air temperature, WTD, PA, and GPP). The selected 25 sites represent wetland types of bog (ombrotrophic), fen (minerotrophic), marsh, wet tundra, and rice paddy in temperate, boreal, and Arctic regions. The mean annual air temperature ranged from -5 °C to 17 °C across the sites, and minimum WTD ranged from -0.62 m to 0.68 m.

We used daily gap-filled CH₄ fluxes and the ancillary biophysical variables at the tower sites. The daily mean values of the gap-filled half-hourly variables were provided in the FLUXNET-CH₄ database (Delwiche et al., 2021). We used two types of daily CH₄ fluxes (i.e., FCH4_F and FCH4_F ANN median) in the database. FCH4_F was gap-filled using a multidimensional scaling (MDS) approach in R EddyProc (Delwiche et al., 2021), but still contained periods of time with long data gaps (<2 months). FCH4_F ANN median was gap-filled based on an artificial neural network method, which fills all data gaps (Knox et al., 2019). As input drivers from the FLUXNET-CH₄ database, daytime-based GPP (GPP_DT) in the database (Lasslop et al., 2010), air temperature (TA_F), barometric pressure (PA_F), soil temperature (TS), and WTD (WTD_F) were used. The gaps in the meteorology (i.e., TA_F, and PA_F) were filled using the ERA-Interim reanalysis data (Vuichard and Papale, 2015), whereas those of WTD and soil temperatures were filled using the MDS method. We used soil temperature at two depths for representing the surface and deep layers in the model. For sites affected by permafrost (RU-Ch2, US-Ics, and US-Uaf), we assumed that the deepest soil temperature measurement was representative of the bottom of the active layer. Data for RU-Ch2, US-Ics, US-Bzf, and US-Bzb sites did not include WTD data in the FLUXNET-CH₄ database, but WTD data were directly provided from principal investigators.
Since WTD for RU-Ch2 was based on discrete manual measurements, we linearly interpolated the data to the daily timescale.

We prepared daily LAI as a model input based on satellite-based LAI smoothed using GPP. First, the four-day LAI data (MCD15A3H; collection 6) was downloaded from MODIS land products subsets. The spatial resolution of the product is 500 m. We used a single grid cell of data centered on the site location. The LAI data were first set to zero for the snow periods, and were then smoothed using a Savitzky-Golay filter (Chen et al., 2004). The snow conditions were determined based on the MODIS reflectance products (MCD43A4; collection 6) from the MODIS land products subsets. Because smoothed LAI often failed to explain seasonal peaks when peak LAI was missing, daily LAI was then modeled using the smoothed LAI and daily GPP normalized with a maximum GPP (nGPP). LAI at day (i) was modeled with a non-centered moving mean of the normalized GPP multiplying a scale factor.

$$LAI_i = L_s \sum_{j=i-D}^{i+D} nGPP_j/(D+1) \quad (1)$$

Two empirical parameters of the scale factor for explaining maximum LAI ($L_s$) and moving window for explaining a lag between GPP and LAI ($D$) were the parameters determined based on a differential evolution method. Since there was no clear relationship between LAI and GPP for NZ-Kop, LAI for NZ-Kop was estimated simply based on 10-day moving mean of the satellite-based LAI. The smoothed LAI well mimicked the satellite-derived LAI, where mean and standard deviation of root mean square error (RMSE) and correlation coefficient (R) were $0.46 \pm 0.24$ and $0.84 \pm 0.11$, respectively, across the sites.

## 2.2 The iPEACE model

Partitioning CH$_4$ emissions from the EC measurements was conducted by optimization of a process-based model with the data. We used the iPEACE model (Ueyama et al., 2022), which was originally proposed to infer CH$_4$ dynamics at a temperate bog in Japan, but has been generalized for the current analysis (Fig. 1).

The iPEACE model consists of two soil layers, a surface layer susceptible to oxic conditions and a deep layer prone to anoxic conditions, and considers CH$_4$ production and oxidation in each layer, as well as three transport pathways: plant-mediated transport, ebullition, and diffusion. The
modeled mechanisms are similar to those used in current process-based models (Ralvonen et al., 2017; Walter and Heimann, 2000; Wania et al., 2010; Riley et al., 2011). The simple formulation of iPEACE allows to effectively fit the model to data at reduced computational costs. The model is driven with GPP for substrate availability, LAI for transport potential through plant stems, soil temperature in the two layers for driving kinetics, oxygen ($O_2$) concentration for redox potential, WTD for diffusivity and hydrostatic pressure that drives ebullition, and PA for ebullitive transport. The $O_2$ concentration was not included in the FLUXNET-$CH_4$ database, and thus was determined from WTD. When the water table position is above or below a soil layer, the layer is assumed to be anoxic or fully oxic, respectively. When WTD is within a soil layer, $O_2$ concentration in that layer is linearly related to that fraction of the layer that is inundated between fully oxic to anoxic conditions.

To explore the underlying processes, the model contains 10 parameters and two initial values of the $CH_4$ pools in each soil layer (mol-$CH_4$ m$^{-3}$), which are calibrated with data (Table 3). For adapting the model to the current analysis, the thickness of the surface layer and root fraction (described below) in the surface layer are calibrated for each site, whereas the previous study (Ueyama et al., 2022) used a fixed value.

Thickness of the surface layer ($z_{surf}$; m) is the parameter constrained by the data. Thickness of the deep layer is calculated as the difference between total soil thickness (1 m, except for permafrost sites) and the thickness of the surface layer. For sites affected by permafrost, total soil thickness is defined as the active layer depth (0.5 m for RU-Ch2, 1.0 m for US-ICs, and 0.6 m for US-Uaf). Seasonal changes in soil thickness associated with soil thaw are not considered in the model for simplicity. Surface root fraction ($f_{root}$) is the parameter explaining how roots are concentrated in the surface layer relative to the total roots. The model assumes that root density is higher in the surface layer than the deep layer.

Methane production is assumed to depend on substrate availability from GPP, kinetics as determined by soil temperature, and anaerobic status as determined by $O_2$ concentration. The fraction of GPP to $CH_4$ substrate ($p_{production}$; mmol-$CH_4$ g$^{-1}$ C) and temperature sensitivity ($Q_{10}$) are both empirical parameters. Modeled $CH_4$ production increases with soil temperature and substrate availability but decreases with increased $O_2$ concentration. The $p_{production}$ parameter is the aggregated parameter for explaining the fraction of root exudates from GPP and the efficiency from exudates to $CH_4$ production and relates to the base production rate in a $Q_{10}$ equation (Chen,
2021). The model does not explicitly consider anaerobic peat decomposition; thus, CH₄ production by decomposition are implicitly incorporated through a decrease in the CH₄ pools. Partitioning of CH₄ substrate in each soil layer is assumed to be a function of the root distribution between the surface and deep soil layers. CH₄ oxidation is calculated with a Michaelis-Menten equation (Wania et al., 2010) with CH₄ concentration and O₂ concentration, where the maximum CH₄ oxidation rate (p_oxidation; mol-CH₄ m⁻³ s⁻¹) is a calibrated parameter.

Plant-mediated transport is calculated by the concentration gradient between a soil layer and the atmosphere, root fraction in each layer, and LAI. The transfer efficiency under a given concentration gradient (p_plant; 10⁻³ d⁻¹) is a calibrated parameter. The model does not consider CH₄ transport by dead plants, which are not accounted for by LAI, with the assumption that collapsed aerenchymatous tissue in senesced leaves has low transport capacity (Korrensalo et al., 2022).

Ebullitive transport is calculated based on a concentration threshold scheme (Peltola et al., 2018), which has two empirical parameters: nondimensional conductivity for bubble transport (p_ebullition) and a parameter for explaining episodic CH₄ bubble transport driven by barometric pressure changes (p pressure; hPa⁻¹). Since the model assumes that CH₄ is not immediately emitted as ebullition but accumulated as bubbles, p_ebullition represents the transport efficiency of bubbles. The p pressure parameter empirically explains the sensitivity to decreasing barometric pressure, i.e., the relative increase in ebullition per 1 hPa decrease in mean PA. In the model, the ebullition flux from each layer is assumed to be directly transported to the atmosphere, when WTD is within the top 10 cm of the soil based on a field study (Stanley et al., 2019). When WTD is deeper than 10 cm, CH₄ transport through ebullition is added to the surface layer CH₄ pool, which is a modification from the original model of Ueyama et al. (2022).

Diffusive flux is calculated using Fick's first law. The diffusion coefficients for gas and water are calculated based on Riley et al. (2011), and then their calibrated correction factors (p diffusion-gas and p diffusion-water) are multiplied to the respective diffusion coefficients.

## 2.3 Model applications

The model parameters, initial conditions, and model error (σ) were determined from the observed data by the Bayesian method as follows:

$$F_{\text{OBS}} \sim \text{Normal}(F_{\text{MODEL}}, \sigma^2) \qquad (2)$$
where the function *Normal* represents the normal distribution, $F_{\text{OBS}}$ is the observed CH$_4$ emission, and $F_{\text{MODEL}}$ is the modeled CH$_4$ emission. The *a priori* distribution of $\sigma$ was assumed to be a log normal distribution with mean of log(0.5) mg CH$_4$ m$^{-2}$ d$^{-1}$ and standard deviation of 0.1 mg CH$_4$ m$^{-2}$ d$^{-1}$, where the hierarchical structure was used to reduce computational costs. Equation 2 assumes that variance for the model-observation mismatch was temporally uniform without incorporating temporal correlation in the observed data.

The *a priori* distributions of the parameters were generally assumed to be uniform (Table 3). The range of uniform distributions were determined by adding plus/minus to the values determined by the differential evolution method for each site (Table S1). The pre-constraint of *a priori* distribution effectively reduces computational costs without decreasing model performance and improves model convergence, based on a preliminary analysis. For constraining the behavior that root density must be higher in the surface layer than the deep layer in the Bayesian optimization, the thickness of the surface layer and root distribution were determined without results from the mathematical optimization. For the parameter optimization, we did not assume the hierarchy in the statistical model.

The posterior distributions of the parameters were estimated using a Markov Chain Monte Carlo (MCMC) method with the No-U-Turn Sampler (NUTS). NUTS is an extension of Hamiltonian Monte Carlo and provides very effective samples without requiring user intervention or costly tuning runs (Hoffman and Gelman, 2014). The efficiency of NUTS was more than 1,000 times that of Metropolis or Gibbs sampling. Posterior distributions of the parameters were estimated using four chains with 1,000 samples after warm-up based on 1,000 sampling. Bayesian inference was performed using the PyStan library (version 2.19.1.1). Owing to a complex and multimodal parameter space, consistent solutions from each chain were not obtained or some chains were not converged for some sites. In this case, we used results from chains that estimated the lowest model errors. The conservative treatment was required because bad chains seem to converge to local minima rather than to mathematically meaningful multimodal distributions and the problem was not fixed using different *a priori*, different initial values or further sampling. The trace plots and probability density functions for all parameters in all sites are shown in Fig. S1, which shows that at least two chains were well converged. Convergence of MCMC was evaluated by the Gelman–Rubin method with the potential scale reduction factor (PSRF), which showed that
all parameters for all sites were well converged (PSRF < 1.05) except slightly high PSRF for two parameters for US-Uaf (PSRF < 1.12; Table S3). Computational costs of the Bayesian inference ranged from 0.35 hours to 2.5 days per site with an average of 6.16 hours (Table S4).

Model parameters were estimated using daily CH₄ fluxes and the ancillary biophysical variables. Specifically, we used daily gap-filled CH₄ flux (FCH4_F), which contained only long data gaps (>2 months), and did not assume embedded functional relationships. In addition, we used FCH4_F_ANNarf when uncertainties in the neural network (FCH4_uncertainty) were less than absolute of FCH4_F_ANNarf. The use of gap-filled fluxes with low uncertainties could prevent propagating uncertainties associated with long-term gap-filling data into the parameter estimation. We also evaluated how the gap-filled data influenced modeled processes, where we eliminated data records where daily CH₄ emission contained more than 80% gaps in half-hourly data, in constraining the model. Apart from this issue, some high-latitude and rice paddy sites provided only growing-season fluxes, which hampered constraining the model for cold non-growing and fallow seasons, respectively. We also found that flux data for the first few days of a model run were important for constraining the initial CH₄ pools (i.e., initial conditions). Without the data, initial conditions were not well converged, and estimated dormant season emissions were unrealistic. Consequently, when FCH4_F was missing, we used the gap-filled CH₄ flux (FCH4_F_ANNarf) during the first six days of a model run and for the winter period (air temperature < -10 °C). The benefits of selectively using gap-filled data could outweigh the propagation of gap-filled errors, where unrealistic CH₄ emissions were not estimated.

The model constraints for each site were evaluated by RMSE normalized by mean, R, and
normalized standard deviation (SD) in daily CH₄ flux. For further interpreting and analyzing
modeled results, we eliminated unconstrained site-data where normalized RMSE was >0.9, R was
<0.6, normalized SD was <0.7, or normalized SD was >1.3.

The sensitivities to the forcing variables were performed using the models successfully constrained for each site. First, we applied perturbations to the inputs of: 1) 1 °C increase to the observed soil temperatures, 2) 10% increase in GPP and LAI, 3) 10 cm increase in WTD, and 4) 10 cm decrease in WTD with all other inputs held at measured conditions. Next, we examined the changes in modeled CH₄ emissions with unperturbed input (control experiment). We conducted the sensitivity analysis for sites spanning at least three years of data because the uncertainties are high in models constrained by short-term data (Ueyama et al., 2022).
To understand the variabilities in the estimated parameters across the sites, we applied principal component analysis (PCA) toward seven parameters: $p_{production}$, $Q_{10production}$, $p_{oxidation}$, $p_{pebullition}$, $p_{plant}$, $p_{diffusion-water}$, and $p_{pressure}$. The parameter for gas diffusion ($p_{diffusion-gas}$) was not included in the PCA because $p_{diffusion-gas}$ did not show a bell-shaped density curve at approximately half of the sites (Fig. S1). The parameters were first standardized with mean and SD and then compressed into two principal components (PC) using the scikit-learn library in python. We chose two principal components because they explained more than 70% of the variance in the parameters across the sites.

# 3. Results

## 3.1. *Model performance*

Across the 25 sites, 19 sites had reliable performance that satisfied the criterion for normalized RMSE, R, and normalized SD (section 2.3). According to the Taylor diagram (Fig. 2), model-data agreement was the best ($R > 0.9$) for RU-Ch2, FI-Lom, SE-Deg, FI-Sii, and CA-SCB. Among the accepted 19 sites, the median of normalized RMSE, R, and normalized SD were 0.59, 0.82, and 0.87, respectively. Except for the five sites with good model fit noted above, the model underestimated the SD of CH$_4$ flux, where the mean and SD of the normalized SD was $0.84 \pm 0.13$ across all sites. For the six sites excluded from subsequent analyses due to low performance (US-Sne, DE-Hte, DE-Zrk, DE-Sfn, US-Bzf, and US-Wpt), the mean seasonality was inconsistent between observations and models (Fig. 3), despite a moderate R and normalized RMSE. The low performance may represent a lack of important processes in the model and insufficient data to constrain the model. For example, US-Sne is a newly restored wetland and has a heterogeneous surface of open water and emergent vegetation, which make it difficult to constrain the processes based only on measured CH$_4$ fluxes for three years. Overall, there was no significant difference in the model performance in terms of wetland type and the number of years used for calibration.

In general, there were no obvious differences in modeled results with the optimized data containing fully gap-filled data or data when excluding days with >80% gaps. However, five sites (US-Sne, DE-Hte, DE-Zrk, FR-LGt, and NZ-Kop) did not meet the standard for a well constrained model with the non-gap-filled data (Fig. S2). The median of normalized RMSE, R, and normalized SD were 0.57, 0.83, and 0.90, respectively, in the model with the data not containing fully gap-filled data. The estimated CH$_4$ transport, production, and oxidation were also consistent among the
two models constrained with two data criteria, except for sites having low record numbers (e.g., RU-Ch2 and JP-Mse) (Fig. S3). Other results, including inter-site differences in CH₄ emission processes and sensitivity to biophysical drivers, were generally consistent among the two models constrained with two data criteria.

## 3.2. Estimated transport processes

Based on model results, plant-mediated transport and ebullition were more important pathways for CH₄ emissions than diffusive transport across sites (Fig. 3, 4; Table 2). In most cases, plant-mediated transport tended to be the major pathway for fen sites (72% ± 10%, n = 8; mean ± SD) and bog sites (55% ± 16%, n = 8; mean ± SD) (Fig. 4). Ebullition accounted for 27% ± 10% of the total emission for the fen sites and 26% ± 10% for the bog sites. In contrast, ebullition was estimated to be the major pathway at the two tundra sites (64% ± 4%) owing to shallow WTD (Fig. 4). Because the modeled plant-mediated transport increased with LAI, relative contribution of ebullition and/or diffusion was found high during periods of low LAI. When LAI was ≥20% of the annual peak, plant-mediated transport was the major pathway (70% ± 14%), except for three sites (RU-Ch2, US-Bzb, and KR-Crk) during the growing season (Fig. 3; Table 2). Diffusion was a minor pathway at most sites, but tended to be high in two marsh sites (US-Myb and US-Tw1) and a bog site (SE-Sto). For the three sites, the model predicted an anoxic surface layer, negligible oxidation, and high CH₄ concentrations in the surface layer at high WTD sites, allowing for surface diffusion. Since US-Myb was a restored wetland, the contribution of diffusion was approximately half of the CH₄ emissions in open water conditions (2010-2011) and then decreased to 31% ± 6% with the expansion of emergent vegetation from 2012 to 2018.

Based on cross-correlation analysis, CH₄ emissions lagged CH₄ production by 1-32 days (Fig. 5). There was more than a 30-day lag between CH₄ production and CH₄ emissions at US-Uaf. Lags tended to be, on average, longer in bogs (13 ± 10 days; n = 7; mean ± SD) than in fens (3 ± 2 days; n = 5), rice paddies (11 ± 3 days; n = 2), or tundra (6 ± 3 days; n = 2). Even in a long-lagged site (> 30 days for US-Uaf), the correlation between CH₄ production and CH₄ emission was good (R > 0.70), indicating that CH₄ production controlled temporal variations in CH₄ emission.

Inter-site variations in CH₄ production explained inter-site variations in CH₄ emissions (R² = 0.72; p = 0.01), except for sites where the ratio of oxidation to production was high (Fig. 6a). For sites with high oxidative fraction to production, CH₄ emissions were relatively low considering
their production (Fig. 6a). These sites with high oxidation generally exhibited low minimum WTD (Fig. 6b). CH₄ production and emission were positively correlated with soil temperature and GPP across the sites having low oxidation (Fig. 6c-f). This result is unexpected because the model was constrained in each site using temporal variations in the variables, as there was no assumption about inter-site variations in constraining the model. Based on the variable importance analysis using random forest regression, soil temperature and GPP almost equally explained the inter-site variations in CH₄ production. In contrast to production and oxidation, inter-site variations in three transport pathways did not correlate with CH₄ emissions.

### 3.3. *Estimated parameters*

Most parameters in our model were well converged (Table S3), but $p_{\text{diffusion-gas}}$ did not show a bell-shaped density curve with a single peak at 8 of the 19 sites (Fig. S1). Substrates for CH₄ production per GPP ($p_{\text{production}}$) were converged on the lower end of a *priori* range (median = 1.1 mmol m⁻² gC⁻¹ m²) over the 19 sites. The median and SD of $Q_{10}$ of CH₄ production was 3.7 ± 1.9, where there was a weak negative correlation between $p_{\text{production}}$ and $Q_{10}$ across the sites ($R^2 = 0.31$; $p < 0.01$). The maximum oxidation parameter was estimated to be in the middle of the prescribed upper and lower range at most sites. Estimated $p_{\text{ebullition}}$ and $p_{\text{plant}}$ were not correlated with contributions from ebullition and plant-mediated transport to CH₄ emission, respectively. Ebullition from 9 sites had a marginal sensitivity to pressure decline (<2% hPa⁻¹), where there was no correlation between $p_{\text{pressure}}$ and contributions of ebullition to the total emission across the sites. There was no significant difference ($p < 0.05$) in all parameters aggregated by aerenchymatous and moss vegetation.

Based on the PCA analysis, 77% of the variance in the parameters among the sites was compressed with two PCs (Fig. 7). The first PC represented a tradeoff of two parameters for CH₄ production between high $p_{\text{production}}$ and low $Q_{10}$ and *vice versa*, representing 61% of the parameter distribution across the sites. The second PC explained 16% of the distribution and represented a tradeoff between CH₄ production and transport through plants and gas diffusion. There were weak clusters for bog sites with relatively high $Q_{10}$, tundra sites with low transport parameters, and rice paddies with high transport parameters. No clusters were apparent for fen and marsh sites.

The thickness of the surface layer, $z_{\text{surf}}$, was the conceptual depth separating surface oxic and deeper anoxic layers, and thus negatively correlated to WTD for sites where minimum WTD was
below -0.1 m ($z_{surf}$ = -1.2 * WTD - 0.05 m; R² = 0.48; p=0.03; n = 10). The regression analysis showed that $z_{surf}$ was close to minimum WTD. In contrast, there was no significant trend in the surface layer thickness for sites with high mean annual WTD (> -0.1 m). For sites with high WTD (i.e., always above the ground surface), the thickness of the soil layers did not control the degree of redox conditions for the two layers because the surface layer was always anaerobic.

### 3.4. Sensitivity to biophysical variables

Based on the sensitivity analysis, CH₄ emissions increased by 9.6% or 3.5 gCH₄ m⁻² yr⁻¹ (median relative increase), with 10% increase in GPP across the sites, with the increases higher in the sites with high annual soil temperatures (Fig. 8a). The sensitivity analysis was performed on sites that had at least three years of data (14 sites) among the 19 sites. The sensitivities aggregated for high or low WTD sites (sites having mean water table position above or below the ground surface) indicated that the relative increases in CH₄ emissions did not differ significantly between the two WTD classes (p = 0.35 in Welch's t test; inset in Fig. 8a).

The 1 °C increases in soil temperatures increased CH₄ emissions by 6.6% or 2.5 gCH₄ m⁻² yr⁻¹ (median relative increase) (Fig. 8b). The increases were similar in magnitude to those from the 10% increase in GPP. Compared with the sensitivity to GPP, the increased magnitudes appeared to not be clearly related to the mean annual soil temperatures and WTD, likely because temperature sensitivity ($Q_{10}$) for CH₄ production differed by site. The increases in CH₄ emissions also did not differ significantly between the two WTD classes (p = 0.80; inset in Fig. 8b).

The increase in CH₄ emissions with 1 °C increases were lower than those estimated based on an empirical $Q_{10}$ relationship between daily mean soil temperature and CH₄ emissions (Fig. 9). Eight of the 14 sites were estimated to have higher CH₄ emission sensitivity using the empirical $Q_{10}$ model than iPEACE. Across all 14 sites, the relative increases in CH₄ emissions tended to be higher in the empirical $Q_{10}$ model (12%) than the iPEACE model (8%) across the sites (p = 0.12) (US-Uaf was not included in relative changes in emission owing to the small magnitude in emission).

Decreased CH₄ emissions associated with a 10 cm decrease in WTD were greater than increased CH₄ emissions with a 10 cm increase in WTD (Fig. 8c, d). A decrease in WTD decreased CH₄ emissions at most sites and *vice versa*, where the median changes by the decrease and increase in WTD were -31% and +6.5%, respectively. A site with a WTD permanently well above the
ground surface (US-Myb) did not exhibit significant responses to changing WTD, as WTD always remained above the surface. The relative changes in CH₄ emissions did not differ significantly between sites with low and high WTD with 10 cm increases in WTD (inset in Fig. 8c; p = 0.34) and 10 cm decrease in WTD (inset in Fig. 8d; p = 0.15).

There were two mechanisms for reduced CH₄ emissions by decreased WTD. The first mechanism is associated with changes in the frequency with which the surface layer becomes oxic conditions. In this mechanism, CH₄ production from the surface layer decreases when the WTD decreases with the perturbed input mostly fluctuating within the surface layer throughout the year. The second mechanism is related to the long-lasting change in redox conditions in the deep layer. We argue that reduced anaerobic conditions in the deep layer, which was rarely affected by oxic conditions with the unperturbed WTD, but was affected by the perturbed decrease in WTD. Owing to the loss of anaerobic conditions, CH₄ in the deep layer was consumed through oxidation; thus, the effects were relatively long-lasting until CH₄ concentrations built-up again. The median decrease in CH₄ production was -6.9 gCH₄ m⁻² yr⁻¹, and median increase in CH₄ oxidation was 12.9 gCH₄ m⁻² yr⁻¹, indicating that the second mechanism was the major process responsible for the reduction in CH₄ emissions. As an exceptional response examined at NZ-Kop, the decreased WTD could change sustained anoxic conditions to oxic conditions in the deep layer, resulting in decreased CH₄ production, reduced CH₄ pool, and finally decreased oxidation.

## 4. Discussion

The estimated processes for CH₄ emissions provide meaningful insights for interpreting observed data and estimating sensitivities to the forcing variables. The current analysis aims to shed light on the relative importance of processes involved in CH₄ production, transport, and oxidation across 25 freshwater wetland sites in temperate, boreal, and Arctic regions. The observed data included in the FLUXNET-CH₄ database were used to constrain a process-based model which has a similar structure used in previous modeling studies (Walter and Heimann, 2000; Wania et al., 2010; Riley et al., 2011). Flux partitioning is typically applied to net CO₂ fluxes for estimating GPP and ecosystem respiration (Reichstein et al., 2005), and has successfully provided deeper insights on their biotic and abiotic controls (Jung et al., 2017; Mahecha et al., 2010). Compared to the partitioning of CO₂ fluxes, more complex models are required to explain wetland CH₄ emissions and partition net CH₄ flux observations (Wania et al., 2010; Riley et al., 2011; Grant et al., 2019;
541 Chen 2021). Partitioned CH₄ fluxes can be useful for evaluating inter-site differences in fluxes
542 (Figs. 3, 4), time lags between surface emissions and production (Fig. 5), different responses of
543 CH₄ processes (e.g., production, oxidation, and transport) to biophysical variables (Figs. 6 and 8),
544 and model parameterizations (Fig. 7). Key processes and parameters estimated in this study need
545 to be better constrained with further long-term observations and different data streams.
546

547 4.1. Inter-site variations in estimated processes
548 The inter-site variations in CH₄ emissions were found to be primarily associated with those in CH₄
549 production rather than those in oxidation and transport (Fig. 6), especially for sites with high WTD
550 and low CH₄ oxidation. These results could explain the correlation of annual CH₄ emissions with
551 mean annual air or soil temperature across global wetlands in the FLUXNET-CH₄ database (Knox
552 et al., 2019; Delwiche et al., 2021), where temperature was found to be an important driver of
553 methanogenesis substrates (Chang et al. 2021) and CH₄ production (Yvon-Durocher et al., 2014).
554 In contrast, oxidation increased with decreasing WTD (Fig. 6b), resulting in oxidation as the
555 second most important process for explaining inter-site variations in CH₄ emissions. These results
556 are also consistent with global syntheses, which showed that a positive correlation between CH₄
557 emissions and WTD was only detected in sites with relatively low WTD (i.e., mean annual WTD
558 was below the soil surface) (Knox et al., 2019, 2021).
559

560 Transport processes were estimated to regulate the time-lag between CH₄ production and
561 emissions (Fig. 5), albeit we found no significant effect on total CH₄ emissions because annual
562 emissions were mainly controlled by CH₄ production (Fig. 6). The lag between production and
563 emission occurred due to the time required to increase the CH₄ concentrations to drive CH₄
564 transport. The lag of CH₄ emissions to soil temperature or GPP was reported in studies using
565 FLUXNET-CH₄ (Chang et al., 2019; Delwiche et al., 2021; Knox et al., 2021; Yuan et al., 2022).
566 For example, Knox et al. (2021) estimated that on average CH₄ emissions lagged soil temperature
567 and GPP by 5.4 days and 20.7 days, respectively, across wetlands globally. The lag between CH₄
568 emission and production (Fig. 5) nonetheless partly explained the lag between emission and
569 biophysical variables, as time is required for building up sufficient CH₄ concentrations driving
570 CH₄ emissions.
571

571 4.2. Sensitivities of CH₄ emissions to biophysical drivers
The estimated sensitivity of CH₄ emissions to GPP (Fig. 8a, b) indicates the importance of substrate availability. A strong relationship between net ecosystem production and CH₄ emissions was previously reported across wetlands extending from subarctic peatlands to subtropical marshes associated with substrate availability (Whiting and Chanton, 1993). The estimated sensitivity occurred because CH₄ production in iPEACE was driven by GPP and soil temperature, reflecting the concept that increased GPP will increase substrate availability and thereby CH₄ emissions. The strong relationship with GPP (Fig. 8a) was unexpected, however, because the sensitivity to GPP (p̂<sub>production</sub>) was calibrated in each site and thus was expected to show high variability among the sites. It is worth noting that the estimated sensitivity to GPP might be caused by model assumptions. The model assumed that substrates for CH₄ were only provided by GPP, but old peat previously fixed is also known to be a substrate for CH₄ production (Chasar et al., 2000). Substrates from recent primary production and peat organic carbon should be incorporated into future modeling with iPEACE.

Based on our sensitivity analyses, CH₄ emissions were sensitive to a decrease in WTD for most sites (Fig. 8). The most important mechanism associated with decreased WTD was increased oxidation at the deep layer. Because the buildup of the CH₄ pool after loss of anaerobic conditions is time consuming, the effects can be long-lasting. This result is consistent with previous studies. Brown et al. (2014) indicated that a long recovery time was required for CH₄ emissions after rewetting following a drop in WTD at a site where the mean WTD was below the surface. They proposed a reason for the long recovery time as breaking the critical zone for CH₄ emissions by low WTD conditions. Simultaneously, when increased WTD resulted in aerobic layers switching to anaerobic conditions, CH₄ emissions increased, but the response was smaller than those to a decreasing WTD. This difference occurs because increased WTD increased the frequency of anoxic conditions at the surface layer, but the surface layer was still susceptible to oxic conditions even with perturbation increase in WTD, resulting in limited increases in CH₄ concentration. When deep soil remained anaerobic owing to shallow WTD, increases in soil temperature and GPP were equally important drivers of CH₄ emissions through kinetics and substrate availability, respectively (Fig. 8).

## 4.3. Comparison of estimated processes to observations from previous studies

Estimated transport flux was compared to EC measurements at various sites (Table 1). The high
ebullition (50% of total emissions) was measured with chamber measurements at JP-Bby (Tokida et al., 2007a, b), which was consistent with the current study. Windham-Myers et al. (2018) measured ebullition with a static chamber during five days in summer at US-Tw1, and ebullition contributions to the total emission (10-30%) were comparable to those by the current study (26%). In contrast, plant-mediated transport estimated with chambers for FI-Sii (31%) and FI-Si2 (21%) was smaller than our model estimates (91% for FI-Sii and 67% for FI-Si2). However, Susiluoto et al. (2018) reported contributions similar to the current study based on process-based models also constrained using EC data (75-95%) for FI-Sii. Kwon et al., (2017) measured lower contributions of plant-mediated transport (25%) and ebullition (2%) in RU-Ch2 than the model estimates. McNicol et al. (2017) measured ebullition by bubble traps (< 1.3%) and diffusion by dissolved CH₄ (< 4.1%) from open water bodies within the flux footprint US-Myb, values which are smaller than the current estimates (18% and 24%, respectively). One reason for the inconsistency might be the spatial heterogeneity at US-Myb. Their study did not consider areas of emergent vegetation where contributions by ebullition can be higher (Villa et al., 2021). Hwang et al. (2020) estimated smaller ebullition (10-17%) than the current study (61%) based on the wavelet analysis of EC data at KR-Crk. For KR-Crk data in the FLUXNET-CH₄ database, WTD under drainage was provided as 0 cm; thus, the model predicted more saturated conditions at the surface than the actual conditions, resulting in higher contributions by ebullition.

Based on the site-scale validation, iPEACE estimates were consistent with production, ebullition, or diffusive flux observations at two sites, but inconsistent with observations from four sites. A comprehensive validation of estimated transport fluxes is challenging at the site scale owing to limited sites with both EC data and process studies available at the same location (Table 1). Furthermore, no study has in-situ measured the three transport fluxes simultaneously, resulting in uncertainties in how transport fluxes by process studies are consistent with CH₄ emissions measured with EC towers. Plant-mediated transport could be the priority for in-situ measured transport fluxes to validate CH₄ emissions, since it was estimated to be a major pathway in most sites (Table 2) and in other modeling studies (Table 1). Differences in spatial representativeness between EC towers and process studies could also contribute to inconsistencies.

Our estimated wetland CH₄ emissions were within the range of those measured or predicted with process-based models regardless of difficulties in direct comparisons at the site scale. Although the contributions of each transport flux were highly variable among previous studies
(Table 1), plant-mediated transport and ebullition tended to be major transport pathways, consistent with our current estimates (Fig. 4). Previous models also estimated plant-mediated transport as the major pathway (Table 1), although the VISIT model predicted ebullition as the major pathway for Arctic wetlands (Ito, 2019). In contrast, iPEACE tended to estimate higher contributions from ebullition and lower contributions from diffusion. This difference could be caused by the assumption that ebullition occurs when WTD is greater than 10 cm below the ground (Stanley et al., 2019). The contribution of plant-mediated transport was similar to previous modeling studies because of similar model structure, but tended to be higher than measurements (Table 1). Rhizospheric oxidation (Bansal et al., 2020; Korrensalo et al., 2022) is a potential reason for low CH₄ emissions through vegetation, which was not considered in the current version of iPEACE.

## 4.4. Toward refined parameterizations

Based on the PCA (Fig. 7), modeling wetland CH₄ emissions could be improved with refined parameterization and representation of CH₄ production, plant-mediated transport, and diffusion through water. The importance of parameterizations for production and plant-mediated transport was also estimated in a study constraining a global CH₄ model with observed CH₄ emissions at 16 wetland sites (Müller et al., 2015). The high explanatory power in the first PC by the production parameters suggests that CH₄ production was important for inter-site variations in CH₄ emissions. Considering the structure of iPEACE, sites with high $p_{production}$ could be more limited by substrate availability, whereas sites with high $Q_{10}$ could be more limited by kinetics. The second PC explained CH₄ emissions that are limited by production and/or transport. A similar trade-off between parameters of production and plant-mediated transport was also inferred in an optimized process-based model (Salmon et al., 2022). These results suggest that a model for explaining variabilities in parameters of production and plant-mediated transport across wetlands is needed for refined simulations rather than determining one set of parameters.

## 4.5. Next steps in modeling wetland CH₄ emissions

The estimated processes were the most likely processes for explaining observed CH₄ emissions under the model structure of iPEACE (section 2.2), suggesting that careful interpretation is required. iPEACE considers important processes to explain CH₄ emissions that have been
incorporated in some previous modeling studies (Walter and Heimann, 2000; Wania et al., 2010; Riley et al., 2011). However, definitions and formulations of CH₄-related processes are often different among models (Melton et al., 2013). For instance, iPEACE does not include processes included in more mechanistic models (e.g., Salmon et al., 2022; Susiluto et al., 2018). We need to better define processes in the model and to validate modeled processes, where the model-data fusion could be useful to bridge model and observation.

To improve our understanding of CH₄ emissions from wetlands, future improvements are possible with increased availability of EC data, additional observations, and by incorporating more processes into the model. First, in-situ observations of transport fluxes and production parameters with incubations would be useful to constrain the model because Bayesian optimization can effectively incorporate the additional constraints from observations. Second, more long-term data are required for better constraining the model. The period of the current study ranged from one to nine years with a median of four years. Ueyama et al. (2022) indicated that long-term data (e.g., >3 years) effectively constrained the partitioned fluxes. Furthermore, we did not focus on tree-dominated wetlands (e.g., swamps) owing to the importance of unaccounted processes, such as CH₄ transport to the atmosphere by tree stems (Pangala et al., 2012), or from O₂ transport to the rhizosphere via aerial roots (Purvaja et al., 2004). In this study, we predicted O₂ concentration in the soil based on WTD, but the relationship between O₂ concentrations and WTD is complex (Ueyama et al., 2020). Thus, measurements of WTD and O₂ concentrations are strongly recommended for evaluating CH₄ emissions in wetlands. The current model considers a 1 m thick soil, but anaerobic peat deeper than 1 m could play a role in CH₄ emissions (Tokida et al., 2007a, b; Peltola et al., 2018). Since flux tower measurements did not continuously monitor the O₂ and CH₄ concentrations in the deep peat, constraining processes at the deep peat were difficult in this study. Finally, refined modeling wetland CH₄ emission will be possible by evaluating how partitioned emissions are consistent across different models constrained with the same data.

The Bayesian inference in this study might be improved after considering the outlined limitations. We did not obtain reliable results for 6 of 23 sites. The inability could be caused by lack of important processes, but might be resolved with improved mathematical techniques. The error distribution was assumed with Gaussian distribution, which lacked the ability to fit long-tail, such as data containing outliers. Use of other error distributions might improve posterior inference (Hamura et al., 2022). For 12 sites, at least one chain was not well converged (Fig. S1), possibly
696 due to a problem of slow convergence associated with complex multimodal parameter distributions.
697 Introducing Extended Ensemble Monte Carlo (Iba, 2001), such as the replica exchange method,
698 could improve convergence. The techniques for complex parameter distributions could improve
699 the parameter optimization, where some parameters in the current study hit the range of prior
700 distributions (Fig. S1) possibly owing to the equifinality problem (Schulz et al., 2001).

701 Acknowledgements

702 We acknowledge support from the John Wesley Powell Center for Analysis and Synthesis of the
703 U.S. Geological Survey (USGS)(“Wetland FLUXNET Synthesis for Methane” working group,
704 https://www.usgs.gov/centers/john-wesley-powell-center-for-analysis-and-
705 synthesis/science/wetland-fluxnet-synthesis) and the USGS Ecosystem Mission Area, Land
706 Change Science programs. We thank Dr. Eric J Ward for providing useful comments on the
707 manuscript. ESE was supported by the grants from the Arctic Observatory Program of the National
708 Science Foundation (grant numbers 1936752, 1503912, 1107892) and by the US Geological
709 Survey, Research Work Order 224 to the University of Alaska Fairbanks, the Bonanza Creek
710 Long-Term Ecological Research Program funded by the National Science Foundation (NSF DEB-
711 1026415, DEB-1636476) and the NSF Long-Term Research in Environmental Biology Program
712 (NSF LTREB 2011276). MU was supported by the Arctic Challenge for Sustainability II (ArCS
713 II; JPMXD1420318865) and JSPS KAKENHI (20K21849). IM was supported by ICOS-Finland,
714 Academy of Finland project N-PERM and Horizon Europe project GreenFeedBack. MK was
715 supported by the Rural Development Administration (PJ014892022022). TS was supported by the
716 Helmholtz Association of German Research Centres (VH-NG-821). DE-Zrk is a Terrestrial
717 Environmental Observatories Network (TERENO) site. SD was supported by the U.S. Geological
718 Survey, Ecosystems Mission Area, Land Change Science Program. SG was supported by the SNO
719 Tourbières, CNRS-INSU. OS was supported by the Canada Research Chairs, Canada Foundation
720 for Innovation Leaders Opportunity Fund, and Natural Sciences and Engineering Research
721 Council Discovery Grant Programs. WJR and KYC were supported by the Reducing
722 Uncertainties in Biogeochemical Interactions through Synthesis and Computation (RUBISCO)
723 Scientific Focus Area, Office of Biological and Environmental Research of the U.S. Department
724 of Energy Office of Science. YR was supported by the Ministry of Environment of Korea
725 (2022003640002). Lawrence Berkeley National Laboratory (LBNL) is managed by the University
726 of California for the U.S. Department of Energy under contract DE-AC02-05CH11231, California
727 Department of Water Resources, CA Fish and Wildlife, and US DOE Ameriflux. Any use of trade,
728 firm, or product names is for descriptive purposes only and does not imply endorsement by the
729 U.S. Government.

730
731 **Data availability statement**
732 The data that support the findings of this study are available in the FLUXNET-CH₄ Community
733 Product, available at https://fluxnet.org/data/fluxnet-ch4-community-product/. DOIs for
734 individual site data are provided in Table 2. The iPEACE source code is available upon request to
735 the authors.

736
737 **References**
738 Baldocchi, D. (2014). Measuring fluxes of trace gases and energy between ecosystems and the atmosphere
739 – the state and future of the eddy covariance method. Global Change Biology, 20, 3600-3609.

740 Bansal, S., Johnson, O. F., Meier, J., & Zhu, X. (2020). Vegetation affects timing and location of wetland
741 methane emissions. Journal of Geophysical Research: Biogeosciences, 125, e2020JG005777.

742 Bohn T. J., Melton, J. R., Ito, A., Kleinen, T., Spahni, R., Stocker, B. D., Zhang, B., Zhu, X., Schroeder,
743 R., Glagolev, M. V., Maksyutov, S., Brovkin, V., Chen, G., Denisov, S. N., Eliseev, A. V., Gallego-
744 Sala, A., McDonald, K. C., Rawlins, M. A., Riley, W. J., ... Kaplan, J. O. (2015). WETCHIMP-WSL:
745 intercomparison of wetland methane emissions models over West Siberia. Biogeosciences, 12, 3321-
746 3349.

747 Bridgham, S. D., Cadillo-Quiroz, H., Keller, J. K., & Zhuang, Q. (2013). Methane emissions from wetlands:
748 biogeochemical, microbial, and modeling perspectives from local to global scales. Global Change
749 Biology, 19, 1325-1346.

750 Brown, M. G., Humphreys, E. R., Moore, T. R., Roulet, N. T., & Lafleur, P. M. (2014). Evidence for a
751 nonmonotonic relationship between ecosystem-scale peatland methane emissions and water table depth.
752 Journal of Geophysical Research Biogeosciences, 119, 826-835.

753 Butterbach-Bahl, K., Papen, H., & Rennenberg, H. (1997). Impact of gas transport through rice cultivars
754 on methane emission from rice paddy fields. Plant, Cell and Environment, 20, 1175-1183.

755 Campbell, D., & Goodrich, J. (2020). FLUXNET-CH₄ NZ-Kop Kopuatai. New Zealand.
756 https://doi.org/10.18140/FLX/1669652

757 Castro-Morales, K., Kleinen, T., Kaiser, S., Zaehle, S., Kittler, F., Kwon, M. J., Beer, C., & Göckede, M.
758 (2018). Year-round simulated methane emissions from a permafrost ecosystem in Northeast Siberia,
759 Biogeosciences, 15, 2691-2722.

760 Chang, K.-Y., Riley, W., Knox, S., Jackson, R., McNicol, G., Poulter, B., Aurela, M., Baldocchi, D., Bansal, S., Bohrer, G., Campbell, D., Cescatti, A., Chu, H., Delwiche, K., Desai, A., Euskirchen, E., Friborg, T., Gockede, M., Helbig, M., ... Zona, D. (2021). Substantial hysteresis in emergent temperature sensitivity of global wetland CH₄ emissions. Nature Communications, 12, 2266.

764 Chang, K.-Y., Riley, Crill, P. M., Grant, R. E., & Saleska, S. R. (2020). Hysteretic temperature sensitivity of wetland CH₄ fluxes explained by substrate availability and microbial activity. Biogeosciences, 17, 5849-5860.

767 Chasar, L. S., Chanton, J. P., Glaser, P. H., Siegel, D. I., & Rivers, J. S. (2000). Radiocarbon and stable carbon isotopic evidence for transport and transformation of dissolved organic carbon, dissolved inorganic carbon, and CH₄ in a northern Minnesota peatland. Global Biogeochemical Cycles, 14, 1095-1108.

771 Chen, J. (2021). *Biophysical Models and Applications in Ecosystem Analysis*. Michigan State University Press.

773 Chen, J., & Chu, H. (2020). FLUXNET-CH₄ US-WPT Winous Point North Marsh. United States. https://doi.org/10.18140/FLX/1669702

775 Chen, J., Jönsson, P., Tamura, M., Gu, Z., Matsushita, B., & Eklundh, L. (2004). A simple method for reconstructing a high-quality NDVI time-series data set based on the Savitzky-Golay filter. Remote Sensing of Environment, 91, 332-344.

778 Chu, H., Chen, J., Gottgens, J. F., Ouyang, Z., John, R., Czajkowski, K., & Becker, R. (2014). Net ecosystem methane and carbon dioxide exchanges in a Lake Erie coastal marsh and a nearby cropland. *Journal of Geophysical Research: Biogeosciences*, 119(5), 722-740.

781 Conrad, R. (2009). The global methane cycle: recent advance in understanding the microbial processes involved. Environmental Microbiology Reports, 1, 285-292.

783 Delwiche, K. B., Knox, S. H., Malhotra, A., Fluet-Chouinard, E., McNicol, G., Feron, S., Ouyang, Z., Papale, D., Trotta, C., Canfora, E., Cheah, Y.-W., Christianson, D., Alberto, M. C. R., Alekseychik, P., Aurela, M., Baldocchi, D., Bansal, S., Billesbach, D. P., Bohrer, G., ... Jackson, R. B. (2021). FLUXNET-CH₄: A global, multi-ecosystem dataset and analysis of methane seasonality from freshwater wetlands. Earth System Science Data, 13, 3607-3689.

788 Desai, A. R. (2020). FLUXNET-CH₄ US-Los Lost Creek. United States. https://doi.org/10.18140/FLX/1669682

790 Dorodnikov, M., Knorr, K.-H., Kuzyakov, Y., & Wilmking, M. (2011). Plant-mediated CH₄ transport and contribution of photosynthates to methanogenesis at a boreal mire: a ¹⁴C pulse-labeling study. Biogeosciences, 8, 2365-2375.
793 Eichelmann, E., Knox, S., Rey-Sanchez, A. C., Valach, A., Sturtevant, C., Szutu, D., Verfaillie, J., & 794 Baldocchi, D. (2020). FLUXNET-CH₄ US-Tw4 twitchell east end wetland. United States. 795 https://doi.org/10.18140/FLX/1669698

796 Euskirchen, E., Bret-Harte, M., & Edgar, C. (2020). FLUXNET-CH₄ US-ICs Imnavait Creek Watershed 797 Wet Sedge Tundra, United States. https://doi.org/10.18140/FLX/1669678

798 Euskirchen, E., & Edgar, C. (2020a). FLUXNET-CH₄ US-BZF Bonanza Creek Rich Fen, United States. 799 https://doi.org/10.18140/FLX/1669669

800 Euskirchen, E., & Edgar, C. (2020b). FLUXNET-CH₄ USBZB Bonanza Creek Thermokarst Bog, United 801 States. https://doi.org/10.18140/FLX/1669668

802 Glaser, P. H., Morin, C. P., Rosenberry, D. O., Siegel, D. I., Ruud, O., Chasar, L. I., & Reeve, A. S. (2004). 803 Surface deformations as indicators of deep ebullition fluxes in a large northern peatland. Global 804 Biogeochemical Cycles 19, GB1003, doi:10.1029/2003GB002069.

805 Goeckede, M. (2020). FLUXNET-CH₄ RU-Ch2 Chersky reference, Russian Federation. 806 https://doi.org/10.18140/FLX/1669654

807 Göckede, M., Kittler, F., & Schaller, C., (2019). Quantifying the impact of emission outbursts and non- 808 stationary flow on eddy covariance CH₄ flux measurements using wavelet techniques. Biogeosciences, 809 16, 3113-3131.

810 Gogo, S., Guimbaud, C., Laggoun-Défarge, F., Catoire, V., & Robert, C. (2011). In situ quantification of 811 CH₄ bubbling events from a peat soil using a new infrared laser spectrometer. Journal of Soils 812 Sediments, 11, 545-551.

813 Grant, R. F., Mekonnen, Z. A., Riley, W. J., Arora, B., & Torn, M. S., (2019). Modeling climate change 814 impacts on an arctic polygonal tundra: 2. changes in CO₂ and CH₄ exchange depend on rates of 815 permafrost thaw as affected by changes in vegetation and drainage. Journal of Geophysical Research: 816 Biogeosciences, 124, 1323-1341.

817 Hamura, Y., Irie, K., Sugasawa, S., (2022). Log-regularly varying scale mixture of normals for robust 818 regression. Computational Statistics & Data Analysis, 173, 107517.

819 Hoffman, M. D., & Gelman, A., (2014). The No-U-Turn sampler: adaptively setting path lengths in 820 Hamiltonian Monte Carlo. Journal of Machine Learning Research, 15, 1593-1623.

821 Hwang, Y., Ryu, Y., Huang, Y., Kim, J., Iwata, H., & Kang, M. (2020). Comprehensive assessments of 822 carbon dynamics in an intermittently-irrigated rice paddy. Agricultural and Forest Meteorology, 285- 823 286, 107933.

824 Iba, Y., (2001). Extended Ensemble Monte Carlo. International Journal of Modern Physics C, 12, 623-656.

825 Ito, A. (2019). Methane emission from pan-Arctic natural wetlands estimated using a process-based model, 826 1901–2016. Polar Science, 21, 26-36.
827 Iwata, H. (2020). FLUXNET-CH₄ JP-Mse Mase rice paddy field. Japan.
https://doi.org/10.18140/FLX/1669647

829 Iwata, H., Hirata, R., Takahashi, Y., Miyabara, Y., Itoh, M., & Iizuka, K. (2018). Partitioning eddy-covariance methane fluxes from a shallow lake into diffusive and ebullitive fluxes. Boundary-Layer Meteorology, 169, 416-428.

832 Iwata, H., Ueyama, M., & Harazono, Y. (2020). FLUXNET-CH₄ US-Uaf University of Alaska, Fairbanks. United States. https://doi.org/10.18140/FLX/1669701

834 Jackson, R. B., Saunois, M., Bousquet, P., Canadell, J. G., Poulter, B., Stavert, A. R., Bergamaschi, P., Niwa, Y., Segers, A., & Tsuruta, A. (2020). Increasing anthropogenic methane emissions arise equally from agricultural and fossil fuel sources. Environmental Research Letters, 15, 071002

837 Jacotot, A., Gogo, S., & Laggoun-Défarge, F. (2020). FLUXNET-CH₄ FR-LGt La Guette, France. https://doi.org/10.18140/FLX/1669641

839 Jung, M., Reichstein, M., Schwalm, C. R., Huntingford, C., Sitch, S., Ahlström, A., Arneth, A., Camps-Valls, G., Ciais, P., Friedlingstein, P., Gans, F., Ichii, K., Jain, A. K., Kato, E., Papale, D., Poulter, B., Raduly, B., Rödenbeck, C., Tramontana, G., ... Zeng, N. (2017). Compensatory water effects link yearly global land CO₂ sink changes to temperature. Nature, 541, 516-520.

843 Kajiura, M., & Tokida, T. (2021). Quantifying bubbling emission (ebullition) of methane from a rice paddy using high-time-resolution concentration data obtained during a closed-chamber measurement. Journal of Agricultural Meteorology, 77, 245-252.

846 Karofeld, E., & Tónisson, H. (2014). Spatio-temporal changes in bog pool bottom topography – temperature effect and its influence on pool development: an example from a raised bog in Estonia. Hydrological Processes 28:958-968.

849 Knox, S. H., Jackson, R. B., Poulter, B., McNicol, G., Fluet-Chouinard, E., Zhang, Z., Hugelius, G., Bousquet, P., Canadell, J. G., Saunois, M., Papale, D., Chu, H., Keenan, T. F., Baldocchi, D., Torn, M. S., Mammarella, I., Trotta, C., Aurela, M., Bohrer, G., ... Zona, D. (2019). FLUXNET-CH₄ synthesis activity: objective, observations, and future directions. Bulletin of the American Meteorological Society, 100, 2607-2632.

854 Knox, S. H., Bansal, S., McNicol, G., Schafer, K., Sturtevant, C., Ueyama, M., Valach, A. C., Baldocchi, D., Delwiche, K., Desai, A. R., Euskirchen, E., Liu, J., Lohila, A., Malhotra, A., Melling, L., Riley, W., Runkle, B. R. K., Turner, J., Vargas, R., ... Jackson, R. B. (2021). Identifying dominant environmental predictors of freshwater wetland methane fluxes across diurnal to seasonal time scales. Global Change Biology, 27, 3582-3604.

859 Koebsch, F., & Jurasinski, G. (2020). FLUXNET-CH₄ DE-Hte Huetelmoor. Germany. https://doi.org/10.18140/FLX/1669634
861 Koebsch, F., Jurasinski, G., Koch, M., Hofmann, J., & Glatzel, S. (2015). Controls for multi-scale temporal variation in ecosystem methane exchange during the growing season of a permanently inundated fen. 862 Agricultural and Forest Meteorology, 204, 94-105. 863

864 Korrensalo, A., Mammarella, I., Alekseychik, P., Vesala, T., & Tuittila, E-S. (2022). Plant mediated methane efflux from a boreal peatland complex. Plant Soil, 471, 375-392. 865

866 Kutzback, L., Wagner, D., & Pfeiffer, E.-M., (2004). Effect of microrelief and vegetation on methane emission from wet polygonal tundra, Lena Delta, Northern Siberia. Biogeochemistry 69:341-362. 867

868 Kwon, M. J., Beulig, F., Ilie, I., Wildner, M., Küsel, K., Merbold, L., Mahecha, M. D., Zimov, N., Zimov, S. A., Heimann, M., Schuur, E. A. G., Kostka, J. E., Kolle, O., Hilke, I., & Göckede, M. (2017). Plants, microorganisms, and soil temperatures contribute to a decrease in methane fluxes on a drained Arctic floodplain, Global Change Biology, 23, 2396-2412. 869 870 871

872 Lasslop, G., Reichstein, M., Papale, D., Richardson, A. D., Arneth, A., Barr, A., Stoy, P., & Wohlfahrt, G. (2010). Separation of net ecosystem exchange into assimilation and respiration using a light response curve approach: critical issues and global evaluation, Global Change Biology, 16, 187-208. 873 874

875 Lohila, A., Aurela, M., Tuovinen, J.-P., Laurila, T., Hatakka, J., Rainne, J., & Mäkelä, T. (2020). FLUXNET-CH₄ FI-Lom Lompolojankka. Finland. https://doi.org/10.18140/FLX/1669638 876

877 Ma, S., Jiang, J., Huang, Y., Shi, Z., Wilson, R. M., Ricciuto, D., Sebestyen, S. D., Hanson, P. J., & Luo, Y. (2017). Data-constrained projections of methane fluxes in a Northern Minnesota peatland in response to elevated CO₂ and warming. Journal of Geophysical Research: Biogeosciences, 122, 2841-2861. 878 879 880

881 Mahecha, M. D., Reichstein, M., Carvlhais, N., Lasslop, G., Lange, H., Seneviratne, S. I., Vargas, R., Ammann, C., Arain, M. A., Cescatti, A., Janssens, I. A., Migliavacca, M., Montagnani, L., & Richardson, A. D. (2010). Global convergence in the temperature sensitivity of respiration at ecosystem level. Science, 329, 838-840. 882 883 884

885 Männistö, E., Korrensalo, A., Alekseychik, P., Mammarella, I., Peltola, O., Vesala, T., & Tuittila, E.-S., (2019). Multi-year methane ebullition measurements from water and bare peat surfaces of a patterned boreal bog. Biogeosciences 16:2409-2421. 886 887

888 Matthes, J. H., Sturtevant, C., Oikawa, P., Chamberlain, S. D., Szutu, D., Ortiz, A. A., Verfaillie, J., & Baldocchi, D. (2020). FLUXNET-CH₄ US-Myb Mayberry Wetland. United States. https://doi.org/10.18140/FLX/1669685 889 890

891 McNicol, G., Sturtevant, C. S., Knox, S. H., Dronova, I., Baldocchi, D. D., & Silver, W. L. (2017). Effects of seasonality, transport pathway, and spatial structure on greenhouse gas fluxes in a restored wetland. 892 Global Change Biology 23:2768-2782. 893

894 Melton, J. R., Wania, R., Hodson, E. L., Poulter, B., Ringeval, B., Spahni, R., Bohn, T., Avis, C. A.,
895 Beerling, D. J., Chen, G., Eliseev, A. V., Denisov, S. N., Hopcroft, P. O., Lettenmaier, D. P., Riley, W.
896 J., Singarayer, J. S., Subin, Z. M., Tian, H., Zürcher, S., ... Kaplan, J. O. (2013). Present state of global
897 wetland extent and wetland methane modelling: conclusions from a model inter-comparison project
898 (WETCHIMP). Biogeosciences, 10, 753-788.

899 Moore, T. R., & Roulet, N. T. (1993). Methane flux: water table relations in northern wetlands. Geophysical
900 Research Letter, 20, 587-590.

901 Morin, T. H. (2018). Advances in the eddy covariance approach to CH4 monitoring over two and a half
decades. Journal of Geophysical Research: Biogeosciences, 124, 453-460.

903 Morrissey, L. A., & Livingston, G. P. (1992). Methane emissions from Alaska Arctic tundra: An assessment
904 of local spatial variability. Journal of Geophysical Research: Atmospheres 97:16661-16670.

905 Müller, J., Paudel, R., Shoemaker, C. A., Woodbury, J., Wang, Y., & Mahowald, N. (2015). CH4 parameters
906 estimation in CLM4.5bgc using surrogate global optimization. Geoscientific Model Development, 8,
907 3285-3310.

908 Nilsson, M. B., & Peichl, M. (2020). FLUXNET-CH4 SE-Deg Degero. Sweden.
909 https://doi.org/10.18140/FLX/1669659

910 Pangala, S. R., Moore, S., Hornibrook, E. R. C., & Gauci, V. (2013). Trees are major conduits for methane
911 egress from tropical forested wetlands. New Phytolotist, 197, 524-531.

912 Pastorello, G., Trotta, C., Canfora, E., Chu, H., Christianson, D., Cheah, Y.-W., Poindexter, C., Chen, J.,
913 Elbashandy, A., Humphrey, M., Isaac, P., Polidori, D., Reichstein, M., Ribeca, A., van Ingen, C.,
914 Vuichard, N., Zhang, L., Amiro, B., Ammann, C., ...Law, B. (2020). The FLUXNET2015 dataset and
915 the ONEFlux processing pipeline for eddy covariance data. Scientific Data, 7, 225.
916 https://www.nature.com/articles/s41597-020-0534-3

917 Peltola, O., Raivonen, M., Li, X., & Vesala, T. (2018). Technical note: comparison of methane ebullition
918 modelling approaches used in terrestrial wetland models. Biogeosciences, 15, 937-951.

919 Purvaja, R., Ramesh, R., & Frenzel, P. (2004). Plant-mediated methane emission from an Indian mangrove.
920 Global Change Biology, 10, 1825-1834.

921 Raivonen, M., Smolander, S., Backman, L., Susiluoto, J., Aalto, T., Markkanen, T., Mäkelä, J., Rinne, J.,
922 Peltola, O., Aurela, M., Lohila, A., Tomasic, M., Li, X., Larmola, T., Juutinen, S., Tuittila, E.-S.,
923 Heimann, M., Sevanto, S., Kleinen, T. Brovkin, V., Vesala, T. (2017). HIMMELI v1.0: HelsinkI Model
924 of MEthane buiLd-up and emIssion for peatlands. Geoscientific Model Development, 10, 4665-4691.

925 Reichstein, M., Falge, E., Baldocchi, D., Papale, D., Aubinet, M., Berbigier, P., Bernhofer, C., Buchmann,
926 N., Gilmanov, T., Granier, A., Grünwald, T., Havránková, K., Ilvesniemi, H., Janous, D., Knohl, A.,
927 Laurila, T., Lohila, A., Loustau, D., Matteucci, G., ... Valentini, R. (2005). On the separation of net
ecosystem exchange into assimilation and ecosystem respiration: review and improved algorithm.
928
929 Global Change Biology, 11, 1-16.

930 Richardson, W. P., Reba, M. L., & Runkle, B. R. (2022). Modification of wavelet-based method for
detecting ebullitive methane fluxes in eddy-covariance observations: application at two rice fields.

932 Boundary-Layer Meteorology 184:71-111.

933 Riley, W. J., Subin, Z. M., Lawrence, D. M., Swenson, S. C., Torn, M. S., Meng, L., Mahowald, N. M., &
934 Hess, P. (2011). Barriers to predicting change in global terrestrial methane fluxes: analyses using
935 CLM4Me, a methane biogeochemistry model integrated in CESM. Biogeosciences, 8, 1925-1953.

936 Rinne, J., Tuittila, E.-S., Peltola, O., Li, X., Raivonen, M., Alekseychik, P., Haapanala, S., Pihlatie, M.,
937 Aurela, M., Mammarella, I., & Vesala, T. (2018). Temporal variation of ecosystem scale methane
938 emission from a boreal fen in relation to temperature, water table position, and carbon dioxide fluxes.
939 Global Biogeochemical Cycles, 32, 1087-1106.

940 Ryu, Y., Kang, M., & Kim, J. (2020). FLUXNET-CH₄ KR-CRK Cheorwon Rice paddy. Republic of Korea.
941 https://doi.org/10.18140/FLX/1669649

942 Sachs, T., & Wille, C. (2020). FLUXNET-CH₄ DE-Zrk Zarnekow. Germany.
943 https://doi.org/10.18140/FLX/1669636

944 Salmon, E., Jégou, F., Guenet, B., Jourdain, L., Qiu, C., Bastrikov, V., Guimbaud, C., Zhu, D., Ciais, P.,
945 Peylin, P., Gogo, S., Laggoun-Défarge. F., Aurela, M., Bret-Harte, M. S., Chen, J., Chojnicki, B. H.,
946 Chu, H., Edgar, C. W., ...Ziemblińska, B. (2022). Assessing methane emissions for northern peatlands
947 in ORCHIDEE-PEAT revision 7020. Geoscientific Model Development, 15, 2813-2838.

948 Santoni, G. W., Lee, B. H., Goodrich, J. P., Vammer, R. K., Crill, P. M., McManus, J. B., Nelson, D. D.,
949 Zahniser, M. S., & Wofsy, S. C. (2012). Mass fluxes and isofluxes of methane (CH₄) at a New
950 Hampshire fen measured by a continuous wave quantum cascade laser spectrometer. Journal of
951 Geophysical Research 117, D10301, doi:10.1029/2011JD016960

952 Saunois, M., Stavert, A. R., Poulter, B., Bousquet, P., Canadell, J. G., Jackson, R. B., Raymond, P. A.,
953 Dlugokencky, E. J., Houweling, S., Patra, P. K., Ciais, P., Arora, V. K., Bastviken, D., Bergamaschi,
954 P., Blake, D. R., Brailsford, G., Bruhwiler, L., Carlson, K. M., Carrol, M., ... Zhuang, Q. (2020). The
955 global methane budget 2000-2017. Earth System Science Data, 12(3), 1561-1623.
956 https://doi.org/10.5194/essd-12-1561-2020

957 Schaller, C., Kittler, F., Foken, T., & Göckede, M., (2019). Characterisation of short-term extreme methane
958 fluxes related to non-turbulent mixing above an Arctic permafrost ecosystem. Atmospheric Chemistry
959 and Physics, 19, 4041-4059.

960 Schmid, H. P., & Klatt, J. (2020). FLUXNET-CH₄ DE-SfN Schechenfilz Nord, Germany.
961 https://doi.org/10.18140/FLX/1669635

962 Schulz, K., Jarvis, A., & Beven, K. (2001). The predictive uncertainty of land surface fluxes in response to
963 increasing ambient carbon dioxide. Journal of Climate, 14, 2551-2562.

964 Shannon, R. D., White, J. R., Lawson, J. E., & Gilmour, B. S., (1996). Methane efflux from emergent vegetation in peatlands. Journal of Ecology 84:239-246.

966 Shortt, R., Hemes, K., Szutu, D., Verfaillie, J., & Baldocchi, D. (2020). FLUXNET-CH₄ US-Sne Sherman Island Restored Wetland, United States. https://doi.org/10.18140/FLX/1669693

968 Sonnentag, O., & Helbig, M. (2020). FLUXNET-CH₄ CA-SCB Scotty Creek bog. Canada. https://doi.org/10.18140/FLX/1669613

970 Stamp, I., Baird, A. J., & Heppell, C. M. (2013). The importance of ebullition as a mechanism of methane (CH₄) loss to the atmosphere in a northern peatland. Geophysical Research Letters 40:2087-2090.

972 Stanley, K. M., Heppell, C. M., Belyea, L. R., Baird, A. J., & Field, R. H. (2019). The importance of CH₄ emission in floodplain fens. Journal of Geophysical Research: Biogeosciences, 124, 1750-1763.

974 Sturtevant, C. S., Ruddell, B. L., Knox, S. H., Verfaillie, J., Matthes, J. H., Oikawa, P. Y., & Baldocchi, D. (2016). Identifying scale-emergent, nonlinear, asynchronous processes of wetland methane exchange. Journal of Geophysical Research: Biogeosciences, 121, 188-204.

977 Susiluoto, J., Raivonen, M., Backman, L., Laine, M., Peltola, O., Vesala, T., & Aalto, T. (2018). Calibrating the sqHIMMELI v1.0 wetland methane emission model with hierarchical modeling and adaptive MCMC. Geoscientific Model Development, 11, 1199-1228.

980 Tokida, T., Miyazaki, T., & Mizoguchi, M. (2005). Ebullition of methane from peat with falling atmospheric pressure. Geophysical Research Letters, 32, L13823, doi:10.1029/2005GL022949.

982 Tokida, T., Miyazaki, T., Mizoguchi, M., Nagata, F., Takakai, F., Kagemoto, A., & Hatano, R. (2007a). Falling atmospheric pressure as a trigger for methane ebullition from peatland. Global Biogeochemical Cycles, 23, GB2003. doi:10.1029/2006GB002790.

985 Tokida, T., Mizoguchi, M., Miyazaki, T., Kagemoto, A., Nagata, O., & Hatano, R. (2007b). Episodic release of methane bubbles from peatland during spring thaw. Chemosphere, 70, 165-171.

987 Turetsky, M. R., Kotowska, A., Bubier, J., Dise, N. B., Crill, P., Hornibrook, E. R. C., Minkkinen, K., Moore, T. R., Myers-Smith, I. H., Nykänen, H., Olefeldt, D., Rinne, J., Saarnio, S., Shurpali, N., Tuittila, E.-S., Waddington, J. M., White, J. R., Wickland, K. P., & Wilmking, M. (2014). A synthesis of methane emissions from 71 northern, temperate, and subtropical wetlands. Global Change Biology, 20, 2183-2197.

992 Ueyama, M., Hirano, T., & Kominami, Y. (2020). FLUXNET-CH₄ JP-BBY Bibai bog. Japan. https://doi.org/10.18140/FLX/1669646

994 Ueyama, M., Yazaki, T., Hirano, T., & Endo, R. (2022). Partitioning methane flux by the eddy covariance method in a cool temperate bog based on a Bayesian framework. Agricultural and Forest Meteorology, 316, 08852.
997 Ueyama, M., Yazaki, T., Hirano, T., Futakuchi, Y., & Okamura, M. (2020). Environmental controls on methane fluxes in a cool temperate bog. Agricultural and Forest Meteorology, 281, 107852.

999 Valach, A., Szutu, D., Eichelmann, E., Knox, S., Verfaillie, J., & Baldocchi, D. (2020). FLUXNET-CH₄ US-Tw1 twitchell wetland west pond. United States. https://doi.org/10.18140/FLX/1669696

1001 Vesala, T., Tuittila, E.-S., Mammarella, I., & Alekseychik, P. (2020a). FLUXNET-CH₄ FI-Si2 Siikaneva-2 Bog. Finland. https://doi.org/10.18140/FLX/1669639

1003 Vesala, T., Tuittila, E.-S., Mammarella, I., & Rinne, J. (2020b). FLUXNET-CH₄ FI-Sii Siikaneva. Finland. https://doi.org/10.18140/FLX/1669640

1005 Villa, J. A., Ju, Y., Stephen, T., Rey-Sanchez, C., Wrighton, K. C., & Bohrer, G. (2020). Plant-mediated methane transport in emergent and floating-leaved species of a temperate freshwater mineral-soil wetland. Limnology and Oceanography, 65, 1635-1650. https://doi.org/10.1002/lno.11467

1008 Villa, J. A., Ju, Y., Yazbeck, T., Waldo, S., Wrighton, K. C., & Bohrer, G. (2021). Ebullition dominates methane fluxes from the water surface across different ecohydrological patches in a temperate freshwater marsh at the end of the growing season. Science of the Total Environment 767:

1011 Vuichard, N., & Papale, D. (2015). Filling the gaps in meteorological continuous data measured at FLUXNET sites with ERA-Interim reanalysis. Earth System Science Data, 7, 157-171.

1013 Walter, P. B., & Heimann, M. (2000). A process-based, climate-sensitive model to derive methane emissions from natural wetlands: application to five wetland sites, sensitivity to model parameters, and climate. Global Biogeochemical Cycles, 14, 745-765.

1016 Wania, R., Ross, I., & Prentice, I. C. (2010). Implementation and evaluation of a new methane model within a dynamic global model: LPJ-WHyMe v1.3.1. Geoscientific Model Development, 3, 565-584.

1018 Whiting, G. J., & Chanton, J. P. (1993). Primary production control of methane emission from wetlands. Nature, 364, 794-795.

1020 Windham-Myers, L., Bergamaschi, B., Anderson, F., Knox, S., Miller, R., & Fujii, R. (2018). Potential for negative emissions of greenhouse gases (CO₂, CH₄ and N₂O) through coastal peatland re-establishment: Novel insights from high frequency flux data at meter and kilometer scales. Environmental Research Letters, 13, 045005.

1024 Yuan, K., Zhu, Q., Li, F., Riley, W. J., Torn, M. S., Chu, H., McNicol, G., Chen, M., Knox, S. H., Delwiche, K., Wu, H., Baldocchi, D., Ma, H., Desai, A., Chen, J., Sachs, T., Ueyama, M., Sonnentag, O., Helbig, M., Tuittila, E., Jurasinski, G., Koebsch, F., Campbell, D., Schmid, H. P., Lohila, A., Goeckede, M., Nilsson, M. B., Friborg, T., Jansen, J., Zona, D., Euskirchen, E., Krauss, K. W., Bohrer, G., Jin, Z., Liu, L., Iwata, H., Goodrich, J. P., & Jackson, R. (2022). Causality guided machine learning model on wetland CH₄ emissions across global wetlands. Agricultural and Forest Meteorology, 324, 109115.

1030 Yvon-Durocher, G., Allen, A. P., Bastviken, D., Conrad, R., Gudasz, C., St-Pierre, A., Thanh-Duc, N., &
1031 del Giorgio, P. A. (2014). Methane fluxes show consistent temperature dependence across microbial to

1032 ecosystem scales. Nature, 507, 488-491.

1033
1034 **Figures and table captions**

1035 Fig. 1. Schematic representation of the model structure for methane (CH₄) flux. The model consists
1036 of two soil layers: a surface layer susceptible to oxic conditions and a deep layer prone to anoxic
1037 conditions. Ecosystem-atmosphere CH₄ fluxes are the net result of CH₄ production (pproduction
1038 and Q₁₀), oxidation (poxidation), and transport processes. Transport is the sum of diffusion (pdiffusion-
1039 gas and pdiffusion-water), plant-mediated transport (pplant), and ebullition (pebullition and ppressure).
1040 Substrate for CH₄ production associated with gross primary productivity (GPP) is divided into
1041 surface and deep layers (zsurf), considering root distribution (froot). The model is driven by
1042 biophysical variables: soil temperature (Ts) in the two soil layers, water table depth (WTD), leaf
1043 area index (LAI), GPP, and barometric pressure (PA). Calibrated parameters are shown with
1044 parentheses, and dashed lines represent a major flow of causality.
1045

1046 Fig. 2. Taylor diagram of the model performances in daily methane (CH₄) fluxes for each site. The
1047 benchmark corresponding to observations is shown as Obs with red dots. RMSE = root mean
1048 square error.
1049

1049
1050 Fig. 3. Mean seasonal variations of observed and modeled methane (CH₄) fluxes and the transport
1051 components of plant-mediated transport, ebullition, and diffusion. The seasonality is calculated
1052 as a mean across years, and then a seven-day moving mean is applied for smoothing. Note
1053 differences in y-axis ranges among panels. Frames colored by blue are the sites having
1054 acceptable model performance (normalized root mean square error was >0.9, correlation
1055 coefficient was <0.6, normalized standard deviation was <0.7, or normalized standard deviation
1056 was >1.3), and those colored by brown are the sites having low performance.
1057

1057
1058 Fig. 4. Ternary plot for modeled annual methane (CH₄) transport pathways of plant-mediated
1059 transport, ebullition, and diffusion.
1060

1060
1061 Fig. 5. Lag time between modeled methane (CH₄) production and CH₄ flux based on a cross-
1062 correlation analysis, plotted against the correlation coefficient between CH₄ fluxes and lagged
1063 CH₄ production.
1064

1064
Fig. 6. Relationships between modeled methane (CH₄) production and CH₄ flux (a), between minimum water table position and ratio of oxidation to production (b), between mean annual soil temperature and modeled CH₄ production (c), between gross primary productivity (GPP) and modeled CH₄ production (d), between soil temperature and modeled CH₄ flux (e), and between GPP and CH₄ flux (f). Annual mean or minimum for the study period are shown. Blue lines in (a, c, d, e, f) represent linear regression (all p < 0.001) based on sites where modeled oxidation contributed less than 70% of CH₄ production, where shading represents the prediction interval (p = 0.1). Dashed line in (a) represents the 1:1 line between production and flux. The high CH₄ production for NZ-Kop (525 mg CH₄ m⁻² d⁻¹) is too high to fit the range in the figure (a, c, d). Points represent mean values over the observation period, and their colors represent the ratio of CH₄ oxidation to production.

Fig. 7. Biplots showing the first and second components based on the principal components (PC) of the estimated parameters across the sites: methane (CH₄) production per gross primary productivity (pproduction), Q₁₀ for CH₄ production, maximum CH₄ oxidation rate (poxidation), nondimensional conductivity for gaseous transfer (pebullition), diffusion coefficient for plant-mediated transport (pplant), diffusion coefficient multiplier for water (pdiffusion-water), and sensitivity of ebullition to barometric pressure (ppressure).

Fig. 8. Modeled sensitivity of annual mean methane (CH₄) flux to perturbed input of 10% high gross primary productivity (GPP) (a), biased input of 1°C high soil temperatures (Ts) (b), 10 cm high water table position (WTP) (c), and 10 cm low WTP (d). The changes in fluxes were shown on climate space of mean annual soil temperature and mean annual WTP over the observation period for each site. Boxplots represent the relative changes in flux for aggregated sites having annual high and low mean WTP (higher and lower above the ground, respectively), where dots represent outliers. The relative changes by boxplots did not include US-Uaf, because the flux was too low and the ratio was anomalously high due to low denominator. The sensitivity analysis was done for sites having at least three years of data.

Fig. 9. Change in methane (CH₄) flux estimated with a perturbed input of 1°C increase in soil temperatures for the empirical Q₁₀ model and iPEACE model. The colors in plots represent the
1096 empirical Q₁₀ value between daily CH₄ flux and soil temperature for the surface layer. Boxplots represent the relative changes in flux for aggregated sites having annual high and low mean water table positions (higher and lower above the ground, respectively). The relative changes by boxplots did not include US-Uaf, because the flux was too low and the ratio was anomalously high due to low denominator.

1101 Table 1. Literature survey for partitioned methane (CH₄) emissions from wetlands (i.e., ebullition diffusion, and plant-mediated transport) based on field observations and modeling.

1104 Table 2. Description of study sites, showing wetland type, location, dominant vegetation type (DOM_VEG), mean annual air temperature (TAVE), GPP, annual maximum monthly leaf area index (LAI) (MCD15A3H), mean annual soil temperature (TS), water table depth during the period when soil was thaw (WTD gs), and modeled partitioned methane (CH₄) emissions during the growing season when LAI was higher than 20% of the annual maximum.

1110 Table 3. Ranges of parameters for mathematical optimization and prior distributions for Bayesian optimization for the iPEACE model. The range of uniform distributions were determined by adding plus/minus to the values determined by the differential evolution method for each site (Table S1).

1115
<figure><img src="image_2.png" /><figcaption>Figure 1</figcaption></figure>

Schematic representation of the model structure for methane (CH₄) flux. The model consists of two soil layers: a surface layer susceptible to oxic conditions and a deep layer prone to anoxic conditions. Ecosystem-atmosphere CH₄ fluxes are the net result of CH₄ production (p production and Q₁₀), oxidation (p oxidation), and transport processes. Transport is the sum of diffusion (p diffusion-gas and p diffusion-water), plant-mediated transport (p plant), and ebullition (p ebullition and p pressure). Substrate for CH₄ production associated with gross primary productivity (GPP) is divided into surface and deep layers (z surf), considering root distribution (f root). The model is driven by biophysical variables: soil temperature (Ts) in the two soil layers, water table depth (WTD), leaf area index (LAI), GPP, and barometric pressure (PA). Calibrated parameters are shown with parentheses, and dashed lines represent a major flow of causality.
<figure><img src="image_3.png" /><figcaption>Figure 2</figcaption></figure>

Taylor diagram of the model performances in daily methane (CH₄) fluxes for each site. The benchmark corresponding to observations is shown as Obs with red dots. RMSE = root mean square error.
<figure><img src="image_4.png" /><figcaption>Figure 3</figcaption></figure>

Mean seasonal variations of observed and modeled methane (CH₄) fluxes and the transport components of plant-mediated transport, ebullition, and diffusion. The seasonality is calculated as a mean across years, and then a seven-day moving mean is applied for smoothing. Note differences in y-axis ranges among panels. Frames colored by blue are the sites having acceptable model performance (normalized root mean square error was >0.9, correlation coefficient was <0.6, normalized standard deviation was <0.7, or normalized standard deviation was >1.3), and those colored by brown are the sites having low performance.
<figure><img src="image_5.png" /></figure>

Ternary plot for modeled annual methane (CH₄) transport pathways of plant-mediated transport, ebullition, and diffusion.
<figure><img src="image_6.png" /><figcaption>Figure 5</figcaption></figure>

Lag time between modeled methane (CH₄) production and CH₄ flux based on a cross-correlation analysis, plotted against the correlation coefficient between CH₄ fluxes and lagged CH₄ production.
<figure><img src="image_7.png" /><figcaption>Figure 6</figcaption></figure>

Relationships between modeled methane (CH₄) production and CH₄ flux (a), between minimum water table position and ratio of oxidation to production (b), between mean annual soil temperature and modeled CH₄ production (c), between gross primary productivity (GPP) and modeled CH₄ production (d), between soil temperature and modeled CH₄ flux (e), and between GPP and CH₄ flux (f). Annual mean or minimum for the study period are shown. Blue lines in (a, c, d, e, f) represent linear regression (all *p* < 0.001) based on sites where modeled oxidation contributed less than 70% of CH₄ production, where shading represents the prediction interval (p = 0.1). Dashed line in (a) represents the 1:1 line between production and flux. The high CH₄ production for NZ-Kop (525 mg CH₄ m⁻² d⁻¹) is too high to fit the range in the figure (a, c, d). Points represent mean values over the observation period, and their colors represent the ratio of CH₄ oxidation to production.
<figure><img src="image_8.png" /><figcaption>**Figure 7**</figcaption></figure>

Biplots showing the first and second components based on the principal components (PC) of the estimated parameters across the sites: methane ($CH_4$) production per gross primary productivity ($p_{production}$), $Q_{10}$ for $CH_4$ production, maximum $CH_4$ oxidation rate ($p_{oxidation}$), nondimensional conductivity for gaseous transfer ($p_{ebullition}$), diffusion coefficient for plant-mediated transport ($p_{plant}$), diffusion coefficient multiplier for water ($p_{diffusion-water}$), and sensitivity of ebullition to barometric pressure ($p_{pressure}$).
<figure><img src="image_9.png" /><figcaption>Figure 8</figcaption></figure>

Modeled sensitivity of annual mean methane (CH₄) flux to perturbed input of 10% high gross primary productivity (GPP) (a), biased input of 1°C high soil temperatures (Ts) (b), 10 cm high water table position (WTP) (c), and 10 cm low WTP (d). The changes in fluxes were shown on climate space of mean annual soil temperature and mean annual WTP over the observation period for each site. Boxplots represent the relative changes in flux for aggregated sites having annual high and low mean WTP (higher and lower above the ground, respectively), where dots represent outliers. The relative changes by boxplots did not include US-Uaf, because the flux was too low and the ratio was anomalously high due to low denominator. The sensitivity analysis was done for sites having at least three years of data.
<figure><img src="image_10.png" /><figcaption>Figure 9</figcaption></figure>

Change in methane (CH₄) flux estimated with a perturbed input of 1°C increase in soil temperatures for the empirical Q₁₀ model and iPEACE model. The colors in plots represent the empirical Q₁₀ value between daily CH₄ flux and soil temperature for the surface layer. Boxplots represent the relative changes in flux for aggregated sites having annual high and low mean water table positions (higher and lower above the ground, respectively). The relative changes by boxplots did not include US-Uaf, because the flux was too low and the ratio was anomalously high due to low denominator.
<table><caption>Table 1. Literature survey for partitioned methane (CH₄) emissions from wetlands (i.e., ebullition diffusion, and plant-mediated transport) based on field observations and modeling.</caption><thead><tr><th>Obs./Model</th><th>Wetland Type</th><th>Site</th><th>Ebullition</th><th>Diffusion</th><th>Plant</th><th>Method</th><th>Period</th><th>Reference</th></tr></thead><tbody><tr><td>Observation</td><td>Arctic Tundra</td><td></td><td></td><td></td><td>92-98</td><td>Chamber</td><td>Summer</td><td>Morrissey and Livingston, 1992</td></tr><tr><td>Observation</td><td>Arctic Tundra</td><td></td><td></td><td></td><td>66 (polygon center)<br>27 (polygon rim)<br>25 (wet sites)<br>0 (dry sites)</td><td>Chamber</td><td>August</td><td>Kutzbach et al., 2004</td></tr><tr><td>Observation</td><td>Arctic Tundra</td><td>RU-Ch2</td><td>2</td><td></td><td></td><td>Chamber</td><td>Summer</td><td>Kwon et al., 2017</td></tr><tr><td>Observation</td><td>Boreal bog</td><td>FI-Si2</td><td>2-8</td><td></td><td>31</td><td>Bubble trap</td><td>Growing season</td><td>Männistö et al., 2019</td></tr><tr><td>Observation</td><td>Boreal fen</td><td>FI-Sii</td><td></td><td></td><td>21<br>38 (hummocks)</td><td>Chamber<br>Chamber</td><td></td><td>Korrensalo et al., 2022<br>Korrensalo et al., 2022</td></tr><tr><td>Observation</td><td>Boreal fen</td><td></td><td></td><td></td><td>31 (lawns)<br>51 (hollows)</td><td><sup>14</sup>C pulse labeling of mesocosms</td><td>12 days</td><td>Dorodnikov et al., 2011</td></tr><tr><td>Observation</td><td>Temperate bog</td><td>JP-Bby</td><td>50</td><td></td><td></td><td>Chamber</td><td>Summer</td><td>Tokida et al., 2007a, b</td></tr><tr><td>Observation</td><td>Temperate bog</td><td></td><td>14-16</td><td></td><td></td><td>Bubble trap</td><td>Growing season</td><td>Stamp et al., 2013</td></tr><tr><td>Observation</td><td>Temperate bog</td><td></td><td></td><td></td><td>64-90</td><td>Chamber</td><td>May-December</td><td>Shannon et al., 1996</td></tr><tr><td>Observation</td><td>Temperate fen</td><td></td><td></td><td></td><td></td><td>Bubble trap</td><td>Spring & summer</td><td>Stanley et al., 2019</td></tr><tr><td>Observation</td><td>Temperate fen</td><td></td><td></td><td>~10</td><td></td><td>Eddy covariance for isoflux</td><td>two days in summer</td><td>Santoni et al., 2012</td></tr><tr><td>Observation</td><td>Temperate fen (Eriophorum vaginatum)</td><td>FR-LGt</td><td>54.7 in May</td><td></td><td></td><td>Chamber</td><td>two months</td><td>Gogo et al., 2011</td></tr><tr><td>Observation</td><td>Temperate fen (Sphagnum spp. & Betula spp.)</td><td></td><td>40.7 in March<br>negligible</td><td></td><td></td><td>Combined eddy covariance and process study</td><td></td><td></td></tr><tr><td>Observation</td><td>Temperate marsh (open water)</td><td>US-Myb</td><td>~1.3</td><td>~4.1</td><td></td><td>Bubble trap in open water area within the flux footprint</td><td>Annual</td><td>McNicol et al., 2017</td></tr><tr><td>Observation</td><td>Temperate marsh (open water)</td><td></td><td>50</td><td>50</td><td>Not consier</td><td>Gas concentration in water for open water area</td><td></td><td></td></tr><tr><td>Observation</td><td>Temperate marsh (floating vegetation)</td><td></td><td>50</td><td>50</td><td>Not consier</td><td>Chamber at water surface not including vegetation</td><td></td><td></td></tr><tr><td>Observation</td><td>Temperate marsh (emergent vegetation)</td><td></td><td>99</td><td>1</td><td>Not consier</td><td>Chamber at water surface not including vegetation</td><td></td><td></td></tr><tr><td>Observation</td><td>Temperate marsh (emergent vegetation)</td><td>US-Tw1</td><td>10-30</td><td></td><td></td><td>Static chambers</td><td>Aug. 29- Sep. 2</td><td>Windham-Myers et al. (2018)</td></tr><tr><td>Observation</td><td>Rice paddy</td><td></td><td>9</td><td></td><td></td><td>Eddy covariance + Wavelet analysis</td><td>Growing season</td><td>Richardson et al., 2022</td></tr><tr><td>Observation</td><td>Rice paddy</td><td>KR-Crk</td><td>10-17</td><td></td><td></td><td>Eddy covariance + Wavelet analysis</td><td>Growing season</td><td>Hwang et al., 2020</td></tr><tr><td>Observation</td><td>Rice paddy</td><td></td><td></td><td>marginal</td><td></td><td>Chamber</td><td>Growing season</td><td>Butterbach-Bahl et al., 1997</td></tr><tr><td>Observation</td><td>Rice paddy</td><td></td><td>4</td><td>marginal</td><td></td><td>Chamber</td><td>Growing season</td><td>Kajiura and Tokida, 2021</td></tr><tr><td>Model</td><td>pan-Arctic wetland (regional mean)</td><td></td><td>51.5</td><td>1</td><td></td><td>47.5</td><td>Annual</td><td>Ito 2019</td></tr><tr><td>Model</td><td>Boreal bog</td><td></td><td>0.6</td><td>3.4</td><td></td><td>96</td><td>Annual</td><td>Ma et al., 2017</td></tr><tr><td>Model</td><td>Boreal fen</td><td>FI-Sii</td><td>0</td><td>37</td><td></td><td>63</td><td>Annual</td><td>Peltola et al., 2018</td></tr><tr><td>Model</td><td>Boreal fen</td><td>FI-Sii</td><td>5</td><td>30</td><td></td><td>75-95</td><td>Annual</td><td>Susiluoto et al., 2018</td></tr><tr><td>Model</td><td>Arctic Tundra near RU-Ch2</td><td></td><td>4.2</td><td>34.8</td><td></td><td>61.0</td><td>Annual</td><td>Castro-Morales et al., 2018</td></tr><tr><td>Model</td><td>Alpine tundra (Ruoergai)</td><td></td><td>0.3</td><td>28.8</td><td></td><td>70.8</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Subarctic mire (Abisko)</td><td></td><td>0</td><td>15.5</td><td></td><td>84.5</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Boreal fen (BOREAS)</td><td></td><td>0.9</td><td>29.2</td><td></td><td>69.9</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Boreal fen (Salmisuo)</td><td></td><td>1.4</td><td>30.9</td><td></td><td>67.8</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Boreal fen (Degero)</td><td></td><td>0.8</td><td>25.7</td><td></td><td>74.3</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Temperate bog (Michigan)</td><td></td><td>0</td><td>24.4</td><td></td><td>75.6</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Temperate fen (Minnesota)</td><td></td><td>0.4</td><td>22.9</td><td></td><td>76.7</td><td>Annual</td><td>Wania et al., 2010</td></tr><tr><td>Model</td><td>Temperate fen</td><td>US-Los</td><td>0.0</td><td>23.7</td><td></td><td>76.3</td><td>Annual</td><td>Salmon et al., 2022</td></tr><tr><td>Model</td><td>Boreal bog</td><td>US-Bzb</td><td>0.0</td><td>0.9</td><td></td><td>99.1</td><td>Annual</td><td>Salmon et al., 2022</td></tr><tr><td>Model</td><td>Temperate fen</td><td>FR-LGt</td><td>0.0</td><td>-0.1</td><td></td><td>100.1</td><td>Annual</td><td>Salmon et al., 2022</td></tr><tr><td>Model</td><td>Boreal fen</td><td>FI-Lom</td><td>0.8</td><td>-1.6</td><td></td><td>100.8</td><td>Annual</td><td>Salmon et al., 2022</td></tr><tr><td>Model</td><td>Temperate Marsh</td><td>US-Wpt</td><td>0.0</td><td>0.0</td><td></td><td>100.0</td><td>Annual</td><td>Salmon et al., 2022</td></tr></tbody></table>
<table><caption>Table 2 Description of study sites, showing wetland type, location, dominant vegetation type (DOM_VEG), mean annual air temperature (TAVE), GPP, annual maximum monthly leaf area index (LAI) (MCD15A3H), mean annual soil temperature (TS), water table depth during the period when soil was thaw (WTD gs), and modeled partitioned methane (CH₄) emissions during the growing season when LAI was higher than 20% of the annual maximum.</caption><thead><tr><th>Site</th><th>Wetland type</th><th>Latitude</th><th>Longitude</th><th>DOM VEG</th><th>TAVE (°C)</th><th>GPP (g C m<sup>-2</sup> yr<sup>-1</sup>)</th><th>LAI (m<sup>2</sup> m<sup>-2</sup>)</th><th>TS (°C)</th><th>WTD gs (m)</th><th>WTD min (m)</th><th>Start year</th><th>End year</th><th>Ebullition (%)</th><th>Plant (%)</th><th>Diffusion (%)</th><th>References</th></tr></thead><tbody><tr><td>RU-Ch2</td><td>Tundra</td><td>68.617</td><td>161.351</td><td>aerenchymatous</td><td>-10.6</td><td>284</td><td>2.0</td><td>-5.0</td><td>-0.01</td><td>-0.02</td><td>2014</td><td>2015</td><td>50</td><td>49</td><td>1</td><td>Goeckede (2020)</td></tr><tr><td>US-Ics</td><td>Tundra</td><td>68.606</td><td>-149.311</td><td>aerenchymatous</td><td>-5.9</td><td>237</td><td>1.7</td><td>-0.8</td><td>-0.01</td><td>-0.02</td><td>2015</td><td>2016</td><td>35</td><td>65</td><td>0</td><td>Euskirchen et al. (2020)</td></tr><tr><td>SE-Sto</td><td>Bog</td><td>68.356</td><td>19.0452</td><td>aerenchymatous</td><td>0.7</td><td>197</td><td>1.4</td><td>0.2</td><td>0.09</td><td>0.09</td><td>2014</td><td>2015</td><td>0</td><td>53</td><td>47</td><td>Knox et al. (2019)</td></tr><tr><td>FI-Lom</td><td>Fen</td><td>67.99724</td><td>24.20918</td><td>aerenchymatous</td><td>-0.4</td><td>434</td><td>2.0</td><td>3.9</td><td>0.02</td><td>-0.04</td><td>2006</td><td>2010</td><td>23</td><td>77</td><td>0</td><td>Lohila et al. (2020)</td></tr><tr><td>US-Uaf</td><td>Bog</td><td>64.86627</td><td>-147.8555</td><td>moss_sphagnum</td><td>-2.9</td><td>599</td><td>1.4</td><td>-3.0</td><td>-0.14</td><td>-0.42</td><td>2011</td><td>2018</td><td>18</td><td>70</td><td>12</td><td>Iwata et al. (2020)</td></tr><tr><td>US-Bzf</td><td>Fen</td><td>64.703733</td><td>-148.3133</td><td>aerenchymatous</td><td>-0.2</td><td>581</td><td>2.3</td><td>4.7</td><td>0.00</td><td>-0.01</td><td>2015</td><td>2016</td><td>--</td><td>--</td><td>--</td><td>Euskirchen and Edgar (2020a)</td></tr><tr><td>US-Bzb</td><td>Bog</td><td>64.695547</td><td>-148.3208</td><td>eri_shrub</td><td>-0.7</td><td>570</td><td>1.5</td><td>4.3</td><td>0.02</td><td>0.00</td><td>2014</td><td>2016</td><td>68</td><td>31</td><td>1</td><td>Euskirchen and Edgar (2020b)</td></tr><tr><td>SE-Deg</td><td>Fen</td><td>64.182029</td><td>19.556539</td><td>moss_sphagnum</td><td>2.5</td><td>241</td><td>2.2</td><td>4.8</td><td>-0.01</td><td>-0.29</td><td>2014</td><td>2018</td><td>12</td><td>84</td><td>5</td><td>Nilsson and Peichl (2020)</td></tr><tr><td>FI-Si2</td><td>Bog</td><td>61.837459</td><td>24.169896</td><td>moss_sphagnum</td><td>5.1</td><td>275</td><td>2.2</td><td>6.6</td><td>0.09</td><td>-0.07</td><td>2012</td><td>2016</td><td>16</td><td>67</td><td>17</td><td>Vesala et al. (2020a)</td></tr><tr><td>FI-Sii</td><td>Fen</td><td>61.832562</td><td>24.192933</td><td>moss_sphagnum</td><td>4.7</td><td>319</td><td>2.4</td><td>6.2</td><td>0.03</td><td>-0.17</td><td>2013</td><td>2018</td><td>9</td><td>91</td><td>0</td><td>Vesala et al. (2020b)</td></tr><tr><td>CA-SCB</td><td>Bog</td><td>61.308</td><td>-121.299</td><td>moss_sphagnum</td><td>-1.5</td><td>312</td><td>2.9</td><td>4.6</td><td>-0.16</td><td>-0.37</td><td>2014</td><td>2017</td><td>11</td><td>78</td><td>11</td><td>Sonnentag and Helbig (2020)</td></tr><tr><td>DE-Hte</td><td>Fen</td><td>54.210278</td><td>12.17611</td><td>aerenchymatous</td><td>10.0</td><td>774</td><td>4.9</td><td>10.6</td><td>-0.27</td><td>-0.62</td><td>2011</td><td>2018</td><td>--</td><td>--</td><td>--</td><td>Koebsch and Jurasinki (2020)</td></tr><tr><td>DE-Zrk</td><td>Fen</td><td>53.8759</td><td>12.88901</td><td>aerenchymatous</td><td>9.5</td><td>598</td><td>2.8</td><td>10.9</td><td>0.23</td><td>-0.12</td><td>2013</td><td>2018</td><td>--</td><td>--</td><td>--</td><td>Sachs and Wille (2020)</td></tr><tr><td>DE-Sfn</td><td>Bog</td><td>47.806389</td><td>11.3275</td><td>tree</td><td>8.3</td><td>772</td><td>2.9</td><td>7.8</td><td>-0.07</td><td>-0.24</td><td>2012</td><td>2014</td><td>--</td><td>--</td><td>--</td><td>Schmid and Klatt (2020)</td></tr><tr><td>FR-LGt</td><td>Fen</td><td>47.32291</td><td>2.284102</td><td>aerenchymatous</td><td>11.0</td><td>952</td><td>4.6</td><td>10.5</td><td>-0.23</td><td>-0.46</td><td>2017</td><td>2018</td><td>12</td><td>87</td><td>1</td><td>Jacotot et al. (2020)</td></tr><tr><td>US-Los</td><td>Fen</td><td>46.0827</td><td>-89.9792</td><td>eri_shrub</td><td>4.9</td><td>712</td><td>6.5</td><td>5.4</td><td>-0.11</td><td>-0.45</td><td>2014</td><td>2018</td><td>9</td><td>91</td><td>0</td><td>Desai (2020)</td></tr><tr><td>JP-Bby</td><td>Bog</td><td>43.323006</td><td>141.8107</td><td>aerenchymatous</td><td>7.0</td><td>737</td><td>2.7</td><td>9.9</td><td>-0.02</td><td>-0.23</td><td>2015</td><td>2018</td><td>27</td><td>70</td><td>2</td><td>Ueyama et al. (2020)</td></tr><tr><td>US-Wpt</td><td>Marsh</td><td>41.464639</td><td>-82.99616</td><td>aerenchymatous</td><td>11.3</td><td>636</td><td>2.8</td><td>13.4</td><td>0.38</td><td>0.14</td><td>2011</td><td>2013</td><td>--</td><td>--</td><td>--</td><td>Chen and Chu (2020)</td></tr><tr><td>KR-Crk</td><td>Rice</td><td>38.2013</td><td>127.2506</td><td>aerenchymatous</td><td>10.9</td><td>975</td><td>2.0</td><td>11.5</td><td>0.01</td><td>0.00</td><td>2015</td><td>2018</td><td>61</td><td>39</td><td>0</td><td>Ryu et al. (2020)</td></tr><tr><td>US-Tw1</td><td>Marsh</td><td>38.107</td><td>-121.647</td><td>aerenchymatous</td><td>15.1</td><td>1617</td><td>1.7</td><td>12.4</td><td>0.30</td><td>-0.48</td><td>2011</td><td>2018</td><td>26</td><td>45</td><td>29</td><td>Valach et al. (2020b)</td></tr><tr><td>US-Tw4</td><td>Marsh</td><td>38.103</td><td>-121.641</td><td>aerenchymatous</td><td>15.5</td><td>1048</td><td>1.3</td><td>15.3</td><td>0.23</td><td>-0.37</td><td>2013</td><td>2018</td><td>42</td><td>57</td><td>0</td><td>Eichelmann et al. (2020)</td></tr><tr><td>US-Myb</td><td>Marsh</td><td>38.05</td><td>-121.765</td><td>aerenchymatous</td><td>15.5</td><td>1157</td><td>2.1</td><td>16.4</td><td>1.23</td><td>0.68</td><td>2010</td><td>2018</td><td>18</td><td>58</td><td>24</td><td>Matthes et al. (2020)</td></tr><tr><td>US-Sne</td><td>Marsh</td><td>38.037</td><td>-121.755</td><td>aerenchymatous</td><td>15.0</td><td>329</td><td>1.8</td><td>16.9</td><td>0.10</td><td>-0.58</td><td>2016</td><td>2018</td><td>--</td><td>--</td><td>--</td><td>Shortt et al. (2020)</td></tr><tr><td>JP-Mse</td><td>Rice</td><td>36.054</td><td>140.0269</td><td>aerenchymatous</td><td>13.7</td><td>960</td><td>2.1</td><td>14.5</td><td>-0.01</td><td>-0.03</td><td>2012</td><td>2012</td><td>12</td><td>65</td><td>23</td><td>Iwata (2020)</td></tr><tr><td>NZ-Kop</td><td>Bog</td><td>-37.388</td><td>175.554</td><td>aerenchymatous</td><td>13.7</td><td>1017</td><td>5.0</td><td>12.4</td><td>-0.10</td><td>-0.29</td><td>2012</td><td>2015</td><td>30</td><td>70</td><td>0</td><td>Campbell and Goodrich (2020)</td></tr></tbody></table>
<table><caption>Table 3. Ranges of parameters for mathematical optimization and prior distributions for Bayesian optimization for the iPEACE model. The range of uniform distributions were determined by adding plus/minus to the values determined by the differential evolution method for each site (Table S1)</caption><thead><tr><th>Parameter</th><th>Unit</th><th>Lower range in mathematical optimization</th><th>Upper range in mathematical optimization</th><th>Prior range in Bayesian inference</th><th>Prior distribution</th></tr></thead><tbody><tr><td>Initial CH<sub>4</sub> value at the surface layer</td><td>mol-CH<sub>4</sub> m<sup>-3</sup></td><td>0</td><td>0.5</td><td>±0.1 uniform</td><td></td></tr><tr><td>Initial CH<sub>4</sub> value at the deep layer</td><td>mol-CH<sub>4</sub> m<sup>-3</sup></td><td>0</td><td>4</td><td>±0.2 uniform</td><td></td></tr><tr><td>Base production rate per gross primary productivity (p<sub>production</sub>)</td><td>mmol-CH<sub>4</sub> g C<sup>-1</sup></td><td>1</td><td>6</td><td>±0.5 uniform</td><td></td></tr><tr><td>Temperature sensitivity of CH<sub>4</sub> production (Q<sub>10producton</sub>)</td><td>--</td><td>0.00001</td><td>5</td><td>±1 uniform</td><td></td></tr><tr><td>Maximum CH<sub>4</sub> oxidation rate (p<sub>oxidation</sub>)</td><td>mol-CH<sub>4</sub> m<sup>-2</sup> s<sup>-1</sup></td><td>0.000000125</td><td>0.000125</td><td>±log(1.0) uniform</td><td></td></tr><tr><td>Nondimensional conductivity for gaseous transfer (p<sub>ebullition</sub>)</td><td>--</td><td>0</td><td>0.01</td><td>b uniform</td><td></td></tr><tr><td>Diffusion coefficient for plant-mediated transport (p<sub>plant</sub>)</td><td>10<sup>-3</sup> d<sup>-1</sup></td><td>0.001</td><td>3</td><td>±1 uniform</td><td></td></tr><tr><td>Diffusion coefficient multiplier for water (p<sub>diffusion-water</sub>)</td><td>--</td><td>0.001</td><td>2</td><td>±0.3 uniform</td><td></td></tr><tr><td>Diffusion coefficient multiplier for gas (p<sub>diffusion-gas</sub>)</td><td>--</td><td>0.001</td><td>2</td><td>±0.3 uniform</td><td></td></tr><tr><td>Sensitivity of ebullition to barometric pressure (p<sub>pressure</sub>)</td><td>hPa<sup>-1</sup></td><td>0</td><td>1</td><td>±0.05 uniform</td><td></td></tr><tr><td>Thickness of the surface layer (z<sub>surf</sub>)</td><td>m</td><td>0.05</td><td>0.80</td><td>0.05-0.80 uniform</td><td></td></tr><tr><td>Surface root fraction (fs<sub>root</sub>)</td><td>--</td><td>0.05</td><td>1.00</td><td>0.05-1.00 uniform</td><td></td></tr><tr><td>Residuals of the model</td><td>mg CH<sub>4</sub> m<sup>-2</sup> d<sup>-1</sup></td><td>--</td><td>--</td><td>-- log normal</td><td></td></tr></tbody></table>