Methanogens and Methanotrophs Show Nutrient-Dependent Community Assemblage Patterns Across Tropical Peatlands of the Pastaza-Marañón Basin, Peruvian Amazonia

Tropical peatlands are globally important carbon reservoirs that play a crucial role in fluxes of atmospheric greenhouse gases. Amazon peatlands are expected to be large source of atmospheric methane (CH4) emissions, however little is understood about the rates of CH4 flux or the microorganisms that mediate it in these environments. Here we studied a mineral nutrient gradient across peatlands in the Pastaza-Marañón Basin, the largest tropical peatland in South America, to describe CH4 fluxes and environmental factors that regulate species assemblages of methanogenic and methanotrophic microorganisms. Peatlands were grouped as minerotrophic, mixed and ombrotrophic categories by their general water source leading to different mineral nutrient content (rich, mixed and poor) quantified by trace elements abundance. Microbial communities clustered dependent on nutrient content (ANOSIM p < 0.001). Higher CH4 flux was associated with minerotrophic communities compared to the other categories. The most dominant methanogens and methanotrophs were represented by Methanobacteriaceae, and Methylocystaceae, respectively. Weighted network analysis demonstrated tight clustering of most methanogen families with minerotrophic-associated microbial families. Populations of Methylocystaceae were present across all peatlands. Null model testing for species assemblage patterns and species rank distributions confirmed non-random aggregations of Methylococcacae methanotroph and methanogen families (p < 0.05). We conclude that in studied amazon peatlands increasing mineral nutrient content provides favorable habitats for Methanobacteriaceae, while Methylocystaceae populations seem to broadly distribute independent of nutrient content.


INTRODUCTION
Peatlands are ecosystems where plant primary productivity exceeds organic matter decomposition, resulting in the accumulation of partially decomposed soil organic matter, acting as globally important organic carbon (OC) reserves Sjogersten et al., 2014;Page and Baird, 2016). While northern, low-temperature peatlands have been extensively characterized in terms of their microbial communities (Andersen et al., 2013;Serkebaeva et al., 2013;Zhou et al., 2017), little is understood about tropical peatlands, particularly for those in the Amazon basin. Amazon peatlands have been heavily understudied due to the lack of reports of their existence prior to 2008 (Lahteenoja et al., 2009), however over 80 550 km 2 containing 7.07 Gt OC have been reported so far, representing 18% of the surface area and 8.3% of the OC content of global tropical peatlands (Rieley and Page, 2016). Due to the high OC content, Amazon peatlands are expected to produce significant amounts of greenhouse gases (GHG) such as carbon dioxide (CO 2 ), methane (CH 4 ) and nitrous oxide (N 2 O). In fact, recent early reports on GHGs fluxes from a few Amazon peatlands show significant yet highly variable levels of CH 4 and N 2 O across sites, across seasons in a site, or across gradients within in a site (Teh et al., 2017;Winton et al., 2017). Regionally, it has been coarsely estimated that the Amazon Basin can emit 31.6 -41.1 Tg CH 4 year −1 , or approximately 7% of global CH 4 emissions (Wilson et al., 2016;Teh et al., 2017;Winton et al., 2017) but such an estimate is based on too few data points. To our knowledge there is a limited availability of GHG flux reports from Amazon peatlands and no study has addressed the composition and potential functional patterns of soil microbial communities and their association with CH 4 flux in Amazon peatlands.
To better understand the role of microbes participating in CH 4 flux it is important to evaluate the environmental controls (i.e., pH, nutrient content, others) as well as microbemicrobe interactions affecting the activity of CH 4 producers or consumers. Methane production is performed by methanogenic Archaea (methanogens) as the final stage of the organic matter decomposition cascade, where a broad range of heterotrophic organisms metabolize complex organic molecules releasing simpler methanogenic substrates (acetate, formate, methanol and CO 2 ) (Garcia et al., 2000). Consequently, methanogens are metabolically dependent on active heterotrophs in peatlands. Methanogens are metabolically and physiologically diverse, and are primarily grouped taxonomically in the Orders Methanopyrales, Methanococcales, Methanobacteriales, Methanomicrobiales, and Methanosarcinales (Garcia et al., 2000;Thauer et al., 2008). Conversely, CH 4 can also act as a carbon and energy source for specific microbial species known as methanotrophs. These predominantly include aerobic Gamma-and Alpha-proteobacteria of the Methylococcaceae and Methylocystaceae families, respectively (Hanson and Hanson, 1996), aerobic Verrucomicrobia (Dunfield et al., 2007), Archaea that anaerobically oxidize CH 4 via coupling with reduction of alternative terminal electron acceptors [like sulfate (Orphan et al., 2002;Orphan et al., 2009) nitrate (Haroon et al., 2013), Fe or Mn (Ettwig et al., 2016)], and NC10 bacteria that is capable to produce its own oxygen from nitrite and oxidize CH 4 in anaerobic environments (Ettwig et al., 2010).
Studies of the ecological diversity of methanogens in tropical peatlands globally are few and in some cases with contrasting results. In Peat Swamp Forest in Thailand, hydrogenotrophic members of the Class Methanomicrobia (Kanokratana et al., 2011) were found dominant, while an equivalent peatland in Malaysia putative methylotrophic Methanomassiliicoccales methanogens were identified as the dominant (Too et al., 2018). The first study is in agreement with well-characterized northern peatlands, in which several studies demonstrated that the Order Methanomicrobiales and, to lesser extent, hydrogenotrophic members of the Methanobacteriales and acetoclastic Methanosarcinales are important for CH 4 production (Juottonen et al., 2005;Cadillo-Quiroz et al., 2006;Steinberg and Regan, 2008;Lin et al., 2012). However, the result of the second study indicates that there is significant methanogen's (and plausible methanotrophs) variation that is yet to be accounted in tropical peatlands. To our knowledge, the diversity of methanotrophs in Amazon peatlands has not been documented. Contrasting reports have been observed in other tropical regions where a tropical peatland in Indonesia showed a dominance of Methylomonas sp. from the Methylococcaceae (Arai et al., 2014), while a Malaysian swamp forest peatland showed Methylocystaceae (Too et al., 2018) as the dominant. In northern peatlands, the dominant methanotrophic bacteria tend to befree living or Sphagnum-moss symbionts of the family Methylocystaceae (Chen et al., 2008;Dedysh, 2009;Wieczorek et al., 2011;Kip et al., 2012).
Naturally existing nutrient gradients have historically been used by ecologists to identify environmental parameters that govern observable species assemblages (Tilman and Wedin, 1991;Keddy, 1992;McGill et al., 2006). Geochemical and nutrient gradients in northern and tropical peatlands are typically associated to the dominant water source described as follows: a) minerotrophic -groundwater fed fens or swamps are relatively nutrient-rich; b) ombrotrophic -rain fed bogs are nutrient-poor; and c) mixed water source poor fens or swamps have intermediate nutrient levels (Juottonen et al., 2005;Cadillo-Quiroz et al., 2006;Martini et al., 2007;Lahteenoja and Page, 2011). The biodiversity of methanogenic species has previously been shown to change across peatland geochemical gradients (Juottonen et al., 2005;Cadillo-Quiroz et al., 2006;Lin et al., 2012), but attempts to describe a set of community-assembly associations for tropical peatland methanogens has yet to be considered.
Here we describe GHG flux variation and the first accounts of prokaryotic community compositions across a tropical peatland nutrient gradient in the Pastaza-Marañón Basin, in the western Amazonia, Peru. The Pastaza-Marañón Basin holds the most extensive tropical peatland distribution in South America (Lahteenoja et al., 2012;Draper et al., 2014). A rich diversity of ecosystem types, based on nutrient content, exists here due to mineral deposition from dynamic hydrological activity of river channels originating in the Andes, or rainfed, or mix water source-fed peat formation (Lahteenoja et al., 2012). The soils included in this study thus originated from peatlands showing minerotrophic (Charo, San Roque, and Buena Vista), mixed (Nueva York and Quistococha) and ombrotrophic (Miraflores and San Jorge) conditions. Sites were assigned to categories as in previous reports (Lahteenoja et al., 2009;Teh et al., 2017) based on soil analysis and site characteristics, plus an additional assessment based on soil mineral content of sodium, magnesium, phosphorous, sulfate, potassium, calcium, manganese, iron, nickel, copper, and zinc ( Table 1) as quantifiable geochemical proxies to the nutrient conditions provided by ground water-fed (minerotrophic) or rain-fed (ombrotrophic) or mixed water sources in peatlands (Whitfield et al., 2009).
We hypothesized the following: 1) nutrient concentration variation, measured as elemental composition via inductively coupled plasma-mass spectrometry (ICP-MS), among peatlands with contrasting geochemistry influenced by their source of water shapes the overall microbial community composition; and 2) such variations in nutrient concentration can be associated with CH 4 fluxes as well as explain discrete (non-random) assemblages of methanogen and methanotroph families in the study sites.

Site Description and Field Sampling
Seven study sites were selected in the Pastaza-Marañon basin, in the Peruvian Western Amazon region (Table 1). Petland information regarding vegetation, geochemistry and peat depths have been described previously for Quistococha , San Jorge and Buena Vista (Kelly et al., 2013), San Roque, Miraflores, and Charo (Lahteenoja et al., 2009) and Nueva York (Lahteenoja and Page, 2011). Annual climatological data was obtained from the database described in Espinoza Villar et al. (2009). Samples were collected in December of 2011 and December of 2012 near the end of the "dry" season or a transitional period to rainy conditions. On-site measurements included pH, air and soil temperature using a HI 99121N meter (Hanna Instruments, RI, United States), and fixed soil volume collection (400 cm 3 ) for bulk density analysis (Chambers et al., 2010). A 100-m long transect near the center of the peatland with north to south direction was done using center location and soil characterization from previous studies (Lahteenoja et al., 2009;Lahteenoja et al., 2012) to better represent peat geochemical conditions. In transect two pair of independent points (∼10 m apart in east to west direction) were placed every 50 m (6 sampling points per site with variable distances from 20 to ∼160 m). Soil samples were collected from the center of sampling points shortly after flux measurements were completed and when relevant at two depths (0-15 and 15-30 cm) (12 samples per site). For each of the six sampling points in the transect, a 15 cm deep and 5 cm diameter mini soil core was collected in sterile bags and mixed, then subsamples were taken for corresponding analysis: triplicated 0.5 gr were aseptically placed in sterile tubes and mixed with 1 ml of MoBio Lifeguard Soil stabilization solution (MoBio, United States), ∼50 g of soil was sealed in a bag for physical or chemical analysis. All samples were stored in refrigerated conditions within 24-48 h of collection during field expedition, and subsequently frozen and stored at −20C (for chemical) or −80C (for DNA) until appropriate analysis. GHG fluxes were measured at each of the six soil sampling points at the soil surface. All flux measurements were conducted in a time window of ±3 h from noon (9 am-3 pm) using metal static vented soil chambers (Keller et al., 2005) enclosing a 0.045 m 2 area, where collars with 0.11 ± 0 0.09 ± 0 0.11 ± 0 0.11 ± 0 0.09 ± 0 0.13 ± 0 0.1 ± 0 pH 5.9 ± 0.2 a 5.6 ± 0.1 a 5.9 ± 0.1 a 3.7 ± 0.1 b 3.2 ± 0 b 2.9 ± 0.1 c 2.5 ± 0 c Alpha Diversity (H') 3.78 ± 0.23 a 3.84 ± 0.04 a 4.07 ± 0.14 a 3.16 ± 0.03 b 3.43 ± 0.13 b 3.52 ± 0.2 c 3.55 ± 0. airproof rubber seals in their base were affix to ground with 30 cm stakes at least 30 min before measurements. After the chamber reached pressure equilibration, air samples were taken with plastic syringes at 10, 20, 30, and 40 min and transferred into nitrogen-flushed, pre-evacuated 6 mL glass vials sealed with butyl stoppers (Bellco, NJ, United States).

Soil and Gas Chemistry
A Quadrupole ICP-MS (Thermo Electron X series with CCT, Fisher Scientific, United States) was used to analyze soil elemental composition in the Metals, Environmental and Terrestrial Analytical Laboratory, part of the Chemical and Environmental Characterization Core Facilities at Arizona State University. Between 1.5 -2.5 g of soil was dried at 80 • C for 72 h, digested in hydrofluoric, nitric and hydrochloric acid (12 and 3 h steps), followed by microwave heat digestion (30 min) and final evaporation. ICP-MS analysis of diluted samples in nitric acid along with standards for all elements evaluated done as described elsewhere (Potts, 1992). GHG were analyzed on a modified SRI Greenhouse Gas GC (SRI Instruments, CA, United States) with dual valco valves for separate gas streams for a Flame Ionization Detector (FID) for CO 2 and CH 4 and an Electron Conductivity Detector (ECD) for N 2 O. The flux rates were calculated from a linear change in trace gas concentration over time with reference to the internal volume of the chamber and the soil area covered (Hutchinson and Livingston, 1993).

Molecular Analyses
MoBio Lifeguard soil stabilization solution was removed by centrifugation and supernatant removal. Then, genomic DNA was extracted from 0.5 g of soil (wet weight) using a MoBio UltraClean Soil DNA isolation kit (MoBio, CA, United States), with the addition of aluminum ammonium sulfate (5 mM) prior to bead beating. Extracted DNA was quantified with a Qubit fluorometer (Invitrogen, MA, United States) assessed via gel electrophoresis. PCR amplifications were carried out using the 515F (5 -GTGYCAGCMGCCGCGGTA) (Baker et al., 2003) and 909R (5 -CCCCGYCAATTCMTTTRAGT) (Wang and Qian, 2009;Tamaki et al., 2011) primers that cover the V4 and V5 regions of the 16S rRNA gene of Bacteria and Archaea. Unique barcodes of 6 -10 base pairs were attached to the 5 end of primers. The PCR was carried out with 1. Paired-end sequences were quality filtered as described elsewhere (Bokulich et al., 2013). Paired-end reads were demultiplexed with an in house script. Downstream analyses for OTU annotation was performed with the Quantitative Insights into Microbial Ecology Two (Qiime2) pipeline (Caporaso et al., 2010). Chimeric sequences were removed with Uchime (Edgar et al., 2011). Closed-reference operational taxonomic unit (OTU) picking at 97% identity was used in relation to the Silva 128 database (Quast et al., 2013). Four 16S rRNA gene analyses performed on Northern peatlands from Sphagnum-dominated bog at Grand Rapids, Minnesota, United States (SRA IDs SRR6514003, SRR6514014, SRR6514002 and SRR6514033), were included in Qiime2 analyses for coarse comparative purposes against Amazon peatlands. Sequences of this study were deposited under SRA Bioproject ID PRJNA501909.

Statistical Analyses
All statistical analyses were performed in R v 3.4.1 (R Development Core Team, 2013). Analysis of variance (ANOVA) was used to compare soil physicochemical variables followed by least-significant difference (LSD) testing (p < 0.05) across defined site categories (minerotrophic, mixed, ombrotrophic). Also, a Permutational multivariate analysis of variance (PERMANOVA) analysis was completed on soil physicochemical variables to test whether the tree site categories will better explained soil conditions than any two-way combination (R2, p < 0.001). ANOVA was also used to compare CO 2 gas fluxes, and Kruskal Wallis for CH 4 and N 2 O gas fluxes. Fluxes were visualized as a stacked barplot, with LSD (p < 0.05) using the defined categories plus a Northern peatland one. Alpha diversity (H ) and evenness, as well as ordination of microbial communities with non-metric multidimensional scaling (NMDS) and fitting of climatological variables, environmental variables and GHG fluxes was performed as described elsewhere (Oksanen et al., 2013). Chloroplast and mitochondrial OTUs were identified and removed from the family level matrix. Weighted network analysis was performed on OTU covariance at the family level, with edge weights filtered at the p = 0.05 level determined by Spearman correlation (Csardi and Nepusz, 2006). The Spinglass technique was used to calculate modularity of the network (Reichardt and Bornholdt, 2006). Networks for minerotrophic, mixed and ombrotrophic sites were generated as above, and the frequency of methanogen and methanotroph connections to other nodes were summed, with corresponding nodes reported at the phylum level. Kruskal-Wallis tests were performed for methanogen and methanotroph family abundance with site classification as a categorical variable. Rank Abundance calculations were performed as described elsewhere (Kindt and Coe, 2005). A lognormal model (Equation 1) and a geometric model (Equation  2) were fit to log family abundances over ranks with nonlinear least squares, and model fits were compared via Akaike Information Criteria (AIC) (R Development Core Team, 2013).
The log-normal model is as follows: Where S(R) is abundance of the R th rank, S o is the abundance of the first rank, a is an inverse measure of the distribution width, and R corresponds to the R th rank (Ludwig and Reynolds, 1988). The geometric model is as follows: Frontiers in Microbiology | www.frontiersin.org Where S(R) is abundance of the R th rank, S 0 is the abundance of the first rank, p is the proportional decrease in abundance and R corresponds to the R th rank (Dolan et al., 2007). Finally, community assembly patterns of methanogen and methanotroph families were conducted using variance ratio (V-ratio) (Gotelli, 2000). Null model testing of 1000 random permutations of the family matrix was performed with fixed sums and equiprobable distribution across sites (Gotelli, 2000).

Soil Properties and Gas Fluxes
In soil properties measured across study sites (Table 1), a general trend was observed across the geochemical and nutrient gradient (minerotrophic to ombrotrophic) of decreasing pH and soil elemental nutrient concentration, with the ombrotrophic sites demonstrating low pH (<3) (LSD, p < 0.05). No trends were noted for mean annual temperature (MAT), mean annual precipitation (MAP), H or evenness. The concentration of elemental soil nutrients consistently supported the separation of the minerotrophic soil category from mixed and ombrotrophic (LSD, p = 005); but not of the last two due to replicate variability despite lower means in ombrotrophic soils (Table 1).
However, multivariate analysis of mineral nutrient, pH and bulk soil measurements (PERMANOVA, p < 0.001), showed the separation of mixed and ombrotrophic categories (in addition to minerotrophic), has the highest R2 fit (33%) than any pairwise (up to 26%) and thus better represent the variation across sites supporting our use of tree site categories for Amazon peatland as in this study. For gas flux comparisons (Figure 1), carbon dioxide fluxes differed among all peatlands (ANOVA, p = 0.04), but these differences were independent of the elemental nutrient gradient (ANOVA, p = 0.09). Similarly, N 2 O flux had high variation among sites but the recorded levels (∼0.5-100 ug N m −2 h −1 ) were highly specific to individual peatlands (Kruskal-Wallis, p = < 0.001), independent of the nutrient gradient (Kruskal-Wallis, p = 0.05). Conversely, CH 4 flux showed a dependence on the elemental nutrient gradient (Kruskal-Wallis, p = 0.03) and independent of individual peatlands (Kruskal-Wallis, p = 0.1). San Roque and Buena Vista (minerotrophic), which had an order of magnitude greater CH 4 flux than the other sites, drove this. Charo, also a minerotrophic site, had an unusually low CH 4 flux compared to the San Roque and Buena Vista sites. CH 4 flux remained relatively low for all mixed and ombrotrophic sites.

Relationships Between Communities Across the Nutrient Gradient
Rarefaction of OTUs at the 97% sequence homology and sequencing depth of 50 000 sequences per sample show sufficient sampling was achieved (Supplementary Figure 1) and from the 32 Bacterial and Archaeal phyla with a relative abundance of greater than 0.005% across the soils over a third (mark with asterisks in Figure 2) differed significantly in their frequency across site's categories (LSD, p < 0.05). Proteobacteria were the dominant phylum fraction FIGURE 1 | Greenhouse gas fluxes of CO 2 (mg C m 2 h −1 ), CH 4 (µg C m 2 h −1 ), and N 2 O (µg N m 2 h −1 ) across the Pastaza-Marañón Basin. Box plots show the mean and quartiles of six independent measurements. All flux measurements were evaluated for statistically significant differences among categories with a Least Significant Difference post hoc test where different letters (a, b) indicate significant differences across site categories (minerotrophic, mixed and ombrotrophic) at the p = 0.05 level.
FIGURE 2 | Stacked barplot of relative abundances (%) of Bacterial and Archaeal phyla across the Pastaza-Marañón Basin. Significant differences between phyla in the minerotrophic, mixed, ombrotrophic and Northern peatland groups are shown as: (**) differences between two groups; (***) differences between three groups; and (****) all groups differed. Phyla above 0.005% are shown. Relative abundances are the mean of six independent measurements for Amazon soils, and four for the Northern peatland soil.
Environmental variables showed pH and elemental nutrient concentration separated all microbial communities along the first NMDS axis. Of the GHG, only CH 4 flux showed an association with minerotrophic communities.

Methanogen and Methanotroph Relationships Across the Nutrient Gradient
The values of relative abundances of methanogen and methanotroph families ( and Methanocellaceae all decreased with decreasing nutrient content (p < 0.005). Members of GOM Arc I (previously termed ANME-2d but not shown to be methanotrophic) were highly abundant in Charo but otherwise not affected by geochemical and nutrient conditions. Only members of the unclassified Methanomicrobia preferred oligotrophic and ombrotrophic peats to rich minerotrophic peats, however these organisms were in very low relative abundance (0.0006%). The OTUs assigned as unclassified Methanomicrobia, Methanomicrobiales or Methanosarcinales are based on environmental, uncultured 16S rRNA genes in the Silva 128 database that lack a close relative to a characterized methanogen family. Covariance network analysis of microbial families across peats showed three separate modules ( Figure 4A). The abundant methanogens in the minerotrophic sites all clustered in Module  Table 1) and gas fluxes (CO 2 , CH 4 , and N 2 O in units as per Figure 1) as loadings. The dotted circles represent 95% confidence intervals of minerotrophic, mixed, and ombrotrophic communities.
One (blue). The Unclassified Methanomicrobia grouped with Module Three (green). Methanotrophic Methylococcacae, only identified in the minerotrophic soils, grouped with Module One. The more evenly dispersed Methylocystaceae grouped with Module Two (red), which was intermediate between Modules One and Three. The analysis of individual networks for minerotrophic, mixed and ombrotrophic communities (Figures 4B-D), showed that the total number of nodes (i.e., taxon families) decreased along the gradient (totals of 495, 327, and 301 in the minerotrophic, mixed and ombrotrophic networks, respectively). The variance to mean ratio (VMR) of node connectivity also decreased along the gradient as the number of highly connected, central nodes dispersed (VMRs of 0.81, 0.66, and 0.54 for minerotrophic, mixed and ombrotrophic networks, respectively). Finally, mean path length also decreased as nodes became more dispersed (lengths of 5.88, 2.95, and 2.81 for minerotrophic, mixed and ombrotrophic networks, respectively). A heatmap analysis of the connection frequency between methanogen/methanotroph nodes to significantly covarying nodes in individual networks (Figure 5), shows that methanogen nodes were most often connected to Actinobacteria, Chloroflexi, Firmicutes, Planctomycetes, Alpha-and Delta-Proteobacteria across all categories, and to Bacteroidetes, Thaumarchaeota and Verrucomicrobia and other Euryarchaeota to a lesser degree. Connections to Delta-Proteobacteria were primarily to sulfate reducers (Desulfaculaceae, Desulfobacteraceae, Desulfurellaceae, Desulfovibrionaceae) and the family Syntrophorhabdaceae (Supplementary Table 2). Connections to Alpha-Proteobacteria were between methanogens and methanotrophic Methylocystaceae, Rhizobiaceae, and Sphingomonads (Supplementary Table 2). Methylocystaceae methanotrophs were predominantly connected to varied Alpha-, Beta-, Delta-and Gamma-Proteobacteria. Methylococcaceae showed sparse connections to other Gamma-Proteobacteria, Saccharibacteria, Spirochaetae, Tectomicrobia, TM6, and Verrucomicrobia in the minerotrophic network (Figure 5).
Rank abundance curves of methanogen and methanotroph families ( Figure 6A) and AIC model comparisons demonstrated geometric distributions (i.e., a poorly structured community) for minerotrophic (blue) (R 2 = 0.95) and mixed (red) (R 2 = 0.93) communities, whereas the ombrotrophic community had a lognormal distribution (green) (R 2 = 0.97) indicative of a highly structured community. Finally, the V-ratio of co-occurrence of methanotroph and methanogen families was 1.6131, which was significantly greater than the upper 95% tail of 1000 permutations of the family matrix (mean of 0.99, standard deviation of 0.04, Figure 6B.).

Variation in CH 4 Flux Across the Nutrient Gradient
Methane fluxes from tropical wetlands can vary greatly, from atmospheric uptake to 30 mg CH 4 -C m 2 h −1 , indicative of the underlying complexities that contribute to CH 4 flux (Bartlett and Harriss, 1993;Sjogersten et al., 2014;Zhu et al., 2015). Both flight-derived measurements above the Amazon Basin and ground-level field-scale measurements agree that regional differences arising from vegetation type, hydrology and soil pH are important contributors to CH 4 flux (Wilson et al., 2016;Teh et al., 2017;Winton et al., 2017). Consistently greater fluxes occur from high-productivity vegetation and nutrient rich soils types that experience seasonal flooding (minerotrophic), compared to rainfall-fed (ombrotrophic) poor and low-productivity vegetation types, which emphasizes the complex interactions between seasonality, hydrology and biological processes in the Amazon Basin (Lahteenoja et al., 2009;Teh et al., 2017).
The measurements recorded here, of 0.3 to 2.5 mg CH 4 -C m 2 h −1 , are in similar range to measurements previously collected in diverse wetlands and flooded areas in the Amazon basin, 0.08 to 6 mg CH 4 -C m 2 h −1 (Bartlett et al., 1988;Richey et al., 1988;Bartlett et al., 1990;Miller et al., 2007;Belger et al., 2011). However, CH 4 values in this study are in a lower range than those observed in the only other multisite measurements completed in peatlands of Pastaza-Marañón Basin for their dry season measurements (1.8-31 mg CH 4 -C m 2 h −1 ) (Teh et al., 2017). Nevertheless, the opposite was observed in N 2 O flux rates in comparison to that same study (0.5-100 vs. ∼0.03-0.87 ug N 2 O-N m 2 h −1 respectively) (Teh et al., 2017). The observed differences between studies are likely explained by the different sampling design, inclusion or not of temporal component in either study, and methodological differences in gas collection indicating further studies, common methodologies and common benchmarks are still needed to better understand CH 4 flux variation in the Amazon and develop more accurate regional estimates. Our study did not  aim to record seasonal variation which can significantly change CH 4 flux in the Pastaza-Marañón peatlands (Teh et al., 2017), but focused on differences across sites in the same season with flux sampling during similar times of day (±3 h around midday) to reduce potential effects of day fluctuations. We found a strong relationship between CH 4 flux, nutrient content and pH (LDS test Table 1, and Kruskal-Wallis, p = 0.03) similar to general findings of other studies in both Northern and tropical peatlands (Bartlett and Harriss, 1993;Lai, 2009;Ye et al., 2012;Teh et al., 2017;Winton et al., 2017). The minerotrophic Charo site proved an exception to this, as CH 4 fluxes from this site mirrored the mixed and ombrotrophic sites. High concentrations of Fe (356 ± 90 g kg −1 dry soil) and Mn (411 ± 73.8 mg kg −1 dry soil) were present at Charo. In fact, supplementary water chemistry assessments also showed that the concentration of sulfate was also high at Charo compared to the other peatland sites (16.4 ± 9.3 mg L −1 versus 0.1 -1.67 mg L −1 , Supplementary Table 1). As microbial Fe (III), Mn and sulfate reduction are more energetically favorable than methanogenesis (Dubey, 2005) the low CH 4 flux observed at Charo is explained by the occurrence of microbial respiration based on alternative terminal electron acceptor for anaerobic respiration (like Fe, Mn, SO 4 ) which are more energetically favorable and outcompete methanogens for several substrates (acetate, formate, etc.) (Bridgham et al., 2013). This was not observed in the other two minerotrophic sites and underscore the importance of evaluating the geochemistry of amazon peatlands as a control of methanogenesis.

Composition of Microbial Phyla Across the Nutrient Gradient
In terms of overall community composition, the dominant Proteobacteria decreased from minerotrophic to ombrotrophic sites, counter mirrored by an increase in the fraction of Acidobacteria in the same gradient. A dominance of Proteobacteria in near-neutral, minerotrophic Northern (Zhou et al., 2017) and poor (pH 3.8) tropical peatlands has been reported in previous community surveys (Espenberg et al., 2018). Similarly, a more dramatic shift in the community favoring Acidobacteria has been found in acidic, nutrientpoor ombrotrophic Northern and tropical peatlands (Jackson et al., 2009;Lin et al., 2012;Lin et al., 2014), and is apparent in the Northern peatland included here for general comparison (pH ∼ 3.5). Bacteroidetes, Firmicutes, Actinobacteria and Chloroflexi represent secondarily dominant phyla in Northern and tropical peatlands, as observed here. Unique to this study was to note relatively high abundances of Bathyarchaeota (up to 8% in the mixed) and Thaumarchaeota (up to 5% in the ombrotrophic peats). As the Bathyarchaeota have only recently been reclassified and incorporated into 16S rRNA gene reference databases, their distribution in terrestrial environments is still not yet fully accounted, but appear to dominate Archaeal communities in some peat, freshwater wetland and mangrove sediments (Xiang et al., 2017). In the comparison to reads from a Northern peatland site incorporated here, Bathyarchaeota were detectable at a mean relative abundance of 0.006%, suggesting a preference for tropical peatlands. However, since Northern peatlands have been shown to hold significant variation in microbial diversity, a more comprehensive comparison against more sites is required to better assess Bathyarchaeota or other group's preferences. Thaumarchaeota have recently become of great interest as certain isolates are involved in aerobic ammonia oxidation (Brochier-Armanet et al., 2012). These Archaea have also been shown to form macroscopic filaments attached to surfaces in tropical mangrove sediments, that support symbiotic growth of sulfur reducing Proteobacteria (Muller et al., 2010). The environmental role of Bathyarchaeota and Thaumarchaeota in these tropical peatlands remains to be determined.

Methanogen and Methanotroph Family Assemblages Across the Gradient
Methane flux from peatlands is a combined effect of in situ production versus consumption due to the activity of methanogens and methanotrophs, respectively. Methanogenic activity is further constrained by the supply of methanogenic substrates from heterotrophic activity or competition for methanogenic substrate consumption when alternative terminal electron acceptors are available in the ecosystem (Garcia et al., 2000). These two factors make methanogens sensitive to nutrient conditions affecting heterotrophs or competitors, while aerobic methanotrophs can be effected by levels of CH 4 or competition for O 2 as a shared electron acceptor with other aerobic heterotrophs.
Here we use the relative abundance of a taxonomic unit as a proxy for an organism's fitness under certain conditions. Relative abundance is a reflection of the capacity of an organism to out-compete others for control of a limiting resource (Pielou, 1975). Such competition gives rise to highly structured, non-uniform species assemblages whereby each species can be considered as occupying a distinct niche (Sugihara, 1980). We sought to investigate how the environment FIGURE 5 | Heatmap of connection frequencies between methanogen and methanotroph nodes to significantly covarying nodes (p < 0.05) from minerotrophic (blue), mixed (red), and ombrotrophic (green) networks. Covarying nodes on the y axis are reported at the Phylum level. Scaled frequencies range from zero (dark blue) to 0.7 (yellow). Units are Hellinger transformed discrete counts of node connections. and potential microbial taxa interactions affect methanogenic and methanotrophic communities and whether changes in community assembly were associated with CH 4 flux. While de novo genome assemblies suggest that members of the Bathyarchaeota may be capable of methanogenesis (Evans et al., 2015), other studies have demonstrated acetogenic (He et al., 2016) and heterotrophic (Maus et al., 2018;Yu et al., 2018) metabolism by this phyla, and so this taxonomic group was not included in ecological analyses as methanogens. Furthermore, while GOM Arc I were initially termed ANME-2d, this group is neither monophyletic with other ANME groups nor has it been shown to consume CH 4 (Lloyd et al., 2006) and so GOM Arc I were considered as methanogens here.
Interestingly, optimal modularity of a network of all peatlands in study agreed with their a priori categorization into three categories. By separating this network based on water geochemistry/nutrient categories and constructing sub-networks, it was clear that connectivity within networks decreased as nodes generally became less central with decreasing mineral levels from minerotrophic to ombrotrophic. From these sub-networks it was possible to gauge the frequency of connections between methanogens/methanotrophs and other phyla. Methanogens were most commonly associated with anaerobic sulfate reducing Delta-Proteobacteria, methanotrophic Alpha-Proteobacterial Methylocystaceae, Actinobacteria, Firmicutes and Chloroflexi. Actinobacteria and Firmicutes include many primary and secondary fermenters (including known syntrophs) sustaining methanogens as shown in northern peatlands Wust et al., 2009). The positive covariance between competitive methanogens and sulfate reducers was surprising. However, sulfate reducing Delta-Proteobacteria have been shown to be predominant members of peatland microflora previously (Morales et al., 2006;Lin et al., 2014;Zhou et al., 2017), while it has also been observed that certain hydrogenotrophic methanogens, such as Methanobacteriaceae and Methanomicrobia that dominated these systems, can coexist in vitro and generate CH 4 in a stable manner while in co-culture with sulfate reducers (Chen et al., 2018). Indeed, as sulfate reducers can utilize propionate, butyrate, lactate and acetate as alternative electron donors to hydrogen, the potential for coexistence between sulfate reducers and methanogens has been noted for some time (Bryant et al., 1977;Muyzer and Stams, 2008). A relationship between methanogens and methanotrophic Methylocystaceae was not surprising, although the absence of an interaction with Methylococcaceae was. This may be a result of differences in methanotroph ecophysiology, explained in greater detail below. Actinobacteria play a vital role in the catabolism of plant-derived organic matter, important for the initial stages of decomposition (Goodfellow and Williams, 1983;Bentley et al., 2002;Berlemont and Martiny, 2013). The phylum Firmicutes includes many well-characterized heterotrophs responsible for producing partially fermented organic compounds suitable for methanogens, and have been noted to positively covary with methanogens when enriched in batch culture (Sierocinski et al., 2018). In fact, primary and secondary fermenters in Actinobacteria and Firmicutes have been shown sustaining methanogenic metabolism in northern peatland soils Wust et al., 2009). Finally, while Chloroflexi are less understood, in silico analyses suggest that this phyla is also involved in the anaerobic degradation of plant-derived organic matter (Hug et al., 2013). Thus, these taxa in consort with methanogens may play an important role in the synergistic degradation of plant-derived organic matter and CH 4 cycling in these peatlands. The methanotrophs primarily demonstrated associations with other Proteobacteria, particularly with the more dominant and evenly dispersed Methylocystaceae. This may imply some level of metabolic codependence between Proteobacterial taxa (Morris et al., 2012;Sachs and Hollowell, 2012), which has been demonstrated in co-cultures of methanotrophs and heterotrophic Proteobacteria (Iguchi et al., 2011;Ho et al., 2014;Jeong and Kim, 2018). Also, methanotrophs produce methanol as part of their metabolism and it has been shown they could accumulate it externally upon inhibition of their methanol dehydrogenase (Patel et al., 2016), and putatively could support methylotrophic Proteobacteria. Alternatively, this effect may simply be a product of strong covariance between highly dominant Proteobacteria.
The modularity visualized with the network were supported further as either geometric (minerotrophic and mixed) or lognormal (ombrotrophic) rank distributions ( Figure 6A). A lognormal rank distribution is reflective of highly structured, hierarchical communities whereby taxonomic units are separated based on niche differentiation (Preston, 1962) whereas a geometric distribution is reflective of poorly even distributions whereby few taxa of the highest ranks are highly competitive and monopolize available resources at the expense of others (Ludwig and Reynolds, 1988). Also, the null hypothesis of an absence of structured assemblages was rejected based on the V-ratio of methanogen and methanotroph family distributions (p < 0.05, Figure 6B). A V-ratio significantly greater than one indicates aggregation of taxa assemblages (Gotelli, 2000), which was driven by aggregations of methanogen families and Type I Methylococcaceae methanotrophs in the minerotrophic systems. Taken together, all these evidences support a structured, niche-dependent ecology of methanogens and Type I Methylococcaceae in the Pastaza-Marañón Basin.

Ecophysiology of Methanogen and Methanotroph Families
It has been shown previously that methanogen species richness is higher in minerotrophic versus ombrotrophic northern peatlands (Juottonen et al., 2005;Lin et al., 2012). There are likely multiple underlying reasons for this, and for the structured community assemblages in the Pastaza-Marañón Basin. Plant net primary productivity is highest in minerotrophic peatlands relative to others, and higher quantities of plant root exudates can provide potential methanogenic substrates to communities (Whiting and Chanton, 1993). The ombrotrophic peats of this study (Miraflores and San Jorge) are dominated by 'pole forest' species (Draper et al., 2014) which are expected to provide more phenolic-like compounds to the peat. Additions of Fe, Ni, Co, and Na have an immediate, positive effect on CH 4 production in peatland soils, as these metals can be limiting cofactors for enzymes involved in microbial fermentation and/or methanogenesis (Basiliko and Yavitt, 2001). Finally, methanogenic activity is highly sensitive to pH and both hydrogenotrophic and acetoclastic methanogenesis decreases if the pH is less than 4.5 (Ye et al., 2012).
While the majority of methanogens favored the minerotrophic peats, the abundances of methanogen families differed. The most dominant families in the minerotrophic peatlands were Methanobacteriaceae and Methanosataceae. Methanobacteriaceae are hydrogenotrophic methanogens that typically dominate minerotrophic northern peats (Juottonen et al., 2005;Steinberg and Regan, 2008;Lin et al., 2012). Hydrogenotrophic methanogenesis is more energetically favorable than acetoclastic methanogenesis (−131 vs. −31 G • ' kJ mol −1 CH 4 , respectively) (Garcia et al., 2000;Thauer et al., 2008). This could explain the dominance of the methanogen population by Methanobacteriaceae relative to acetoclastic Methanosataceae. Simultaneously, the greater abundance of acetoclastic Methanosataceae relative to other hydrogenotrophs could be as it occupies a separate, substrate-dependent niche to the most successful hydrogenotroph, Methanobacteriaceae.
The distributions of Methanobacteriaceae and Methanosataceae were primarily constrained to the minerotrophic peatlands. In the mixed and ombrotrophic environments, the dominant hydrogenotrophic and acetoclastic methanogens were Unclassified Methanomicrobia and Unclassified Methanosarcinales, respectively. The Class Methanomicrobia includes a wide range of Orders, such as Methanocellales, Methanomicrobiales and Methanosarcinales. That these OTUs could not be classified to the Family level suggests there are novel, uncultured methanogens in these mixed and ombrotrophic tropical peatlands. A number of traits in these groups may be favored in such environments. For example, Methanocellales may be uniquely equipped for nutrient-poor, low productivity environments, as they effectively utilize hydrogen under partial pressures reaching the thermodynamic limit of methanogenesis, and dinitrogenfixation pathways have been identified in several genomes (Lyu and Lu, 2015). Methanosarcinales are metabolically diverse and demonstrate greater growth yield efficiency than hydrogenotrophic methanogens (7.2 vs. up to 3 g biomass mol −1 CH 4 , respectively) (Thauer et al., 2008).
Similarly to methanogens, methanotroph families demonstrated non-uniform distributions. Methylocystaceae dominance is consistent with observations in northern peatlands (Wieczorek et al., 2011;Andersen et al., 2013) and a forested swamp peatland in Malaysia (Too et al., 2018), and its abundance was unaffected by nutrient content or the relatively low CH 4 fluxes from oligotrophic and ombrotrophic peats. Unlike Methylococcacae, Methylocystaceae can utilize carbon sources other than CH 4 such as acetate, oxidize CH 4 at relatively low oxygen concentrations, and can survive periods of anoxia via storing carbon intracellularly as poly beta-hydroxybutyrate (Vecherskaya et al., 2001;Semrau et al., 2011). As with Methanocellales, these traits allow Methylocystaceae to dominate nutrient-poor environments. Conversely, Methylococcaceae methanotrophs tend to occupy nutrient-rich environmental niches and are dependent on relatively higher concentrations of CH 4 (Whittenbury et al., 1970). ANME or NC10 methanotrophs were not identified in this study at analyzed peat depths (0 -30 cm), which could be explained by absence or low density at this water-saturated shallow depths. One study targeting deeper layers (80 -135 cm) identified anaerobic methanotrophs in a northern peatland (Zhu et al., 2012), while others identified similar groups in shallow layers (top 15 cm) in wetlands (Segarra et al., 2015;Valenzuela et al., 2017). Since this study, did not target deeper layers to account for anaerobic methanotrophs then the role of nutrient concentration on below-surface anaerobic methanotroph communities, if any, is not addressed by the field sampling design used here.

Relevance for Regional CH 4 Fluxes in the Amazon Basin
As has been shown for northern peatlands, our results show a putative coupling between peat geochemical and nutrient content and CH 4 flux (Whiting and Chanton, 1993). This has implications both for improvement of global CH 4 biogeochemical models and can inform environmental regulations concerning peatlands. Although tropical peatlands are known to be fundamental players in global CH 4 cycling, there is a concerning lack of empirical data to support GHG modeling (Sjogersten et al., 2014). Targeted measurements based on historical knowledge of nutrient status and hydrological patterns may assist in continent scale modeling of CH 4 flux. Furthermore, anthropological disturbance of peatlands has consequences for ecosystem scale C and nutrient cycling (Andersen et al., 2013). As nutrient status could be key for regulating microbial communities involved in CH 4 flux, changes to hydrological regimes and agricultural nutrient deposition has the potential to increase atmospheric CH 4 emissions from tropical peatlands.

CONCLUSION
In conclusion, nutrient concentration shapes microbial communities in the largest tropical peatland in South America, the Pastaza-Marañón Basin. Highly structured, niche-dependent distributions of methanogen and methanotroph families were found to exist, dependent on nutrient concentration. Methanogen communities, dominated by hydrogenotrophic Methanobacteriaceae, aggregate strongly with minerotrophic systems that display higher CH 4 flux. This suggests that hydrogenotrophic methanogenesis is the primary source of CH 4 in these minerotrophic systems. Methanotroph communities, dominated by Methylocystaceae, are dispersed across all systems regardless of CH 4 flux, pH and nutrient content. The putative association between CH 4 flux and peatland geochemical and nutrient content can inform modeling of the CH 4 cycle in the tropics and suggest increased atmospheric CH 4 emissions as a consequence of hydrology disturbance and agricultural nutrient deposition.

DATA AVAILABILITY STATEMENT
The 16S rRNA sequence datasets generated for this study can be found in the NCBI repository under the SRA Bioproject ID PRJNA501909 (https://www.ncbi.nlm.nih.gov/sra/SRX4961727).

AUTHOR CONTRIBUTIONS
HC-Q designed the study. HC-Q, DF, MZ-E, JP, JH, JA-P, and JU-M conducted the sample collection, gas, chemical, and molecular analyses. DF and MZ-E equally performed the data analyses, modeling and synthesis, and drafted the manuscript. All authors contributed to the final version of the manuscript.