Sinan Xu^1,2,3

Zhilei Sun^2,3*

Wei Geng^2,3

Hong Cao^2,3

Xilin Zhang^2,3

Bin Zhai^2,3

Zijun Wu¹*

¹State Key Laboratory of Marine Geology, School of Ocean and Earth Science, Tongji University, Shanghai, China
²Laboratory for Marine Geology, Pilot National Laboratory for Marine Science and Technology, Qingdao, China
³The Key Laboratory of Gas Hydrate, Ministry of Natural Resources, Qingdao Institute of Marine Geology, Qingdao, China

Understanding the modern marine methane processes, which can profoundly affect global climate and have far-reaching impacts on human living environments, is critical for research on the global carbon cycle. Thus, modeling of marine methane processes has attracted increasing attention due to models can accurately simulate and predict the environmental effects of methane on marine and atmospheric ecosystems. In this study, we review the applications of modeling works to marine methane processes, including methanogenesis in sediments, transport and reaction of methane in sediments and seawater, and marine methane emissions to the atmosphere. Compiled a large database of global methanogenesis rates and methane fluxes to the sulfate-methane transition zone, we estimate that the global methanogenesis budget in marine sediments is ∼0.87 Tmol yr⁻¹ and global sedimentary dissolved inorganic carbon produced by anaerobic oxidation of methane is ∼8.9 Tmol yr⁻¹. In addition, although anaerobic oxidation of methane in sediments and aerobic oxidation of methane in seawater act as primary filters to prevent methane leakage from sediments to the hydrosphere as well as the atmosphere, large masses of methane in extreme seafloor environments (e.g., mud volcanic eruptions and hydrate leakage) can still escape microbial oxidation and leakage to seawater or the atmosphere. There is still a lack of models that simulate methane in these extreme marine environments. Therefore, more modeling works are needed to assess the efficiency of marine ecosystems, including sediments and hydrosphere, in filtering methane in the event of large-scale methane leakage from the seafloor. This study provides an interdisciplinary view of methane processes in marine systems and helps identify future directions in the modeling of methane processes in marine system.

1 Introduction

Methane (CH₄), the most important greenhouse gas after carbon dioxide (CO₂), plays a key role in the Earth’s carbon cycle over geological time and ongoing global warming (Solomon et al., 2009; Etminan et al., 2016; Allen et al., 2018; Nisbet et al., 2019; Akam et al., 2020). Oceans constitute the largest carbon reservoir on Earth and cover about 70% of the world’s surface area. Approximately 1,000–5,000 Gt of methane is stored in marine sediments, mainly in the form of methane hydrates (Buffett and Archer, 2004; Burwicz et al., 2011). Even 1% of the methane seeping from sediment can increase the amount of methane in seawater and the atmosphere tenfold (Boetius and Wenzhöfer, 2013). Consequently, numerous modeling studies have been conducted to quantify methane processes in marine sediments and seawater. In this study, we review how models describe methane processes in marine systems, including sediments and seawater.

Compared with CO₂, methane has a relatively short mean lifetime (∼9 years) and is very close to a steady state (Dlugokencky et al., 2011). The reduction of methane with hydroxyl radicals (OH⁻) in the troposphere and stratosphere is a vital process of removing atmospheric methane (Dlugokencky et al., 2011). As a result, methane emissions are rapidly fed back into the global climate. The concentration and increased rate of methane in the atmosphere have reached their highest level since the preindustrial period (Saunois et al., 2016b). According to the mole fraction of methane in the geologic-atmospheric records, approximately 215 Tg year⁻¹ of methane was emitted to the atmosphere in the preindustrial era (Lelieveld et al., 2002). Methane emissions are mainly classified as natural (e.g., those formed in wetlands, termites and oceans) and anthropogenic (e.g., those from rice agriculture, biomass burning, domestic sewage, and coal extraction), with the natural source accounting for most methane emissions. In particular, wetlands (e.g., swamps and tundra) are the largest methane emitting regions (∼150 Tg year⁻¹) (Bousquet et al., 2006). Methane emissions from human activities are also not negligible. From the preindustrial era to the current period, anthropogenic methane emissions have increased significantly, reaching ∼503–884 Tg year⁻¹ (Bousquet et al., 2006; Dlugokencky et al., 2011; Saunois et al., 2016a). Coal extraction, crushing, and processing produce ∼120 Tg year⁻¹ of methane emissions (Bousquet et al., 2006). Rice agriculture and biomass burning also contribute ∼30 and ∼50 Tg year⁻¹, respectively (Bousquet et al., 2006). Moreover, human activities can affect natural emissions. For example, human exploitation of marine hydrates might destabilize the methane hydrates in sediments, resulting in substantial methane emissions.

Owing to global warming concerns and the necessity to extract the vast amounts of methane hydrate resources, marine methane processes are drawing growing interest. The scale of leakage from methane hydrate breakdown in marine sediments is large (∼20 Tg year⁻¹) due to changes in numerous geological and environmental conditions (Dlugokencky et al., 2011; Ruppel and Kessler, 2017). The transport of methane-rich fluids released by methane hydrate destabilization involves numerous biochemical reactions and complex environmental geological changes that profoundly affect the entire marine ecosystem (Regnier et al., 2011). Thus, the research on marine methane processes is significant for the following reasons: 1) Study of the mechanisms of global warming. Although the total amount of methane released into the atmosphere from the oceans is limited compared to emissions from wetlands and other factors (Dlugokencky et al., 2011; Ruppel and Kessler, 2017), marine sediments are the largest methane reservoir on Earth, most of which is in the form of methane hydrates (Burwicz et al., 2011). Therefore, the stability of methane hydrate reservoirs and their impact on global warming are particularly interesting. 2) The demand for exploration and development of methane hydrate resources. In marine system, the migration of methane fluids is closely related to the accumulation of mineral resources, and is directly involved in the formation of minerals (Schulz and Zabel, 2006). For example, the cold seeping fluids on the seafloor are mainly composed of methane, and cold vents are paramount exploration proxies for methane hydrate resources, especially for shallow-surface hydrate resources. Therefore, the search for cold vents formed by methane hydrate decomposition leakage fluids should be identified for exploring shallow-surface hydrate resources. Further, environmental considerations for the later development of hydrate resources include observation and study of leakage mechanism (Schulz and Zabel, 2006). 3) Research on critical scientific questions in the life and Earth sciences. The transport and reaction of fluids from methane hydrate decomposition involving essential biochemical reactions are a key link in the study of carbon cycle (Ruppel and Kessler, 2017). Moreover, a series of chemical interactions based on methane as the energy base is a prerequisite of life on many celestial bodies, including planets and moons, and a relevant direction for astrobiological research (Miller and Smythe, 1970; Jakosky et al., 1995; Pavlov et al., 2003). Therefore, marine methane processes should be observed and studied to better understand critical scientific questions about the oceanic evolution in geological history and the origin and evolution of life, including those on other planets. In addition, marine methane processes are interesting research issues at the intersection of Earth system sciences and life sciences.

Due to differences in the Gibbs free energy of electron acceptors, the methanogenesis process related to organic matter (OM) degradation occurs in deep sediments with anaerobic environments. The produced methane is transported to upper sediment by molecular diffusion and advection (Schulz and Zabel, 2006; Dale et al., 2008c). Although most seawater is oxygenated (Repeta et al., 2016), the methanogenesis process can occur in the microbial digestive tract and release fresh fecal particles into the ocean. In seawater, methane is mainly produced in the mixed layer, where methane concentration can reach ∼5 nM (Bastviken et al., 2004) and the maximum methane concentration leads to the “marine methane paradox” (Lenhart et al., 2015). Owing to the difference in methane chemical potential between the atmosphere and seawater (the net particle flux always flows from the higher chemical energy state system to the lower chemical energy state system), methane in supersaturated surface seawater can leak into the atmosphere (Holmes et al., 2000; Mrazovac et al., 2012). Modeling works are essential for studying the marine methane processes. 1) Modeling can overcome the spatial and temporal limitations of laboratory conditions. 2) Modeling can quantitatively analyze the variation of methane fluxes during each process. 3) Modeling can reasonably verify the accuracy of experimental data (Boudreau, 1997). In addition, the combination of model and measured data provides a more accurate quantification of methane production and consumption in each marine process. Considerable work has been done by scholars to simulate the abovementioned marine methane processes. Reaction-transport model (RTM) is commonly used to simulate the transport of methane in sediments coupled with methanogenesis and anaerobic oxidation of methane (AOM). The transport of methane-rich fluids in seawater involves the effects of ocean currents, water pressure, and temperature, so more complex mathematical models have been constructed to simulate it (Yamazaki et al., 2006).

In this review, we summarize the modeling research on marine methane processes, including methanogenesis and transport and reaction of methane in sediments, the flux of methane leakage at the sediment-seawater interface (SWI), methane processes in seawater, and methane flux from seawater to the atmosphere.

2 Methane Processes in Sediments

The main processes involved with methane in sediments are methanogenesis, methane phase transition, and methane transport-reaction. In this section, we review the application of models to these processes and the key factors influencing these methane processes.

2.1 Methanogenesis in Deep Sediment

Deep anaerobic marine sediments are the world’s largest reservoir of methane (∼4.55 × 10⁵ Tg C), with 90% of the produced methane stored in the continental margins, which are approximately four to eight times larger than the terrestrial surface biosphere and soils (Reeburgh, 2007; Wadham et al., 2012; Chen et al., 2017). Methanogenesis is the major source of methane in sediments and occurs when porewater sulfate (SO₄²⁻) is depleted (Jørgensen and Kasten, 2006). There are two main methanogenesis processes (Megonigal et al., 2004; Fenchel et al., 2012; Burdige and Komada, 2015) the autotrophic pathway (using carbon dioxide reduction), and the acetoclastic pathway (using acetate). The latter process usually occurs in fresh inland rivers (Blair, 1998; Conrad, 2005). According to carbon isotopic data combined with model simulations, the main methanogenesis process in marine sediments is carbon dioxide reduction (autotrophic pathway) (Burdige et al., 2016), which is expressed as follows:

2 C H_{2} O \overset{+ 2 H_{2} O}{\to} 2 C O_{2} + 4 H_{2} \overset{- 2 H_{2} O}{\to} C O_{2} + C H_{4} (1)

The reduction of carbon dioxide can be divided into two steps. First, OM is sequentially broken down into smaller molecules, followed by fermentation, which produces hydrogen gas (H₂). Then, the produced H₂ acts as a driving force for methanogenesis via carbon dioxide reduction (Whiticar, 1999). Modeling has been applied to methanogenesis that relies mainly on OM degradation, which is directly or indirectly involved in almost all biochemical reactions in marine sediments (Arndt et al., 2013). Hence, there are two important factors in methanogenesis modeling: 1) OM degradation rate and 2) switching of OM degradation from sulfate reduction to methanogenesis. Three types of models are commonly used to describe OM degradation (Arndt et al., 2013). 1) Discrete model (G-model). This type of model divides OM into a finite number of reactivity fractions, each of which is degraded according to its single degradation constants (Jørgensen, 1978). 2) Reactive continuum model (RCM). This type of model describes OM degradation by assuming a continuous distribution of OM reactivity at the SWI. A commonly used RCM is based on the Gamma distribution function (γ-RCM) (Boudreau and Ruddick, 1991). 3) Power model. This model is an empirical model based on a large dataset of global sedimentary OM data (Middelburg, 1989). According to the OM degradation model, the OM degradation rate can be expressed as follows:

R_{O M} (t) = \frac{d G_{O M} (t)}{d t} (2)

where R_OM denotes the OM degradation rate, t denotes time, and G_OM(t) denotes OM content at time t. The sulfate concentration in porewater is an important indicator of the methanogenesis, as it occurs at a lower sulfate concentration (Burdige et al., 2016). A threshold value for a given sulfate concentration (∼0.5 mM) is typically used to describe the change in OM degradation switched from sulfate reduction to methanogenesis (Boudreau, 1996; Arndt et al., 2009; Chuang et al., 2019; Dale et al., 2019). The methanogenesis rate (R_ME) is expressed as follows:

R_{M e} = \frac{1}{2} \cdot f (S O_{4}^{2 -}) \cdot R_{O M} (3)

where 1/2 denotes the ratio between OM consumed and methane produced (Eq. 1), and f(SO₄^2–) is a function representing to the threshold of sulfate concentration, which is expressed as (Boudreau, 1996; Chuang et al., 2019):

f (S O_{4}^{2 -}) = 1 - \frac{[S O_{4}^{2 -}]}{[S O_{4}^{2 -}] + K_{S}} (4)

or by the complementary error function (Martens et al., 1998; Dale et al., 2019):

f (S O_{4}^{2 -}) = 0.5 \cdot erfc (\frac{[S O_{4}^{2 -}] - K_{S}}{k_{i n}}) (5)

where [SO₄^2–] denotes the sulfate concentration, Ks represents the threshold of the sulfate concentration, and k_in in Eq. 5 is a parameter controlling the steepness of f(SO₄^2–). Typically, sulfate concentration decreases with depth. Hence, Eq. 4 is a monotonically increasing function with depth. However, the OM degradation rate Eq. 2 is a monotonically decreasing function with depth. As a result, substituting Eq. 4 into Eq. 3 cannot strictly guarantee that the methanogenesis rate Eq. 3 decreases monotonically with depth, especially in the surface sediments where the change in OM degradation rate is extremely high (Middelburg, 1989). Compared with Eqs 4, 5 can solve this problem, as the complementary error function is equal to either 1 or 0 (f(SO₄^2–) = 1, if [SO₄^2–] < Ks or f(SO₄^2–) = 0, if [SO₄^2–] > Ks).

From the expression of the methanogenesis rate (Eqs 2, 3), the intensity of methanogenesis in sediments is directly related to the OM degradation rate. Hence, factors affecting OM degradation also influence methanogenesis processes. The content of OM and its reactivity are the main factors that affect OM degradation (Arndt et al., 2013). The higher the OM content at the SWI, the higher is the OM content that can reach in the deep sediments and the higher the methanogenesis rate (Burwicz et al., 2011). Interestingly, the higher the OM reactivity, the lower the methanogenesis rate for the same OM content at the SWI. This is because under higher OM reactivity, more active OM is consumed in the upper sediment, and less OM is transferred to the deep sediment for methanogenesis (Meister et al., 2013). The sedimentation rate is also a paramount factor affecting OM degradation and further methanogenesis (Buffett and Archer, 2004; Burwicz et al., 2011). The sedimentation rate is a critical proxy of depositional environments, as its quantity in shelf regions (∼0.05 cm year⁻¹) is usually greater than that in deep sea regions (<0.001 cm year⁻¹) (Tromp et al., 1995; Burwicz et al., 2011). In addition, a higher sedimentation rate is generally coupled with higher OM fluxes and faster burial processes, promoting methanogenesis (Seiter et al., 2004; Meister et al., 2013).

The above shows the main factors affecting methanogenesis in sediments. We collected methanogenesis data from 48 sites around the global ocean and calculated the total methanogenesis rate in sediments by depth integration. We found that the sedimentation rate (w) and depth of sulfate-methane transition zone (SMTZ) have a well linear regression result with the depth-integrated methanogenesis rate, where methanogenesis is positively correlated with the sedimentation rate (R² = 0.71) and negatively correlated with SMTZ depth (R² = 0.79) (Table 1 and Figure 1). The positive relationship between depth-integrated methanogenesis rate and OM reactivity is weak (R² = 0.53). The distribution of OM reactivity in global marine sediments is complex, and no model can constrain global OM reactivity well (Arndt et al., 2013). For example, shelf regions are generally considered to have high OM reactivity due to high OM flux from inland and marine biosphere, but there is also a large amount of terrestrial refractory OM (LaRowe et al., 2020). Given the deeper water depth, OM fluxes and OM reactivity are generally lower in abyssal regions. However, high OM reactivity was found in the eastern costal of equatorial Pacific regions. The main reasons are warmer water temperature, which enhances net primary production in surface seawater, and lower water column oxygen concentration, which inhibits OM degradation in seawater and promotes more active OM that can reach sediments (Arndt et al., 2013). Based on the empirical formula we obtained (Figure 1) and the mean depth of SMTZ in different regions (Egger et al., 2018), we estimated that approximately 0.87 Tmol year⁻¹ methane is produced in global deep marine sediments and 78% in shelf regions (water depth <200 m) (Table 3).

TABLE 1

TABLE 1. Summary of methanogenesis model studies. The table reports data on sedimentation rate (w), the depth of SMTZ, OM content at the SWI (G_OM), OM apparent reactivity (<k>) and depth-integrated methanogenesis rate (R_ME) at each site.

FIGURE 1

FIGURE 1. Global empirical relationships. (A). Log–log plots of sedimentation rate and depth-integrated methanogenesis rate. (B). Log–log plots of SMTZ depth and depth-integrated methanogenesis rate. (C). Log–log plots of OM apparent reactivity and depth-integrated methanogenesis rate.

2.2 Methane Phase Transition in Sediments

Methane produced in bottom sediments has three phase states: gaseous, liquid (dissolved in porewater), and solid (methane hydrates) (Jørgensen and Kasten, 2006). The phase of methane in sediments varies with environmental factors, typically the ambient temperature, pressure, and salinity (Buffett and Archer, 2004; Regnier et al., 2011). The saturation concentration of methane (C_{Sa_Me}) in porewater has been extensively studied and can be expressed by a polynomial in temperature, pressure and salinity (Eq. 6) (Duan and Weare, 1992).

C_{S a_M e} = 1.437 \times 10^{- 7} S T P - 4.412 \times 10^{- 5} T P - 4.6842 \times 10^{- 5} S P + 4.129 \times 10^{- 9} S T + 1.43465 \times 10^{- 2} P - 1.6027 \times 10^{- 6} T - 1.2676 \times 10^{- 6} S - 4.9581 \times 10^{- 4} (6)

where S, T and p denote ambient salinity (−), temperaure (K), and pressure (atm), respectively. This empirical formula is applied to S of 1–35, T of 273.15–290.15 K, and p of 1–30 atm. If the concentration of dissolved methane in porewater is greater than the methane saturation and the ambient pressure and temperature conforms to the conditions for hydrate formation, methane hydrate will form in the sediment. The zone in the sediment where hydrate stability can be preserved is called the gas hydrate stability zone (GHSZ) (Kvenvolden, 1993). The GHSZ depth is mainly determined by the ambient temperature and pressure. The methane phase transitions from liquid + hydrate to liquid + gas below the GHSZ as geothermal temperature increases with depth.

The saturation of free methane in porewater is a key parameter for assessing the state of methane in sediments (Burwicz et al., 2011; Tishchenko et al., 2005). The transition of the methane phase in the GHSZ involves four main processes: methane hydrate formation (gaseous → solid, Eq. 7), methane hydrate dissolution (solid → gaseous, Eq. 8), free methane gas formation (dissolved → gaseous, Eq. 9), and free methane gas dissolution (gaseous → dissolved, Eq. 10). The simplest formulation for the above transition is based on a linear dependence (Wallmann et al., 2006a; Burwicz et al., 2011; Regnier et al., 2011), expressed as follows:

R_{G H} = k_{G H} \cdot (\frac{[C H_{4}]}{C_{S}^{d i s s}} - 1), i f [C H_{4}] \geq C_{S}^{d i s s} (7)

R_{D G H} = k_{D G H} \cdot (1 - \frac{[C H_{4}]}{C_{S}^{d i s s}}), i f [C H_{4}] < C_{S}^{d i s s} (8)

R_{F G} = k_{F G} \cdot (\frac{[C H_{4}]}{C_{S}^{f r e e}} - 1), i f [C H_{4}] \geq C_{S}^{f r e e} (9)

R_{D F G} = k_{D F G} \cdot (1 - \frac{[C H_{4}]}{C_{S}^{f r e e}}), i f [C H_{4}] < C_{S}^{f r e e} (10)

where [CH₄] denotes the methane concentration, C_S^diss denotes the solubility of dissolved methane, C_S^free denotes the solubility of gaseous methane, and k_GH, k_DGH, k_FG, and k_DFG denote the kinetic constants of gas hydrate formation, free methane gas formation, gas hydrate dissolution, and free methane gas dissolution, respectively.

The thickness of the GHSZ is an essential parameter for assessing hydrate resources (Burwicz et al., 2011). In addition, a sufficient methane source (e.g., deeper hydrate reservoirs and methanogenesis (Burdige et al., 2016)), low temperature, and high pressure are keys to methane hydrate reservoir formation. The minimum content of surface OM content for methane hydrate formation is ∼1 wt% (Buffett and Archer, 2004). It is difficult to form methane hydrates in regions within water depth of 600 m because of the high ambient temperature in sediments (Buffett and Archer, 2004), and methane hydrates are commonly found within water depth of 1,000–3,000 m with a bottom water temperature of ∼2°C (Kvenvolden, 1993). Methane hydrate formation regions can be divided into passive and active margins (Dale et al., 2008c). Active regions have sufficient overlying methane-rich fluids (methane from deeper hydrate reservoirs), and the abundance of methane hydrates within the GHSZ reaches 30%–50%. By contrast, diffusion dominates the transport of methane in passive margins and hydrate abundances in these regions are smaller (Dale et al., 2008c). The amount of methane stored in marine sediments is within 500–57,000 Gt C (Dickens, 2001; Buffett and Archer, 2004; Milkov, 2004; Klauda and Sandler, 2005; Archer et al., 2009; Burwicz et al., 2011; Pinero et al., 2013; Kretschmer et al., 2015). Buffett and Archer. (2004) first estimated methane hydrate reservoirs in global ocean sediments using the 2-G OM degradation model coupled with sulfate reduction, methanogenesis and AOM. They found that the rain rate of particulate OM is key in the global inventory of methane hydrate reservoirs (Archer et al., 2002). In addition, they estimated a total methane hydrate inventory of 3,000 Gt C, but it decreased to 600 Gt C when the methane in the overlying fluid was disregarded (Buffett and Archer, 2004). Klauda and Sandler, (2005) slightly modified Buffett’s model and obtained a global methane hydrate inventory of ∼57,000 Gt C (Klauda and Sandler, 2005), which was almost two orders of magnitude higher than that of Buffett and Archer. (2004). The main reason for the large error was that the OM degradation model was a 1-G model that the entire particulate organic carbon (POC) pool was assumed to be one group with a single degradation rate constant (4.7 × 10^–7 year⁻¹). Moreover, biogeochemical reactions associated with the consumption of methane in porewater were ignored (e.g., sulfate reduction and AOM). OM reactivity decreased with sediment depth, whereas this feature cannot be reflected by the G model (Middelburg, 1989; Arndt et al., 2013). Burwicz et al. (2011) estimated the global methane hydrate inventory in the range of 4.18–995 Gt C, where the first value refers to present-day conditions (using the relatively low Holocene sedimentation rate) and the second value corresponds to a higher Quaternary sedimentation rate. An RCM (γ-RCM) was used to describe OM degradation in sediments coupled with sulfate reduction and AOM. In addition, an improved power model was used to simulate the OM degradation rate (R_OM) (Eq. 11) (Wallmann et al., 2006a; Wallmann et al., 2012; Kretschmer et al., 2015).

R_{O M} (t) = \frac{K_{C}}{[D I C] + [C H_{4}] + K_{C}} \cdot k_{x} \cdot G_{O M} (t) (11)

where [DIC] and [CH₄] denote the concentration of ambit dissolved inorganic carbon (DIC) and methane, respectively, K_c denotes the Monod inhibition constant (mM), k_x denotes the age-dependent kinetic constant (year⁻¹), and G_OM(t) denotes the OM content. These simulations found that the main factors affecting hydrate storage in the GHSZ were the ambient temperature, oxygen content, sedimentation rate, and upward fluid especially in the active margins (Buffett and Archer, 2004; Burwicz et al., 2011). In particular, ambient temperature is a vital factor influencing the stability of methane hydrate reservoirs, as global methane storage may be reduced by 85% if sedimentary temperature increases by even 3°C (Buffett and Archer, 2004).

According to the global exploration of methane hydrate reservoirs, many sites with hydrates have been discovered globally. Figure 2 shows the main regions in the world where methane hydrates have been found, such as the Gulf of Mexico (Davidson et al., 1986), Caribbean Sea (Reed et al., 1990), the eastern margin of South America (Jahren et al., 2001), western margin of Africa (Ben-Avraham et al., 2002), Bering Sea (Scholl and Hart, 1993), Sea of Okhotsk (Shoji et al., 2005), Okinawa Trough (Sakai et al., 1990), Sea of Japan (Yun et al., 2011), Shikoku Trough (Saito and Suzuki, 2007), South China Sea (Li et al., 2018), eastern Pacific Ocean (Inagaki et al., 2006), California Coast (Dickens and Quinby-Hunt, 1994), Peruvian margin (Kvenvolden and Kastner, 1990), the Gulf of Oman (White, 1979), Ross Sea and Weddell Sea in the Antarctic (Stoll and Bains, 2003; Giustiniani et al., 2018), and Barents Sea (Andreassen et al., 1990). In addition, the methanogenesis model and methane phase transition models were employed to predict the global distribution of methane hydrate reservoirs (Burwicz et al., 2011; Kretschmer et al., 2015). Both exploration and modeling works reveal that methane hydrate reservoirs occur in shallow marine sediments, mainly distributed along the coastal regions of continents (Figure 2).

FIGURE 2

FIGURE 2. Global distribution of inferred sites, drilling/sampling sites and real production tests for methane hydrate. The underlying map shows global distribution of methane hydrates in marine sediments estimated by model (Kretschmer et al., 2015).

2.3 Transport-Reaction of Methane in Marine Sediments

The methane transported in marine sediments can be divided into dissolved methane (in porewater) and bubbling methane (Regnier et al., 2011). The latter is generated when methane hydrates in the deeper layers break down or when methanogenesis rates are extremely high (Martens and Klump, 1980). Methane transport in the sediment is described by RTM (Meister et al., 2013; Berner, 2020), as follows:

φ \cdot \frac{\partial [C H_{4}] (x, t)}{\partial t} = \frac{\partial (φ \cdot \frac{D_{m}}{τ^{2}} \cdot \frac{\partial [C H_{4}] (x, t)}{\partial x})}{\partial x} - \frac{\partial (φ \cdot w \cdot [C H_{4}] (x, t))}{\partial x} + v \cdot λ \cdot \frac{\partial ([C H_{4}] (x, t) - C_{S a_M e})}{\partial x} + φ \cdot \sum R (x, t) (12)

where x denotes the depth below the seafloor, t denotes time, φ denotes porosity, D_m denotes the molecular diffusion coefficient of methane, τ denotes tortuosity calculated as τ² = 1 − ln(φ²) (Boudreau, 1997), w denotes the burial velocity, [CH₄] denotes the methane concentration, and R denotes methane related biogeochemical reactions. The methane-rich upward fluid is expressed by the third term on the right side of Eq. 11, v denotes the upward fluid velocity, and C_{Sa_Me} denotes the saturation concentration of methane (Eq. 6). The parameter λ = 1, if [CH₄] > C_{Sa_Me} and λ = 0, if [CH₄] < C_{Sa_Me}. Methane transport distinctly differs under passive and active conditions. Under passive conditions, the diffusion process dominates the transport of methane on a long-time scale, and methane cannot reach the SWI. By contrast, the advection process under active conditions rapidly pushes the methane from the deeper to the upper sediments within a short time scale (∼12 years), and methane can reach the seawater (Dale et al., 2008c).

The most common reaction involved in methane transport is the AOM (Eq. 13), which consumes 90% of the methane in the overlying fluids and is thus the biggest methane sink (Reeburgh, 2007).

C H_{4} + S O_{4}^{2 -} \to H C O_{3}^{-} + H S^{-} + H_{2} O (13)

From Eq. 13, the AOM rate is related to the concentration of methane and sulfate. The simplest model describing the AOM rate is follows (Regnier et al., 2011; Van Cappellen and Gaillard, 1996):

R_{A O M} = k_{A O M} \cdot [C H_{4}] \cdot [S O_{4}^{2 -}] (14)

where k_AOM denotes the first-order rate constant for the AOM, [CH₄] denotes the methane concentrations, and [SO₄^2-] denotes the sulfate concentrations. This model can be easily used to fit measured data and describe the distributions of AOM rates as well as the location and thicknesses of the SMTZs in sediments. In addition, various factors can affect the AOM rate that cannot be reflected by Eq. 14, such as temperature, substrate concentration, enzyme reactivity, and bioenergy (Regnier et al., 2005; Dale et al., 2006; Regnier et al., 2011). The concentrations of microorganisms and enzymes involved in the AOM are related to substrate concentration. Specifically, they are limited at low substrate concentrations. In light of this, the AOM rate can be expressed as follows (Regnier et al., 2005):

R_{A O M} = v_{m a x} \cdot \frac{[C H_{4}]}{[C H_{4}] + K_{m}} \cdot \frac{[S O_{4}^{2 -}]}{[S O_{4}^{2 -}] + K_{s}} (15)

where v_max denotes the maximum AOM rate, and K_m and K_s represent the half-saturation constants of methane and sulfate, respectively. Furthermore, the minimum energy supply that can facilitate the AOM is ∼11 kJ mol⁻¹ (Dale et al., 2006), whereas bioenergy available within the SMTZ is limited. Considering the limitation of bioenergy, the AOM rate can be expressed as follows:

R_{A O M} = v_{m a x} \cdot F_{K} \cdot F_{T} (16)

where F_K and F_T denote the kinetic and thermodynamic driving forces for AOM, respectively. The F_K is the abbreviation of the latter two terms of Eq. 15, and F_T is expressed as follows:

F_{T} = 1 - exp (\frac{Δ G_{r} + Δ G_{B Q}}{χ \cdot R \cdot T}) (17)

where ΔG_r denotes the Gibbs energy of the AOM reaction, ΔG_BQ denotes the minimum energy required to sustain the AOM, χ denotes the number of protons translocated across the cell membrane, R denotes the gas constant, and T denotes temperature.

When the methane hydrate in sediments is destabilized or the rate of methanogenesis is high, methane can exist as bubbles in supersaturated porewater (Martens and Klump, 1980). Bubbling methane leakages are widely observed and investigated in coastal sediments (Chanton et al., 1989; Anderson et al., 1998; Veloso‐Alarcón et al., 2019). The rate of bubble growth is described by the first-order model (Davie and Buffett, 2001):

Φ = R_{b} \cdot ([C H_{4}] - C_{S a_M e}) (18)

where R_b denotes the rate constant, [CH₄] denotes the methane concentration, and C_{Sa_Me} denotes the saturation concentration of methane in porewater. A continuous three-phase system (solid: hydrate, liquids: porewater, and gas: bubble) has been applied to describe gaseous methane transport in sediments. This approach has been widely applied to soils, aquifers, petroleum and shale gas extraction based on the Darcian flow theory (Schowalter, 1979; Molins and Mayer, 2007; Reagan and Moridis, 2008; Molins et al., 2010). Unlike that of dissolved methane in the porewater, the process of bubble ascent is also influenced by eddy diffusion, which can be described as follows (Haeckel et al., 2007):

K_{e d d y} \approx 0.928 \cdot \sqrt{g \cdot r_{b u b b l e}^{3}} (19)

where K_eddy denotes the first-order eddy diffusivity constant, g denotes the acceleration due to gravity, and r_bubble denotes the bubble radius. The velocity of eddy diffusion (K_eddy > 1 × 10⁵ cm2 year⁻¹) is several orders of magnitude greater than that of molecular diffusion, resulting in a sevenfold increase in the maximum AOM rate in the sediments while methane is released into the overlying water column as bubbles (Haeckel et al., 2007). However, the continuous three-phase system assumes that the fluid is transported through a static solid matrix. This limits the application of the system to marine sediments, which are affected by compaction, burial, and particle mixing (Regnier et al., 2011).

Numerous studies have been conducted to simulate AOM reaction in global marine sediments based on Eqs 12–19 (Regnier et al., 2011). AOM reaction mainly occurs in SMTZ, which is commonly defined as a segment close to the maximum AOM rate. The SMTZ depth varies from 2 cm in the Black Sea (Wallmann et al., 2006b) to 400 m in the equatorial Atlantic (Arndt et al., 2006). Three factors impact the SMTZ depth: 1) methane flux into the SMTZ, 2) upward fluid velocity, and 3) OM content and reactivity (Regnier et al., 2011; Meister et al., 2013; Chuang et al., 2019). Higher methane flux and upward fluid velocity can generate a shallow and narrow SMTZ (Regnier et al., 2011). For example, the SMTZ depth in active margins (∼10 m) is much shallower than that in passive margins (∼100 m) (Jørgensen and Kasten, 2006; Dale et al., 2008c). Moreover, the higher the content and reactivity of OM in surface sediment, the shallower is the SMTZ, attributable to the consumption of sulfate by OM degradation, which promotes the upward transport of methane-bearing fluids (Meister et al., 2013). The SMTZ is also related to the sedimentation rate. The higher the sedimentation rate, the shallower is the SMTZ (Meister et al., 2013).

The SMTZ acts as the main methane barrier where most electron acceptors are depleted, and OM degradation initiates to dominate methanogenesis (Egger et al., 2018). Approximately ∼191 Tmol C year⁻¹ of OM is transported to the global seafloor, where 3%–4% of OM is converted into methane (5.7–34.4 Tmol year⁻¹) in continental margin sediments (Muller-Karger et al., 2005; Dunne et al., 2007). Egger et al. (2018) compiled methane and sulfate data from 740 sites of global marine sediments and suggested ∼3.8 Tmol C annual flux of methane into the SMTZ. According to this methane budget, approximately 3.8 Tmol year⁻¹ DIC is produced by the AOM because 1 mol of methane consumed produces 1 mol of DIC (Akam et al., 2020). However, there are some errors in this budget that calculated only by stoichiometric relationship of the AOM and methane flux data in Egger et al. (2018). The main reason is that the methane flux calculated by Egger et al. (2018) is based on the gradient of methane at locations where porewater methane and sulfate concentrations are equal. However, the flux calculated by Egger et al. (2018) cannot indicate the methane flux involved in the AOM throughout the sediment, and we found that the ratio between this methane flux and depth integrated AOM rate is not 1:1(Table 2). Mathematically, the AOM exists throughout sediments, but the concentrations of methane/sulfate in the surface/bottom sediments are too low to detect. Therefore, it is difficult to quantify the total methane flux in the sediments involved in the AOM. Compiled methane flux into the SMTZ data calculated by Egger et al. (2018) and the depth-integrated AOM rate data in the global marine sediments, we found a good regression in the log-log coordinate between the two (Figure 3). According to the methane flux distribution in different marine regions and the empirical formula (Figure 3), we estimated that ∼8.9 Tmol of AOM-related DIC is produced annually, with ∼82% being in the continental shelf sediments (water depth <200 m) (Table 3).

TABLE 2

TABLE 2. Summary of AOM rate studies. The table reports data on the depth of SMTZ, methane flux into SMTZ (J_CH) and depth integrated AOM rate (R_AOM) at each site. The check marks indicate the presence of methane hydrates in the bottom sediments.

FIGURE 3

FIGURE 3. Log–log plots of methane flux into the SMTZ and integrated AOM rate. The gray shaded area indicates the error range of the linear regression, where the upper bound is log(R_AOM) = 0.91 × log(J_CH4) + 0.95, and the lower bound is log(R_AOM) = 0.91 × log(J_CH4) −0.40.

TABLE 3

TABLE 3. Global budgets of methanogenesis rate and DIC flux related to AOM.

3 Methane Transport and Reaction in Seawater

Although sediment AOM consumes nearly all upward methane, particularly in passive margins, methane leakages are commonly observed in active margins (Dale et al., 2008d), such as the northern Barents Sea (Andreassen et al., 2017), East Siberian Arctic Shelf (Shakhova et al., 2010), and Gulf of Mexico (MacDonald et al., 1994). The main pathway of methane leakage from sediments to seawater in active margins (e.g., cold seep and pockmarks) is methane plume (Schulz and Zabel, 2006). The formation of authigenic carbonates promoted by AOM processes is an essential factor in the structure of cold seeps and pockmarks (Bohrmann et al., 1998; Luff and Wallmann, 2003; Aloisi et al., 2004; Nöthen and Kasten, 2011).

Cold seeps are often called the “windows to the deep geosphere” (Boetius and Wenzhöfer, 2013), and they form a central link of material and energy exchange between the lithosphere and hydrosphere together with hydrothermal vents systems. Cold seeps formed by hydrate decomposition leakages have globally been found in the entire bathymetric range of the continental slopes from high to low latitudes. For example, more than 1,200 leakages have been detected in the Hornsund Fault Zone of the Barents Sea (Waage et al., 2019), more than 600 leakages have been discovered in the polar North Atlantic (Bear Island Trough, northern Bjørnøyrenna) (Andreassen et al., 2017), approximately 5,000 leakages have been identified in the northern Gulf of Mexico at water depths of >200 m (Solomon et al., 2009), and 27,000 leakages have been detected in the shallow waters of eastern Siberia (Shakhova et al., 2014). Preliminary estimates suggest that there may be over hundreds of thousands of cold-spring leakage sites globally, erupting simultaneously due to hydrate leakage. Large-scale methane leakages on the seafloor are usually caused by the destabilization of hydrate reservoirs located at shallow water depths or buried at shallow depths and connected to fractures (Freire et al., 2011). The following main factors influence these large leaks. 1) Abnormal overpressure. When hydrocarbons accumulate in the pores of marine sediments and the pressure reaches a high enough level, these hydrocarbons will migrate upward through the GHSZ (Tréhu et al., 2004). The decomposition of methane hydrates at the edge of this stability zone can also form a considerable number of hydraulic fractures in its upper part, which then become an ideal channel for gas migration (Xu and Germanovich, 2006). 2) Regional geological environment fluctuation. These changes are observed in seismic activity (Fischer et al., 2013), tidal cycles (Boles et al., 2001), glacial melting (Andreassen et al., 2017), and bottom water temperature fluctuations (Ferré et al., 2020). 3) Erosion of geological bodies. The erosion of the submarine canyons and destabilization of the canyon sidewall sediments lead to the erosion of the hydrate-bearing caps or reservoirs, resulting in the destabilization and decomposition of the hydrate-bearing layers (Paull et al., 2005). 4) Dramatic changes in the global environment. Dramatic global sea level fluctuations and rapid climate warming during geological history could have triggered the catastrophic release of hydrates by disintegration. For example, the amount of methane hydrates that decomposed and escaped to the atmosphere during the Paleocene/Eocene thermal maximum (PETM) was estimated to be ∼2,100 Gt C (Dickens et al., 1997).

In the cold seeps, mud volcanoes, or pockmark regions, methane generally seeps into the hydrosphere in the form of gas flares, also called bubble plumes (Boetius and Wenzhöfer, 2013). Only a few in situ measurements approaches can be used to calculate methane flux from sediments to the hydrosphere, and they commonly use remotely operated vehicle (ROV). Combined with ROV observation, methane flux (F_methane) at the SWI can be estimated as follows (Blomberg et al., 2016; Lohrberg et al., 2020; Mau et al., 2020):

F_{m e t h a n e} = N \cdot n_{m e t h a n e_b u b b l e} \cdot M_{m e t h a n e} \cdot f_{m e t h a n e} (20)

where N denotes the number of seep locations, n_{methane_bubble} denotes the amount of methane per gas bubble, M_methane denotes the molar mass of methane, and f_methane denotes the gas bubble emission frequency, which can be obtained by ROV observation. Another method of estimating methane efflux is the box model that the box is defined by the hydrocasts model (Mau et al., 2020), as follows:

F_{m e t h a n e} = \frac{I \cdot u (z) / l_{p a t h}}{A} (21)

where I denotes the estimated methane inventory within the box, u(z) denotes the current velocity, l_path denotes the migration of the entire box away from its original position, and A denotes the surface area of the grid. In addition, according to porewater methane profiles, methane flux at the SWI can be calculated by Fick’s first law (Eq. 22) (Haese et al., 2003; Chen et al., 2017), as follows:

F_{m e t h a n e} = φ \cdot D_{m} \cdot \frac{d [C H_{4}]}{d z} (22)

where φ denotes porosity, D_m denotes the molecular diffusion coefficient of methane, and the last term represents the gradient of methane concentration at SWI.

We collected methane flux data in cold seep and pockmarks (Table 4). The largest methane flux (1,169–1,175 mmol m⁻² d⁻¹) is at the center of the giant pockmark REGAB located in the West African margin at a water depth of 3,160 m (Pop-Ristova et al., 2012; Boetius and Wenzhöfer, 2013). The regional methane flux can be obtained on the basis of gas flares observed on the seabed and the estimation of the methane flux at each flare. As shown in Table 4, 2.5–169.9 Mg year⁻¹ of methane flux is identified in the Gulf of Mexico regions (∼6,041.25 km²), of which more than 90% is occupied by the northern Gulf of Mexico (Weber et al., 2014; Römer et al., 2019; Wei et al., 2021).

TABLE 4

TABLE 4. Synthesis of methane flux into hydrosphere and their detection methods.

The transport of methane in seawater, commonly in the form of gas flares or plumes, is an extremely complex process, where bubble plumes are rapidly dissolved and subsequently diluted by mixing with overlying seawater and then dispersed by ocean currents (Graves et al., 2015). The Navier-Stokes equations are commonly used to describe hydrothermal plumes in the seawater where momentum, mass, heat, methane saturation, and microbially-mediated chemical reactions are considered (Yamazaki et al., 2006; Jiang and Breier, 2014). The radius of methane bubbles (R_bubble) typically ranges from 0.001 to 0.015 cm (Shakhova et al., 2014; Higgs et al., 2019). The bubbles less than 10 mm in diameter are expected to dissolve before they reach the surface mixing layer (Gentz et al., 2014). The rate of bubble dissolution depends mainly on the initial bubble size, water temperature, salinity, pressure, and bubble rise rate (Leifer and Patro, 2002; Rehder et al., 2009; Shagapov et al., 2017). The bubble dissolution rate can be written as follows (Fu et al., 2021):

\frac{d M}{d z} = \frac{- (4 π R_{b u b b l e}^{2}) \cdot K \cdot (C_{S} - C_{0})}{V_{b u b b l e}} (23)

where dM/dz denotes the change in a bubble’s methane content over the rise interval dz, V_bubble denotes the bubble rise velocity, and K denotes the mass transfer rate.

Based on large ROV observations of methane bubble radius (R_bubble) and bubble rising velocity (V_bubble), empirical formulae were used to describe bubble rise velocity, as follows (Clift et al., 2005; Leifer et al., 2006):

V_{b u b b l e} = 276 R_{b u b b l e} - 1648 R_{b u b b l e}^{2} + 4882 R_{b u b b l e}^{3} - 7429 R_{b u b b l e}^{4} + 5618 R_{b u b b l e}^{5} (24)

The well-known theory of Morton, Taylor and Turner (MTT model) has been used to estimate the maximum rising height of a single plume (Morton et al., 1956), and the scaling of the maximum plume rise height (Z_max) is given as follows:

Z_{m a x} = C_{e} \cdot {(\frac{B_{e x i t}}{N^{3}})}^{\frac{1}{4}} (25)

where B_exit denotes the source flux, N denotes the ambient frequency, and C_e denotes the scaling coefficient, which was estimated to be 3.76 by analyzing literature data from laboratory plume experiments (Briggs, 1969).

Compared with methane reaction in sediments with sulfate (AOM) in anoxic environments, bottom waters in the continental margin are generally completely oxygenated, with oxygen concentrations ranging from 250 to 350 μM (Boetius and Wenzhöfer, 2013). The consumption of methane in the hydrosphere is related to oxygen, as aerobic oxidation of methane (AeOM):

C H_{4} + O_{2} \to C O_{2} + H_{2} O (26)

The AeOM rate can be calculated as follows (Reeburgh et al., 1991; Valentine et al., 2010; Mau et al., 2020):

R_{A e O M} = k \cdot [C H_{4}] (27)

where k denotes the first-order rate constant, and [CH₄] denotes the methane concentration. Many factors that influence AeOM, such as (1) water depth; (2) temporal and spatial distributions of fluid release flux; (3) bubble characteristics (e.g., size, oil film and hydrate coating) (Veloso-Alarcón et al., 2019); (4) ocean euphotic layer; (5) dissolved oxygen concentration, temperature and salinity of seawater (Crespo-Medina et al., 2014); and (6) ocean currents (Steinle et al., 2015). According to the data compiled by Boetius and Wenzhöfer, (2013), the contribution of AeOM to methane consumption in seeps regions can exceed the contribution of AOM, particularly at seeps with low microbial abundance. In regions with low methane fluxes, the AOM consumes 90% of the methane in the overlying fluid (Reeburgh, 2007), whereas most of the methane transported to the seafloor in the seepage regions is consumed by the AeOM (Boetius and Wenzhöfer, 2013).

In cold seep regions, the efficiency of the methane filter (AOM and AeOM) decreases from ∼80% for low fluid flow systems (methane flux: ∼35 mmol m⁻² d⁻¹) to ∼20% for moderate flow systems (methane flux: ∼100 mmol m⁻² d⁻¹). With microorganisms and bacterial mats absent and with intense gas eruptions, the filter efficiency becomes lower than 10% in high fluid flow systems (e.g., giant pockmark systems) (Pop-Ristova et al., 2012; Boetius and Wenzhöfer, 2013). When such a large amount of methane escapes from microbial oxidation on the seafloor, it will be consumed aerobically in the seawater or transferred to the upper mixed layer and then to the atmosphere (Reeburgh, 2007; Boetius and Wenzhöfer, 2013).

4 Methane Flux to the Atmosphere

Owing to the human activities, the concentration of greenhouse gases in the atmosphere has been increasing since the 20th century, causing global warming (Sommer et al., 2009; Anderson et al., 2016). The greenhouse effect causes some ecological problems such as glacier retreat, sea-level rise, and the northward shift of the climate zone, which will cause great harm to the natural environment (Larcombe et al., 1995). As the main greenhouse gas, it is significant to calculate the amount of methane transported from seawater to the atmosphere.

The methane flux released from seawater to the atmosphere has been mainly predicted by digital simulation and remote sensing technology (Bovensmann et al., 2010). These methods usually use scientific ship positioning measurements and reverse atmospheric models to estimate methane release. For example, during the 2012 North Sea Elgin blowout accident, such a method was used to assess methane emissions to the atmosphere (Gerilowski et al., 2015). Considering this technology is relatively clumsy and inflexible, passive remote sensing techniques have been developed to collect atmospheric methane concentrations around the study areas, such as wireless remote sensing to collect atmospheric methane concentrations around the study region (Somov et al., 2013; van Kessel et al., 2018). However, when surveying with remote sensing instruments using short-wave infrared radiation, the weak reflectivity of water affects the information (Seelig et al., 2008).

The model calculates the methane flux at the sea-air interface mainly based on the differences in methane chemical potential between seawater and the atmosphere (Seelig et al., 2008). When methane is supersaturated in seawater, the methane in seawater can be emitted to the atmosphere driven by chemical potential. The methane fluxes to the atmosphere are calculated by the diffusive exchange equation (Solomon et al., 2009; Michel et al., 2021):

F l u x_{M e t h a n e} = k_{a v g} \cdot (C_{p l u m e} - C_{e q}) (28)

where k_avg denotes the gas transfer coefficient at the average wind speed, and C_eq denotes the seawater methane concentration in equilibrium with air under ambient conditions (Yamamoto et al., 1976). The gas transfer coefficients are calculated using the empirical formula (Wanninkhof, 1992), as follows:

k_{a v g} = 0.31 \cdot u_{a v g} \cdot {(\frac{S_{C}}{600})}^{- 0.5} (29)

where u_avg denotes the average wind speed at 10 m above the sea surface, and S_c denotes the Schmidt number (a function of salinity and temperature).

The combined cycle model for CO₂ and CH₄ (CMCDMC) was used to assess the role of different environmental parameters in the natural and anthropogenic components of methane leakage and climate change (Krapivin et al., 2017). Given that low-rate methane leaks are nearly completely consumed by the AOM, the climatic impact of methane may be overlooked. However, the high-rate leakage enables methane to enter the seawater and atmosphere directly, thereby causing global warming (Buffett and Archer, 2004). Although there is no conclusive evidence that hydrate-derived methane presently enters the atmosphere, more observational data and improved numerical models will help better describe the climate-hydrate synergy in the future (Ruppel and Kessler, 2017).

5 Challenges and Outlook

The application of the model to the marine methane processes is summarized as shown in Figure 4, including methanogenesis, methane transport and reaction in sediments and seawater, and methane flux from seawater to the atmosphere. The global ocean methane cycle has been increasingly studied, and considerable progress has been made in understanding the marine methane cycle. However, many environmental factors involved in methanogenesis, methane reaction and transport have not been described well by the models. Human understanding of the global methane cycle and the capability to assess the contribution of methane leakage to past and future global changes depend heavily on the accuracy of model construction and calculations. Future progress will rely heavily on additional observational data from different marine environments in the global ocean and on linking models to the observed complexities. Below, we list the remaining problems concerning methane cycle processes to help direct future research.

FIGURE 4

FIGURE 4. Schematic diagram of the application of model in the marine methane processes.

5.1 Metal-dependent AOM (Metal-AOM)

Sulfate–driven AOM is widely found in global marine sediments, but scholars have found that AOM driven by some active metals (e.g., manganese and iron) is also quite common (Beal et al., 2009). The reduction of the other metal ions with methane is often neglected because the concentration of sulfate in porewater is several orders of magnitude higher than those of other electron acceptors and the almost complete consumption of methane (Reeburgh, 2007). If the manganese flux (∼19 Tg year⁻¹) and iron flux (∼730 Tg year⁻¹) of the whole world were used to oxidize methane, the result might account for around a fourth of today’s AOM consumption. Even if only a small portion of manganese and iron fluxes is used for the AOM, the process can be a massive methane sink, as manganese and iron may be oxidized and reduced by 100–300 times before being buried (Canfield et al., 1993; Beal et al., 2009). Although several experimental studies have been conducted on metal-AOM (Sivan et al., 2011; Egger et al., 2015; Ettwig et al., 2016), there is no single study that assesses metal-AOM via a modelling approach. Thus, more modelling work about metal–AOM needs to be done to help quantify global metal-AOM and enhance the understanding of the changes of δ¹³C-DIC in porewater.

5.2 Model Construction, Boundary Conditions, and Marine Biological Environments

As methane plumes mainly originate from the decomposition of methane hydrates, transport to the seawater, and eventually reach the atmosphere, a continuum model should be established to describe methane transport in different media (sediment, seawater and atmosphere). At present, the boundary conditions, the initial conditions, and parameters of the mathematical models used for marine geological investigations are ideal (Boudreau, 1997). The transport of porewater species and the initial and boundary conditions are always set as a constant or closed boundary. Further, the pressure and temperature boundary conditions given for simulating the decomposition and release of methane from seabed natural gas hydrates are within a limited predicted range. For complex geological processes, the selection of boundary conditions should be closely combined with field monitoring or laboratory experimental data rather than a simple boundary prediction range or constant. The boundary conditions for methane leakage from sediments and methane migration through different media should be a set of complexes, dynamic boundary conditions that biochemical and geological factors should be considered.

Methane transport in the seabed and seawater is affected by marine and geological organisms (Keppler et al., 2009). Therefore, model accuracy and sensitivity are strongly related to the consideration of biological factors. At present, most mathematical models are lacking in accounting for the effects of marine organisms on methane transport, or the effects of these marine biological factors are treated as a constant (Chuang et al., 2019). By contrast, the activity of marine organisms is influenced by multiple factors. Thus, further work is needed to describe the impact of marine organisms on the methane cycle. The results of field investigation and experimental analyses of recent biochemical reactions of methane should be employed to characterize the models.

5.3 Methane Leakage on the Seafloor

Hydrate can be formed when the seabed temperature and pressure reach certain conditions and if there is a leak that adds alkane fluids to the seabed. Therefore, the stable boundary of hydrates may be at a certain depth in the water body, and the reservoir location may not necessarily be in the sediment. The Gulf of Mexico (Brooks et al., 1994; Boswell, 2009), Joetsu basin in Japan, and the South China Sea have exposed hydrate outcrops or shallow hydrates on the seafloor (Hiruta et al., 2009; Zhang et al., 2017). When methane bubbles ascend on the seafloor in some areas, their surfaces can rapidly form hydrate crusts during upward migration. Therefore, the bubble plumes at the leakage point of cold seeps can be classified as “clean” or “dirty” (Barnes and Goldberg, 1976). There is no numerical simulation for bubble leakage under these conditions. For some anoxic “euxinic basins,” such as Cariaco Basin and Black Sea (Reeburgh et al., 1991; Van Rensbergen et al., 2002), research on the processes of methane transport also lacks simulation. Therefore, there is a large room for improvement in accurately simulating modern seabed methane processes, especially methane leaks associated with gas hydrates.

5.4 Predictive Capability of Accidental Large-Scale Methane Leakage

The major advantage of the modeling approach is that it can reasonably predict the trend of methane transport under predetermined conditions (Boudreau, 1997). Owing to the AOM process, methane in sediments hardly seeps into the seawater, or it seeps into the seawater only for a short period (<60 years) and is then quickly re-limited in sediments by the AOM process (Dale et al., 2008c; Regnier et al., 2011). Many ROV observations reveal that strong gas ebullition in the cold seeps or mud volcano regions allows a massive amount of methane to escape microbial oxidation in sediments and be ultimately transported to seawater or even the atmosphere (Boetius and Wenzhöfer, 2013; Andreassen et al., 2017). Nevertheless, there has been no systematic analysis of modeling work combined with ROV observation data to simulate methane processes in cold seeps and mud volcano regions. Considering that diffusion-dominated methane processes in sediments are extremely slow and insignificantly impact on human life (Jørgensen and Kasten, 2006; Regnier et al., 2011), it would be more meaningful to simulate large-scale methane leakage in cold seeps and mud volcano regions. Better modeling work on methane dynamics in extreme environments and under changing environmental conditions can help improve the methane cycle predictive capability.

Author Contributions

ZS and SX jointly put forward the research objectives of this task, formulated the research route, wrote and organized this research together with ZW, WG, and HC: Data curation, Formal analysis. XZ and BZ: Review and editing.

Funding

This study was supported by the Natural Science Foundation of China (42176057, 41376077, and 41976057), the Shandong Special Fund of Pilot National Laboratory for Marine Science and Technology (Qingdao) (2021QNLM020002), and the Marine Geological Survey Program of China Geological Survey (DD20221707).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Adler, M., Hensen, C., Kasten, S., and Schulz, H. D. (2000). Computer Simulation of Deep Sulfate Reduction in Sediments of the Amazon Fan. Int. J. earth Sci. 88 (4), 641–654. doi:10.1007/s005310050294