Geochemical and Microbial Signatures of Siboglinid Tubeworm Habitats at an Active Mud Volcano in the Canadian Beaufort Sea

During the ARA08C expedition in 2017, sediment push cores were collected at an active mud volcano (420 m water depth) in the Canadian Beaufort Sea from two visually discriminative siboglinid tubeworm (ST) habitats that were colonized densely and less densely (ST1 and ST2, respectively). In this study, we investigated the biogeochemical and microbial community characteristics at ST1 by analyzing the geochemical properties, microbial lipids, and nucleic acid signatures, and comparing them with the data previously reported from ST2. The two ST sites showed distinct differences in vertical geochemical gradients [methane, sulfate, dissolved inorganic carbon (DIC), total organic carbon, and total sulfur], with a higher methane flux recorded at ST1 (0.05 mmol cm–2 y–1) than at ST2 (0.01 mmol cm–2 y–1). Notably, the δ13C values of DIC were more depleted at ST1 than at ST2, resulting in a higher proportion of DIC derived from the anaerobic oxidation of methane (AOM) at ST1 than at ST2. Moreover, both the ST1 and ST2 sites revealed the dominance of AOM-related lipid biomarkers (especially sn-2-hydroxyarchaeol), showing highly 13C-depleted values. The 16S rRNA analyses showed the presence of AOM-related archaea, predominantly anaerobic methanotrophic archaea (ANME)-3 at ST1 and ST2. Our results suggest that AOM-related byproducts (sulfide and DIC) potentially derived from ANME-3 were more abundant at ST1 than at ST2. This variation was attributed to the intensity and persistence of ascending methane. Therefore, our study suggests that AOM-derived byproducts are possibly an essential energy source for tubeworms during chemosynthetic metabolism, shaping different colony types on the seafloor.


INTRODUCTION
Cold seeps can host diverse faunal and microbial communities (e.g., siboglinid tubeworms (STs), mytilid and vesicomyid bivalves, and giant sulfide-oxidizing bacteria; Levin, 2005;Jørgensen and Boetius, 2007;Vanreusel et al., 2009;Niemann and Boetius, 2010). STs (referred to as typical chemosynthetic organisms) lack a digestive organ and thus rely entirely on their endosymbionts for nutrient and energy supply (Cavanaugh et al., 1981;Felbeck et al., 1981;Southward et al., 1986 and references therein). In most cases, bacterial symbionts are embedded within host cells in the interior of the trophosomes (highly vascularized tissue) (Bright and Giere, 2005). For instance, the bottom sections (i.e., roots) of tubeworms have a distinct function for taking up the reduced sulfur species (Julian et al., 1999;Freytag et al., 2001;Cordes et al., 2005), and transport various substrates to the endosymbionts mainly inhabiting trophosomes Duperron et al., 2014). Furthermore, the seawater-exposed tube of tubeworms may act as a direct conduit of sulfate from seawater to sediment layers (Duperron et al., 2014). Therefore, the growth of these tubeworms inhabiting cold seeps may provide the basis for chemosynthetic ecosystems that are decoupled from photosynthesis (Sibuet and Olu, 1998).
In methane-rich seeping systems, one of the most dominant microbial processes in sediments is the anaerobic oxidation of methane (AOM) associated with sulfate reduction (SR), which produces byproducts, such as sulfide and dissolved inorganic carbon (DIC) (Reeburgh, 1996;Hinrichs and Boetius, 2002). Tubeworms take up sulfide compounds, oxidize these compounds into sulfate, and then release sulfate by the posterior end to sediments, which may enhance the AOM (Felden et al., 2010;Hilário et al., 2011). Therefore, the long-term development of tubeworm colonization may play an important ecological role in the reduction of methane venting (Sibuet and Olu, 1998;Cordes et al., 2005;de Beer et al., 2006;Boetius and Wenzhöfer, 2013). The aerobic oxidation of methane (MOx) is another biological sink for ascending methane. The oxidized surface sediments may support specific MOx-related bacterial groups, which depend on the availability of oxygen (Omoregie et al., 2009;Boetius and Wenzhöfer, 2013). Such communities not only form dense mats at the sediment-water interface, but also coexist as endosymbiotic associations with benthic fauna such as bivalves and tubeworms (Fisher, 1990;Duperron et al., 2007;Petersen and Dubilier, 2009;Rodrigues et al., 2011). Therefore, the chemosynthetic interaction between benthic organisms and methanotrophs may play an essential role in efficiently reducing methane emissions from the seafloor.
Since methane assimilation performed by AOM-and MOxrelated microbial communities is influenced by varying gasrich fluid fluxes (Niemann et al., 2006), the ecological function of chemosynthetic ecosystems appears to be driven by local environmental factors, such as morphologic setting and physicochemical conditions (Rossel et al., 2011). In this context, active mud volcanoes (MVs) discovered in the Canadian Beaufort Sea can be considered as a natural laboratory for studying the biogeochemical dynamics of methane in sediments displaying various chemosynthetic habitats (Paull et al., 2015). In such settings, the activities of benthic organisms potentially control methane emissions from the seafloor (Boetius and Wenzhöfer, 2013). In this regard, Lee et al. (2019b) recently investigated vertical variations in AOM-and MOx-related methanotrophic communities in the mud volcano MV420 (420 m water depth) in the Canadian Beaufort Sea using sediment push cores collected from visually discriminative habitats devoid of megafauna and microbial mats to the naked eye, covered with bacterial mats, or colonized by STs. They showed that a niche diversification in MV420 shaped distinct methanotrophic communities due to the availability of electron acceptors in association with varying degrees of methane flux. However, our knowledge of the primary processes controlling the biological succession of STs in this region remains limited.
In this study, we investigated the biogeochemical and microbial community characteristics in a densely colonized ST field (ST1) in the MV420 of the Canadian Beaufort Sea by analyzing geochemical properties, microbial lipids, and nucleic acid signatures. The data from ST1 were then compared with those previously obtained from a less densely colonized ST field (ST2, Lee et al., 2019b) to constrain the primary biological process shaping the visually discriminative ST habitats in an active MV in the Canadian Beaufort Sea.

Sample Collection
In this study, one MV located at a water depth of 420 m (MV420) was mapped in detail (Figure 1) using an autonomous underwater vehicle during the ARA08C expedition of the Korean ice-breaker RV Araon in September 2017. Based on the observation of tubeworm patches dominated by Oligobrachia sp. at this MV (Paull et al., 2011(Paull et al., , 2015Lee et al., 2019a), we collected two push cores (ARA08C-DIVE104-1 and -10) using the remotely operated vehicle (ROV) "Mini" from a densely colonized ST site for investigation (ST1; 70.792 • N, 135.556 • W) (Figure 1). Similar to a previous study conducted at a less densely colonized ST site (ST2; 70.790 • N, 135.565 • W; Lee et al., 2019b), the first sediment core was subsampled at 2-4 cm intervals and the subsamples were transferred into serum vials containing saturated NaCl solution for gas analysis. After subsampling, these vials were immediately sealed with butyl rubber stoppers to prevent gas exchange. The second sediment core was used to extract pore waters at 4-cm intervals using the Rhizon soil moisture sampler. The extracted pore water was filtered through a 0.20 µm disposable polytetrafluorethylene filter before collection for onboard and post-cruise analyses. For the sulfate, DIC, and nutrient (NH 4 + and PO 4 3− ) analyses, an aliquot of pore water was stored at −20 • C until ion chromatography analysis. Another aliquot was collected in a 2-mL septum screw-lid glass vial and preserved with 20 µL of saturated mercuric chloride (HgCl 2 ) solution for DIC carbon isotope analysis (δ 13 C DIC ). Afterward, the same core was sliced at 1-cm intervals for the analyses of bulk elements, lipid biomarkers, and 16S rRNA gene sequences. Sediment samples were transferred into aluminum sheets and stored at −80 • C until analysis.

Gas and Porewater Analyses
The hydrocarbon composition in the headspace gas was measured using a gas chromatograph (GC) 7890A (Agilent Technologies, CA, United States) with a flame ionization detector (FID) at the Korea Institute of Geoscience and Mineral Resources (KIGAM). Precision from the repeated standard analysis was better than 5%. The δ 13 C values of methane (δ 13 C CH4 ) were obtained using an isotope ratio mass spectrometer (IRMS; Finnigan MAT 252, Thermo Fisher Scientific, Waltham, MA, United States) with a Combustion III interface (Thermo Fisher Scientific) at Nagoya University. The δ 13 C values were calibrated using methane standards with certified isotopic values (−35.2 ) obtained from the National Institute of Standards and Technology. The analytical errors were ± 0.2 for δ 13 C CH4 . Sulfate concentrations were measured using ion chromatography (761 compact ion chromatograph; Metrohm, Ionenstrasse, Switzerland). DIC concentrations were analyzed by titration using HCl at KIGAM. Reproducibility was monitored by repeated titrations of < 2% standard seawater supplied by the International Association of Physical Sciences of the Oceans for DIC. The δ 13 C DIC values were obtained using an IRMS (Finnigan DELTAplusXL, Thermo Fisher Scientific, Waltham, MA, United States) with a Finnigan GasBench-II headspace autosampler at Oregon State University. The δ 13 C values were calibrated using analyses of the Wiley CaCO 3 DIC standard with certified isotopic values (−0.4 ). The analytical error was ± 0.02 for δ 13 C DIC . Phosphate and ammonium concentrations were measured at KIGAM using colorimetric methods with a spectrophotometer (UV-2450; Shimadzu, Kyoto, Japan) at 885 and 640 nm (Gieskes et al., 1991). Data on methane and sulfate fluxes were calculated as reported in a previous study (Lee et al., 2019b), which includes a detailed description of the corresponding calculation.

Bulk Element Analysis
Sediment samples were freeze-dried and homogenized using an agate mortar prior to the bulk element analyses. For total organic carbon (TOC), sediment samples (∼1 g) were treated with 8 mL of 1 N HCl to remove carbonates. The TOC contents and their isotopic compositions were measured using an elemental analyzer (vario PYRO cube, Elementar, Germany) connected to an IRMS (Isoprime visION, Isoprime, Manchester, United Kingdom). The stable carbon isotopic values of TOC (δ 13 C TOC ) were reported in standard delta (δ) notation relative to the Vienna Pee-Dee Belemnite (VPDB) scale per mil. The δ 13 C values were calibrated using the International Atomic Energy Agency (IAEA) CH 3 standard with a certified isotopic value of −24.7 . The analytical precisions were ± 0.2 wt.% for TOC and ± 0.1 for δ 13 C TOC . The total sulfur (TS) content and its isotopic composition were measured using an elemental analyzer (EA1110, Thermo) connected to an IRMS (Dual-pumped 20-20-S, Secon) in organic elemental analysis labs in the United Kingdom. The stable sulfur isotopic values of TS (δ 34 S TS ) were reported in standard delta (δ) notation relative to the Vienna Canyon Diablo Troilite (VCDT) scale per mil. The δ 34 C TS values were calibrated using the IAEA S-2 and S-3 standards, with certified isotopic values of 22.7 and −32.3 , respectively. The analytical precisions were ± 0.3 wt.% for TS and ± 0.4 for δ 34 C TS .

Lipid Biomarker Analysis
Lipid analytical procedures were performed as described by Lee et al. (2018). Briefly, total lipids were extracted using an accelerated solvent extractor (Dionex ASE 200, Dionex Corporation, Sunnyvale, CA, United States) and separated into apolar and polar fractions over an Al 2 O 3 column (activated for 2 h at 150 • C). After column separation, 40 µL of 5α-androstane (10 µg mL −1 ) was added to the apolar fraction as an internal standard. The polar fraction was divided into two aliquots to which either C22 7, 16-diol (10 µg mL −1 ), or C46 glycerol dialkyl glycerol tetraethers (GDGTs; 10 µg mL −1 ) were added as internal standards. Half of the polar fraction containing C22 7,16-diol was dried and silylated with 25 µL of N,O-bis(trimethylsilyl) trifluoroacetamide (BSTFA) and 25 µL of pyridine before heating to 60 • C for 20 min to form trimethylsilyl derivatives. The second half of the polar faction containing C 46 GDGT was redissolved by sonication (5 min) in hexane:isopropanol (99:1, v:v) and filtered with a 0.45 µm polytetrafluoroethylene filter. Subsequently, an aliquot of the filtered fraction was treated with iodic acid following the procedure described by Kaneko et al. (2011) to cleave ether bonds from the GDGTs. All the apolar and polar fractions were analyzed using a GC equipped with an FID, GC connected to a mass spectrometer (MS), and high-performance liquid chromatography-atmospheric pressure positive-ion chemical ionization-MS (HPLC-APCI-MS). The δ 13 C values of selected compounds were determined by an IRMS connected to a GC via a combustion interface [glass tube packed with copper oxide (CuO), operated at 850 • C]. Isotopic values were expressed as δ 13 C values per mil relative to the VPDB. The δ 13 C values were corrected for the introduction of additional carbon atoms during salivation. The analytical errors were less than ± 0.4 .

Nucleic Acid Analysis
Genomic DNA was extracted using an Exgene Soil SV mini kit (Cambio, United Kingdom) at GeneAll (Seoul, South Korea). The V4-V5 regions (coverage of 85-95%) of bacterial 16S rRNA genes and the V6-V8 regions (coverage of 66-88%) of archaeal 16S rRNA genes were amplified by polymerase chain reaction using primer pairs 515F/926R and A956F/A1401R, respectively, at the Integrated Microbiome Resource (IMR) at Dalhousie University, Canada 1 (Walters et al., 2016). Sequencing of the amplicons was conducted at IMR using the paired-end (2 × 300 bp) Illumina MiSeq system (Illumina, United States). Paired-end sequences were trimmed based on quality scores using Sickle (Schirmer et al., 2015) followed by error corrections with BayesHammer (Nikolenko et al., 2013). The resultant quality-trimmed and error-corrected paired-end sequences were assembled using PANDAseq (Masella et al., 2012). Sequenceprocessing steps were performed by following the Miseq SOP in mothur. The assembled sequences were aligned against a SILVA alignment (version 132 2 ; Pruesse et al., 2007) and subsequently denoised using the "pre.cluster" command. Chimeric sequences were removed using the "chimera.uchime" command in de novo mode (Edgar et al., 2011). The sequences were further clustered to operational taxonomic units (OTUs) at a 97% sequence similarity level using the OptiClust clustering algorithm. Sequences were normalized by rarefying sequences per sample to the smallest library size (17,667 for bacteria and 6,768 for archaea). Taxonomic assignments of each OTU (2,286 bacterial OTUs and 2,758 archaeal OTUs) were determined by sequence similarity searches against the EzBioCloud database (Yoon et al., 2017). All sequence data used in this study were deposited in the Sequence Read Archive at the National Center for Biotechnology Information under the accession numbers PRJNA559003 and PRJNA643655. Phylogenetic trees of major archaeal OTUs of Methanomicrobia and bacterial OTUs of Gammaproteobacteria and Deltaproteobacteria with greater than 3% relative abundance were constructed using the maximum-likelihood algorithm (Felsenstein, 1981) and the General Time Reversible evolutionary model with MEGA X (Kumar et al., 2018). The robustness of the tree topologies was assessed by bootstrap analyses based on 1,000 replications of the sequences. To infer the putative sulfate-reducing partners of ANMEs, co-occurrence between two dominant archaeal OTUs of ANMEs and 8 bacterial OTUs of Deltaproteobacteria with greater than 3% relative abundances was investigated by calculating Pearson correlation coefficients using the "cor" function of the base R package.

Geochemical Profiles
The gas and porewater properties varied along the sediment depths of ST1 (Figure 2 and Supplementary Table 1). The methane concentrations were below 5.4 mM, and the δ 13 C CH4 values ranged from −42.8 to −15.0 . The sulfate concentrations varied between 2.9 and 27.8 mM. The measured DIC concentrations and δ 13 C DIC values ranged from 2.6 to16.1 mM and from −25.9 to −7.8 , respectively. The phosphate concentrations were in the range of 0.3-0.8 mM, while the ammonium concentrations were below 0.1 mM. For the bulk elements, the TOC contents and δ 13 C TOC values ranged from 1.2 to.4 wt.% and from −27.0 to −26.4 , respectively (Figure 2). The TS contents and δ 34 S TS values varied between 0.3 and 0.5 wt.% and ranged from −18.3 to −2.4 , respectively (Figure 2).

Geochemical Signatures of Methane Oxidation
Overall, the methane concentrations were lower in the upper 5 cm of the sediment core at ST1 (Figure 2), which was comparable to those previously reported at ST2 (Lee et al., 2019b, see Supplementary Figure 1). This depleted profile of methane concentration at ST1 was also similar to those at other cold seepages (e.g., Haakon Mosby MV and Japan Trench) which were characterized by tubeworm patches and clam colonies Felden et al., 2014). Moreover, the δ 13 C CH4 profile investigated at ST1 was enriched only to a sediment depth of 5 cm, in contrast to the relatively consistent δ 13 C CH4 values observed at ST2 (Figure 2 and see also Supplementary  Figure 1). The enriched degree of δ 13 C CH4 values appears to be related to methane oxidation in connection with the penetrated sulfate depths, due to selective isotopic fractionation (Whiticar et al., 1986;Boetius et al., 2000;Knittel and Boetius, 2009). In such cases, most of the methane is reduced via AOM, which occurs at sulfate-penetrated depths within the sulfate-methane transition zone (SMTZ) (Reeburgh, 1980), indicating that the remaining methane can approach relatively enriched isotope values (Pohlman et al., 2009). Notably, the signature of 13 Cenriched methane could be in part attributed to the contribution of methylotrophic methanogenesis, as shown by the metabolism of Methanococcoides that was cultivated from seeping sediments (Wegener et al., 2016). It is also worth noting that tubeworms live in symbiosis with thiotrophic bacteria that utilize sulfide as energy source and release sulfate into sediments (Lösekann et al., 2007;Felden et al., 2010). Hence, the 13 C-depleted value of DIC at 5 cm sediment depth might be associated with the enhanced AOM. The depth and thickness of SMTZ within tubeworm habitats can be thus controlled by various factors, such as the depth of the methane production zone, bioirrigation, the flux of methane and sulfate, and their consumption rates (Jørgensen and Kasten, 2006;Knittel and Boetius, 2009;Felden et al., 2010;Meister et al., 2013). Considering the different methane fluxes calculated at both ST sites (Table 1), the SMTZ may have formed below a sediment depth of approximately 5 cm at ST1, whereas it may occur at deeper sediment depths (i.e., below ca. 9 cm) at ST2. With respect to the incomplete depletion of methane and sulfate in the deeper sediment layers, we could not determine the exact SMTZ layers at both sites. Nevertheless, together with 13 C-enriched methane, the most 13 C-depleted DIC (∼-25.9 ) occurred at a 5-cm sediment depth at ST1 (Figure 2), reflecting the potential increase of AOM-derived byproducts within the sediments Felden et al., 2010). The relative contribution of the AOM-derived DIC along the sediment depths was approximately estimated based on the isotopic mass balance (cf. Coffin et al., 2013), assuming that the δ 13 C signatures of the parent methane and seawater DIC in the Canadian Beaufort Sea were −64 and −0.4 , respectively (Paull et al., 2015). Our results revealed that 11.7-40.2% of the total DIC was contributed by the AOM-derived DIC source at both ST sites (see also Supplementary Figure 2 and Table 1). In particular, the highest contribution of AOM-derived DIC in the depth profiles at ST1 was consistent with the most depleted δ 13 C DIC value at a 5-cm sediment depth, compared to the relatively lower contributions (below 30%) at ST2. Furthermore, the increasing nutrient profiles in the deeper sedimentary layers at both ST sites were typically indicative of the presence of byproducts derived from the degradation of organic matter (Figure 2).
The bulk element (C and S) contents and their isotopic values showed different vertical variations at both ST sites (Figure 2 and see also Supplementary Figure 1). The origin of 13 C-depleted isotopic signatures is commonly regarded as terrestrial plants (Peters et al., 1978;Fry, 1999;Perdue and Koprivnjak, 2007). However, 13 C-depleted TOC has also been found in seeping sediments, likely resulting from the assimilation of methane into bacterial biomasses (Haese et al., 2003;Elvert et al., 2005). In this regard, the bacterial biomass can be incorporated and preserved in solid-phase sediment, and through time can constitute a significant portion of the organic carbon pool (Canuel and Martens, 1993;Gong and Hollander, 1997). As AOM occurs at deeper sediment layers under low methane fluxes (Niemann et al., 2006;Lee et al., 2019b), bacterial biomass involved in MOx are incorporated to some degree into the sedimentary carbon pool at the surface of ST2. Under this setting, the isotopic shift in the TABLE 1 | Geographic coordinates (latitude and longitude), water temperature at the sea bottom ( • C), sediment core length, and methane and sulfate fluxes (mmol cm −2 y −1 ) calculated from concentration profiles and penetration depths (cm) of sulfate and AOM-related DIC production (%).
Siboglinid tubeworms a Siboglinid tubeworms field 1 (ST1) (ARA08C-field 2 (ST2)(ARA08C-DIVE104-1 and 10) DIVE105-12) and 14) Latitude (  TS profiles may also provide indirect evidence for the supply of electron acceptors (sulfate and oxygen) derived from seawater (Canfield and Thamdrup, 1994). Assuming that seawater sulfate with a δ 34 S value of + 21 was supplied to the sediments (Rees et al., 1978), the significant isotopic shift of δ 34 S TS at ST2 can potentially be regarded as the active penetration of downward sulfate with oxygen from seawater. Conversely, with respect to the significant AOM signals at ST1, the depleted δ 34 S signatures (−12.0 ± 2.0 ) along depth profiles may indicate that sulfide or sulfur species are likely produced via microbial SR coupled with methane oxidation (Boetius et al., 2000;Jørgensen et al., 2004;Reeburgh, 2007). Indeed, the δ 34 S values of sulfide may be depleted by as much as 46 relative to those of seawater sulfate during SR (Hoefs, 2007). With respect to the input of seawaterderived sulfate and the oxidation of sulfide to intermediates (e.g., elemental sulfur and thiosulfate), the sulfide pool at ST1 appears to be, at least in part, associated with AOM-derived SR. Although more complicated sulfur cycles involved in alternative pathways may occur in sediments, the range of δ 34 S TS values at both sites suggests that SR coupled with AOM was more active under higher methane fluxes. Therefore, we argue that the different methane fluxes may be potentially associated with the discriminative biogeochemical signals observed at both ST sites.

Molecular Biogeochemical Signatures of Methane Oxidation
For AOM-related biomarkers, we identified sn-2hydroxyarchaeol at ST1 and ST2 (Figure 3 and see also Supplementary Figure 3). The sn-2-hydroxyarchaeol concentrations were generally higher in the SMTZ than in the surface layers at both sites (Figures 6A,C), with 13 Cdepleted values compared to the δ 13 C CH 4 values (δ 13 C CH4 ; −27.1 ± 10.8 (n = 6) at ST1 and −59.2 ± 4.5 (n = 4) at ST2). The difference in the δ 13 C values between sn-2hydroxyarchaeol and methane (i.e., δ 13 C) in the SMTZ was, on average, −69.5 ± 7.2 (n = 3) at ST1 and −48.7 ± 2.4 (n = 2) at ST2 (Figure 6B). The methane-derived microbial biomass was generally depleted in 13 C compared to the source methane as a result of isotopic fractionation during methane assimilation (Whiticar, 1999). The negative δ 13 C values of sn-2-hydroxyarchaeol in the SMTZ at both sites were thus indicative of the AOM occurrence in sediments where sulfate was present (Pancost et al., 2000;Elvert et al., 2005;Niemann et al., 2005;Bradley et al., 2009). Furthermore, the δ 13 C values of predominant archaeol at ST1 (−27.7 ± 2.7 , n = 3) were also negative (Figures 6A,B), providing additional evidence for an active AOM in the SMTZ at ST1. It is worth noting that the δ 13 C values of both sn-2-hydroxyarchaeol and archaeol were more depleted at ST1 than at ST2, suggesting that the AOM activities were stronger at ST1 than at ST2. Interestingly, the SMTZ δ 13 C values of archaeol were more enriched in the Håkon Mosby Mud Volcano (HMMV) setting, which was also covered by beggiatoa organisms, than at the Beaufort sites. Enriched 13 C-archaeol has typically been associated with the lipid production of methanogens that use CO 2 as their carbon source (Belyaev et al., 1983). Accordingly, the isotopic compositions of methanogenic biomarkers are approximately equal to those of methane. Thus, considering the isotopic fractionation involved in methanogenesis, the 13 C-enriched archaeol at the HMMV was attributed to the contribution of methane-producing, rather than methane-consuming organisms (Whiticar, 1999;Niemann et al., 2006). Isoprenoidal GDGTs were found in all the sediments investigated (see Supplementary Figure 3 and Supplementary Table 2). The concentrations of GDGT-1, -2, and -3 remained constant along the core depth; however, both GDGT-0 and crenarchaeol were more abundant near the surface layers at both sites. Moreover, the isoprenoidal GDGT distributions did not show a clear dominance of GDGT-1, −2, and −3 over GDGT-0, which is an indicator of the contribution of methanotrophic Euryarchaeota (Schouten et al., 2013). The values of other methanotroph-related indices such as GDGT-0 to crenarchaeol  and GDGT-2 to crenarchaeol (Weijers et al., 2011), and the methane index values (Zhang et al., 2011) were also low, with ranges of 0.4-1.6, 0.0-0.2, and 0.1-0.4, respectively. Accordingly, it appears that marine pelagic Thaumarchaeota (e.g., Pearson et al., 2001;Wuchter et al., 2005;Schouten et al., 2013) were the major contributors to the total GDGT pool in the Canadian Beaufort Sea MV, in contrast to other cold seeps such as the Hydrate Ridge  and the Marmara Sea (Chevalier et al., 2013). With respect to carbon metabolism at the molecular level, BP-1 has one cyclopentane ring that may originate from GDGT-1, -2, or -3 (Schouten et al., 2013). If methanotrophs are responsible for synthesizing these three GDGT compounds, then the isotopic signatures derived from BP-1 are subject to the presence of AOM (Niemann et al., 2005). The δ 13 C values of BP-1 in the SMTZ at ST1 were FIGURE 5 | Depth profiles of the relative abundances of major OTUs (Methanomicrobia; A_OTU001 and A_OTU008, and Methylococcales; B_OTU004 and B_OTU010) obtained from the (A) ST1 and (B) ST2 sites. The ST2 data were obtained from Lee et al. (2019b). Gray shading indicates potential SMTZ layers. similar to those of archaeol (Figure 6B), suggesting the influence of AOM-related methanotrophs. However, the contribution of the AOM-related methanotrophs should be minor, considering the low methanotroph-related indices mentioned above. In contrast, BP-3 is produced only by crenarchaeol, and thus contains the pelagic Thaumarchaeota signature (Schouten et al., 2013). The 13 C-enriched isotopic signatures of BP-3 relative to other compounds in the SMTZ at both ST sites ( Figure 6B) were comparable to those previously reported in non-seeping marine environments (Hoefs et al., 1997). This indicates that marine pelagic Thaumarchaeota using seawater derived-DIC as a carbon source were the primary contributors to the GDGT pool at both sites (e.g., Hoefs et al., 1997;Könneke et al., 2005;Wuchter et al., 2005).
For MOx-related biomarkers, we identified 4α-methyl sterols and C16:1ω8 at ST1 and ST2 (Figure 3 and see also Supplementary Figure 3). Overall, the concentrations of 4α-methyl sterols increased toward the upper sediment layers at both sites (Figure 3 and see also Supplementary Figure 3). Their occurrences have been associated with aerobic methanotrophs (e.g., Methylococcales belonging to alphaproteobacteria) near the surface where oxygen as an electron acceptor is likely available . The δ 13 C values of the 4α-methyl sterols were more depleted near the surface at both sites (Figure 3 and see also Supplementary Figure 3). Thus, it appears that the MOx process predominantly occurred in the upper sediment layers at both sites. The biosynthesis of 4α-methyl sterols might occur during the oxidative removal of methyl carbons at C-4 and C-14 of lanosterol or of another, yet undescribed, protosterol (Pearson et al., 2003). Furthermore, a previous study suggested that Methylococcus capsulatus uses an alternative enzymatic pathway for the C-4 demethylase reaction (Summons et al., 2006). Therefore, the predominance of these MOx-related methanotrophs near the surface at both STs seems to be influenced by more effective use of nutrients and substrates (e.g., oxygen/nitrate and methane), in contrast to those of AOM-related archaeal and sulfur-oxidizing bacterial communities. Notably, the average δ 13 C value of 4α-methyl sterols investigated at ST1 (−19.4 ± 3.0 , n = 3, Figure 6D) was more depleted than that of ST2 (14.1 ± 8.3 , n = 4, Figure 6D), and most likely resulted from the stronger MOx activities at ST1 than at ST2. As supporting evidence, we also found higher abundances and depleted δ 13 C values of C16:1ω8 in the surface layers at ST1 and ST2 (Figure 3). However, the δ 13 C values of C16:1ω8 were lower at ST1 (−16.4 ± 1.5 , n = 2) than at ST2 (−19.4 ± 3.0 , n = 3), similar to those of the 4α-methyl sterols ( Figure 6D). With respect to the synthesis of C16:1 fatty acids with ω8 double bond positions (Bowman et al., 1991(Bowman et al., , 1993Fang et al., 2000), the 13 C-depleted C16:1ω8 values in the nearsurface layers at ST1 suggest the presence of MOx-related bacteria derived from type I methanotrophs (Makula, 1978;Nichols et al., 1985). In particular, the negative δ 13 C values of C16:1ω8 are indicative of the involvement of type I methanotrophs in the RuMP pathway for methane oxidation (Jahnke et al., 1999).
As evidence for other bacterial communities, we found C16:1ω5 at both sites, which were also more abundant in the upper sediment layers where the 4α-methyl sterols were predominant (Figure 3 and see also Supplementary Figure 3). The δ 13 C values of C16:1ω5, associated with sulfate-reducing bacteria (Elvert et al., 2003;Blumenberg et al., 2004;Niemann et al., 2006), were lower in the upper sediment layers at ST1 (−12.9 , n = 1, Figure 6D) than at ST2 (-1.4 ± 2.8 , n = 3, Figure 6D). The δ 13 C values of C16:1ω5 at both sites were more depleted than those ( δ 13 C values: 30-40 ) in nonseeping Arctic environments (Birgel et al., 2004). Given that the MOx-related methanotrophs relative to the AOM-related methanotrophs were more abundant near the surface layers, as discussed above, the isotopic signature of C16:1ω5 at ST1 might be related to the assimilation of MOx-derived DIC. Hence, the δ 13 C pattern of bacterial fatty acids (i.e., C16:1ω8 and C16:1ω5) at ST1 and ST2 was in line with that of the MOxrelated biomarker, that is, the 4α-methyl sterols, which suggests that MOx-related processes were stronger at ST1 than at ST2.
To further clarify the phylogenetic position within the class Methanomicrobia, phylogenies of the two dominant Methanomicrobia OTUs (A_OTU001 and A_OTU008) were inferred from 16S rRNA gene sequences (Supplementary Figure 4 and Supplementary Table 3). ANME-3 (A_OTU001) accounted for more than 50% of the archaeal sequences at 5-17 cm sediment depths at both sites ( Figure 5). ANME-2 (A_OTU008) was mostly detected at 4-9 cm sediment depths at ST2, accounting for 2.4-8.7% of total OTUs (Figure 5). Together with our lipid data, ANME-2 and ANME-3 groups were identified at 4-11 cm sediment depths at ST2, in contrast to the predominant occurrence of ANME-3 at ST1. To investigate the relationship between ANMEs and SRB, the presence of the sulfate-reducing Deltaproteobacteria was determined by analyzing deltaproteobacteria OTUs (Desulfobacteraceae and Desulfobulbaceae) at both ST sites. Several OTUs of Desulfobacteraceae showed positive Pearson correlations with dominant ANME-2 and ANME-3. However, deltaproteobacterial OTUs obtained from the Beaufort Sea which have positive correlations with ANME-2 and ANME-3 are phylogenetically distant from previously known SRB of Desulfosarcina-Desulfococcus or Desulfobulbus (Knittel et al., 2005;Niemann et al., 2006). Hence, to better assess ANME-SRB consortia, further studies in the Beaufort Sea are required by applying fluorescence in situ hybridization. Notably, both ST sites FIGURE 6 | Comparison of the concentrations and isotopic compositions ( δ 13 C) of selected lipid biomarkers (archaeol, sn-2-hydroxyarchaeol, biphytanes derived from GDGTs, 4α-methyl sterols, and bacterial fatty acids) obtained from ST1, ST2, and the Haakon Mosby mud volcano (HMMV): sites (A,B) for the sediment layers within the SMTZ, and (C,D) for the sediment layers above the SMTZ. Note that the ST2 and HMMV data were obtained from Niemann et al. (2006) and Lee et al. (2019b), respectively. δ 13 C values were calculated by subtracting the δ 13 C values of methane from those of selected lipid biomarkers.
contained the 4α-methyl sterols, which are putatively specific for bacterial groups involved in MOx. The 13 C-depleted 4α-methyl sterols were more abundant in the upper sediment layers at both ST sites in conjunction with high proportions of OTUs belonging to the order Methylococcales. The vertical distributions of Methylococcales varied from those of ANME-2 and ANME-3, showing higher proportions in the shallower sediment layers at both ST sites. Thus, the vertically discriminative trends between the AOM-and MOx-related methanotrophs may be related to the supply of available electron acceptors (i.e., sulfate and oxygen) under the bioirrigation activities of STs.
At MV420, the distinct development of ST patches investigated at ST1 appeared to correlate with the active AOM process associated with mainly ANME-3, while the co-occurrence of ANME-2 and ANME-3 was observed at ST2 which was characterized by lower methane flux (Figure 5 and see also Supplementary Figure 4). In this regard, ANME-3 has been found only in cold seepages (including MVs) with high methane partial pressures and temperatures of < 20 • C (Niemann et al., 2006;Lösekann et al., 2007;Ruff et al., 2015;Gründger et al., 2019). For instance, ANME-3 was reported at the HMMV in the Barents Sea (at a water depth of 1,250 m), Kazan MV in the Eastern Mediterranean (at a water depth of 1,700 m), and surface sediments from the Guaymas Basin in the Gulf of California (at a water depth of 1,500 m) (Niemann et al., 2006;Heijs et al., 2007;Pachiadaki et al., 2010;Vigneron et al., 2013). The bioirrigation of STs sustains the supply of available electron acceptors (especially sulfate) in sediment layers, resulting in increasing AOM rates within the SMTZ (Treude et al., 2003;Niemann et al., 2006;Ruff et al., 2013Ruff et al., , 2018. Therefore, ANME-3 might have been more strongly forced to tolerate high sulfide concentrations produced during the AOM at ST1 than at ST2, thereby supporting the metabolism of sulfur-oxidizing endosymbionts inhabiting STs Duperron et al., 2014;Lee et al., 2019a).
Notably, ANME-2 (mainly ANME-2c) often occurs in sediments bioirrigated by benthic organisms, such as at the Hydrate Ridge in the Northeast Pacific and the REGAB pockmark in the Congo Basin Knittel et al., 2005;Pop Ristova et al., 2012), and in low fluid flux regimes Wegener et al., 2008). Therefore, ANME-2 detected at ST2 may have a preferential niche in relatively low-methane fluxes. With respect to the relatively small contribution of AOM-derived DIC at ST2, AOM activities by ANME-2 and -3 may be weakened for the production of sulfide and DIC, causing the low density of tubeworm patches. Furthermore, the commonly increased abundance of Methylococcales near the surface of both sites may additionally contribute to the DIC supply related to tubeworm metabolism.
A previous MV420 study in the Canadian Beaufort Sea suggested that visually discriminative habitats that were devoid of megafauna and microbial mats (DM) to the naked eye, covered with bacterial mats (BM), or colonized by STs were shaped in association with varying degrees of methane flux (Lee et al., 2019b). Notably, the BM (0.04 mmol cm −2 y −1 ) and ST1 (0.05 mmol cm −2 y −1 ) sites showed similar methane fluxes. Therefore, it appears that the methane flux is not a unique environmental factor for niche diversification of chemosynthetic organisms within an MV system. In this regard, the discriminative colonization of chemosynthetic organisms might be related to substrate availability, according to the mass transfer limitation of porewater fluids in stagnant sediments (Boetius and Suess, 2004;de Beer et al., 2006). Under slower up-flow fluids that permit the diffusional penetration of sulfate, AOM communities can be developed, producing sulfide, which is then aerobically oxidized by chemosynthetic bacteria (e.g., Beggiatoa and tubeworms). In contrast to DM site which was characterized by the recent, active methane outgassing, the sediments at the BM and ST1 sites were influenced by upflowing methane for a longer period, providing sufficient time for a microbial and faunal colonization associated with the AOM. However, ST2 may have had insufficient time for its formation to be colonized by slower growing STs. Therefore, we infer that the differences in colonization between ST1 and ST2 may be related to the varying persistence of ascending methane, together with varied methane fluxes. Our data thus provide a new clue about a succession of potentially competing tubeworm communities over methane exposure time on the same MV.

CONCLUSION
Integrated biogeochemical and microbial investigations were performed on the sediment push cores retrieved from densely (ST1) and less densely (ST2) colonized tubeworm patches at MV420 in the Canadian Beaufort Sea. Compared to ST2, which recorded a low methane flux (0.01 mmol cm −2 y −1 ), the significantly depleted δ 13 C DIC values measured in the depth profiles of ST1 can be regarded as active AOM evidence, indicating the relatively highest contribution (40.2%) of AOMderived DIC within the SMTZ of this site. Similarly, the mostly depleted δ 34 S TS values recorded at ST1 also appear to reflect active SR involved in AOM. More specifically, we found isotopically 13 C-depleted lipid biomarkers (sn-2hydroxyarchaeol and 4α-methyl sterols) and related OTUs, belonging to ANME-2/-3 and Methylococcales, mediating AOM and MOx in the sediments of both STs. Notably, the prevalence of ST patches observed at ST1 may be associated with the supply of AOM-related byproducts (sulfide and DIC) performed mostly by ANME-3 as a result of active AOM. Conversely, at ST2, sulfide and DIC production via AOM may not be well facilitated by ANME-2 and -3, resulting in the low density of ST patches. Taking into consideration the common biogeochemical signatures investigated at another habitat under a similar methane flux as the Canadian Beaufort MV, a possible scenario explaining the biological succession of benthic communities is likely the change in the up-flow persistence of methane-rich and sulfate-depleted pore water. Therefore, we infer that the various byproducts produced by methanotrophic activities play an important role in the life strategy of STs involved in chemosynthetic metabolism. Consequently, our study reveals that the interaction between methanotrophs and benthic organisms is critical for understanding the ecology, diversity, and biogeochemistry of MVs and other types of cold seep ecosystems.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

AUTHOR CONTRIBUTIONS
DHL, YML, and JHK designed and coordinated the study and led the writing and interpreted the data. YKJ, CP, and KHS contributed to sample acquisition. DHL conducted experiments for gas, porewater, and lipid biomarkers, and YML analyzed microbial diversity. J-HK provided the methane and DIC concentration data. All authors commented on the manuscript and contributed to the writing.

ACKNOWLEDGMENTS
We sincerely appreciated the two reviewers for their valuable comments. We would like to thank the captain and crews of the R/V Araon for their help at sea. We also thank S. Kang, D. Kim, and S. Kim for their analytical assistance in the laboratory at Hanyang University and D. Graves, F. Flores, and R. Gwiazda at the Monterey Bay Aquarium Research Institute for remotely operated vehicle operation.