Appl Microbiol Biotechnol DOI 10.1007/s00253-014-5648-0

ENVIRONMENTAL BIOTECHNOLOGY

Metagenomic analysis of sludge from full-scale anaerobic digesters operated in municipal wastewater treatment plants Ying Yang & Ke Yu & Yu Xia & Frankie T. K. Lau & Daniel T. W. Tang & Wing Cheong Fung & Herbert H. P. Fang & Tong Zhang

Received: 13 January 2014 / Revised: 24 February 2014 / Accepted: 26 February 2014 # Springer-Verlag Berlin Heidelberg 2014

Abstract This study applied Illumina high-throughput sequencing to explore the microbial communities and functions in anaerobic digestion sludge (ADS) from two wastewater treatment plants based on a metagenomic view. Taxonomic analysis using SILVA SSU database indicated that Proteobacteria (9.52–13.50 %), Bacteroidetes (7.18 %– 10.65 %) and Firmicutes (7.53 %–9.46 %) were the most abundant phyla in the ADS. Differences of microbial communities between the two types of ADS were identified. Genera of Methanosaeta and Methanosarcina were the major methanogens. Functional analysis by SEED subsystems showed that the basic metabolic functions of metagenomes in the four ADS samples had no significant difference among them, but they were different from other microbial communities from activated sludge, human faeces, ocean and soil. Abundances of genes in methanogenesis pathway were also quantified using a methanogenesis genes database extracted from KEGG. Results showed that acetotrophic was the major methanogenic pathway in the anaerobic sludge digestion. Keywords Anaerobic digestion sludge . Taxonomicanalysis . Functional SEED subsystems . Methanogenesis pathway

Electronic supplementary material The online version of this article (doi:10.1007/s00253-014-5648-0) contains supplementary material, which is available to authorized users. Y. Yang : K. Yu : Y. Xia : H. H. P. Fang : T. Zhang (*) Environmental Biotechnology Laboratory, Department of Civil Engineering, University of Hong Kong, Hong Kong, China e-mail: [email protected] F. T. K. Lau : D. T. W. Tang : W. C. Fung Drainage Services Department, Government of the Hong Kong Special Administrative Region, Hong Kong, China

Introduction One of the major challenges in wastewater treatment to date is the management of excessive sludge generated in the treatment process. The cost of the sludge treatment could be higher than the cost of wastewater treatment in many cases (PérezElvira et al. 2006). Several treatment processes, including anaerobic digestion, aerobic digestion and composting, are usually applied to reduce the amount of sludge as well as the number of pathogenic microorganisms. Among them, anaerobic sludge digestion is most widely applied to reduce the amount of sludge (Novak et al. 2007; Sun et al. 2010), eliminate pathogen (Sahlstrom 2003) and generate biofuel in form of methane (Kim et al. 2003) as well. Various methods had been applied to investigate the microbial communities in anaerobic digesters, including clone library of 16S rRNA genes (Rincón et al. 2008), denaturing gradient gel electrophoresis (DGGE) analysis (Palatsi et al. 2010, 2011; Shin et al. 2010; Supaphol et al. 2011) and fluorescence in situ hybridization (FISH) (Braguglia et al. 2012). However, these methods only gave limited information compared with the emerging metagenomic approaches based on high-throughput sequencing (HTS), including 454 pyrosequencing, Illumina, Ion Torrent, Solid and PacBio. Using 16S rRNA gene fragments obtained by pyrosequencing, Garcia et al. (2011) studied the transition of microbial communities in the carrot waste digesters, and Garcia-Peña et al. (2011) investigated the microbial communities in anaerobic digestion and co-digestion processes of vegetable and fruit residues. 454 pyrosequencing of 16S rRNA genes provides huge information on the microbial communities, but limited information regarding to the functions of the microorganisms. However, using 454 pyrosequencing to study the functions of the environmental microorganisms using metagenomic analysis would be very costly, since large data sizes are usually required due to the complexity of the environmental samples. Li et al. (2013) used 454 pyrosequencing to

Appl Microbiol Biotechnol

perform metagenomic analysis on a biogas reactor. About 647 Mb of data was generated but still not yet close to the natural status for microorganism in the reactors established in laboratory. On the other hand, Illumina sequencing offers an alternative with lower cost (Mardis 2008; Metzker 2009) and had been widely used in the study of metagenome of human gut microbes (Qin et al. 2010), permafrost soil cores (Mackelprang et al. 2011), cow rumen (Hess et al. 2011), enhanced biological phosphorus removal reactors (Albertsen et al. 2011), activated sludge (Fang et al. 2013) and drinking water (Shi et al. 2013) as well. However, information on metagenomic analysis for fullscale anaerobic digesters in municipal wastewater treatment plants (WWTPs) is still very limited. The present study aimed to examine the broad spectrum of microbial population and functions in anaerobic digestion sludge (ADS) from a metagenomic view based on HTS, and learn the differences and common points of two types of fullscale anaerobic digesters. Four anaerobic sludge samples were collected from digesters in two local WWTPs and a total of 11 Gb DNA sequences were obtained. The metagenome data was used to investigate the microbial communities in anaerobic sludge, characterize the functional profiles of the microbial communities and screen possible genes associated with methanogenesis pathway.

Materials and methods Sample collection The ADS samples were collected from Shatin Sewage Treatment Works (ST STW) and Shek Wu Hui Sewage Treatment Works (SWH STW) in September 2011 and March 2012. The samples were designated as ST1109, ST1203, SWH1109 and SWH1203, respectively. Technical replicates of ADS were derived from SWH1203 and designated as SWH1203-A and SWH1203-B. Samples were mixed with 100 % ethanol at a ratio of 1:1 (v/v) immediately after being taken from the digesters and stored at −20 °C before DNA extraction. ST STW is the largest secondary wastewater treatment works in Hong Kong which serves a population of 600,000 and process 340,000 m3 of wastewater per day. While SWH STW processes 93,000 m 3 wastewater per day. Wastewater in ST STW contains ~1.1 % salinity as TDS because seawater is used in the toilet flushing system in that area of Hong Kong. Operation conditions and performance of the digesters are summarized in Table S1. DNA extraction and Illumina sequencing For DNA extraction, 2 ml of each sample was centrifuged at 4,000 rpm for 5 min to collect approximately 200 mg of the pellet. Genomic DNA of microorganisms was extracted using

the FastDNA@ Spin kit for Soil (MP Biomedicals, Illkirch, France) following the manufacturer’s instruction. DNA concentration and purity were determined by NanoDrop microspectrophotometry (ND-1000, Thermo Fisher Scientific, US). Approximately 6 μg for each DNA sample was used for library construction. Sequence libraries of ~180-bp DNA fragments were prepared and then sequenced by Illumina HiSeq 2000 at the Beijing Genomics Institute (BGI, Shenzhen, China). About 1 Gb (giga base pairs) of sequencing data were generated for ST1109 and SWH1109 each, 3 Gb of data were generated for ST1203, SWH1203-A and SWH 1203-B each (please refer to Supporting Information for details of Illumina sequencing). Data from SWH1203-B was used to investigate the profiles of antibiotic resistance genes in ADS previously (Ma et al. 2014). Data processing The platforms of high-throughput sequencing often produce artificial replicates that are nearly identical. Failure to remove these replicated sequences could lead to incorrect conclusions (Gomez-Alvarez et al. 2009). Therefore, the sequencing results were firstly filtered by a custom-made script to remove duplicates which first 50 bp were identical according to the recommendation in MG-RAST QC pipeline. In order to increase the length of the sequences from Illumina paired-end sequencing, the paired-end reads after deduplication were merged to generate longer DNA fragments by a customized script, which were defined as tags (Zhou et al. 2010). In detail, one of the paired-end reads was converted into its reverse-complement counterpart. Then if the reversecomplement counterpart and the other corresponding pairedend read had overlap longer than 10 bp, these two reads were merged into a longer tag as illustrated in Fig. S1 (python script are available at http://web.hku.hk/~zhangt/ZhangT.htm). Bioinformatics analysis In this study, the sample was computationally characterized without a prior assembly step, providing quantitative information into the microbial composition and gene content. Microbial community structures in the four anaerobic sludge samples were analyzed through the server of MG-RAST (version 3) (Meyer et al. 2008). MG-RAST IDs for the datasets are 4485273.3, 4485272.3, 4485327.3, 4485326.3, 4504598.3, 4504600.3, 4504628.3, 4504630.3, 4504632.3 and 4504638.3 for paired-end reads in ST1109, SWH1109, ST1203, SWH1203-A and SWH1203-B, respectively. MGRAST IDs for the merged tags are 4489147.3, 4489184.3, 4505688.3, 4505689.3 and 4505690.3, respectively. Taxa were assigned by Best Hit classification at the E-value cutoff of 10−5 with minimum alignment length of 50 bp (Pfister et al. 2010; Tringe 2005) using three rRNA gene

Appl Microbiol Biotechnol

databases in MG-RAST, i.e., SILVA Small Subunit (SSU) rRNA database (version 10.4) (Pruesse et al. 2007), Greengenes database (version 2011) (DeSantis et al. 2006) and Ribosomal Data Project database (RDP, version 10.22) (Cole et al. 2009). SILVA SSU database was also used to characterize the microbial community structure by local BLAST (Altschul et al. 1997) with E-value cutoff of 10−20 (Mackelprang et al. 2011). The sequences from the BLAST results were then assigned to NCBI taxonomies with MEGAN (version 4.67.5) (Huson et al. 2007) using the Lowest Common Ancestor (LCA) algorithm, and the absolute cutoff was BLAST bitscore 50, while the relative cutoff was 10 % of the top 50 hits. Functions of ADS were studied through MG-RAST by assigning sequences based on SEED subsystems database (Overbeek et al. 2005). Alignment of Illumina reads and tags to the SEED subsystems database was conducted using Hierarchical classification at the E-value cutoff of 10−5 and minimum alignment length of 17 amino acids (aa) (Jung et al. 2011; Pfister et al. 2010). The number of raw reads and tags assigned to different taxa or functions were normalized against the total number of annotated sequences for each database to avoid bias introduced by the absolute number. The results could be therefore used to compare in parallel, as well as to explore the microbial compositions and functions in the ADS samples. Methane production is a noteworthy advantage of anaerobic digestion. To characterize the abundances of genes encoding enzymes in methanogenesis pathway, a database of these genes was constructed using the method described previously (Yu and Zhang 2013). Sequences in the database were extracted from KEGG according to the KO number of genes in the methanogenesis pathway. Metagenomic datasets of the ADS samples were aligned against this methanogenesis database using BLAST with E-value cutoff of 10−5. Confirmation of methanogenesis genes was made by aligning the matched sequences against NCBI non-redundant database using BLAST with E-value cutoff of 10−10.

Results The reproducibility of metagenomic analysis Tags were used in the analysis since the percentages of annotated tags were higher in both taxonomic and functional annotations (Table S3; for details of the comparison of reads and tags, please refer to Supporting Information). Many researches had demonstrated the applications of high throughput sequencing in the study of microbial distributions and functional genes in different microbial communities (Hess et al. 2011; Jung et al. 2011). However, most of them did not have replicates for statistical analysis to check the technical

reproducibility of the applied methods (Knight et al. 2012; Prosser 2010). In the present study, technical replicates of ADS sample (from DNA extraction to final Illumina sequencing), i.e., SWH1203-A and SWH1203-B from the Shek Wu Hui STW (SWH ADS) were used to investigate the reproducibility of metagenomic analysis combined with Illumina high throughput sequencing. Tags in the two replicate datasets were annotated using three rRNA gene databases (SSU, RDP and Greengenes) and one protein database (SEED subsystems) through MGRAST. They were annotated by local BLAST using SSU database and LCA algorithm and the default parameters in MEGAN as well. The reproducibility of the annotated results was shown in Figs. S2 and S3, respectively. The results presented a very satisfactory reproducibility of the replicates (R2 >0.99), based on the linear fitting of taxonomic distribution percentage at different levels (phylum and genus using Best Hit classification in MG-RAST; phylum, class, order, family, and genus using LCA algorithm in BLAST-MEGAN approach) using various 16S rRNA gene databases, and functions distribution percentages using SEED subsystems. Microbial communities in ADS As stated above, the results of SWH1203-A and SWH1203-B suggested satisfactory reproducibility of the technical replicates. The sequencing data from these technical datasets were combined as SWH1203 for the following analyses. For each STW, samples were collected at two time points, i.e., September 2011 and March 2012, respectively. Among all the tags of the four samples, 0.08 % to 0.11 % of the sequences were identified as 16S rRNA genes at the E-value cutoff of 10−20 using the SILVA SSU database by BLAST. SILVA SSU database provide comprehensive information of rRNA sequences (Quast et al. 2012) and BLAST was applied because it is more sensitive than BLAT, the algorithm used in MG-RAST. Using this method, most of these sequences were annotated as Bacteria (83.91–94.77 %), followed by Archaea (2.81–4.64 %) and Eukaryota (0.31–1.62 %) (Table S4). In the domain of Bacteria, the phyla of Proteobacteria (9.52–13.50 % of the identified 16S rRNA gene fragments), Bacteroidetes (7.18–10.65 %) and Firmicutes (7.53–9.46 %) were dominant in all ADS samples (Fig. 1 and Table S5). Beta- and Delta-proteobacteria were the major classes in the phylum of Proteobacteria (Table S6A), having the average percentages of 25.06 % and 20.63 %, respectively. The percentages of Alpha- and Gamma-proteobacteria were higher in the September samples while more Delta-proteobacteria was found in the March samples for both STWs. The major classes in Bacteroidetes phylum were Bacteroidia and Bacteroidetes (class), Shingobacteriia and Flavobacteriia (Table S6B), taking 45.24 %, 25.29 %, 17.01 % and 10.86 % of all the Bacteroidetes (phylum) sequences in the four samples,

Appl Microbiol Biotechnol

Fig. 1 Relative abundances of detected phyla in the four ADS samples based on the taxonomy annotation from SSU using MEGAN. The data were visualized via Circos software (Krzywinski et al. 2009). The width of bars from each phylum indicate the relative abundances of that phylum in the ADS samples

respectively. The percentages of Bacteroidetes (class) were higher in the ST samples while the percentages of Sphingobacteriia were higher in the SWH samples. The family of Bacteroidaceae in Bacteroidetes (class) has been known as a group of fermentative bacteria in the acidogenic phase of the digestion process (Traversi et al. 2012). Firmicutes are well-known fermenters and syntrophic bacteria that can degrade various substrates (Garcia et al. 2011). Within the phylum of Firmicutes, the major classes were Bacilli and Clostridia with average percentages of 6.15 % and 75.35 %, respectively (Table S6C). Percentages of Baccilli were about two times higher in ST samples than that of SWH samples while percentages of Clostridia were higher in SWH samples by approximately one third than those of ST samples.

Other abundant phyla of Bacteria included Actinobacteria (5.05–7.89 %), Thermotogae (2.42–6.18 %), Chloroflexi (1.75–4.54 %), Spirochaetes (1.27–3.42 %), Tenericutes (0.04–3.93 %) and Planctomycetes (0.80–1.89 %) (Fig. 1 and Table S5). Euryarchaeota was the only phylum found in Archaea (Fig. 1), accounting for 2.66 % to 4.05 % of tags from the four ADS samples. 64.37 % to 73.20 % of Euryarchaeota were class of Methanomicrobia in ST and SWH, respectively (Table S6E). Moreover, Methanobacteria was only found in SWH samples with the percentages of 7.19 % and 5.48 % while Halobacteria was only found in ST samples with the percentages of 3.48 % and 3.24 %. The genera of Methanosaeta and Methanosarcina in the class of Methanomicrobia were major methanogens in the

Appl Microbiol Biotechnol

samples. The other methanogens, including Methanolinea, Methanobrevibacter, Methanospirillum and Methanobacterium, were only detected in one or two of the ADS samples with low abundances (Table S7). Abundances of SEED subsystems in ADS The annotation of functional genes was conducted using SEED subsystems in MG-RAST with E-value cutoff of 10−5 and minimum alignment length of 17 aa (Jung et al. 2011; Pfister et al. 2010). The distributions of the 28 Level 1 subsystems in four ADS samples were highly alike (Table S8). For the 462 Level 2 subsystems, the abundances of major (more than 1 %) Level 2 subsystems also showed no significant difference (Table S9), suggesting the consistent major function profiles in different anaerobic digesters. Among the major Level 1 subsystems besides clusteringbased subsystem, the subsystem of carbohydrates was the most abundant, followed by protein metabolism, amino acids and derivatives. Therefore, the subsystem of carbohydrate was subjected for further analysis. As illustrated by Fig. S4, the major Level 2 subsystems in carbohydrate subsystem were central carbohydrate metabolism and one-carbon metabolism. Central carbohydrate metabolism is the integration of pathways of transporting and oxidation of main carbon sources into the cell (Papagianni 2012). It uses a complex series of enzymatic steps to convert sugars into metabolic precursors. These precursors are then used to generate the entire biomass of the cell (Noor et al. 2010). Abundances of genes in methanogenesis pathways The alignment result of the methanogenesis genes was showed in Fig. 2. The three types of methanogenic pathways are (1) hydrogenotrophic, (2) acetotrophic, (3) methylotrophic (Whiticar 1999). Acetotrophic pathway is known to be the major one for more than 70 % of methane production in most Fig. 2 Percentages of genes in methanogenesis pathways from the four ADS metagenomic datasets

engineering anaerobic digestion process (Yu et al. 2005). Generally, the genes involved in acetotrophic pathway have higher abundances (Fig. 2). However, the most abundant gene was in the hydrogenotrophic pathway, transferring formate to CO2. The abundances of genes in methylotrophic pathway were the least.

Discussion Annotation using tags Merging short reads into tags increased annotated percentage in both taxonomic and functional annotation because of the increased length of the sequences (Huson et al. 2007), which were close to 200 bp. As indicated by the results, the increase of functional annotation percentage was remarkable. The average annotation percentage in functional annotation increased from 2.87 % in reads to 38.24 % in tags. Sequences with longer length can be annotated into finer taxonomic levels in taxonomic analysis as well. One drawback of analysis using tags is that some reads cannot be merged and therefore the throughput of the dataset was compromised. The merging percentage of the reads largely depended on the strategy in the library construction. In the present study, 61.43 % to 84.47 % of the reads can be merged into longer tags (Table S2), suggesting a good representation. Taxonomic analysis of anaerobic sludge The phyla of Proteobacteria, Bacteroidetes and Firmicutes were abundant in both ST and SWH samples (Fig. 1), followed by Bacteroidetes, Actinobacteria, Thermotogae and Chloroflexi. Rivière et al. (2009) suggested that Proteobacteria, Bacteroidetes, Firmicutes and Chloroflexi were the dominant groups in seven anaerobic sludge digesters using PCR approach. In a study of anaerobic digestion of

Appl Microbiol Biotechnol

carrot waste (Garcia et al. 2011), the major phyla were Firmicutes (45.4 %) and Bacteroidetes (31.3 %) while Proteobacteria took only 4.9 %. In a laboratory-scale biogas reactor being fed with kitchen waste, pig manure and sludge, Firmicutes is the most prevalent phylum, followed by Proteobacteria and Bacteroidetes (Li et al. 2013). This indicated that the abundant Proteobacteria in ADS from ST and SWH was possibly resulted from the feed sludge since Proteobacteria was quite rich in activated sludge (Zhang et al. 2011). Proteobacteria and Firmicutes are also important microbes in anaerobic digesters because Beta-proteobacteria and Gamma-proteobacteria were the propionate, butyrate, and acetate-utilizing microbial communities and Firmicutes was the butyrate-utilizing microbial community (Ariesyady et al. 2007). These two classes of Proteobacteria were abundant in the present study (Table S6A) as well. The phylum of Actinobacteria is known to express different ligninases (Jaenicke et al. 2011). The high relative abundance of Actinobacteria than many other phyla in the ADS samples also suggested that they had played an important role in the anaerobic digestion process. The abundances of Thermotoga were similar in ST1109, ST1203 and SWH1203, but relatively lower in SWH1109. Further analysis showed that the major genus in Thermotoga was Kosmotaga in ST and Fervidobacterium in SWH (Fig. 3). A previous study reported that nearly all of the Thermotoga were affiliated with Kosmotoga in ST ADS (Ye and Zhang 2012). Kosmotoga was a new genus proposed recently, which contains members isolated from marine environments, including oil production fluid in the North Sea (DiPippo et al. 2009), and a shallow hydrothermal system occurring within a coral reef (Nunoura et al. 2010). The salinity of the wastewater in ST (1.1 % salinity) (Zhang et al. 2011) might be one of the reasons why Kosmotoga was the dominant genus in Thermotoga in ST other than Fervidobacterium in SWH. Core genera in the ADS samples In the previous studies, microbial communities were commonly investigated by PCR using specific primers designed for Bacteria or Archaea (Rivière et al. 2009) and cannot be compared directly. However, using metagenomic approach, the abundances of Bacteria and Archaea can be revealed on one united platform without PCR bias, although the sensitivity of metagenomic approach was compromised by the sequencing depth. Even though the ADS samples were collected from two digesters with similar operation and performance, the shared genera in ADS samples were limited under the current sequencing depth since microbial community composition is very complex in the full-scale anaerobic digesters. Only 20 genera were shared by the two ST samples and 16 genera were shared by the two SWH samples (Fig. S5 and Table S10).

Only seven of those genera were shared by the four ADS samples, including Bacteroides, Clostridium, Mycobacterium, Planctomyces, Spirochaeta, Syntrophomonas in the domain of Bacteria, and Methanosaeta in the domain of Archaea. The total abundances of the core genera in the samples were relatively low, ranging from 1.94 % (ST1203) to 3.36 % (ST1109). These core genera belonged to the phyla of Bacteroidetes, Firmicutes, Actinobacteria, Planctomycetes and Spirochaetes in the Bacteria domain. However, a study using PCR approach suggested that the core microorganisms in Bacteria domain were affiliated with Chloroflexi, Betaproteobacteria, Bacteroidetes and Synergistetes in ADS (Rivière et al. 2009), which was not consistent with the metagenomic results. Considering both the potential PCR bias and the low sensitivity of metagenomic approach, more works are needed to study the actual core microorganism compositions in ADS samples, which is critical to understand the relationship among biodiversity, operational conditions and digester performance. Functional analysis of the ADS samples The major Level 1 subsystems in ADS were carbohydrates, protein metabolism, amino acids and derivatives, besides clustering-based subsystem (Table S8). Several other studies also suggested similar major subsystems in the various ecosystems, including grassland (Delmont et al. 2012), kimchi (Jung et al. 2011), marine, fresh water and microbialite (Breitbart et al. 2009). Principal component analysis (PCA) (Fig. 4) was carried out based on the Level 1 subsystems distribution to compare five different environmental samples, including ADS, AS, oceans, human faeces and soil (Table S11). The functional characteristic of ADS was close to human faeces and distant from AS and ocean ecosystems, which was reasonable since ADS and human faeces are both in anaerobic environment. PCA is able to identify which variables are driving the separation between different kinds of samples. Functions of clustering-based subsystems, cell wall and capsule, virulence, disease and defense subsystems appeared to be more significant related to ADS, human faeces and soil than in other ecosystems. In contrast, sequences related to protein metabolism, respiration, and amino acid and derivatives were more represented in AS, indicating the significance of functional difference between aerobic environment and anaerobic environment. Heterogeneous microenvironment may exist within one system. Previous studies showed the presence of heterogeneous microenvironment within sediment biofilm (Nguyen et al. 2012) and AS (Han et al. 2012). The dissolved oxygen (DO) concentrations within the center of the flocs from activated sludge were less than at the surface (Han et al. 2012), which make the center of flocs an anoxic microenvironment.

Appl Microbiol Biotechnol Fig. 3 Heat map of percentages of to 30 most abundant genera in ADS samples. The neighborjoining tree was drawn with MEGA 5.2

However, samples of AS were still distant from other anaerobic environment according to the PCA result based on the functional distribution in different environment, which suggested the importance of integrating metagneomic analysis with biogeochemical measurements for microbial ecology studies.

Fig. 4 Principle components analysis (PCA) of the five ecosystems using the abundances of Level 1 subsystems. Data included metagenomic data of human faeces, ocean and soil, activated sludge (AS) and ADS samples in this study. Detailed information on the ecosystems is showed in Table S11

Methane metabolism and methanogens The abundances of the genes in methane metabolism were similar in all of the ADS samples (Fig. 2). Most of sequences matched to two genes in hydrogenotrophic pathway and two genes in acetotrophic pathway. The two genes in

Appl Microbiol Biotechnol

hydrogenotrophic pathway were genes encoding formate dehydrogenase (EC: 1.2.1.2) and formylmethanofuran dehydrogenase (EC: 1.2.99.5). Both of which are involved in the initial step of the hydrogenotrophic pathway (Li et al. 2013). The two abundant genes in acetotrophic pathway included genes encoding acetyl-CoA synthetase (EC: 6.2.1.1) and acetyl-CoA decarbonylase/synthase complex (ACDS), which are essential in the synthesis of acetyl-CoA from acetate (Li et al. 2013). This result supported that these genes play a key role in generation of methane in anaerobic digesters. As indicated in Fig. 2, genes in acetotrophic pathway were generally more abundant, which was consistent with the previous study that acetotrophic pathway was the major pathway in methanogenesis (Yu et al. 2005). However, it should be kept in mind that the abundance of genes in methanogenesis pathway was discussed on the DNA level, which cannot show the expression level of the genes or the activity of the relevant enzymes. Metatranscriptomics or metaproteomics are needed to explore the expression of genes or activity of the enzymes involved in the methanogenesis pathway in future study. Among methanogens, the order of Methanosarcinales is aceticlastic. Two genera in this order were detected in the ADS samples, Methanosaeta and Methanosarcina. Methanosaeta was abundant in the samples of ST1109, ST1203 and SWH1109, but very low in SWH1203. However, Methanosarcina was abundant in SWH1203 and it was not detectable in the other three samples (Fig. 3). These two genera are different in the physiology and biokinetics. The genus of Methanosaeta has a high affinity for acetate but a relatively low maximum acetate utilization rate. It would tend to live in environment with low acetate concentration. On the other hand, the genus of Methanosarcina has a much lower affinity for acetate with a very higher substrate utilization rate. It would be favored in the environment with high acetate concentration (Yu et al. 2005). The variation of abundances of Methanosaeta and Methanosarcina was highly possible due to the acetate concentration in the anaerobic digester, which was affected by the treatment efficiency. In this study, we had characterized the microbial communities, performed analysis of the functions from the metagenomic data of the of four ADS samples. To the authors’ knowledge, this is the first study of metagenomic analysis of full-scale anaerobic sludge digester from municipal wastewater treatment plants. Based on the analyzed results, the following conclusions were drawn. (1) The phyla of Proteobacteria, Bacteroidetes, Firmicutes and Actinobacteria were the predominant phyla in ADS, playing important roles in the anaerobic sludge digestion process. (2) Functional analysis showed that the abundances of the Level 1 subsystems defined in SEED had no significant difference among different environmental samples, although PCA showed that three Level 1 subsystems: clustering-based subsystems, cell wall and capsule and virulence, disease and defense subsystems were more represented in ADS. (3) The

results of genes in methanogenesis pathway shown that acetotrophic was the major methanogenic pathway. It was also consistent with the taxonomic results that Methanosaeta and Methanosarcina, the aceticlastic methanogens, were the major genus in the detectable methanogens in ADS samples. Acknowledgments This work is supported by Hong Kong General Research Fund (7198/10E). Ying Yang, Ke Yu, and Yu Xia thank The University of Hong Kong for the postgraduate studentship. We thank the Drainage Services Department of the Government of the Hong Kong Special Administrative Region for their support for this study. Conflict of interest The authors declare that they have no conflict of interest.

References Albertsen M, Hansen LBS, Saunders AM, Nielsen PH, Nielsen KL (2011) A metagenome of a full-scale microbial community carrying out enhanced biological phosphorus removal. ISME J 6(6):1094– 1106. doi:10.1038/ismej.2011.176 Altschul SF, Madden TL, Schäffer AA, Zhang J, Zhang Z, Miller W, Lipman DJ (1997) Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res 25(17):3389–3402 Ariesyady H, Ito T, Okabe S (2007) Functional bacterial and archaeal community structures of major trophic groups in a full-scale anaerobic sludge digester. Water Res 41(7):1554–1568. doi:10.1016/j. watres.2006.12.036 Braguglia CM, Gagliano MC, Rossetti S (2012) High frequency ultrasound pretreatment for sludge anaerobic digestion: effect on floc structure and microbial population. Bioresour Technol 110:43–49. doi:10.1016/j.biortech.2012.01.074 Breitbart M, Hoare A, Nitti A, Siefert J, Haynes M, Dinsdale E, Edwards R, Souza V, Rohwer F, Hollander D (2009) Metagenomic and stable isotopic analyses of modern freshwater microbialites in Cuatro Ciénegas, Mexico. Environ Microbiol 11(1):16–34. doi:10.1111/j. 1462-2920.2008.01725.x Cole JR, Wang Q, Cardenas E, Fish J, Chai B, Farris RJ, Kulam-SyedMohideen AS, McGarrell DM, Marsh T, Garrity GM, Tiedje JM (2009) The Ribosomal Database Project: improved alignments and new tools for rRNA analysis. Nucleic Acids Res 37(Database issue): D141–D145. doi:10.1093/nar/gkn879 Delmont TO, Prestat E, Keegan KP, Faubladier M, Robe P, Clark IM, Pelletier E, Hirsch PR, Meyer F, Gilbert JA, Le Paslier D, Simonet P, Vogel TM (2012) Structure, fluctuation and magnitude of a natural grassland soil metagenome. ISME J 6:1677–1687. doi:10.1038/ ismej.2011.197 DeSantis TZ, Hugenholtz P, Larsen N, Rojas M, Brodie EL, Keller K, Huber T, Dalevi D, Hu P, Andersen GL (2006) Greengenes, a chimera-checked 16S rRNA gene database and workbench compatible with ARB. Appl Environ Microbiol 72(7):5069–5072. doi:10. 1128/AEM.03006-05 DiPippo JL, Nesbo CL, Dahle H, Doolittle WF, Birkland NK, Noll KM (2009) Kosmotoga olearia gen. nov., sp. nov., a thermophilic, anaerobic heterotroph isolated from an oil production fluid. Int J Syst Evol Micr 59(12):2991–3000. doi:10.1099/ijs.0.008045-0 Fang H, Cai L, Yu Y, Zhang T (2013) Metagenomic analysis reveals the prevalence of biodegradation genes for organic pollutants in activated sludge. Bioresour Technol 129:209–218. doi:10.1016/j.biortech. 2012.11.054

Appl Microbiol Biotechnol Garcia SL, Jangid K, Whitman WB, Das KC (2011) Transition of microbial communities during the adaption to anaerobic digestion of carrot waste. Bioresour Technol 102(15):7249–7256. doi:10. 1016/j.biortech.2011.04.098 Garcia-Peña EI, Parameswaran P, Kang DW, Canul-Chan M, KrajmalnikBrown R (2011) Anaerobic digestion and co-digestion processes of vegetable and fruit residues: process and microbial ecology. Bioresour Technol 102(20):9447–9455. doi:10.1016/j.biortech. 2011.07.068 Gomez-Alvarez V, Teal TK, Schmidt TM (2009) Systematic artifacts in metagenomes from complex microbial communities. ISME J 3(11): 1314–1317. doi:10.1038/ismej.2009.72 Han Y, Liu J, Guo X, Li L (2012) Micro-environment characteristics and microbial communities in activated sludge flocs of different particle size. Bioresour Technol 124:252–258. doi:10.1016/j.biortech.2012. 08.008 Hess M, Sczyrba A, Egan R, Kim TW, Chokhawala H, Schroth G, Luo S, Clark DS, Chen F, Zhang T, Mackie RI, Pennacchio LA, Tringe SG, Visel A, Woyke T, Wang Z, Rubin EM (2011) Metagenomic discovery of biomass-degrading genes and genomes from cow rumen. Science 331(6016):463–467. doi:10.1126/science.1200387 Huson DH, Auch AF, Qi J, Schuster SC (2007) MEGAN analysis of metagenomic data. Genome Res 17(3):377–386. doi:10.1101/gr. 5969107 Jaenicke S, Ander C, Bekel T, Bisdorf R, Dröge M, Gartemann K-H, Jünemann S, Kaiser O, Krause L, Tille F, Zakrzewski M, Pühler A, Schlüter A, Goesmann A (2011) Comparative and joint analysis of two metagenomic datasets from a biogas fermenter obtained by 454pyrosequencing. PLoS ONE 6(1):e14519. doi:10.1371/journal. pone.0014519 Jung JY, Lee SH, Kim JM, Park MS, Bae JW, Hahn Y, Madsen EL, Jeon CO (2011) Metagenomic analysis of kimchi, a traditional Korean fermented food. Appl Environ Microb 77(7):2264–2274. doi:10. 1128/aem.02157-10 Kim J, Park C, Kim T-H, Lee M, Kim S, Kim E-W, Lee J (2003) Effects of various pretreatments for enhanced anaerobic digestion with waste activated sludge. J Biosci Bioeng 95(3):271–275 Knight R, Jansson J, Field D, Fierer N, Desai N, Fuhrman JA, Hugenholtz P, van der Lelie D, Meyer F, Stevens R, Bailey MJ, Gordon JI, Kowalchuk GA, Gilbert JA (2012) Unlocking the potential of metagenomics through replicated experimental design. Nat Biotechnol 30(6):513–520. doi:10.1038/nbt.2235 Krzywinski M, Schein J, Birol I, Connors J, Gascoyne R, Horsman D, Jones SJ, Marra MA (2009) Circos: an information aesthetic for comparative genomics. Genome Res 19(9):1639–1645. doi:10. 1101/gr.092759.109 Li A, Yn C, Wang X, Ren L, Yu J, Liu X, Yan J, Zhang L, Wu S, Li S (2013) A pyrosequencing-based metagenomic study of methaneproducing microbial community in solid-state biogas reactor. Biotechnol Biofuels. doi:10.1186/1754-6834-6-3 Ma LP, Li B, Zhang T (2014) Abundant rifampin resistance genes and significant correlations of antibiotic resistance genes and plasmids in various environments revealed by metagenomic analysis. Appl Microbiol Biot. doi:10.1007/s00253-014-5511-3 Mackelprang R, Waldrop MP, DeAngelis KM, David MM, Chavarria KL, Blazewicz SJ, Rubin EM, Jansson JK (2011) Metagenomic analysis of a permafrost microbial community reveals a rapid response to thaw. Nature 480(7377):368–371. doi:10.1038/ nature10576 Mardis ER (2008) The impact of next-generation sequencing technology on genetics. Trends Genet 24(3):133–141. doi:10.1016/j.tig.2007. 12.007 Metzker ML (2009) Sequencing technologies — the next generation. Nat Rev Genet 11(1):31–46. doi:10.1038/nrg2626 Meyer F, Paarmann D, D’Souza M, Olson R, Glass EM, Kubal M, Paczian T, Rodriguez A, Stevens R, Wilke A, Wilkening J,

Edwards RA (2008) The metagenomics RAST server – a public resource for the automatic phylogenetic and functional analysis of metagenomes. BMC Bioinforma 9(1):386. doi:10.1186/1471-21059-386 Nguyen HD, Cao B, Mishra B, Boyanov MI, Kemner KM, Fredrickson JK, Beyenal H (2012) Microscale geochemical gradients in Hanford 300 Area sediment biofilms and influence of uranium. Water Res 46(1):227–234. doi:10.1016/j.watres.2011.10.054 Noor E, Eden E, Milo R, Alon U (2010) Central carbon metabolism as a minimal biochemical walk between precursors for biomass and energy. Mol Cell 39(5):809–820. doi:10.1016/j.molcel.2010. 08.031 Novak JT, Chon DH, Curtis B-A, Doyle M (2007) Biological solids reduction using the cannibal process. Water Environ Res 79(12): 2380–2386. doi:10.2175/106143007x183862 Nunoura T, Hirai M, Imachi H, Miyazaki M, Makita H, Hirayama H, Furushima Y, Yamamoto H, Takai K (2010) Kosmotoga arenicorallina sp. nov. a thermophilic and obligately anaerobic heterotroph isolated from a shallow hydrothermal system occurring within a coral reef, southern part of the Yaeyama Archipelago, Japan, reclassification of Thermococcoides shengliensis as Kosmotoga shengliensis comb. nov., and emended description of the genus Kosmotoga. Arch Microbiol 192(10):811–819. doi:10. 1007/s00203-010-0611-7 Overbeek R, Begley T, Butler RM, Choudhuri JV, Chuang HY, Cohoon M, de Crecy-Lagard V, Diaz N, Disz T, Edwards R, Fonstein M, Frank ED, Gerdes S, Glass EM, Goesmann A, Hanson A, IwataReuyl D, Jensen R, Jamshidi N, Krause L, Kubal M, Larsen N, Linke B, McHardy AC, Meyer F, Neuweger H, Olsen G, Olson R, Osterman A, Portnoy V, Pusch GD, Rodionov DA, Ruckert C, Steiner J, Stevens R, Thiele I, Vassieva O, Ye Y, Zagnitko O, Vonstein V (2005) The subsystems approach to genome annotation and its use in the project to annotate 1000 genomes. Nucleic Acids Res 33(17):5691–5702. doi:10.1093/nar/gki866 Palatsi J, Illa J, Prenafeta-Boldú FX, Laureni M, Fernandez B, Angelidaki I, Flotats X (2010) Long-chain fatty acids inhibition and adaptation process in anaerobic thermophilic digestion: batch tests, microbial community structure and mathematical modelling. Bioresour Technol 101(7):2243–2251. doi:10.1016/j.biortech.2009.11.069 Palatsi J, Viñas M, Guivernau M, Fernandez B, Flotats X (2011) Anaerobic digestion of slaughterhouse waste: main process limitations and microbial community interactions. Bioresour Technol 102(3):2219–2227. doi:10.1016/j.biortech.2010.09.121 Papagianni M (2012) Recent advances in engineering the central carbon metabolism of industrially important bacteria. Microb Cell Fact 11(50) Pérez-Elvira SI, Nieto Diez P, Fdz-Polanco F (2006) Sludge minimisation technologies. Rev Environ Sci Bio Technol 5(4):375–398. doi:10. 1007/s11157-005-5728-9 Pfister CA, Meyer F, Antonopoulos DA (2010) Metagenomic profiling of a microbial assemblage associated with the California Mussel: a node in networks of carbon and nitrogen cycling. PLoS ONE 5(5): e10518. doi:10.1371/journal.pone.0010518 Prosser JI (2010) Replicate or lie. Environ Microbiol 12(7):1806–1810. doi:10.1111/j.1462-2920.2010.02201.x Pruesse E, Quast C, Knittel K, Fuchs BM, Ludwig W, Peplies J, Glockner FO (2007) SILVA: a comprehensive online resource for quality checked and aligned ribosomal RNA sequence data compatible with ARB. Nucleic Acids Res 35(21):7188–7196. doi:10.1093/nar/ gkm864 Qin J, Li R, Raes J, Arumugam M, Burgdorf KS, Manichanh C, Nielsen T, Pons N, Levenez F, Yamada T, Mende DR, Li J, Xu J, Li S, Li D, Cao J, Wang B, Liang H, Zheng H, Xie Y, Tap J, Lepage P, Bertalan M, Batto J-M, Hansen T, Le Paslier D, Linneberg A, Nielsen HB, Pelletier E, Renault P, Sicheritz-Ponten T, Turner K, Zhu H, Yu C, Li S, Jian M, Zhou Y, Li Y, Zhang X, Li S, Qin N, Yang H, Wang J,

Appl Microbiol Biotechnol Brunak S, Doré J, Guarner F, Kristiansen K, Pedersen O, Parkhill J, Weissenbach J, Antolin M, Artiguenave F, Blottiere H, Borruel N, Bruls T, Casellas F, Chervaux C, Cultrone A, Delorme C, Denariaz G, Dervyn R, Forte M, Friss C, van de Guchte M, Guedon E, Haimet F, Jamet A, Juste C, Kaci G, Kleerebezem M, Knol J, Kristensen M, Layec S, Le Roux K, Leclerc M, Maguin E, Melo Minardi R, Oozeer R, Rescigno M, Sanchez N, Tims S, Torrejon T, Varela E, de Vos W, Winogradsky Y, Zoetendal E, Bork P, Ehrlich SD, Wang J (2010) A human gut microbial gene catalogue established by metagenomic sequencing. Nature 464(7285):59–65. doi:10.1038/ nature08821 Quast C, Pruesse E, Yilmaz P, Gerken J, Schweer T, Yarza P, Peplies J, Glockner FO (2012) The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res 41(D1):D590–D596. doi:10.1093/nar/gks1219 Rincón B, Borja R, González JM, Portillo MC, Sáiz-Jiménez C (2008) Influence of organic loading rate and hydraulic retention time on the performance, stability and microbial communities of one-stage anaerobic digestion of two-phase olive mill solid residue. Biochem Eng J 40(2):253–261. doi:10.1016/j.bej.2007.12.019 Rivière D, Desvignes V, Pelletier E, Chaussonnerie S, Guermazi S, Weissenbach J, Li T, Camacho P, Sghir A (2009) Towards the definition of a core of microorganisms involved in anaerobic digestion of sludge. ISME J 3(6):700–714. doi:10.1038/ismej.2009.2 Sahlstrom L (2003) A review of survival of pathogenic bacteria in organic waste used in biogas plants. Bioresour Technol 87:161–166 Shi P, Jia S, Zhang X-X, Zhang T, Cheng S, Li A (2013) Metagenomic insights into chlorination effects on microbial antibiotic resistance in drinking water. Water Res 47(1):111–120. doi:10.1016/j.watres. 2012.09.046 Shin SG, Lee S, Lee C, Hwang K, Hwang S (2010) Qualitative and quantitative assessment of microbial community in batch anaerobic digestion of secondary sludge. Bioresour Technol 101(24):9461– 9470. doi:10.1016/j.biortech.2010.07.081

Sun L, Randall CW, Novak JT (2010) The influence of sludge interchange times on the oxic–settling–anoxic process. Water Environ Res 82(6):519–523. doi:10.2175/106143009x12487095236711 Supaphol S, Jenkins SN, Intomo P, Waite IS, O’Donnell AG (2011) Microbial community dynamics in mesophilic anaerobic codigestion of mixed waste. Bioresour Technol 102(5):4021–4027. doi:10.1016/j.biortech.2010.11.124 Traversi D, Villa S, Lorenzi E, Degan R, Gilli G (2012) Application of a real-time qPCR method to measure the methanogen concentration during anaerobic digestion as an indicator of biogas production capacity. J Environ Manag 111:173–177. doi:10.1016/j.jenvman. 2012.07.021 Tringe SG (2005) Comparative metagenomics of microbial communities. Science 308(5721):554–557. doi:10.1126/science.1107851 Whiticar MJ (1999) Carbon and hydrogen isotope systematics of bacterial formation and oxidation of methane. Chem Geol 161(1–3):291–314 Ye L, Zhang T (2012) Bacterial communities in different sections of a municipal wastewater treatment plant revealed by 16S rDNA 454 pyrosequencing. Appl Microbiol Biot 97(6):2381–2690. doi:10. 1007/s00253-012-4082-4 Yu K, Zhang T (2013) Construction of customized sub-databases from NCBI-nr database for rapid annotation of huge metagenomic datasets using a combined BLAST and MEGAN approach. PLoS ONE 8(4):e59831. doi:10.1371/journal.pone.0059831 Yu Y, Lee C, Hwang S (2005) Analysis of community structures in anaerobic processes using a quantitative real-time PCR method. Water Sci Technol 52(1–2):85–91 Zhang T, Shao M-F, Ye L (2011) 454 Pyrosequencing reveals bacterial diversity of activated sludge from 14 sewage treatment plants. ISME J 6(6):1137–1147. doi:10.1038/ismej.2011.188 Zhou H-W, Li D-F, Tam NF-Y, Jiang X-T, Zhang H, Sheng H-F, Qin J, Liu X, Zou F (2010) BIPES, a cost-effective high-throughput method for assessing microbial diversity. ISME J 5(4):741–749. doi:10. 1038/ismej.2010.160