Representing methane emissions from wet tropical forest soils using microbial functional groups constrained by soil diffusivity

Tropical ecosystems contribute significantly to global emissions of methane (CH₄), and landscape topography influences the rate of CH₄ emissions from wet tropical forest soils. However, extreme events such as drought can alter normal topographic patterns of emissions. Here we explain the dynamics of CH₄ emissions during normal and drought conditions across a catena in the Luquillo Experimental Forest, Puerto Rico. Valley soils served as the major source of CH₄ emissions in a normal precipitation year (2016), but drought recovery in 2015 resulted in dramatic pulses in CH₄ emissions from all topographic positions. Geochemical parameters including (i) dissolved organic carbon (C), acetate, and soil pH and (ii) hydrological parameters like soil moisture and oxygen (O₂) concentrations varied across the catena. During the drought, soil moisture decreased in the slope and ridge, and O₂ concentrations increased in the valley. We simulated the dynamics of CH₄ emissions with the Microbial Model for Methane Dynamics-Dual Arrhenius and Michaelis–Menten (M3D-DAMM), which couples a microbial functional group CH₄ model with a diffusivity module for solute and gas transport within soil microsites. Contrasting patterns of soil moisture, O₂, acetate, and associated changes in soil pH with topography regulated simulated CH₄ emissions, but emissions were also altered by rate-limited diffusion in soil microsites. Changes in simulated available substrate for CH₄ production (acetate, CO₂, and H₂) and oxidation (O₂ and CH₄) increased the predicted biomass of methanotrophs during the drought event and methanogens during drought recovery, which in turn affected net emissions of CH₄. A variance-based sensitivity analysis suggested that parameters related to aceticlastic methanogenesis and methanotrophy were most critical to simulate net CH₄ emissions. This study enhanced the predictive capability for CH₄ emissions associated with complex topography and drought in wet tropical forest soils.

1 Introduction

Wet tropical forest soils contribute significantly to global emissions of methane (CH₄; Pachauri et al., 2014). Although net emissions of CH₄ from upland soils are infrequent in temperate climates, studies show that CH₄ emissions are common in wet tropical forests, even in upland soils (Cattânio et al., 2002; Keller and Matson, 1994; Silver et al., 1999; Teh et al., 2005; Verchot et al., 2000). Landscape topography can strongly influence the proportions of CH₄ production and oxidation in mountainous tropical regions, affecting net emissions (Silver et al., 1999; O'Connell et al., 2018). Climate, and specifically patterns in rainfall, also affects emissions from tropical forests. Climate change may increase the frequency and severity of extreme rainfall and drought events, altering the spatial and temporal dynamics of CH₄ emissions through changes in redox dynamics and substrate availability (Silver et al., 1999; Chadwick et al., 2016; Neelin et al., 2006). Thus, accurately estimating CH₄ emissions under a variety of climatic and topographic conditions is important for predicting soil carbon–climate feedbacks in the humid tropical biome.

Several studies have reported the effect of drought events on biogenic CH₄ emissions across different wet tropical forest soils. For example, Aronson et al. (2019) demonstrated that the lower soil moisture conditions during the 2015–2016 El Niño event increased consumption of atmospheric CH₄ in a wet tropical forest Oxisol of Costa Rica. Similarly, a large-scale, 5-year throughfall exclusion experiment in a moist tropical forest Oxisol in Brazil also reported increased consumption of atmospheric CH₄ under the drought treatment, followed by a recovery of CH₄ emissions to pre-treatment values after the experiment ceased (Davidson et al., 2004, 2008). Using rainout shelters, Wood and Silver (2012) found spatial variability in CH₄ oxidation rates, with an increase of 480 % uptake in valleys in an Ultisol in Puerto Rico. More recently, in a similar Puerto Rico Ultisol, O'Connell et al. (2018) reported increasing consumption of atmospheric CH₄ during a Caribbean drought event, followed by increased production of CH₄ after the drought was over. The post-drought net CH₄ emission rates were higher than the pre-drought emissions, such that the benefits to atmospheric radiation imparted by the lowered emissions during the drought were eliminated. The sharp differences between pre- and post-drought emissions suggested that drought affected the balance of methanogenesis and methanotrophy in the soils, but the study lacked analysis of the microbial community's contributions to these two separate processes.

The concept of “microsites” inside soil aggregates or within soil micropores can help explain the coexistence of oxidative and reductive processes in soils (Silver et al., 1999; Teh and Silver, 2006), which may have occurred in the post-drought period in the O'Connell et al. (2018) study. Oxygen can remain inside micropores during saturated conditions and thereby maintain aerobic microbial respiration; likewise, hypoxic conditions can persist in microsites under extended droughts and thereby maintain anaerobiosis. Additionally, liquid substrates for methanogenesis such as acetate can accumulate in microsites under dry conditions because their diffusion to hungry microbial communities may be restricted. Conversely, gaseous substrates such as CO₂ and H₂ may accumulate in microsites under saturated conditions because gaseous diffusion can be limited. The observed rapid flush of CH₄ in response to a post-drought wetting event (O'Connell et al., 2018) suggests methanogenesis continued during the drought in the Ultisol's microsites, despite low soil moisture and high O₂ supply (Andersen et al., 1998; Bosse and Frenzel, 1998; Teh et al., 2005; von Fischer and Hedin, 2002). Finely textured soils common to the humid tropics can facilitate the co-existence of reduced solute and gas species with O₂ because the rate of solute and gaseous exchanges is controlled by diffusion into and out of microaggregates (Hall and Silver, 2013; Liptzin et al., 2011; Silver et al., 2013). In particular, hematite precipitation on clay minerals, found in both Oxisols and Ultisols, can enhance formation of soil aggregates because of their high surface area and charge properties (Hall et al., 2016). Soil organic matter can also enhance aggregation and at the same time consume O₂ (Six et al., 2004). However, few if any measurements of microsites exist in real field soils.

To explain the diverse observations of CH₄ emissions during and after drought across a wet tropical forest catena, we hypothesized that explicit representations of diffusion into and out of microsites for gas and solute transport would be required. To account for the balance of methanotrophy and methanogenesis, separate microbial functional groups for CH₄ production and oxidation would need to be defined. Therefore, a microbial functional group model for CH₄ production and consumption (Xu et al., 2015) was merged with a soil diffusivity module (Davidson et al., 2012; Sihi et al., 2018) to simulate the dynamics of net in situ CH₄ emissions from soil microsites (Sihi et al., 2020a). This module considers three key mechanisms for CH₄ production and consumption: aceticlastic methanogenesis (production from acetate), hydrogenotrophic methanogenesis (production from H₂ and CO₂), and aerobic methanotrophy (oxidation of CH₄ and reduction of O₂) (Fig. 1). Here we report a modeling experiment to explain contrasting patterns of observed CH₄ emissions following a severe drought in 2015, and we provide new data to describe CH₄ emissions under non-drought conditions in 2016. We explicitly account for changes in soil moisture, O₂, acetate, and microbial functional group dynamics within soil microsites in the model.

Figure 1 Conceptual figure of the modeling approach. Panel (a) shows the model representation of soil microsite distribution (modified from Sihi et al., 2020a; also see Eq. 14). The cylinder refers to the volume beneath the soil chambers. The intensity of different cylinder colors in the figure refers to rate of a process or the intensity of a concentration inside microsites in each theoretical cylinder; e.g., a dark color means a higher rate/intensity, and a light color means a lower rate/intensity for a given process. The 2D graph on the right refers to the probability density function of the rate of the process or intensity of the concentration in the bulk soil. A wide distribution skewed to the right (dark line) implies higher bulk rates of the process or higher concentrations, and a narrow distribution skewed to the left (light line) implies lower bulk rates of the process or lower concentrations, of any of the following: solute concentration [S_i], gas concentration [G_i], soil moisture (SoilM_i), gas and solute diffusion (Diff_i), methane production (Prod_i), and methane oxidation (Ox_i). Panel (b) is the schematic of the microbial-functional-group-based model for simulating soil methane (CH₄) dynamics in field soils (modified from Xu et al., 2015). The schematic represents the decomposition of soil organic matter (SOM) and plant litter into carbon dioxide (CO₂) and dissolved organic matter (DOC); the production of acetate and hydronium ion (H⁺) from decomposition and fermentation of DOC, which also decreases pH; the production of acetate and hydronium ion (H⁺) from homoacetogenesis, which decreases pH; and the production of dihydrogen ion (H₂) and CO₂ from decomposition of DOC. The intermediary products then have three possible non-mutually exclusive pathways: (1) aceticlastic methanogenesis, which is the production of methane from aqueous acetate found in soil solutions; (2) hydrogenotrophic methanogenesis, which is the production of methane from hydrogen; and (3) methanotrophy, which is the oxidation of methane into carbon dioxide.

2 Materials and methods

2.1 Study site

The study was conducted across a tropical forest catena near the El Verde Field Station in the Luquillo Experimental Forest (LEF) in northeastern Puerto Rico in the United States (latitude 18^∘19^′16.83^′′ N, longitude 65^∘49^′10.13^′′ W). The site is part of a National Science Foundation Long-Term Ecological Research (LTER) and Critical Zone Observatory (CZO) site and is also part of the US Department of Energy's Next Generation Ecosystem Experiment–Tropics. The mean annual temperature at the site is 23 ^∘C, and the long-term mean rainfall is approximately 3500 mm yr⁻¹ with low seasonality (Scatena, 1989). Inter-annual variability of rainfall ranges between 2600 and 5800 mm yr⁻¹, sometimes associated with extreme rainfall events (approximately 100 mm d⁻¹) from Caribbean storm systems (Heartsill-Scalley et al., 2007). The LEF is classified as a wet tropical forest according to the Holdridge life zone system, which considers rainfall, elevation, latitude, humidity, and evapotranspiration (Harris et al., 2012).

The landscape at the field site is highly dissected with short catenas, characterized by a land surface distance of <30 m from ridgetop to valley (O'Connell et al., 2018). This study partitioned sampling along a catena from ridgetop, slope, and valley topographic positions (Fig. S1). The soils are clay-rich Ultisols, which were derived from basaltic and andesitic volcanoclastic parent materials. Soils are acidic (average pH is 4.3 and 5.1 in ridge and valley topographic positions, respectively, Fig. 2). The valley soils have approximately 30 % clay and approximately 15 % sand, while the ridge soils have approximately 22 % clay and approximately 30 % sand (Brenner et al., 2019). The soils contain high concentrations of iron (Fe) and aluminum (Al) (oxy)hydroxides where their relative concentrations vary along the catena, and differences in Fe speciation are associated with variable redox conditions (Hall and Silver, 2013, 2015). A detailed soil survey in the immediate vicinity lists three soil types (Zarzal, Cristal, and Prieto) with a minimal litter layer due to rapid decomposition (Parton et al., 2007; Cusack et al., 2009), surface (A) horizons 5 cm thick, and B horizons 130 to 150 cm thick (Soil Survey Staff, 1995). The surface soil bulk density ranged from 0.5 to 0.7 g cm⁻³, and by 25 cm depth it ranged from 0.8 to 1.1 g cm⁻³ (Cabugao et al., 2021), similar to previous observations (Johnson et al., 2014; Silver et al., 1999). The forest composition is relatively diverse, with the mature tabonuco (Dacryodes excelsa Vahl) and sierra palm (Prestoea montana) trees being most dominant (Scatena and Lugo, 1995; Wadsworth, 1951).

Figure 2 Soil and porewater chemistry (dissolved organic carbon [DOC], a; acetate, b; and pH, c) along the ridge–slope–valley topographic gradient.

2.2 Soil and porewater sampling

Previous CH₄ measurements in the LEF at the soil surface, 10 cm depth, and 35 cm depth found the highest CH₄ concentrations at 10 cm depth (Silver et al., 1999), while 30 cm depth was the location of maximum soil organic carbon (SOC) concentrations (Johnson et al., 2014). To initialize the model, soil and soil water samples were collected from depths ranging from 0–30 cm in accordance with these previous studies. Soil samples were collected in triplicate from a depth of 0–10 cm and on a quarterly timeframe from the ridgetop, slope, and valley positions for over 2 years. The soil pH was determined using a 1 : 2 ratio of soil : solution using a glass electrode with 0.005 M CaCl₂ as the equilibrated soil solution (Thomas, 1996; Sihi et al., 2020b). Porewater samples were collected approximately weekly for over 2 years using MacroRhizon soil water samplers (length = 5 cm) (Rhizosphere Research Products B.V., Wageningen, the Netherlands) installed at both 5–10 and 25–30 cm depth in triplicate in the ridge, slope, and valley topographic positions (Sihi et al., 2020c). The soil water samples were analyzed for organic acid concentrations (acetate) using high-performance liquid chromatography (Dionex ICS-5000+ Thermo-Fisher Waltham, MA, USA) with the Dionex IonPac AS11-HC column using a potassium hydroxide eluent and gradient elution. The samples were analyzed for total dissolved organic carbon (DOC) using a Shimadzu total organic C analyzer (Shimadzu TOC-L CSH/CSN analyzer, Baltimore, MD, USA). The soil and porewater measurements were conducted in 2017–2018 (the number of samples n ranged between 20 and 35, Fig. 2) to initialize different model parameters for the catena, because measurements were not available for 2015–2016. To that end, the chemical data were used as the reference characteristics of the bulk soil, and the temporal evolution of DOC, acetate, and soil pH at the microsites was calculated using probability distributions of soil moisture and O₂ across soil microsites over the 2-year measurement window. Soil bulk density and particle density values were taken from O'Connell et al. (2018).

2.3 In situ methane flux and soil driver measurements

Campbell Scientific CS 655 soil moisture and temperature sensors and Apogee SO-110 O₂ sensors were co-located with soil gas flux chambers at 15 cm soil depth along the catena, each with five replications along five transects (Fig. S1) (O'Connell et al., 2018). Following Liptzin et al. (2011), soil O₂ sensors were installed in gas-permeable soil equilibration chambers (295 cm³). Data from these sensors were collected hourly using Campbell Scientific CR10000 data loggers and AM16/32B multiplexers (Campbell Scientific, Logan, UT, USA), which were processed using site-based calibration equations.

Soil flux chambers were placed on the top of the soil surface. Soil CH₄ emissions along the catena were measured during 2015 (26 February to 23 December, O'Connell et al., 2018; Silver, 2019) and 2016 (5 April to 18 July) (Sihi et al., 2020d) using a cavity ring-down spectroscopy gas analyzer (Picarro G2508, Santa Clara, CA, USA) connected to 12 automated eosAC closed dynamic soil chambers (Pumpanen et al., 2004) using a multiplexer (Eosense Inc., Dartmouth, Nova Scotia, Canada). Data for soil CH₄ emissions were processed using eosAnalyze-AC (v3.5.0) software followed by a series of quality control protocols (O'Connell et al., 2018). We used daily average values of drivers (soil temperature, soil moisture, and O₂ concentrations) and CH₄ emissions in the modeling exercise. See O'Connell et al. (2018) for more information on the soil sensor, chamber arrays, and the data analysis pipeline.

The data from the 2015 Caribbean drought was partitioned into four distinct periods (O'Connell et al., 2018): (1) pre-drought from day of year (DOY) 57 to 115 (dark gray on Fig. 3), (2) the drought from DOY 116 to 236 (medium gray on Fig. 3), (3) drought recovery from DOY 237 to 328 (light gray on Fig. 3), and (4) post-drought from DOY 329 to 354 (white on Fig. 3). Total precipitation during the drought period was 700 mm in 2015 and 1088 mm during the same time frame in 2016 (Zimmerman, 2000).

Figure 3 Temporal dynamics of observed meteorological drivers (soil temperature, a–c; soil moisture, d–f; soil oxygen, g–i) and net methane emissions (j–l) for 2015 (Data are taken from O'Connell et al., 2018). For methane emissions, symbols represent observed data and lines represent model simulations. Dark gray, medium gray, light gray, and white shading represent pre-drought, drought, drought recovery, and post-drought events (O'Connell et al., 2018).

2.4 Modeling approach

2.4.1 Microbial functional group model for methane production and oxidation

An existing microbial-functional-group-based model for CH₄ production and consumption (Xu et al., 2015) was adopted for this research (Sihi, 2020). As shown in Fig. 1, acetate and H₂ $/$ CO₂ represent substrate $\left[{\mathrm{Substrate}}_{{\mathrm{func}}_i}\right]$ (nM cm⁻³) for aceticlastic and hydrogenotrophic methanogenesis reactions, respectively. On the other hand, CH₄ and O₂ concentrations represent substrate for the methanotrophy reaction. Acetate and CO₂ are inputs based on measurements of soil water and pH described in Sect. 2.2. In the model, acetate is formed by fermentation and by homoacetogenesis (but not by syntrophic acetate oxidation) as defined in Xu et al. (2015) in their Appendix in Eqs. (A15) and (A16) (Fig. 1b). Methylotrophic methanogenesis (Narrowe et al., 2019) is neglected in the model. The overall reaction rates are represented as

$$\begin{array}{}\text{(1)}& \begin{array}{rl}{\mathrm{Reaction}}_{{\mathrm{rate}}_i}& ={\mathrm{Biomass}}_{{\mathrm{func}}_i}\times {\displaystyle \frac{{\mathrm{GrowR}}_{{\mathrm{func}}_i}}{{\mathrm{Efficiency}}_{{\mathrm{func}}_i}}}\\ & \times {\displaystyle \frac{\left[{\mathrm{Substrate}}_{{\mathrm{func}}_{\mathrm{1}\mathrm{\dots }n}}\right]}{\left[{\mathrm{Substrate}}_{{\mathrm{func}}_{\mathrm{1}\mathrm{\dots }n}}\right]+{\mathrm{KM}}_{{\mathrm{func}}_{\mathrm{1}\mathrm{\dots }n}}}}\\ & \times f\left(T\right)\times f\left(\mathrm{pH}\right),\end{array}\end{array}$$

where ${\mathrm{Reaction}}_{{\mathrm{rate}}_i}$ (in nM cm⁻³ h⁻¹) is the rate of CH₄ production and/or consumption under variable substrate concentrations. ${\mathrm{Biomass}}_{{\mathrm{func}}_i}$ (nM cm⁻³) represents microbial functional groups: aceticlastic methanogens, hydrogenotrophic methanogens, and aerobic methanotrophs, respectively. Growth rates and substrate use efficiencies of microbial functional groups are represented as ${\mathrm{GrowR}}_{{\mathrm{func}}_i}$ (h⁻¹) and ${\mathrm{Efficiency}}_{{\mathrm{func}}_i}$ (unitless), respectively (Table 1). The substrate limitation on CH₄ production is imposed by assuming a Michaelis–Menten relationship with the half-saturation constants for CH₄ production and oxidation being ${\mathrm{KM}}_{{\mathrm{func}}_{\mathrm{1}\mathrm{\dots }n}}$ (nM cm⁻³). Although minor contributions of iron-dependent anaerobic CH₄ oxidation to net CH₄ emissions can be expected in our study site (Ettwig et al., 2016), we did not represent this process here as anaerobic oxidation of CH₄ is still not fully understood, and it is generally low in most ecosystems.

Table 1 Fitted values of M3D-DAMM model parameters.

Initial values of model parameters were collected from literature (“Source”). Also see Xu et al. (2015) for detailed information on model parameters.

The extent of change in ${\mathrm{Biomass}}_{{\mathrm{func}}_i}\left({\mathrm{dBiomass}}_{{\mathrm{func}}_i}\right)$ is controlled by the balance between ${\mathrm{Growth}}_{{\mathrm{func}}_i}$ and ${\mathrm{Death}}_{{\mathrm{func}}_i}$ following

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{{\mathrm{dBiomass}}_{{\mathrm{func}}_i}}{{\mathrm{dt}}_{{\mathrm{func}}_i}}}={\mathrm{Growth}}_{{\mathrm{func}}_i}-{\mathrm{Death}}_{{\mathrm{func}}_i},\text{(3)}& {\displaystyle }{\mathrm{Growth}}_{\mathrm{func}}={\mathrm{Efficiency}}_{{\mathrm{func}}_i}\times {\mathrm{Reaction}}_{{\mathrm{rate}}_i},\end{array}$$

where ${\mathrm{Growth}}_{{\mathrm{func}}_i}$ is calculated as a multiplicative function of ${\mathrm{Efficiency}}_{{\mathrm{func}}_i}$ and the ${\mathrm{Reaction}}_{{\mathrm{rate}}_i}$,

$$\begin{array}{}\text{(4)}& {\mathrm{Death}}_{{\mathrm{func}}_i}={\mathrm{DeadR}}_{{\mathrm{func}}_i}\times {\mathrm{Biomass}}_{{\mathrm{func}}_i},\end{array}$$

and ${\mathrm{Death}}_{{\mathrm{func}}_i}$ is a function of ${\mathrm{DeadR}}_{{\mathrm{func}}_i}$ (death rate, Table 1) and ${\mathrm{Biomass}}_{{\mathrm{func}}_i}$ (microbial biomass).

All rate equations were modified by the scalers for temperature, f(T) and pH, f(pH) functions, described below. We represented the temperature effect, f(T), using a classic Q₁₀ function:

$$\begin{array}{}\text{(5)}& f\left(T\right)=Q_{\mathrm{10},i}^{\frac{{\mathrm{Temperature}}_{\mathrm{soil}}-{\mathrm{Temperature}}_{\mathrm{reference}}}{\mathrm{10}}}.\end{array}$$

We represented the pH effect, f(pH), based on Cao et al. (1995):

$$\begin{array}{}\text{(6)}& f\left(\mathrm{pH}\right)={\displaystyle \frac{\left(\mathrm{pH}-{\mathrm{pH}}_{\mathrm{minimum}}\right)\cdot \left(\mathrm{pH}-{\mathrm{pH}}_{\mathrm{maximum}}\right)}{{\begin{array}{c}\left(\mathrm{pH}-{\mathrm{pH}}_{\mathrm{minimum}}\right)\cdot \left(\mathrm{pH}-{\mathrm{pH}}_{\mathrm{maximum}}\right)\\ -\left(\mathrm{pH}-{\mathrm{pH}}_{\mathrm{optimum}}\right)\end{array}}^{\mathrm{2}}}},\end{array}$$

where we set the minimum, optimum, and maximum soil pH values to 4, 7, and 10, respectively. Following Xu et al. (2015), we considered the contribution of acetate to pH as follows:

$$\begin{array}{}\text{(7)}& \mathrm{pH}=-\mathrm{1}\cdot \log \left({\mathrm{10}}^{{\mathrm{pH}}_{\mathrm{initial}}}+\mathrm{4.2}E-\mathrm{9}\cdot \mathrm{Acetate}\right).\end{array}$$

Although other mechanisms to alter soil pH are present at the site, e.g., Fe reduction and oxidation (Teh et al., 2005; Hall and Silver, 2013), these are not considered in the model at this time. Calibrated values of ${\mathrm{GrowR}}_{{\mathrm{func}}_i}$ ${\mathrm{DeadR}}_{{\mathrm{func}}_i}$, ${\mathrm{Efficiency}}_{{\mathrm{func}}_i},{\mathrm{KM}}_{{\mathrm{func}}_i}$, and Q_10,i are presented in Table 1.

2.4.2 Diffusion module for gaseous and solute transport in soil profile and across soil–air boundary

In order to account for the diffusion of gases across the soil–air boundary and solutes (e.g., acetate) through soil water films (Fig. 1), we added the diffusion module of the Dual Arrhenius and Michaelis–Menten (DAMM) model (Davidson et al., 2012; Sihi, 2020; Sihi et al., 2018, 2020a) to the existing microbial functional group model, which we refer to as M3D-DAMM. We calculated initial concentration of gases like O_2, H₂, CO_2, and CH₄, [Gas_conc] (unit: V V⁻¹), as a function of a unitless diffusion coefficient of gas in air (D_gas), volume fraction of gas in air (V V⁻¹), and gas diffusivity ($a^{\mathrm{4}/\mathrm{3}}$ ) as follows:

$$\begin{array}{}\text{(8)}& \left[{\mathrm{Gas}}_{\mathrm{conc}}\right]={\mathrm{D}}_{\mathrm{gas}}\times \text{atmospheric concentration}\times a^{\mathrm{4}/\mathrm{3}},\end{array}$$

where $a^{\mathrm{4}/\mathrm{3}}$ represents the tortuosity of diffusion pathway for gases as a function of soil water (SoilM) and temperature (SoilT):

$$\begin{array}{}\text{(9)}& \begin{array}{rl}{a}^{\mathrm{4}/\mathrm{3}}& ={\left(\mathrm{Porosity}-{\displaystyle \frac{\mathrm{SoilM}}{\mathrm{100}}}\right)}^{\mathrm{4}/\mathrm{3}}\\ & \times {\left({\displaystyle \frac{\mathrm{SoilT}+\mathrm{273.15}}{\mathrm{293.15}}}\right)}^{\mathrm{1.75}},\end{array}\end{array}$$

where the air-filled porosity (a) was calculated by subtracting the volume fraction of soil moisture (V V⁻¹) from total porosity. Porosity was calculated as

$$\begin{array}{}\text{(10)}& \left(\mathrm{1}-{\displaystyle \frac{\mathrm{Bulk}\mathrm{density}}{\mathrm{Particle}\mathrm{density}}}\right).\end{array}$$

The exponent of $\mathrm{4}/\mathrm{3}$ accounts for diffusivity of gases through porous media (Davidson and Trumbore, 1995). The exponent of 1.75 represents the temperature response of gaseous diffusion (Massman, 1998; Davidson et al., 2006). Following Davidson et al. (2012), the value used for gaseous diffusivity coefficient (D_gas) was calculated based on an assumed boundary condition such that the concentration of gaseous substrates in the soil pore space would be equivalent to the volume fraction of gases in air under completely dry conditions.

We assumed another boundary condition to determine the value of the aqueous diffusion coefficient, D_liq, such that soluble substrates like acetate would be available at the enzymatic reaction site under conditions with saturating soil water content (Davidson et al., 2012):

$$\begin{array}{}\text{(11)}& D_{\mathrm{liq}}={\displaystyle \frac{\mathrm{1}}{{\mathrm{Porosity}}^{\mathrm{3}}}}.\end{array}$$

We represented soluble substrates (acetate) diffused through a soil water film as Aqueous-substrate_av (µM L⁻¹), which we calculated as follows:

$$\begin{array}{}\text{(12)}& \begin{array}{rl}{\text{Aqueous-substrate}}_{\mathrm{av}}& =\text{Aqueous-substrate}\\ & \times D_{\mathrm{liq}}\times {\left({\displaystyle \frac{\mathrm{SoilM}}{\mathrm{100}}}\right)}^{\mathrm{3}},\end{array}\end{array}$$

where the ${\left(\frac{\mathrm{SoilM}}{\mathrm{100}}\right)}^{\mathrm{3}}$ term represents the diffusion rate of aqueous substrates to the enzymatic active site (Papendick and Campbell, 1981). Concentrations of acetate in the aqueous phase (µM L⁻¹) were obtained from the measurements across the catena averaged by depths (10 and 30 cm) of Rhizon samplers.

We calculated CH₄ emissions, CH_4emission (unit: µmole m⁻² h⁻¹), as a function of concentration ([CH_4conc]), production (CH_4prod), and oxidation (CH_4ox) of CH₄, multiplied by the equivalent depth (set to 15 cm) (for cm⁻³ volume to cm⁻² area conversion) and 10⁴ (for m² to cm² conversion) as follows:

$$\begin{array}{}\text{(13)}& \begin{array}{rl}{\mathrm{CH}}_{\mathrm{4}\mathrm{emission}}& =\left[{\mathrm{CH}}_{\mathrm{4}\mathrm{conc}}\right]+\left({\mathrm{CH}}_{\mathrm{4}\mathrm{prod}}-{\mathrm{CH}}_{\mathrm{4}\mathrm{ox}}\right)\\ & \times {\mathrm{10}}^{\mathrm{4}}\times \mathrm{depth}.\end{array}\end{array}$$

2.4.3 Soil microsites

The importance of diverse microsite conditions was inferred based on many previous observations in the field and the lab of co-occurrences of oxic soil concentrations and reduced redox-active species (Silver et al., 1999; Teh et al., 2005; Megonigal and Geunther 2008; Hall et al. 2013, 2016; Sihi et al., 2020a). The high clay content, abundant Fe oxides, and visible redox mottling, particularly in the valley and slope soils, facilitate a diversity of soil microenvironments where O₂ and CH₄ can seemingly co-occur, albeit in different microsite locations (Silver et al., 1999; Teh and Silver, 2006). Microsite diversity was also invoked to help explain the rapid CH₄ emissions following drought at the field site (O'Connell et al., 2018). Techniques for accurately measuring in situ microsite activities remain very limited to date, here or elsewhere. Therefore, we simulated production, consumption, and diffusion processes within soil microsites using a log-normal probability distribution function of soil moisture and available C based on these previously observed relationships (Fig. 1). The average values of individual processes across simulated microsites (represented by “i”) represent the reaction in the bulk soil, which we constrained using the net measured CH₄ emissions:

$$\begin{array}{}\text{(14)}& \mathrm{Bulk}{\mathrm{soil}}_{\mathrm{average}}={\displaystyle \frac{\sum {\mathrm{Frequency}}_i\times [\mathrm{microsite}{]}_i}{\mathrm{Total}\mathrm{microsites}}}.\end{array}$$

We directly adopted the probability distribution function of soil moisture and C from Sihi et al. (2020a), which constrained values of Frequency_i of soil microsites. We a priori assigned the size of the microsites to be at least an order magnitude smaller than the diameter used for bulk measurements of CH₄ fluxes. Thus, the mean diameter of microsites was assumed to be at the millimeter scale (the size class of small stable aggregates in these soils), as the diameter of soil chambers was 15.24 cm. Thus, the resultant number of total microsites below each soil flux chamber was 10 000.

2.4.4 Sensitivity analysis

We evaluated the sensitivity of model parameters with a global variance-based sensitivity analysis using the R-multisensi package. This method uses a global sensitivity index (0 < GSI < 1) to determine the sensitivity of CH₄ emissions to model parameter values (Bidot et al., 2018). We conducted a multivariate technique to estimate GSI values in sequential steps. First, we implemented a factorial design on the uncertain model parameters, which is followed by a principal component analysis on model outputs. Then, we extracted GSI values by an ANOVA-based sensitivity analysis on the first principal component. To that end, parameters with high GSI values may explain high temporal variations of the observed CH₄ emissions, and those with low GSI values are insignificant to reproduce the temporal dynamics of CH₄ emissions.

2.4.5 Statistical analysis

We used R (version 3.5.1) for statistical analyses, modeling, and visualization purposes (R Core Team, 2018). Statistical analyses and figures were produced using R-ggstatsplot (Patil, 2018) and R-ggplot2 (Wickham, 2016) packages. Differences in soil and porewater chemistry across the catena were compared using robust t test. Correlograms for soil temperature, soil moisture, O₂, and soil CH₄ emissions were created using adjusted Holm correlation coefficients. All statistical analyses were conducted at the 5 % significance level. We implemented the M3D-DAMM model using R-FME package (Soetaert, 2016).

3 Results

3.1 Observational dynamics of soil biogeochemistry

Soil and porewater chemistry varied along the catena (Fig. 2). Dissolved organic carbon (DOC) values followed the trend of ridge > slope > valley (p≤0.001). Soil DOC concentrations (mean ± SE) were 0.55 ± 0.10, 0.30 ± 0.03, and 0.18 ± 0.03 mg g⁻¹ in ridge, slope, and valley soils, respectively. Organic acid (acetate) concentrations were significantly higher in the ridge (6.57 ± 1.48 µM L⁻¹) and slope (6.42 ± 2.19 µM L⁻¹) than in the valley (1.80 ± 0.20 µM L⁻¹) (p=0.003). Soil pH followed the trend of valley > slope > ridge (p<0.001). Average soil pH ranged from 4.25 ± 0.11 in the ridge to 4.49 ± 0.08 in the slope and to 5.05 ± 0.09 in the valley.

Soil moisture and soil O₂ concentrations were distinctly different in the drought year (2015) compared to 2016. The drought in 2015 decreased soil moisture in the slope and ridge soils and increased O₂ concentrations in the valley soils (Fig. 3) (also see O'Connell et al., 2018). Generally, average soil moisture was higher in the valley (0.47 ± 0.05 in 2015 and 0.51 ± 0.01 % in 2016) as compared to the ridge (0.31 ± 0.12 in 2015 and 0.39 ± 0.03 % in 2016) and slope (0.30 ± 0.16 in 2015 and 0.41 ± 0.04 % in 2016). Average O₂ concentrations were generally lower in the valley (11.54 ± 5.94 in 2015 and 6.30 ± 2.96 % in 2016) as compared to the ridge (18.37 ± 0.72 in 2015 and 17.52 ± 0.42 % in 2016) and slope (18.09 ± 1.22 in 2015 and 16.89 ± 0.58 % in 2016). After the drought ended, the recovery of soil moisture in the ridge and slope soils proceeded more quickly than the recovery of O₂ concentrations in the valley soils (Fig. 3). Soil temperature ranges were averaged across the topographic gradient and were similar in both years (average was 21.58 ± 1.88 in 2015 and 22.97 ± 1.04 ^∘C in 2016).

In 2016, net CH₄ emissions were generally positive in the valley and were marginally negative in the ridge and slope (Fig. 4). The dynamics of CH₄ were very different following the 2015 drought, resulting in net positive CH₄ emissions in the post-drought period for all topographic positions (Fig. 3) (as described in more detail in O'Connell et al., 2018). The magnitude of CH₄ emissions was greater in the valley, followed by the slope and then the ridge.

Figure 4 Temporal dynamics of observed meteorological drivers (soil temperature, a–c; soil moisture, d–f; soil oxygen, g–i) and net methane emissions (j–l) for 2016. For methane emissions, symbols represent observed data and lines represent model simulations.

The strength of the relationships between net CH₄ emissions and soil temperature, moisture, and O₂ concentrations was contingent on both topographic position and year (2015 vs. 2016) (Fig. 5). For example, the relation between CH₄ emissions and soil moisture was stronger in 2016 (normal year) than in 2015 (drought year). The correlation between CH₄ emissions and O₂ concentrations was stronger and more negative in 2015 than 2016. Correlations between soil moisture and O₂ concentrations were negative and stronger in 2016 than in 2015. Correlation coefficients between soil O₂ concentrations and CH₄ emissions were negative and strongest for valley soils and lowest for ridge soils in 2015, but they were uncorrelated in 2016 for ridge and slope soils (Fig. S2).

Figure 5 Relation between soil meteorology and methane emissions for 2015 (a) and 2016 (b). SoilM, SoilT, O₂, and CH₄ represent soil moisture, soil temperature, oxygen, and methane, respectively. Numbers represent adjusted Holm correlation coefficients, and numbers with “X” indicate a non-significant correlation at p<0.05.

3.2 Model simulations of methanogenesis and methanotrophy

In general, there was little bias in the relationships between the observed and simulated CH₄ emissions (Fig. 6). The model explained 72 % and 67 % of the variation in soil CH₄ emissions for 2015 and 2016, respectively, although the model performance varied across the catena (Figs. 6, S3, S4). Overall, simulated CH₄ emissions captured the trend of valley > slope ≥ ridge for 2016. The model also captured the dramatically different dynamics of field CH₄ emissions as a function of topography during and after the 2015 drought. Net positive CH₄ emissions were simulated in the drought recovery and post-drought periods in the ridge and slope in 2015, while net negative emissions were simulated in the other times for these landscape positions. Additionally, simulated net CH₄ emissions were decreased during the drought and drought recovery in the valley soils, as well as the strong net CH₄ emissions in the valley soils in the post-drought period.

Figure 6 Observed versus simulated methane (CH₄) emissions and model residuals for 2015 (a, b) and 2016 (c, d).

The ridge and slope positions were more similar to each other than to the valley soils. Simulated decreased production of acetate and hydrogen during the 2015 drought in the ridge and slope positions resulted in decreased biomass of aceticlastic methanogens and hydrogenotrophic methanogens (Figs. S5, S6). Gross CH₄ production therefore decreased during these time periods (Fig. S7). Simultaneously, as soil moisture decreased, simulated methanotrophic biomass increased during the drought (Fig. S5). The simulated biomass of both aceticlastic methanogens and hydrogenotrophic methanogens increased dramatically in the ridge and slope soils during the drought recovery (aceticlastic methanogens: 3.3 and 5.3 times higher than drought period for ridge and slope, respectively; hydrogenotrophic methanogens: 6.1 and 12 times higher than drought period for ridge and slope, respectively) and post-drought (aceticlastic methanogens: 5.2 and 8.8 times higher than drought period for ridge and slope, respectively; hydrogenotrophic methanogens: 12 and 24 times higher than drought period for ridge and slope, respectively) period. Concomitantly, production of acetate and H₂ was much higher in the ridge and slope soils during the drought recovery (acetate: 1.8 and 2.4 times higher than the drought period for ridge and slope soils, respectively; H₂: 3.5 and 6.0 times higher than the drought period for ridge and slope soils, respectively) and the post-drought (acetate: 2.3 and 3.2 times higher than the drought period for ridge and slope, respectively; H₂: 5.6 and 10 times higher than the drought period for ridge and slope, respectively) period. Together, gross CH₄ production in the ridge and slope soils was significantly higher during the drought recovery (1.9 and 2.5 times higher than the drought period for ridge and slope, respectively) and post-drought periods (3.4 and 4.6 times higher than the drought period for ridge and slope, respectively) compared to the drought (Fig. S7). Simulated production of acetate was increased, which also lowered soil pH values during drought recovery (Fig. S6), with a more pronounced effect in the ridge and slope soils. Additionally, simulated methanotrophic biomass and CH₄ oxidation decreased during the post-drought period (Figs. S5, S7), which is the same time period during which net CH₄ production increased strongly.

For the valley soils, simulated values of aceticlastic methanogens and concomitant acetate production increased during the 2015 drought (Figs. S5, S6). During the drought recovery and post-drought period, both aceticlastic methanogens and acetate production decreased in the valley, while hydrogenotrophic methanogens and H₂ production were stable. Gross CH₄ production, however, remained relatively flat during the drought event in the valley and only increased during the post-drought period (Fig. S7). Simulated CH₄ oxidation and methanotrophic biomass, on the other hand, increased dramatically during the drought and drought recovery period (Figs. S5, S7) and then decreased strongly during the post-drought period. However, simulated methanotrophic biomass was smaller in the valley soils compared to the ridge and slope soils. Methane oxidation by methanotrophs exerted strong controls on simulated net CH₄ emissions, not only in the valley but in all the topographic positions.

3.3 The influence of microsites on net methane emissions

Concomitant with decreased soil moisture, the simulated diffusion of gases (O₂, H₂) was enhanced during the drought event in 2015, while diffusion of the solute (acetate) was dramatically decreased, particularly for the ridge and slope soils (Fig. S8). However, reduction in soil moisture and increase in O₂ can inhibit fermentative hydrogen production (Cabrol et al., 2017). Consequently, simulated gross CH₄ production through hydrogenotrophic and aceticlastic pathways both decreased during the drought event for the ridge and slope positions (Figs. S7, S9). As soil moisture increased during the drought recovery and post-drought periods, the diffusion of gases decreased, and diffusion of acetate increased in the ridge and slope soils (Fig. S8). Consequently, simulated values of gross CH₄ production increased and gross CH₄ oxidation decreased during drought recovery and the post-drought period (Fig. S7). These factors likely contribute to the large pulses of net CH₄ emissions during the post-drought period for ridge and slope positions (Fig. 3).

Overall, the valley soils were relatively insensitive to changes in the diffusion rate of either gases or solutes (Fig. S8), most likely because soil moisture remained relatively stable, regardless of drought conditions (Fig. 3). The lower sand and higher clay contents in the valley soils (Brenner et al., 2019), as well as the lower topographic position, likely caused the valley soils to remain wetter than the slope and ridge soils. Therefore, simulated values of gross CH₄ production were fairly stable in the valley soils (Fig. S7) during the drought and drought recovery period.

Simulated production, oxidation, and net flux of CH₄ was further modified by reactions occurring within soil microsites. For example, during the drought (∼ DOY 200 in 2015), gross CH₄ production was more frequent in soil microsites in the valley compared to the slope and ridge (Fig. 7). Simulated values of CH₄ oxidation were much greater in microsites in the slope and ridge positions, so the net CH₄ emissions were positive in the valley soils and negative in the ridge and slope positions. During the 2015 post-drought period (DOY 345), the frequency of CH₄ production was much greater in all topographic positions compared to the drought period (DOY 200), and it was also more enhanced in the valley soils compared to the slope and ridge. Thus, net positive CH₄ emissions were observed in all topographic positions in the post-drought period (Fig. 3). Methane oxidation at DOY 345 was much greater in the ridge and slope compared to the valley, similar to predictions at DOY 200. Therefore, the prominent CH₄ emissions from all three topographic positions were primarily due to increased production (CH₄ production on DOY 345 was 150 %, 248 %, and 80 % higher than DOY 200 in ridge, slope, and valley, respectively) rather than decreased oxidation (CH₄ oxidation was 32 %, 31 %, and 43 % lower on DOY 345 than DOY 200 in ridge, slope, and valley, respectively), which agrees with previous studies in our site (Teh et al., 2005, 2008; von Fischer and Hedin, 2002).

Figure 7 Rates of gross methane (CH₄) production (a, b), oxidation (c, d), and net flux (e, f) across simulated soil microsites. Day of year 200 and 345 represent drought and post-drought recovery, respectively (see medium gray and white shading in Fig. 3).

Diffusion into microsites strongly affected the concentrations of gases and solutes experienced by microbes, and differences as a function of topographic position were again predicted. Acetate production and diffusion were enhanced in valley soils during the drought, when compared to the slope and ridge soils (Fig. S10). The H₂ production was also enhanced in the valley soils during the drought, but the wetter valley soils experienced lower rates of H₂ diffusion compared to the ridge and slope soils. Increases in O₂ diffusion were also apparent in the ridge and slope soils during the drought, and those increases were greater than in the valley soils. During the post-drought period, however, the frequency of H₂ and O₂ diffusion was much greater for the ridge soils compared to the valley soils (Fig. S10).

Of all parameters, the most sensitive ones were those that controlled CH₄ production through the aceticlastic pathway, followed by the parameters related to CH₄ oxidation (Fig. 8). The GSI values for parameters related to aceticlastic methanogenesis and methanotrophy ranged between 0.25–0.75, whereas the corresponding GSI values for hydrogenotrophic methanogenesis were always <0.1.

Figure 8 Global sensitivity indices of M3D-DAMM model parameters (defined in Table 1). Gray, yellow, and blue colors represent parameters for aceticlastic methanogenesis, hydrogenotrophic methanogenesis, and methanotrophy, respectively.

4 Discussion

4.1 Mechanisms governing net methane emissions

Although the initial concentrations of available C for fermentation (i.e., DOC) and substrate for aceticlastic methanogenesis (i.e., acetate) in the bulk soil followed the trend of ridge > slope > valley (Fig. 2), the pattern of net CH₄ emissions across the catena was opposite (valley > slope ≥ ridge), especially in 2016 (Fig. 4). The seemingly counterintuitive relations of substrate concentrations in the bulk soil versus net CH₄ emissions can be explained by modeling the differing redox conditions across soil microsites. Diffusion promoted the availability of the acetate substrate through more connected soil water films in the wetter valley soils and caused higher gross CH₄ production in 2016, as compared to the relatively drier slope and ridge soils (Figs. S7, S8). In contrast, diffusion of gaseous methanotrophic substrates (CH₄ and O₂) was promoted in the air-filled pore spaces in the drier ridge and slope soils (Fig. S8), resulting in reduced net CH₄ emissions for these two topographic positions in 2016 (Fig. 4). Further, reduced diffusion of O₂ in the wetter valley soils decreased gross methanotrophy compared to the slope and ridge soils (Figs. S7, S8). Consequently, in 2016, net CH₄ emissions dominated the valley soils but were minimal in the ridge and slope soils.

On the other hand, the drought event in 2015 decreased the simulated CH₄ emission in the slope and ridge soils by decreasing H₂ production and both production (Fig. S6) and diffusion of acetate (Fig. S8). The drought increased the CH₄ sink strength of both ridge and slope soils as the observed net CH₄ emissions became more negative during the drought compared to the pre-drought period (Fig. 3). Contributing factors predicted by the model include enhanced O₂ diffusion into the drier ridge and valley soils (Fig. S8), as well as enhanced methanotrophic biomass (Fig. S5). In the valley, the primary impact of the drought appeared to be due to increased methanotrophy (Fig. S7), since acetate, H₂, and gross CH₄ production were predicted to continue unabated (Figs. S6, S7). This suggests that drought enhanced consumption of atmospheric CH₄ in our site, which is consistent with findings from natural droughts and throughfall exclusion experiments in other wet tropical forest soils (Aronson et al., 2019; Davidson et al., 2004, 2008; Wood and Silver, 2012).

However, simulation of observed CH₄ emission during drought recovery in 2015 required explicit representations of the complex interaction of the diffusive supply of solute and gases, dynamics of the microbial functional groups, and the associated acetate–pH feedback loop across the distribution of soil microsites (Fig. 3). The drought recovery increased soil moisture, which likely prompted anaerobiosis across all topographic locations by significantly reducing gas diffusivity in a fraction of the simulated microsites (11, 17, and 21 % in ridge, slope, and valley, respectively) (McNicol and Silver, 2014; Sihi et al., 2020a; Teh et al., 2005). The return to dominantly reducing conditions also was predicted to stimulate fermentation and the production of acetate through homoacetogenesis (Fig. S6). Enhanced production and diffusion of acetate during recovery (Fig.  S8) triggered growth in the predicted biomass of aceticlastic methanogens (Fig. S5), which in turn increased rates of aceticlastic methanogenesis (Fig. S9).

Simulated rates of hydrogenotrophic methanogenesis also increased in anaerobic microsites (Figs. S9, S10), mediated by increased production of H₂ and subsequent stimulation of the biomass of hydrogenotrophic methanogens during the drought recovery in 2015 (Fig. S5). Overall, the absolute values of simulated gross CH₄ production through hydrogenotrophic and aceticlastic pathways (Fig. S9) outweighed the simulated gross CH₄ oxidation rates (Fig. S7), resulting in net soil CH₄ emissions across the catena during the post-drought period (Fig. 3).

Acetate-driven CH₄ increases, decreases in methanotrophy due to decreasing O₂, and increasing hydrogenotrophic methanogenesis all contributed to the post-drought pulses of CH₄ (Fig. 8). Both kinds of methanogens increase during drought recovery and post-drought, but aceticlastic methanogens were 2 orders of magnitude more abundant than hydrogenotrophic methanogens. Additionally, acetate may accumulate in microsites during drought and then become more available with drought recovery due to enhanced solute diffusion (Fig. S8). The model simulations suggest that hydrogen diffusion was lessened under the drought recovery, which is consistent with decreasing rates of gas diffusion through saturated soils (Fig. S8). Further, H₂ has a faster turnover rate compared to acetate (Xu et al., 2015), and therefore accumulation in soils, especially shallow soils which are the subject of this study, is minimized. So, acetate versus hydrogen substrate availability in microsites better explains the observations of higher CH₄ production under the post-drought conditions.

Additionally, acetate is a source of proton and should reduce soil pH (Amaral et al., 1998; Conrad and Klose, 1999; Jones et al., 2003). Previous studies (Xu et al., 2015; Xu et al., 2010) demonstrated that acetate-driven soil pH reduction can reduce net CH₄ production by as much as 30 %, especially in systems with low initial soil pH like our study site. Given that optimal pH for biological activities peaks near neutral pH, the relatively higher soil pH in the valley versus ridge and slope soil further enhanced the topographic patterns of CH₄ emissions (Conrad, 1996; also see Figs. 2, 3, and 4). Note that the initial soil pH across the landscape was already in the acidic range (Fig. 2); consequently, the simulated acetate production and concomitant decrease in soil pH during the 2015 drought recovery further suppressed gross CH₄ production in ridge soils in comparison to the valley soils (Figs. S6 and S7). Iron-reducing bacteria can also suppress CH₄ production either by competing with aceticlastic methanogens for acetate substrate or controlling the flow of acetate to both hydrogenotrophic and aceticlastic methanogens by dissimilatory iron reduction (Teh et al., 2008). Additionally, Fe reduction can increase soil pH by proton consumption and colloid dispersion, whereas Fe oxidation can lead to more acidic conditions (Hall and Silver, 2013; Thompson et al., 2006). None of the Fe-associated mechanisms are currently represented in the M3D-DAMM model.

Hence, high-temporal-resolution field-scale measurements of CH₄ emissions and soil and porewater chemistry facilitated evaluation of the combined effects of soil redox conditions (moisture and O₂ concentrations) and associated pH feedbacks on underlying processes occurring across soil microsites while accounting for variation along the catena as a result of changing climatic drivers over time. The M3D-DAMM model captured the Birch-type effect by quantifying the pulses in soil CH₄ emissions as a function of increases in soil moisture following a strong drought (Birch, 1958). Specifically, the model coupled with microsite diffusivity explained CH₄ emissions common to wet valley soils and rare in comparatively drier ridge and slope soils and predicted the net release of CH₄ emissions from all topographic positions following a strong drought.

4.2 Sensitivity analysis

The variance-based sensitivity analysis confirmed the importance of microbial functional groups and their complex interactions with the surrounding biophysical and chemical environments in controlling CH₄ production and oxidation. For example, the growth and death of aceticlastic methanogens and the relative efficiency of aceticlastic methanogenesis were the most sensitive parameters (Fig. 8), which is consistent with another modeling effort on CH₄ fluxes across the Arctic landscape (Wang et al., 2019). Although from completely different ecosystem types, Wang et al. (2019) and the present study confirmed the importance of simulating soil topographies and microbial mechanisms when evaluating the heterogeneities in CH₄ fluxes. Representations of both direct (methanogenic substrate) and indirect (soil pH feedback) effects of acetate may have contributed to higher GSI values for parameters representing aceticlastic methanogenesis, which is similar to a previous study (Xu et al., 2015). The sensitivity of CH₄ emissions to the parameters representing methanotrophy was secondary to that of those representing aceticlastic methanogenesis, which is consistent with the increase in methanotrophic biomass during the drought. Our predicted changes in microbial biomass might be unacceptably large for the entire soil microbial community, which may only double or perhaps quadruple in response to changes in conditions, but individuals can grow exponentially (Goberna et al. 2010; Pavlov and Ehrenberg 2013; Roussel et al. 2015; Buan 2018).

4.3 Other processes

We did not completely reproduce the net emissions of soil CH₄ during the 2015 post-drought period across the catena with the M3D-DAMM model. To capture the full potential of net emissions of CH₄ (white shading in Fig. 3) from sesquioxide-rich soils, future modeling efforts may need to explicitly include the dynamics of redox-sensitive elements such as Fe and associated pH feedback under contrasting redox conditions (Barcellos et al., 2018; Bhattacharyya et al., 2018; Hall and Silver, 2013, 2015, 2016; O'Connell et al., 2018; Parfitt et al., 1975; and Silver et al., 1999). Wetting events can lower soil redox potential and reduce electron acceptors like Fe(III) to Fe(II). This concomitant reduction of Fe may increase soil pH, especially in anaerobic microsites, which could further increase net emissions of soil CH₄ (Tang et al., 2016; Zheng et al., 2019). Accounting for these effects may allow model simulations to better match the highest observed net CH₄ emissions in the post-drought period (Fig. 3).

Additionally, the reduction of Fe(III) to Fe(II) has supported anaerobic CH₄ oxidation in other ecosystems (Ettwig et al., 2016). Within this context, a measurable amount of anaerobic oxidation of CH₄ has previously been reported at our study site (Blazewicz et al., 2012). Additionally, Fe-reducing microorganisms can utilize acetate as a substrate and thereby compete with methanogens and reduce net methane emissions (Teh et al., 2008). Given the gradient of Fe in our study site, it is likely that biogeochemical cycling of Fe and CH₄ is coupled (O'Connell et al., 2018), which should be accounted for in future modeling efforts. For example, a modeling study supported the importance of Fe in simulating CH₄ cycling in an Arctic soil (Tang et al., 2016). To that end, building a comprehensive framework that also includes Fe biogeochemistry will afford greater confidence in projected CH₄ emissions from wet tropical forests under future climatic conditions (Bonan, 2008; Pachauri et al., 2014; Xu et al., 2016).

5 Conclusions

High-frequency CH₄ emission measurements coupled with real-time soil chemical measurements identified spatial and temporal variations affecting CH₄ production and oxidation in wet tropical forest soils of Puerto Rico. Overall, contrasting patterns of soil moisture between ridge and valley soils played an instrumental role in governing net CH₄ emissions. For example, consistently greater soil moisture likely favored methanogenesis by lowering the availability of O₂ in valley soils compared to ridgetop soils, especially in microsites with high soil moisture and soil C content. However, soil porewater chemistry, particularly the concentrations of acetate and associated soil pH, influenced the pattern of net emissions of CH₄ across the catena (valley > slope > ridge) during wetting after the 2015 drought. Thus, our results provide compelling evidence of the importance of both hot spots and hot moments in generating and mediating CH₄ emissions in wet tropical forest soils. A microbial-functional-group-based model coupled with a diffusivity module and consideration of soil microsites adequately reproduced both the spatial and temporal dynamics of soil CH₄ emissions, although mechanisms involving Fe biogeochemistry were neglected.

This study suggests that representing the microbial mechanisms and the interactions of microbial functional groups with the soil biophysical and chemical environment across soil microsites is critical for modeling CH₄ production and consumption. To that end, explicit consideration of these underlying mechanisms improved predictions of CH₄ dynamics in response to regional climatic events and provided insight into differential dynamics of solute and gas diffusion, different microbial functions, and gross CH₄ production and oxidation as a function of topography. Hence, we contribute to the ongoing development and improvements of Earth system and process models to better simulate microbial roles in CH₄ cycling at regional and global scales. However, observational data concerning the activities of different soil microbial functional groups are still needed to confirm the mechanisms proposed here. Future studies should integrate geochemical and microbiological information relevant for oscillatory redox conditions in wet tropical forests, especially those related to the redox-sensitive elements to build a comprehensive framework for modeling tropical soil CH₄ emissions.