Novel molecular markers for the detection of methanogens and phylogenetic analyses of methanogenic communities
- 1Department of Bacterial Genetics, Institute of Microbiology, Faculty of Biology, University of Warsaw, Warsaw, Poland
- 2Laboratory of RNA Metabolism and Functional Genomics, Institute of Biochemistry and Biophysics, Polish Academy of Sciences, Warsaw, Poland
- 3Institute of Genetics and Biotechnology, Faculty of Biology, University of Warsaw, Warsaw, Poland
- 4Department of Bioinformatics, Institute of Biochemistry and Biophysics, Polish Academy of Sciences, Warsaw, Poland
- 5Department of Systems Biology, Institute of Plant Experimental Biology and Biotechnology, Faculty of Biology, University of Warsaw, Warsaw, Poland
- 6Laboratory of Environmental Pollution Analysis, Faculty of Biology, University of Warsaw, Warsaw, Poland
Methanogenic Archaea produce approximately one billion tons of methane annually, but their biology remains largely unknown. This is partially due to the large phylogenetic and phenotypic diversity of this group of organisms, which inhabit various anoxic environments including peatlands, freshwater sediments, landfills, anaerobic digesters and the intestinal tracts of ruminants. Research is also hampered by the inability to cultivate methanogenic Archaea. Therefore, biodiversity studies have relied on the use of 16S rRNA and mcrA [encoding the α subunit of the methyl coenzyme M (methyl-CoM) reductase] genes as molecular markers for the detection and phylogenetic analysis of methanogens. Here, we describe four novel molecular markers that should prove useful in the detailed analysis of methanogenic consortia, with a special focus on methylotrophic methanogens. We have developed and validated sets of degenerate PCR primers for the amplification of genes encoding key enzymes involved in methanogenesis: mcrB and mcrG (encoding β and γ subunits of the methyl-CoM reductase, involved in the conversion of methyl-CoM to methane), mtaB (encoding methanol-5-hydroxybenzimidazolylcobamide Co-methyltransferase, catalyzing the conversion of methanol to methyl-CoM) and mtbA (encoding methylated [methylamine-specific corrinoid protein]:coenzyme M methyltransferase, involved in the conversion of mono-, di- and trimethylamine into methyl-CoM). The sensitivity of these primers was verified by high-throughput sequencing of PCR products amplified from DNA isolated from microorganisms present in anaerobic digesters. The selectivity of the markers was analyzed using phylogenetic methods. Our results indicate that the selected markers and the PCR primer sets can be used as specific tools for in-depth diversity analyses of methanogenic consortia.
Methanogenesis is a metabolic process driven by obligate anaerobic Archaea. It is responsible for the production of over 90% of methane on Earth (Costa and Leigh, 2014). There are three main methanogenic pathways: (i) hydrogenotrophic methanogenesis using H2/CO2 for methane synthesis, (ii) acetoclastic methanogenesis, in which the methyl group from acetate is transferred to tetrahydrosarcinapterin and then to coenzyme M (CoM), and (iii) methylotrophic methanogenesis, using methyl groups from methanol and methylamines (mono-, di-, and trimethylamine) for the production of methyl-coenzyme M (Figure 1). The final step in all these pathways is common and involves the conversion of methyl-CoM into methane by methyl-coenzyme M reductase, an enzymatic complex that is present in all methanogens (Borrel et al., 2013) (Figure 1).
Figure 1. Schematic diagram of the superpathway of methanogenesis. E.C. numbers for particular enzymes are shown in parentheses. The red star indicates the mcrA gene encoding subunit α of a methyl-coenzyme M reductase I, which is commonly used as a molecular marker for the detection of methanogens. The yellow stars denote molecular markers developed in this study, for which sets of PCR primers were designed.
Methanogenesis is of great importance for biotechnology (e.g., fuel production) and environmental protection (methane emissions contribute to global warming) (Escamilla-Alvaradoa et al., 2012). Therefore, the process has been extensively studied (Gao and Gupta, 2007; Ferry, 2010; Yoon et al., 2013). Consequently, novel species representing particular groups of methanogens are regularly reported (e.g., Dridi et al., 2012; Garcia-Maldonado et al., 2015), and various tools for the genetic and bioinformatic analysis of methanogenic Archaea are being developed (e.g., Farkas et al., 2013; Zakrzewski et al., 2013).
Methanogenic Archaea form complex consortia which remain largely uncharacterized. Methanogens form close relationships with their syntrophic partners and require very specific environmental conditions for growth, so they have proven very difficult to cultivate in the laboratory (Sekiguchi, 2006; Sakai et al., 2009). Therefore, a number of culture-independent methods have been applied to examine methanogenic consortia: (i) community fingerprinting by denaturing gradient gel electrophoresis—DGGE (Watanabe et al., 2004), (ii) single strand conformation polymorphism—SSCP (Delbes et al., 2001), (iii) terminal restriction fragment length polymorphism—T-RFLP (Akuzawa et al., 2011), (iv) fluorescence in situ hybridization—FISH (Diaz et al., 2006), and (v) real-time quantitative PCR—qPCR (Sawayama et al., 2006). However, the most reliable approach for the characterization of methanogenic communities is high-throughput sequencing using either 454 pyrosequencing (e.g., Schlüter et al., 2008; Rademacher et al., 2012; Stolze et al., 2015) or Illumina sequencing technologies (e.g., Caporaso et al., 2011; Zhou et al., 2011; Kuroda et al., 2014; Li et al., 2014).
The most frequently used molecular marker for phylogenetic analyses in metagenomic studies, of methanogenic communities is the 16S rRNA gene. However, low specificity of the oligonucleotide primers employed means that they generate 16S rDNA amplicons for all Archaea (not only methanogens) whose DNA is present in the analyzed sample. In the search for a more specific molecular marker for methanogens, the gene encoding the α subunit of the methyl-CoM reductase (mcrA) was identified and primers were developed for its amplification (Springer et al., 1995; Lueders et al., 2001; Luton et al., 2002; Friedrich et al., 2005; Yu et al., 2005; Denman et al., 2007; Steinberg and Regan, 2009). Of these, primers designed by Luton et al. (2002), are probably the most extensively used in ecological studies, since they produce the lowest bias in amplifying mcrA gene fragments from a wide range of phylogenetically diverse methanogens (e.g., Juottonen et al., 2006).
Several studies have demonstrated that the phylogeny of methanogens based on 16S rDNA and mcrA markers is consistent, although greater richness is usually observed using the latter (Luton et al., 2002; Hallam et al., 2003; Bapteste et al., 2005; Nettmann et al., 2008; Borrel et al., 2013). Interestingly, Wilkins and coworkers showed that these two genes produce different taxonomic profiles for samples taken from anaerobic digesters, i.e., environments extremely rich in methanogens (Wilkins et al., 2015). Clearly, the characterization of methanogenic communities requires a systematic approach using reliable molecular markers.
In this study, we have developed a set of degenerate PCR primers for the amplification of genes encoding key enzymes involved in methanogenesis. Some of these represent an alternative to mcrA primers commonly used for metagenomic analyses of methanogens. These novel primers amplify fragments of other genes of the mcr cluster, i.e., mcrB and mcrG encoding subunits β and γ of methyl-CoM reductase, respectively. Moreover, we have identified appropriate molecular markers for methylotrophic methanogens, which are probably the least explored group of methanogenic Archaea. These primers amplify fragments of the genes mtaB (encoding methanol-5-hydroxybenzimidazolylcobamide Co-methyltransferase, which is responsible for the conversion of methanol to methyl-CoM) and mtbA (encoding methylated [methylamine-specific corrinoid protein]:coenzyme M methyltransferase involved in the conversion of methylated amines into methyl-CoM). The extended panel of molecular markers provided by these novel primer sets should permit a deeper insight into the complex phylogeny, biology, and evolution of methanogens.
Materials and Methods
Standard Genetic Manipulations
PCR was performed in a Mastercycler (Eppendorf) using Taq DNA polymerase (Qiagen; with supplied buffer), dNTP mixture and appropriate primer pairs [Table 1 and additionally primer pairs S-D-Arch-0349-a-S-17/S-D-Arch-0786-a-A-20 for amplification of the variable region (V3V4) of the archaeal 16S rRNA gene (Klindworth et al., 2013), and MLf/MLr for mcrA gene amplification (Luton et al., 2002)]. PCR products of the methanogenesis-linked genes were purified by gel extraction, cloned using the pGEM®-T Easy Vector System (Promega) and transformed into E. coli TG1 (Stratagene) according to a standard procedure (Kushner, 1978). Standard methods were used for the isolation of recombinant plasmid DNA and for common DNA manipulation techniques (Sambrook and Russell, 2001).
Samples of microbial consortia involved in biogas production were collected from (i) the fermenter tank of an agricultural two-stage biogas plant anaerobic digester (AD) in Miedzyrzec Podlaski (Poland) and (ii) an effluent sludge tank from a one-stage wastewater treatment plant anaerobic digester (WD) at MPWIK Pulawy (Poland). In both cases, the samples were centrifuged (8000 g, 4°C, 15 min) and the pellets immediately stored in dry ice prior to DNA extraction.
DNA Extraction and Purification
DNA was isolated from anaerobic digester samples using a modified bead beating protocol. 1 g of pellet material (containing solids and microorganisms) was resuspended in 2 ml of lysis buffer [100 mM Tris-HCl (pH 8.0), 100 mM sodium EDTA (pH 8.0), 100 mM sodium phosphate (pH 8.0), 1.5 M NaCl, 1% (w/v) CTAB] (Zhou et al., 1996). The cells were then disrupted by a 5-step bead beating protocol performed at 1800 rpm (4 × 15 s) and 3200 rpm (1 × 15 s) (MiniBeadBeater 8) using 0.8 g of zirconia/silica beads (ø 0.5 mm, BioSpec). After each round of bead beating the sample was centrifuged (8000 g, 5 min, 4°C), the supernatant retained, and the pellet resuspended in fresh lysis buffer. In addition, after the third round of bead beating, the samples were freeze/thawed five times. The supernatant from each round was extracted with phenol-chloroform-isoamyl alcohol [25:24:1 (vol)]. DNA was then precipitated with one volume of isopropanol, 0.1 volume of 3 M sodium acetate (pH 5.2), recovered by centrifugation at 13,000 g for 20 min, and the pellets washed twice with 70% (v/v) ethanol before resuspending in TE buffer.
The prepared DNA was purified to remove proteins, humic substances, and other impurities by cesium chloride density gradient centrifugation. The concentration and quality of the purified DNA were estimated using a NanoDrop 2000c spectrophotometer (NanoDrop Technologies) and by agarose gel electrophoresis. The applied method yielded highly pure DNA (A260/A280 = 1.8; A260/A230 = 1.9) suitable for metagenomic analysis.
Library Preparation and Amplicon Sequencing
PCR products were analyzed by electrophoresis on 2% agarose gels (1x TAE buffer) with ethidium bromide staining. The amplified DNA fragments from replicate PCRs were pooled and then purified using Agencourt AMPure XP beads (Beckman Coulter). Approximately 250 ng of each amplicon was used for library preparation with an Illumina TruSeq DNA Sample Preparation Kit according to the manufacturer's protocol, except that the final library amplification was omitted. The libraries were verified using a 2100 Bioanalyzer (Agilent) High-Sensitivity DNA Assay and KAPA Library Quantification Kit for the Illumina. Sequencing of amplicon DNA was performed using an Illumina MiSeq (MiSeq Reagent Kit v2, 500 cycles) with a read length of 250 bp.
Designing Oligonucleotide Primers Specific for mcrB, mcrG, mtaB, and mtbA Genes
Data from the NCBI database were used to design degenerate primers to amplify mcrB, mcrG, mtaB, and mtbA gene fragments. A two-stage design strategy was employed. First, nucleotide sequences of genes annotated as mcrB, mcrG, mtaB, and mtbA were retrieved from the NCBI database. These sequences were then used as a query to recover additional gene sequences that were not annotated or were incorrectly annotated. Nucleotide sequences of particular genes were retrieved from genome sequences (completed and drafts) of methanogenic Archaea available on Jan 10th 2014. For each gene, multiple sequence alignments were prepared using ClustalW (Chenna et al., 2003) and MEGA6 (Tamura et al., 2013). Conserved regions within the obtained alignments were identified and used in the design of appropriate degenerate primers. Primer pairs with the lowest degree of degeneracy and producing amplicons not exceeding 500 bp were chosen. This size limit was imposed so that both 454 pyrosequencing and Illumina platforms could be used for amplicon sequencing.
In silico PCR with iPCRess (Slater and Birney, 2005) was done on dataset consisting of complete microbial genomes (5274 in total) obtained from NCBI database. We allowed for two mismatches per primer and required that both primers match and the product length is similar (±50 nucleotides) to expected length. The only exception was the set of mcrG-specific primers, that required allowance of 4 mismatches, due to bigger length of their sequences.
Bioinformatic Analysis of High-throughput Amplicon Sequencing Data
For each selected protein family a reference set of sequences was assembled from the results of searches of the NCBI NR database with BLAST software (Altschul et al., 1997), using known archaeal members of each family as query sequences and an E-value of 0.001 as the threshold. These reference sets were not specifically curated to allow the presence of false positives such as proteins from Bacteria or Eukarya. We consider them false positives, as the process of methanogenesis is limited only to Archaea. A presence of the sequences more similar to bacterial homologs of marker proteins than to archaeal ones would indicate low specificity of the designed primers. We specifically screened for such a cases after phylogenetic placement of reads.
Paired-end reads were merged with FLASH (Magoc and Salzberg, 2011) and then mapped to reference sets using BLASTX, again with an E-value of 0.001 as the threshold. Translated sequences were extracted from the BLAST high scoring pairs (HSPs), and reads with no hits, containing stop codons (presumably generated by frameshifts) or sequences shorter than 30 amino acids were discarded. Therefore, for each primer pair, a reference set of known protein sequences was obtained, as well as a set of protein sequences derived from sequenced amplicons. The latter are referred as “inferred peptides” as they correspond to fragments of target proteins. The ratio of number of inferred peptides to number of all merged reads is the measure of primer sensitivity.
Sequences from the reference sets were aligned with MAFFT (Katoh and Standley, 2013) using default options. Based on these alignments, a maximum likelihood phylogenetic tree was constructed for each protein family using FastTree software (Price et al., 2009) with the Gamma20 model. Sequences inferred from reads were then merged with sequences from reference sets for each protein family and aligned with MAFFT as described above. The resulting alignment and the phylogenetic tree of reference sequences were used as the input to the Evolutionary Placement Algorithm, part of the RAxML package (Stamatakis, 2014). The reads were placed on the reference phylogenetic tree using the PROTGAMMAWAG substitution model. Placements were subsequently trimmed with guppy software (Matsen et al., 2010) using 0.01 as the minimal threshold mass from the leaf to the root. Results underwent guppy “fat” conversion to the PhyloXML file format (Han and Zmasek, 2009) and were then visualized using Archeopteryx software (Han and Zmasek, 2009). The visualization resulted in coloring branches that point to a node or a leaf to which reads were assigned in red. All other branches were colored in black.
Amplicons from 16S rDNA were processed differently. All sequence reads were processed via the NGS analysis pipeline of the SILVA rRNA gene database project (SILVAngs 1.2) (Quast et al., 2013). Using the SILVA Incremental Aligner [SINA SINA v1.2.10 for ARB SVN (revision 21008)] (Pruesse et al., 2012), each read was aligned against the SILVA SSU rRNA SEED and quality controlled (Quast et al., 2013). Reads shorter than 50 aligned nucleotides and reads with more than 2% ambiguities or 2% homopolymers, were excluded from further processing. In addition, putative contaminants and artifacts, and reads with low alignment quality (50 alignment identity, 40 alignment score reported by the SINA), were identified and excluded from downstream analysis.
The classification of each operational taxonomic unit (OTU) reference read was mapped onto all reads that were assigned to the respective OTU. This yielded quantitative information (number of individual reads per taxonomic path), within the limitations of PCR and sequencing technique biases, as well as multiple rRNA operons. Reads without any BLAST hits or those with weak BLAST hits, where the function “(% sequence identity + % alignment coverage)/2” did not exceed the value of 93, remained unclassified.
Raw sequences obtained in this study have been deposited in the SRA (NCBI) database with the accession number PRJNA284604.
Results and Discussion
General Diversity of Archaea in Anaerobic Digesters—16S rRNA and mcrA Molecular Marker Analyses
In the analyses performed in this study metagenomic DNA was extracted from two samples of microbial consortia involved in biogas production (and therefore reach in methanogens). For the description of the overall diversity of Archaea in the analyzed samples, 16S rDNA-specific primers were used (Klindworth et al., 2013). This analysis revealed that methanogens are dominant microorganisms in the studied anaerobic digesters (74% for AD and 95% for WD) and include representatives of four of the seven methanogenic orders (i.e., Methanosarcinales, Methanomicrobiales, Methanobacteriales, Methanomassiliicoccales). The most abundant methanogens in both digesters were Methanosarcinales, represented by the families Methanosaetaceae (~38%) and Methanosarcinaceae (~18%), followed by Methanomicrobiaceae (~20%) of the Methanomicrobiales order (Figure 2).
Figure 2. Relative abundance of archaeal OTUs defined using the 16S rRNA gene hyper-variable region V3V4. The bar chart shows the diversity of Archaea at the lowest reliable taxonomic level (where possible the default family is denoted in the key). AD, agricultural biogas plant anaerobic digester; WD, wastewater treatment plant anaerobic digester.
Abundant non-methanogenic Archaea such as Miscellaneous Crenarchaeotic Group (MCG) (11%) and Halobacteria (7%) represented by Deep Sea Euryarchaeotic Group (DSEG) and Deep Sea Hydrothermal Vent Gp 6 (DHVEG-6) were also detected in the AD sample (Figure 2). These groups are phylogenetically diverse and there is a little knowledge of their ecology and metabolism, however it seems that MCG archaeons are able to ferment wide variety of recalcitrant substrates (Kubo et al., 2012) and DSEG are positively correlated with putative ammonia-oxidizing Thaumarchaeota (Restrepo-Ortiz and Casamayor, 2013).
In addition to the 16S rRNA marker, the mcrA gene was used for taxonomic profiling of methanogenic communities in both digesters. The mcrA gene fragments amplified using primers MLf/MLr (Luton et al., 2002) were sequenced and analyzed. More than half of the sequences (57%) amplified from the AD sample were assigned to uncultured Archaea, belonging to the Methanomassiliicoccales (23%), Methanomicrobiales (13%), Methanobacteriales (11%) and Methanosarcinales (10%) orders (Figure 3), suggesting dominance of hydrogenotrophic methanogens over acetoclastic Archaea. The most abundant genera in AD were Methanobacterium sp. Maddingley MBC34 (11%) followed by Methanosaeta concilli (9%) and Methanoculleus spp. (4%) (Figure 3). Similarly in WD, the majority of the mcrA amplicons were classified as uncultured Archaea belonging to orders Methanomicrobiales (27%) and Methanomassiliicoccales (7%) (Figure 3), while at the genus level most of the methanogens were identified as Methanometylovorans hollandica (21%), Methanosaeta concilli (16%), Methonoculleus spp. (12%), or Methanoplanus petrolearius (3%) (Figure 3).
Figure 3. Phylogenetic placement of mcrA amplicons from AD (A) and WD (B) samples. The width of the red branches corresponds to the number of unique mcrA amplicon sequence reads in that particular branch (this can be either a leaf or node). The collapse of some branches (mapped to uncultured archaeons) to increase the clarity of the trees is indicated by a black triangle. Numbers next to the entries “uncultured archeon” indicate the same microorganism on both trees.
The results obtained for both marker genes (16S rRNA and mcrA) only partially overlapped, probably due to differences in primer affinities and variation in the gene copy numbers. This observation is in agreement with a previous report showing that these two marker genes generate different taxonomic profiles (Wilkins et al., 2015). Therefore, for a greater insight into the structure of methanogenic communities and to verify the obtained results, novel molecular markers specific for other methanogenesis-linked genes were developed.
Development of mcrB-, mcrG-, mtaB-, and mtbA-specific Primers
For the design of degenerate primers specific for the mcrB, mcrG, mtaB, and mtbA genes, sequences were retrieved from the NCBI database [36 nucleotide sequences for mcrB (Figure S1), 61 for mcrG (Figure S2), 26 for mtaB (Figure S3) and 13 for mtbA (Figure S4)]. The mcrG gene turned out to be highly variable, which hampered primer design. Therefore, phylogenetic analysis was performed to distinguish conserved clusters among the analyzed mcrG genes. Two groups of mcrG sequences were distinguished: (i) MCR_G1 (grouping 35 mcrG genes of Methanobacterium spp., Methanobrevibacter spp., Methanocaldococcus spp., Methanococcus spp., Methanothermobacter spp., Methanothermococcus spp., Methanothermus spp., Methanotorris spp., Methanosphaera spp.), and (ii) MCR_G2 (grouping 26 mcrG genes of Methanocella spp., Methanococcoides spp., Methanocorpusculum spp., Methanoculleus spp., Methanohalobium spp., Methanohalophilus spp., Methanolobus spp., Methanoplanus spp., Methanopyrus spp., Methanoregula spp., Methanosalsum spp., Methanosarcina spp., Methanospirillum spp., Methanosphaerula spp.) (Figure S5). The nucleotide sequences of mcrG genes from particular groups were then used to design specific primer pairs.
For the subsequent functional analyses, 28 primers were selected for synthesis, including 6 for mcrB, 9 for mcrG and mtaB, and 4 for mtbA. The initial PCRs were performed with all primer pairs and DNA samples from the AD and WD fermenters as templates. The primer pairs giving the strongest amplification products of the expected size were selected for further analysis. The PCR products were cloned in vector pGEM-T Easy and then inserts of five random clones from each experimental set were sequenced using the sequencing primer M13 Reverse. The BLAST analysis of the resulting sequences revealed the specificity of each primer pair. At this stage, all primers designed for amplification of the mcrG genes of MCR_G2 group methanogens were rejected due to low specificity. Based on those analyses and taking into account the amplification yield, four primer pairs were selected and the optimal PCR conditions were determined (Table 1). Primer pairs specificity was also initially confirmed by in silico PCR analysis using 5274 complete microbial genomes (Table S1).
Since the panel of primers developed in this study was designed to be used in the high-throughput amplicon sequencing analysis of methanogenic communities, their selectivity was tested in the high-throughput sequencing experiments.
Analysis of the Selectivity of the mcrB- and mcrG-specific Primers
DNA fragments were amplified using the developed primer pairs with template DNAs isolated from the anaerobic reactors AD and WD. The raw sequence data obtained from Illumina sequencing were processed and analyzed (Table 2).
This analysis revealed that LMCRB/RMCRB primers, designed to the mcrB gene, amplified DNA fragments comprising sequences representing four methanogenic orders: Methanobacteriales, Methanomassiliicoccales, Methanomicrobiales, and Methanosarcinales (Figure 4). The dominant genus in both digesters was Methanoculleus spp. (48% for AD and 67% for WD), with M. marisnigri as the most abundant species (37 and 53%, respectively). This finding remains in good agreement with previous observations showing that the predominant order in biogas-producing microbial communities in anaerobic digesters is usually Methanomicrobiales, and the most abundant species is hydrogenotrophic M. marisnigri (Wirth et al., 2012). Moreover, in AD, 27% of sequences were classified as uncultured Methanomassiliicoccales (with 4% described as Candidatus Methanoplasma termitum) and 17% as Methanosaeta concilli. The second and third most abundant methanogens in WD were Methanomethylovorans hollandica (19%) and Methanosaeta concilli (6%), respectively (Figure 4).
Figure 4. Phylogenetic placement of mcrB amplicons from AD (A) and WD (B) samples. The width of the red branches corresponds to the number of unique mcrB amplicon sequence reads in that particular branch (this can be either a leaf or node). The leaf for Methanoculleus sp. MH98A was shortened, as indicated by two slashes.
The mcrG gene fragments (amplified with primers LMCRG1/RMCRG1) comprised sequences representing five methanogenic orders: Methanobacteriales, Methanococcales, Methanomicrobiales, Methanomassiliicoccales, and Methanosarcinales. However, representatives of hydrogenotrophic Methanobacteriales were absolutely dominant in both digesters (Figure 5). The most abundant OTUmcrG in AD was assigned to Methanobacterium spp. (97%) (with 7% mapped to M. formicicum), while WD was dominated by Methanosphaera stadtmanae (54%) and Methanobacterium spp. (39%) (with 28% mapped to M. formicicum) and Methanobrevibacter spp. (5%) (Figure 5).
Figure 5. Phylogenetic placement of mcrG amplicons from AD (A) and WD (B) samples. The width of the red branches corresponds to the number of unique mcrG amplicon sequence reads in that particular branch (this can be either a leaf or node).
The above analysis revealed that primers LMCRB/RMCRB are highly specific for mcrB genes of methanogens. Therefore, similarly to the commonly employed mcrA-specific primers, they may be used for an overall characterization of the taxonomic structure of methanogenic communities. The application of both mcrA and mcrB molecular markers permits cross-checking and should give a deeper and more detailed insight into the taxonomic structure of various methanogenic communities. It is worth mentioning that the results obtained using the newly developed primers for mcrB were partially consistent with those obtained by mcrA analysis, and confirmed that the hydrogenotrophic pathway of methane synthesis is employed in the analyzed environments. Moreover, these results demonstrated the importance of the newly described seventh order of methanogenic Methanomassiliicoccales (Iino et al., 2013; Borrel et al., 2014) in the analyzed biogas digesters (Figure 6, Table S2).
Figure 6. Relative abundance of archaeal OTUs defined using the mcrA and mcrB and mcrG (of MCR_G1 group) gene fragments. The bar chart shows the diversity of Archaea at the lowest reliable taxonomic level (mostly genus). AD, agricultural biogas plant anaerobic digester; WD, wastewater treatment plant anaerobic digester.
The mcrG primers LMCRG1/RMCRG1 permitted the analysis of the minority of methanogenic Archaea that were not dominant in mcrA/mcrB analysis (except Methanobacterium for the mcrA marker). Therefore, the obtained results were not consistent with those obtained by mcrA and mcrB analyses. This is the consequence of the fact that the primers LMCRG1/RMCRG1 are specific only for the previously described MCR_G1 group of sequences (Figure S5) and their use could generate programmed bias (Figure 6, Table S2).
Analysis of the Selectivity of the mtaB- and mtbA-specific Primers
In the course of this study, two other marker genes (mtaB and mtbA) specific for methylotrophic methanogens were selected and primer pairs developed. High-throughput sequencing of amplicons obtained using mtaB primers LMTAB/RMTAB detected sequences representing only two orders, Methanosarcinales and Methanobacteriales. In AD, 76% of sequences were assigned to Methanosarcina spp. [including M. barkeri (69%) and M. mazei (7%)] and 23% to Methanosphaera stadtmanae. Reactor WD was dominated by M. hollandica (94%), followed by M. stadtmanae (6%). In comparison, use of mtbA-specific primers LMTBA/RMTBA detected sequences mostly belonging to the Methanosarcinales, with two dominating species: M. barkeri (99%) in AD and M. hollandica (99%) in WD. Single sequences in WD and AD were assigned to Halobacteriales and Methanomassiliicoccales, respectively.
Sequencing of the mtaB and mtbA amplicons clearly indicated that in the analyzed digesters, Methanosarcinales are mainly responsible for the utilization of methylamines, while the conversion of methanol to methane is additionally performed by M. stadtmanae (of Methanobacteriales), which is consistent with previous studies (Fricke et al., 2006; Liu and Whitman, 2008).
Four novel molecular markers were designed and tested for the detection and taxonomic analyses of methanogenic communities. Primers specific to the mcrB and mcrG genes (present in all methanogens), as well as the mtaB and mtbA genes, characteristic for methylotrophic methanogens, were developed. High-throughput sequencing of the amplicons obtained using these primers revealed their high specificity and indicated that these marker genes could be used for taxonomic profiling of methanogenic consortia.
The mcrB and mcrG molecular markers increased the resolution of high-throughput amplicon sequencing analyses of methanogenic communities that until now have only been investigated using the mcrA gene. The use of mcrA, mcrB, and mcrG, together with the 16S rRNA gene marker, should give a much broader overview of the taxonomic diversity of complex methanogenic communities. In addition, the analysis of two other marker genes (mtaB and mtbA) can provide an insight into the metabolic potential of the analyzed methanogens, since they permit the detection and analysis of an enigmatic group of methylotrophic methanogens, which are able to produce methane from methanol or methylamines.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
This work was supported by the National Centre for Research and Development (Poland) grant no. 177481, as well as by the EU European Regional Development Fund, the Operational Program Innovative Economy 2007–2013, agreement POIG.01.01.02-14-054/09-00. Some experiments were carried out with the use of CePT infrastructure financed by the European Union—the European Regional Development Fund [Innovative economy 2007–13, Agreement POIG.02.02.00-14-024/08-00].
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/article/10.3389/fmicb.2015.00694
Figure S1. Alignment of the conserved fragments of the mcrB genes of 36 methanogens used in the design of primers LMCRB and RMCRB.
Figure S2. Alignment of the conserved fragments of the mcrG genes (of MCR_G1 group) of 35 methanogens used in the design of primers LMCRG1 and RMCRG1.
Figure S3. Alignment of the conserved fragments of the mtaB genes of 26 methanogens used in the design of primers LMTAB and RMTAB.
Figure S4. Alignment of the conserved fragments of the mtbA genes of 13 methanogens used in the design of primers LMTBA and RMTBA.
Figure S5. Phylogenetic tree for mcrG nucleotide sequences (from NCBI database). The tree was constructed using the maximum-likelihood algorithm. Statistical support for the internal nodes was determined by 1000 bootstrap replicates and values of >50% are shown.
Akuzawa, M., Hori, T., Haruta, S., Ueno, Y., Ishii, M., and Igarashi, Y. (2011). Distinctive responses of metabolically active microbiota to acidification in a thermophilic anaerobic digester. Microb. Ecol. 61, 595–605. doi: 10.1007/s00248-010-9788-1
Altschul, S. F., Madden, T. L., Schaffer, A. A., Zhang, J., Zhang, Z., Miller, W., et al. (1997). Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res. 25, 3389–3402. doi: 10.1093/nar/25.17.3389
Borrel, G., O'Toole, P. W., Harris, H. M., Peyret, P., Brugere, J. F., and Gribaldo, S. (2013). Phylogenomic data support a seventh order of methylotrophic methanogens and provide insights into the evolution of methanogenesis. Genome Biol. Evol. 5, 1769–1780. doi: 10.1093/gbe/evt128
Borrel, G., Parisot, N., Harris, H. M., Peyretaillade, E., Gaci, N., Tottey, W., et al. (2014). Comparative genomics highlights the unique biology of Methanomassiliicoccales, a Thermoplasmatales-related seventh order of methanogenic archaea that encodes pyrrolysine. BMC Genomics 15:679. doi: 10.1186/1471-2164-15-679
Caporaso, J. G., Lauber, C. L., Walters, W. A., Berg-Lyons, D., Lozupone, C. A., Turnbaugh, P. J., et al. (2011). Global patterns of 16S rRNA diversity at a depth of millions of sequences per sample. Proc. Natl. Acad. Sci. U.S.A. 108(Suppl. 1), 4516–4522. doi: 10.1073/pnas.1000080107
Chenna, R., Sugawara, H., Koike, T., Lopez, R., Gibson, T. J., Higgins, D. G., et al. (2003). Multiple sequence alignment with the Clustal series of programs. Nucleic Acids Res. 31, 3497–3500. doi: 10.1093/nar/gkg500
Delbes, C., Moletta, R., and Godon, J. (2001). Bacterial and archaeal 16S rDNA and 16S rRNA dynamics during an acetate crisis in an anaerobic digestor ecosystem. FEMS Microbiol. Ecol. 35, 19–26. doi: 10.1016/S0168-6496(00)00107-0
Denman, S. E., Tomkins, N. W., and McSweeney, C. S. (2007). Quantitation and diversity analysis of ruminal methanogenic populations in response to the antimethanogenic compound bromochloromethane. FEMS Microbiol. Ecol. 62, 313–322. doi: 10.1111/j.1574-6941.2007.00394.x
Diaz, E. E., Stams, A. J., Amils, R., and Sanz, J. L. (2006). Phenotypic properties and microbial diversity of methanogenic granules from a full-scale upflow anaerobic sludge bed reactor treating brewery wastewater. Appl. Environ. Microbiol. 72, 4942–4949. doi: 10.1128/AEM.02985-05
Dridi, B., Fardeau, M. L., Ollivier, B., Raoult, D., and Drancourt, M. (2012). Methanomassiliicoccus luminyensis gen. nov., sp. nov., a methanogenic archaeon isolated from human faeces. Int. J. Syst. Evol. Microbiol. 62, 1902–1907. doi: 10.1099/ijs.0.033712-0
Escamilla-Alvaradoa, C., Ríos-Lealb, E., Ponce-Noyolac, M. T., and Poggi-Varaldo, H. M. (2012). Gas biofuels from solid substrate hydrogenogenic–methanogenic fermentation of the organic fraction of municipal solid waste. Process Biochem. 47, 1572–1587. doi: 10.1016/j.procbio.2011.12.006
Fricke, W. F., Seedorf, H., Henne, A., Kruer, M., Liesegang, H., Hedderich, R., et al. (2006). The genome sequence of Methanosphaera stadtmanae reveals why this human intestinal archaeon is restricted to methanol and H2 for methane formation and ATP synthesis. J. Bacteriol. 188, 642–658. doi: 10.1128/JB.188.2.642-658.2006
Gao, B., and Gupta, R. S. (2007). Phylogenomic analysis of proteins that are distinctive of Archaea and its main subgroups and the origin of methanogenesis. BMC Genomics 8:86. doi: 10.1186/1471-2164-8-86
Garcia-Maldonado, J. Q., Bebout, B. M., Everroad, R. C., and Lopez-Cortes, A. (2015). Evidence of novel phylogenetic lineages of methanogenic archaea from hypersaline microbial mats. Microb. Ecol. 69, 106–117. doi: 10.1007/s00248-014-0473-7
Hallam, S. J., Girguis, P. R., Preston, C. M., Richardson, P. M., and Delong, E. F. (2003). Identification of methyl coenzyme M reductase A (mcrA) genes associated with methane-oxidizing archaea. Appl. Environ. Microbiol. 69, 5483–5491. doi: 10.1128/AEM.69.9.5483-5491.2003
Iino, T., Tamaki, H., Tamazawa, S., Ueno, Y., Ohkuma, M., Suzuki, K., et al. (2013). Candidatus Methanogranum caenicola: a novel methanogen from the anaerobic digested sludge, and proposal of Methanomassiliicoccaceae fam. nov. and Methanomassiliicoccales ord. nov., for a methanogenic lineage of the class Thermoplasmata. Microbes Environ. 28, 244–250. doi: 10.1264/jsme2.ME12189
Juottonen, H., Galand, P. E., and Yrjala, K. (2006). Detection of methanogenic Archaea in peat: comparison of PCR primers targeting the mcrA gene. Res. Microbiol. 157, 914–921. doi: 10.1016/j.resmic.2006.08.006
Klindworth, A., Pruesse, E., Schweer, T., Peplies, J., Quast, C., Horn, M., et al. (2013). Evaluation of general 16S ribosomal RNA gene PCR primers for classical and next-generation sequencing-based diversity studies. Nucleic Acids Res. 41:e1. doi: 10.1093/nar/gks808
Kubo, K., Lloyd, K. G., Biddle, J. F., Amann, R., Teske, A., and Knittel, K. (2012). Archaea of the Miscellaneous Crenarchaeotal Group are abundant, diverse and widespread in marine sediments. ISME J. 6, 1949–1965. doi: 10.1038/ismej.2012.37
Kuroda, K., Hatamoto, M., Nakahara, N., Abe, K., Takahashi, M., Araki, N., et al. (2014). Community composition of known and uncultured archaeal lineages in anaerobic or anoxic wastewater treatment sludge. Microb. Ecol. 69, 586–596. doi: 10.1007/s00248-014-0525-z
Kushner, S. R. (1978). “An improved method for transformation of E. coli with ColE1 derived plasmids,” in Genetic Engineering, eds H. B. Boyer and S. Nicosia (Amsterdam: Elsevier/North-Holland), 17–23.
Li, Y. F., Nelson, M. C., Chen, P. H., Graf, J., Li, Y., and Yu, Z. (2014). Comparison of the microbial communities in solid-state anaerobic digestion (SS-AD) reactors operated at mesophilic and thermophilic temperatures. Appl. Microbiol. Biotechnol. 99, 969–980. doi: 10.1007/s00253-014-6036-5
Lueders, T., Chin, K. J., Conrad, R., and Friedrich, M. (2001). Molecular analyses of methyl-coenzyme M reductase alpha-subunit (mcrA) genes in rice field soil and enrichment cultures reveal the methanogenic phenotype of a novel archaeal lineage. Environ. Microbiol. 3, 194–204. doi: 10.1046/j.1462-2920.2001.00179.x
Luton, P. E., Wayne, J. M., Sharp, R. J., and Riley, P. W. (2002). The mcrA gene as an alternative to 16S rRNA in the phylogenetic analysis of methanogen populations in landfill. Microbiology 148, 3521–3530. doi: 10.1099/00221287-148-11-3521
Matsen, F. A., Kodner, R. B., and Armbrust, E. V. (2010). pplacer: linear time maximum-likelihood and Bayesian phylogenetic placement of sequences onto a fixed reference tree. BMC Bioinformatics 11:538. doi: 10.1186/1471-2105-11-538
Nettmann, E., Bergmann, I., Mundt, K., Linke, B., and Klocke, M. (2008). Archaea diversity within a commercial biogas plant utilizing herbal biomass determined by 16S rDNA and mcrA analysis. J. Appl. Microbiol. 105, 1835–1850. doi: 10.1111/j.1365-2672.2008.03949.x
Price, M. N., Dehal, P. S., and Arkin, A. P. (2009). FastTree: computing large minimum evolution trees with profiles instead of a distance matrix. Mol. Biol. Evol. 26, 1641–1650. doi: 10.1093/molbev/msp077
Pruesse, E., Peplies, J., and Glockner, F. O. (2012). SINA: accurate high-throughput multiple sequence alignment of ribosomal RNA genes. Bioinformatics 28, 1823–1829. doi: 10.1093/bioinformatics/bts252
Quast, C., Pruesse, E., Yilmaz, P., Gerken, J., Schweer, T., Yarza, P., et al. (2013). The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41, D590–D596. doi: 10.1093/nar/gks1219
Rademacher, A., Zakrzewski, M., Schluter, A., Schonberg, M., Szczepanowski, R., Goesmann, A., et al. (2012). Characterization of microbial biofilms in a thermophilic biogas system by high-throughput metagenome sequencing. FEMS Microbiol. Ecol. 79, 785–799. doi: 10.1111/j.1574-6941.2011.01265.x
Restrepo-Ortiz, C. X., and Casamayor, E. O. (2013). Environmental distribution of two widespread uncultured freshwater Euryarchaeota clades unveiled by specific primers and quantitative PCR. Environ. Microbiol. Rep. 5, 861–867. doi: 10.1111/1758-2229.12088
Sakai, S., Imachi, H., Sekiguchi, Y., Tseng, I. C., Ohashi, A., Harada, H., et al. (2009). Cultivation of methanogens under low-hydrogen conditions by using the coculture method. Appl. Environ. Microbiol. 75, 4892–4896. doi: 10.1128/AEM.02835-08
Sawayama, S., Tsukahara, K., and Yagishita, T. (2006). Phylogenetic description of immobilized methanogenic community using real-time PCR in a fixed-bed anaerobic digester. Bioresour. Technol. 97, 69–76. doi: 10.1016/j.biortech.2005.02.011
Schlüter, A., Bekel, T., Diaz, N. N., Dondrup, M., Eichenlaub, R., Gartemann, K. H., et al. (2008). The metagenome of a biogas-producing microbial community of a production-scale biogas plant fermenter analysed by the 454-pyrosequencing technology. J. Biotechnol. 136, 77–90. doi: 10.1016/j.jbiotec.2008.05.008
Springer, E., Sachs, M. S., Woese, C. R., and Boone, D. R. (1995). Partial gene sequences for the A subunit of methyl-coenzyme M reductase (mcrI) as a phylogenetic tool for the family Methanosarcinaceae. Int. J. Syst. Bacteriol. 45, 554–559. doi: 10.1099/00207713-45-3-554
Stolze, Y., Zakrzewski, M., Maus, I., Eikmeyer, F., Jaenicke, S., Rottmann, N., et al. (2015). Comparative metagenomics of biogas-producing microbial communities from production-scale biogas plants operating under wet or dry fermentation conditions. Biotechnol. Biofuels 8, 14. doi: 10.1186/s13068-014-0193-8
Watanabe, T., Asakawa, S., Nakamura, A., Nagaoka, K., and Kimura, M. (2004). DGGE method for analyzing 16S rDNA of methanogenic archaeal community in paddy field soil. FEMS Microbiol. Lett. 232, 153–163. doi: 10.1016/S0378-1097(04)00045-X
Wilkins, D., Lu, X. Y., Shen, Z., Chen, J., and Lee, P. K. (2015). Pyrosequencing of mcrA and archaeal 16S rRNA genes reveals diversity and substrate preferences of methanogen communities in anaerobic digesters. Appl. Environ. Microbiol. 81, 604–613. doi: 10.1128/AEM.02566-14
Wirth, R., Kovacs, E., Maroti, G., Bagi, Z., Rakhely, G., and Kovacs, K. L. (2012). Characterization of a biogas-producing microbial community by short-read next generation DNA sequencing. Biotechnol. Biofuels 5:41. doi: 10.1186/1754-6834-5-41
Yoon, S. H., Turkarslan, S., Reiss, D. J., Pan, M., Burn, J. A., Costa, K. C., et al. (2013). A systems level predictive model for global gene regulation of methanogenesis in a hydrogenotrophic methanogen. Genome Res. 23, 1839–1851. doi: 10.1101/gr.153916.112
Yu, Y., Lee, C., Kim, J., and Hwang, S. (2005). Group-specific primer and probe sets to detect methanogenic communities using quantitative real-time polymerase chain reaction. Biotechnol. Bioeng. 89, 670–679. doi: 10.1002/bit.20347
Zakrzewski, M., Bekel, T., Ander, C., Pühler, A., Rupp, O., Stoye, J., et al. (2013). MetaSAMS—a novel software platform for taxonomic classification, functional annotation and comparative analysis of metagenome datasets. J. Biotechnol. 167, 156–165. doi: 10.1016/j.jbiotec.2012.09.013
Zhou, H. W., Li, D. F., Tam, N. F., Jiang, X. T., Zhang, H., Sheng, H. F., et al. (2011). BIPES, a cost-effective high-throughput method for assessing microbial diversity. ISME J. 5, 741–749. doi: 10.1038/ismej.2010.160
Keywords: methanogenesis, metagenomics, methanogenic consortia, mcrB, mcrG, mtaB, mtbA
Citation: Dziewit L, Pyzik A, Romaniuk K, Sobczak A, Szczesny P, Lipinski L, Bartosik D and Drewniak L (2015) Novel molecular markers for the detection of methanogens and phylogenetic analyses of methanogenic communities. Front. Microbiol. 6:694. doi: 10.3389/fmicb.2015.00694
Received: 07 May 2015; Accepted: 22 June 2015;
Published: 07 July 2015.
Edited by:Eamonn P. Culligan, University College Cork, Ireland
Copyright © 2015 Dziewit, Pyzik, Romaniuk, Sobczak, Szczesny, Lipinski, Bartosik and Drewniak. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Lukasz Dziewit, Department of Bacterial Genetics, Institute of Microbiology, Faculty of Biology, University of Warsaw, Miecznikowa 1, Warsaw 02-096, Poland, firstname.lastname@example.org
†These authors have contributed equally to this work.