The soil was collected from a raised bog, Touxi (42°17.144 N, 127°50.277E, as indicated by the triangle; Figure 1a), located northwest of the dormant volcano found in the Changbai Mountain range situated at the boundary between northeastern China and North Korea. Topographically, the area is miscellaneous, comprising hills, steep slopes, valleys, and hills. The average annual temperature is around −7 to +3°C. The mean annual precipitation of 700–1,400 mm, and the peatland vegetation was mainly dominated by Sphagnum sp. along with a mixture of Sanguisorba sp., Carex, Thelypteris sp., Iris sp., Aulacomnium palustre, Orchis sp., Eriophorum angustifolium, Oxycoccus palustris, Trichophorum sp. and Euphorbia sp. The sample was collected up to 50 cm in depth, sealed in polyethylene plastic bags, and stored at 0°C before further analysis at the laboratory.

The sample was thawed to room temperature in the laboratory, and all of the peat was milled into small particles (diameter < 1 cm) manually. Then, it was homogenized to ensure an equal amount of organic matter in the whole sample, and plant roots and debris were removed. A subsample was collected at the start and was frozen at −20°C for DNA analysis within a week.

The mesocosms experimental units were polyvinyl chloride cylinders that were 60 cm high and 30 cm wide (Figure 1b) and were placed in the greenhouse to maintain a constant temperature (25°C). Each column was filled with clean quartz at the bottom with a thickness of ~5 cm prior to peat addition, then fine white mesh was plated to isolate the next added peat and quartz. Subsequently, native soil (peat) was added with a total depth of ~50 cm, i.e., about 5 cm below the top of the column. The mesocosm was then filled with water, and the water table was kept 3 cm above the surface of the peat in order to avoid O₂ diffusion from air. Additionally, small pieces of PE shading foil were added to surface water to prevent algal growth. After about 65 days, the concentration of CO₂ and CH₄ have reached a steady state. Therefore, each of the samples, except control, was amended for the following treatments (1) Control (just with sand and peat-filled column and water); (2) Addition of 0.03 mol SO₄²⁻ (Na₂SO₄); (3) Addition of 0.03 mol SO₄²⁻ + 0.5 g HA; (4) Addition of generally 0.5 g HA and (5) Addition of 0.05 mol (generally) Goethite (as a source of Fe). Each mesocosm was also equipped with water and gas samplers. Pore water and silicon gas samplers were installed at 5, 10, 15, 25, 35, and 45 cm depth, as shown in Figure 1b. The mesocosms were allowed to settle for the first 60 days and then incubated in the laboratory for ~223 days (in total) in a climate chamber at 15°C (16 h light/8 h dark cycles, 660 μmol s⁻¹ photosynthetic photon flux). Along with measuring CO₂ and CH₄ at different depths, the flux of these gasses was also measured on similar days (68, 75, 85, 96, 103, 116, 129, 143, 157, 181, 195, 207, and 223 days).

Gas samples were collected, and dark chamber flux measurements were conducted on days 68, 75, 85, 96, 103, 116, 129, 143, 157, 181, 195, 207, and 223 days. Compared to the values obtained at different depths of mesocosm designated as the control, the percentage values have been calculated and given along with actual values in the results. The samples were analyzed to measure the amount of CH₄ and CO₂ were measured using a gas chromatograph (Agilent 7890A, Japan) equipped with FID and a CO₂ methanizer.

Additionally, a vented stated chamber was used to collect surface gas samples from each mesocosm. CH₄ and CO₂ fluxes were measured on similar dates for each mesocosm. A dark chamber (20 cm x 21 cm) was used to determine the flux of these gasses through dynamic chamber measurement. Briefly, three identical floating chambers were placed on the top of mesocosms to estimate a better representative measurement as GHG emission on the water surface is heterogeneous. Fluxes from the headspace were calculated by linear regression of at least six consecutive readings within 30 min (one reading every 5 min) and stored for further analysis. A long measurement period was used to detect even the most minor changes in GHG fluxes. The concentration was later measured by the same method in GC as described above.

For DNA analysis, samples were collected from two places, the subsurface layer, and the bottom layer. The samples were named accordingly (Table 1). E.Z.N.A.^® DNA Kit for Soil (Omega Bio-tek, Norcross, GA, United States) was used to extract the microbial DNA from all mesocosms soil samples, following the given manual. 2% agarose gel was used to evaluate the quality of the DNA. The V4-V5 region of the bacteria 16S ribosomal RNA gene was amplified by PCR using forward primer 515F (5′-barcode-GTGCCAGCMGCCGCGG-3′) and 907R (5′-CCGTCAATTCMTTTRAGTTT-3′). The barcode is a unique eight base sequence for each sample. PCR mixture included of 20 μL mixture comprising four μL of 5 × FastPfu Buffer, 2 μL of 2.5 mM dNTPs, 0.8 μL of each primer (5 μM), 0.4 μL of FastPfu Polymerase, and 10 ng of template DNA. Agarose gel (2%) was used to extract amplicons that were purified by using AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, United States) according to the manufacturer’s instructions and later Quantus^™ Fluorometer (Promega, United States) was used to quantify them.

The PCR products were quantified using Qubit^®3.0 (Life Invitrogen), and amplicons (every 24) with different barcodes were equally mixed. In addition, Illumina’s genomic DNA library preparation procedure was used to construct the Illumina Pair-End library from pooled DNA at Shanghai BIOZERON Biotechnology Co., Ltd. (Shanghai, China) on the NovaSeq PE250 platform.

In order to analyze ASV files obtained from QIIME 2, a microbiome analyst was used, which compiles and analyzes microbiome data online comprehensively and visually (Chong et al., 2020). Briefly, the alpha diversity (diversity of species composition) was determined by Shannon, ACE, Simpson, Chao1 index, and Fisher curves. In contrast, An analysis of similarities (ANOSIM) test was used to compare beta-diversity determined using Bray-Curtis distances and principal coordinate analysis plots (PCoA) based on Bray-Curtis distances. A heatmap of the most abundant classes was generated based on complete hierarchical clustering using Euclidian distances. Bacteria differentially represented in different mesocosms were evaluated by using linear discriminant analysis (LDA) together with effect size (LEfSe). A taxonomic cladogram illustrates differences among genera based on significant taxa (Segata et al., 2011).

In the present study, four different treatments with different TEAs, (1) SO₄²⁻; (2) SO₄²⁻ + HA; (3) HA, and (4) Goethite, along with control, were studied to evaluate the impact of TEAs on CO₂ and CH₄ production as well as the density and diversity of microbial communities in peatland soils. In the first and second mesocosm, SO₄²⁻ was added individually and with HA, respectively. The addition of SO₄²⁻ resulted in higher depletion of CH₄ in the surface layers of mesocosm, it was decreased up to 124.7 (↓56.2%) and 108 (↓ 71.2%) μmol L⁻¹ in the 5 and 10 cm sample, respectively. In contrast, increased concentrations of sulfate were less detrimental to the production of CH₄ (↓34%) in the deeper layers of the sediment (Figure 2Ab). The addition of SO₄²⁻, however, led to an increase in the production of CO₂ (Figure 2Bb). The addition of SO₄²⁻ had increased CO₂ by almost twice the amount of CO₂ produced in control, i.e., 6.75–6.668 mmol L⁻¹ (↑70–86.3%). After the addition of SO₄²⁻ and HA, the production of CH₄ was strongly inhibited, up to 36.06 μmol L⁻¹ (↓90.49%) in the surface layers and ↓47.53% in the bottom layer. On the other hand, only a slight increase in CO₂ production was observed between 3.55 and 5.29 mmol L-1 (↑18.18–24.7%) (Figures 2Ac,Bc). The CH₄ production decreased significantly after HA was added, down to 56.06 μmol L⁻¹ (↓55.9%) in the top layer and 127.7 μmol L⁻¹ (↓41.42%) in the bottom layer. As with other TEAs, the addition of HA increased the production of CO₂ following its addition from 29.5–46.6% (Figure 2Bd).

A drastic reduction in CH₄ production was observed upon the addition of Goethite to 12.38 μmol L⁻¹ (↓90.25%) in the top layer and 109.1 μmol L⁻¹ (↓71.6%) in the bottom layer (Figure 2Ae), which was the highest reduction among all mesocosms. Interestingly as compared to CH₄, goethite addition did not result in a significant increase in CO₂ production, where its value was approximately 3.8 mmol L⁻¹ (↑ ~ 19%) (Figure 2Be).

There were similar trends in surface-emissions at the top as in deeper soils, though the changes were more evident in top layers compared to the bottom layers of mesocosms. Briefly, the addition of TEAs inhibited and enhanced the emission of CH₄ and CO₂, respectively, is similar as explained for mesocosm 1–5 in deeper layers (Supplementary Figure S1).

During lab-scale mesocosm experiments, gas emission varies from top to bottom, therefore, for DNA analysis the samples were taken at two different depths, from the top and near the bottom. A comparison was performed between different treatments based on the microbial diversity and community differences in mesocosms treated with different TEAs. After quality trimming and subsampling, the total number of high-quality reads generated by all the mesocosm samples was 116,231 (7,105–15,603 reads per sample, Supplementary Figure S2). As a result of taxonomic categorization at various levels, 28 phyla, 62 classes, 107 orders, 137 families, and 585 genera were identified.

Among all the 28 bacterial and archaeal phyla found in all the depths, bacterial phyla comprised ~90–97% of the total sequences. Among bacterial phyla, Proteobacteria, Acidobacteria, Chloroflexi, and Bacteroidetes were the most dominant phyla compared to others (Figure 3a). The most abundant phyla in the control group were Acidobacteria (11.76–33.13%), as well as in the HA treated group. Proteobacteria (18.3–44.3%) was the most abundant phyla in all the other treatments, followed by Acidobacteria. Soil-dwelling Proteobacteria play an essential role in the biogeochemical cycle. Chloroflexi (11.50–21.90%) and Bacteroidetes (5.05–20.11%) followed the same abundance pattern as Proteobacteria.

Other dominant phyla included Planctomycetes, Actinobacteria, Spirochetes, Nitrospirae, and Firmicutes, etc., whereas the rest of the phyla were less than 1% in their abundance. Additionally, the abundance pattern of different phyla varies after the addition of different external sources as well as with depth. For example, the abundance of Acidobacteria was decreased in all treatments after the addition of external sources, whereas that of Proteobacteria and other phyla showed relatively opposite trends.

Furthermore, the abundance of Acidobacteria was slightly higher in the top layers of mesocosms compared to the bottom layers. However, in HA treatment, it was higher in the bottom compared to the top layer. The abundance of other phyla also followed the same pattern as Acidobacteria except in HA-treated mesocosm. Our analysis also indicated the presence of 62 classes. The comparison of the dataset among these classes showed that Acidobacteriia was again the most abundant class in control as well as treatment groups, followed by Alphaproteobacteria and Ktenobacteria, Gammaproteobacteria, Bacteroidia, Deltaproteobacteria, Igvanibacteria, Planctomycetacia, Thermosulfovibrionia, and Thermoleophilia, etc., however, the dominance varies individually in every treatment (Figure 3b).

As shown in Figure 4a, the alpha diversity indices, including Ace1, Fisher, Shannon, Simpson, and Chao showed variable trends. However, the HA treatment in all the treatments showed a decrease in the richness and diversity of the microbial community. Statistical analysis showed that alpha diversity was significantly lower in the treatment with HA alone than in the control group (p < 0.05). Additionally, the depth also affects the diversity indices, and the diversity decreased in all the mesocosms with the depth except after SO₄²⁻ + 0.5 g HA addition.

Oxygen (O₂) has a higher oxidation capacity throughout natural ecosystems than redox ions such as SO₄²⁻, Fe³⁺, and quinones of humic acid, giving precedence to redox reactions. Consequently, peat exposed to oxygen undergoes a faster redox reaction and decomposes more quickly. O₂ is scarce in the lower layers of often flooded peatlands, and CO₂ production is minimal (Osaki et al., 2021). The redox process creates CO₂ as well as CH₄ when the peat layer is strictly anaerobic. The response of different microbial metabolic activities in peatlands usually determines the ratio of CO₂ and CH₄ (CO₂:CH₄) production, which should be stoichiometrically 1:1 in the absence of any other TEAs (Wilson et al., 2017). However, a couple of laboratory scale and in-situ experimental studies have previously reported that this ratio is not the same, and CO₂ production is enhanced by different environmental factors (Glaser et al., 2016). Therefore, It’s crucial to evaluate how different terminal electron acceptors (TEAs) influence CO2 and CH4 emissions, as studies exploring their combined effects are limited.

In the present study, the pilot scale mesocosms were used in the present study in order to simulate the natural peatland system as far as possible, with a larger incubator size close to natural status. It was observed that after the addition of SO₄²⁻, a Pungent smell like H₂S was felt in the first two mesocosms. It is important to note that this gas was not detected in deeper layers (results not shown here). The emission of H₂S might be due to organoclastic sulfate reduction. It was observed that the addition of different TEAs has mostly inhibited the emission of CH₄. However, the emission of CO₂ increased. After the addition of SO₄², the reduction in CH₄ production might be due to sulfur-reducing bacteria, which outcompete methanogens for peatland’s available carbon (Liu et al., 2020). Whereas the production of CO₂ was doubled in this mesocosm, these results were in accordance with previous studies, where sulfate’s addition decreased methanogenesis and increased CO₂ production (Estop-Aragonés et al., 2016; Minderlein and Blodau, 2010). The emission of CO₂ and CH₄ followed the same pattern in the third mesocosm, where sulfate was added in combination with HA.

The addition of HA had also impeded CH₄ emission and enhanced CO₂ production. One of the reasons for CH₄ emission reduction might be the inhibitory effect of HA on methanogens, and previously it was reported that the addition of HA could decrease CH₄ production by up to 89% (Khadem et al., 2017). Moreover, a significant part of humic acid’s function acts as an electron acceptor and constrains CH₄ emissions. HA has a negative charge and electron shuttling properties, so a high concentration of HA inside methanogens might also alter their electron transport system. In comparison, HA serves as an electron acceptor for the reducing equivalents transported out of cells through cell membranes of methanogens due to alteration in the electron transport system suppressing microbial growth (Khadem et al., 2017). The peatland soils are generally rich in organic matter and are deficient in inorganic electron acceptors. As a result, the presence of HA would augment the process of microbial respiration, leading to an increase in CO₂ production and deter CH₄ production by suppressing methanogenic activity (Keller and Takagi, 2013; Klüpfel et al., 2014).

Goethite is basically hydrated iron oxide, and its addition showed similar trends to other TEAs in the present experiment. The presence of alternative TEAs such as Fe³⁺ and SO₄²⁻ plays a vital role in deciding how the organic matter would be processed in a peatland. They are also reported to suppress methanogenesis, which is the primary form of anoxic mineralization in peatlands (Smemo and Yavitt, 2011). During water-logged conditions, the peatlands are O₂ deficient, and the iron-reducing microorganisms utilize Fe as an electron acceptor and organic matter as an electron donor, resulting in the production of more CO₂ (Kappler et al., 2021). Briefly, the addition of all the different TEAs increased CO₂ production to some extent, which might have a constrained impact on the microbial communities involved in the respiration or methanogenesis process resulting in CO₂ and CH₄ production. Additionally, the minor changes in CO₂ emission in the deeper layers can be related to the fact that microbial communities in the deeper layers are comparatively more resistant to external factors or nutrient deprivation (Bai et al., 2017a), the resistance against these changes results in the minor emission of CO₂ in deeper layers (Liu et al., 2016) of mesocosm. Briefly, the addition of TEAs had significantly affected the amount CO₂ and CH₄ emission. Moreover, after TEAs addition, the stoichiometric 1:1 ratio between CH₄ and CO₂ was changed completely, sometimes doubling CO₂ and/or suppressing CH₄ more than 80%.

The microbial community structure differences were based on adding different TEAs. Nevertheless, the leading bacterial phyla were the same in all the mesocosm, but their abundance varied in the different mesocosm, Proteobacteria was more abundant in treatment 4a compared to control or 5a treatment. Similarly, Bacteroidetes abundance was much less in 4a compared to 3a (Figure 3a). Similar results were reported by Potter et al. (2017) that the proportion of microbial communities, especially, Acidobacteria and Proteobacteria, is significantly affected by the depth of peatland soils. The variation in the microbial community structure along depth is also due to oxygen deprivation in deeper layers. When oxygen is present in ample amounts, the aerobic microorganisms capable of methane oxidation are active (as in surface layers). Moreover, the aerobic environment, along with substrate addition, profoundly affects the community structure (Tian et al., 2022). The effect can also be seen in the previous section, where the addition of substrate has significantly affected gas emission in surface layers compared to deeper layers.

Additionally, the abundance results are also familiar with previous studies with the same bacterial phyla in different types of peatlands (Luo et al., 2021). Three microbial groups, Bacteroidetes, Actinobacteria, and Verrucomicrobia, were previously reported to contribute most to soil organic carbon degradation in northern peatlands (Tveit et al., 2013). Compared to this study, Actinobacteria and Bacteroidetes comprised a higher percentage compared to Verrucomicrobia in the present study. According to the present study, in peatland bacterial communities, the effects of TEAs appear to be rather subtle. Further investigations are needed to determine how these TEAs affect bacterial populations in specific ways. Due to their reliance on primary nutritional substrates for growth and metabolism, the microbial communities on the upper layers of the soil are more diverse.

The diversity indices also showed that the addition of TEAs has a profound impact on the density and diversity of microbial community structure (Figure 4a). Shannon diversity was comparatively higher in the top layer compared to the bottom layer, except for goethite addition which showed the opposite trends and vice versa. This might be due to physico-chemical changes in the properties of soil that are not feasible for the survival of the microbial communities, and only some bacterial strains survived that were able to grow in harsh climatic conditions (Zhang et al., 2019).

According to the results, the main methanotrophs found in the present study were Methylomirabilales, Methylococcales, and uncultured_methanogenic_archaeon (Figure 4b). As the results showed, the most abundant of these were Methylomirabilales (Figure 4b). The microbial species in the order Methylomirabilales order can oxidize methane anaerobically as well as reduce nitrogen to dinitrogen (Ettwig et al., 2010). There is evidence that particulate methane monooxygenase (pMMO) encoding genes are present in the genomes of Methylomirabilales (Versantvoort et al., 2018). Methylococcales bacteria, primarily aerobic methane-oxidizing bacteria, are reported to be involved in AMO (anaerobic methane oxidation) (Cabrol et al., 2020). Briefly, an increase in microbial CH₄ consumption is most likely behind the observed reduction in CH₄ emissions under oxygen and other terminal electron acceptors.

The LEfSe analysis revealed the most abundant groups in 5 mesocosms were different from each other, While there is a limited amount of knowledge concerning the specific ecology of these groups in peatlands (Figure 5b). It is well documented that Thaumarchaeota species are found in soil samples worldwide, and their contribution to global biogeochemical cycles is crucial (Lin et al., 2015). In addition to having mixotrophic, chemoautotrophic, and heterotrophic lifestyles, Thaumarchaeota is capable of producing energy through diverse mechanisms. Peat biotopes are home to Thaumarchaeota where they produce living compounds by oxidizing organic fatty acids, amino acids, and inorganic ammonia. However, there is still uncertainty regarding the environmental role of Thaumarchaeota in peatlands. Proteobacteria, Bacteroidetes, and Chloroflexi, the most abundant groups in mesocosm 3, were also found in Northern peat lands. It has already been reported that a large proportion of peatlands throughout northern and tropical regions were dominated by Proteobacteria while Bacteroidetes and Chloroflexi constituted secondary dominant phyla (Finn et al., 2020). Rhodocyclaceaebacterium, Betaproteobacterials, and Holophagales were found in mesocosm 4. There are few available reports about detecting Rhodocyclaceaebacterium, Betaproteobacterials, and Holophagales, especially in peatlands. However, there has been evidence that Betaproteobacterials methylotrophs can utilize and assimilate ¹³C-methanol (Osaka et al., 2006). Whereas TRA3_20 (uncultured Betaproteobacteria) was abundant in mesocosm 5, and the information available about TRA3_20 is also limited.