Microbial Community Analysis Using MiSeq Sequencing and Pathway of Methane Production in Tehran WWTP: A Full-Scale Anaerobic Digester

1 Research Center for Environmental Health Technology, Iran University of Medical Sciences, Tehran, Iran. 2 Department of Environmental Health Engineering, School of Public Health, Iran University of Medical Sciences, Tehran, Iran. 3 Department of Civil Engineering, Amirkabir University of Medical Sciences, Tehran, Iran. 4 Managing Director of Tehran Water and Wastewater Co., Tehran, Iran. 5 Research Fellow, the Experimental Forest, College of Bio-Resources and Agriculture, National Taiwan University, Nantou 55750, Taiwan. 6 Student Research Committee, Saveh University of Medical Sciences, Saveh, Iran.


Introduction
Activated sludge process uses a biological process to treat municipal and industrial wastewater 1 . This process is noticed by process designers due to its high efficiency in removing/decreasing most pollutants established as the discharge standards to acceptable levels 2,3 . The wastewater treatment plant operation reports show that the sludge generation through biological wastewater treatment is the main activated sludge side effects because this process needs more than 50-60% of the integrated wastewater management budget 4 . The activated sludge process should be moderated to decrease sludge generation. Anaerobic digestion is highly efficient in sludge treating and bioenergy generating simultaneously [5][6][7] . By applying this process, the treated sludge had acceptable organic matter, pathogens, and nutrient content for use in land application. In addition, the sludge quantity decreased and led to easy final management and generated methane as a valuable biogas. Generally, the anaerobic digestion process is capable of converting 50% bio-solids (in average) to methane 8 . The microbial communities in the sludge are the heart of efficient methane generation through anaerobic digestion and sludge reduction. In other words, we need to have enough information about the sludge microbial community to have a selective pathway to operate and generate methane. Furthermore, if some microbial community is predominant, the pathways change that leads to unpleased results. Generally, assessing the microbial community of sludge is one of the necessary steps to have an efficient and cheap anaerobic sludge digestion process 9,10 . Based on the available reports, some molecular approaches exist to find predominant microbial communities in anaerobic systems. Among them, Illumina sequencing and 454 pyrosequencing technologies at high-throughout sequencing method category are the most efficient and user-friendly methods to microbial community analysis anaerobic digested sludge [11][12][13] . By exploring the literature, we found the Illumina sequencing as the most interesting method to assess microbial communities among soil, ocean, human gut microbes, and activated sludge 14,15 . The use of Illumina sequencing to assess anaerobic sludge community is rare, but its application is rising. By in-depth analysis of the microbial community structure, the microbial functions can be identified in anaerobic sludge digestion and methane production pathway. The present study aimed to explore the microbial community arrangement in a full-scale anaerobic digester of Tehran wastewater treatment plant. To this end, DNA extraction was done from a fullscale anaerobic digestion sludge. The dominant microorganisms, functional profiles, and methane production pathway from AD of sludge were investigated. Especially, significant microorganisms involved in sludge digestion and methane production were analyzed comprehensively according to the acquired molecular data. This study aimed to assess the microbial community structure of a full-sized wastewater treatment having AD for sludge managing and methane production, simultaneously.

Sampling of AD sludge from a full-scale WWTP
In this study a full-scale WWTP, located at south of Tehran, Iran selected. There is used AD to manage activated sludge. The south wastewater treatment plant of Tehran is located in a region of about 110 hectares in the south of Rey. The plan is envisioned for the population of 2,100,000 people. The total capacity of the plant is 5.2 m 3 /s having four modules. Wastewater treatment method is activated sludge conventional system. The extra sludge from the primary sedimentation goes into gravity sludge thickening. Later, the excess sludge from the secondary sedimentation was thickened using a mechanical thicker and went into the sludge digester. The studied fullscale anaerobic reactor digesting activated sludge including thickening tanks, anaerobic mesophilic digestion and dewatering. Figure 1 presents wastewater treatment process diagram.

DNA extraction, library construction and sequencing
At first, the sampled sludge mixed with ethanol (with purity of ~100%) at 1:1 ratio (v/v %). Afterward, it was transferred to lab in standard condition (using a cool box with -20 0 C temperature). All the qualified DNAs (that were checked by quality test) were applied to construct libraries. Gel electrophoresis (1% agarose gel, Voltage: 150 V; Electrophoresis Time: 40 min) was applied to evaluate the quality of used to specify the concentrations of DNA. The DNA concentration of the sample was 52.6 ng/μL. The microbial community sequencing was done with Illumina Hi Seq/Mi Seq platform (by NIGEB Genome Center, China). Meta rDNA Amplicon sequencing includes Meta rDNA V3, V6, V4, V1-V3, V3-V4, ITS1, and ITS2 regions amplicon library construction. The paired-end reads with overlap were merged to tags and the tags were clustered to the Operational Taxonomic Unit (OUT) at 97% sequence similarity by scripts of USEARCH software (v7.0.1090). The library was sequenced on an Illumina Hi Seq/Mi Seq platform. If the two paired-end reads overlapped, the consensus sequence would be generated by FLASH (Fast Length Adjustment of Short reads, v1.2.11). The parameters adopted for overlapping were the minimal overlaping length of 15 bp and the allowed mismatching ratio of overlapped region was 0.1.

Microbial Community analysis
USEARCH software (v. 7.0.1090) scripts were used to cluster the tags into OTU. The tags were clustered into OTU with a 97% threshold using UPARSE and the OTU unique representative sequences were obtained. UCHIME (v4.2.40) was used to filter out the chimeras. By matching the gold database (V20110519, UNITE (V20140703)) with 16S rDNA and ITS sequences distinctly, de novo chimera detection was done for 18S rDNA sequences. All tags were mapped to each OTU representative sequences using USEARCH GLOBAL and then the number of tags for each OTU in each sample was summarized to OTU

Species phylogenetic and Diversity analysis
The Silva core set (Silva_108_core_ aligned_ seqs) using PyNAST by 'align_seqs.py' used to represent lined up sequences. The QIIME (v1.80) built-in scripts were used to construct the representative OTU phylogenetic tree. The tags with a maximum abundance of each Genus were selected as the matching Genus-representative sequences. In addition, the Genus level phylogenetic tree was obtained the same as OTU phylogenetic tree. R (v3.0.3) software was used to draw the phylogeny tree. The rarefaction curve was drown using the indices calculated with extracted tags of the samples. Mothur (v1.31.2) was applied to calculate the indices and the corresponding rarefaction curve was drawn by R (v3.0.3) rarefaction curve. Each index was formulated and calculated by the following address: http://www. mothur.org/wiki/ Calculators. Furthermore, the microbial diversity indices were Chao-the Chao1 estimator, ACE-the ACE estimator, Shannon-the Shannon index, Simpson-the Simpson index, and goods coverage index that can be calculated using the following addressws, respectively (http://www.mothur. org/ wiki/Chao), (http://www. mothur.org/wiki/Ace), (http:// www. mothur.org/ wiki/Shannon),and (http://www. mothur. org/ wiki/ Simpson, (http:// www. mothur. org/ wiki/ Coverage). To analyze the complexity of the species' diversity, the Alpha diversity was used (http://www.mothur.org/wiki/ Schloss_SOP#Alpha_diversity).

Microbial Richness and Diversity Evolution
Microbial community analysis performed by MiSeq sequencing method. The diversity estimators of Chao1, ACE, and Shannon index, Simpson index, and P-value of the sequencing for each sample are shown in Table 1.
A representative OTU phylogenetic tree was constructed using the QIIME (v1.80) built-in scripts including the fast tree method for tree construction. The tags with the highest abundance of each Genus were chosen as the corresponding Genus representative sequences and the Genus level phylogenetic tree was obtained similar to OTU phylogenetic tree. The phylogeny tree was imaged by R (v3.0.3) software at last.

Discussion
In this study, the microbial community within a full-Scale AD sludge determined using the MiSeq sequencing method. Meanwhile, the pathways of methane production in the studied WWTP were proposed.

Microbial Richness and Diversity Evolution
Logged species value, chao1 value, and ACE value can reflect the species' richness of the community. It is possible to use the rarefaction curve according to the three values to assess if produced data is sufficient for covering all species in the community 16,17 . The community richness of the sludge had diverse values based on the community diversity of sludge. Rationally, the increased loading or changing the sludge property that may be as a result of alternation in wastewater biodegradable organic matter content, can affect the diversity and richness in the sludge, considerably 18 . The communities of sludge showed a minor decrease in richness and diversity for sample two. Accordong to table 1, the highest richness and diversity was logged for sample three. According to our results, the species richness and phylogenetic diversity increased for the bulk sludge; then, it reduced by increase of loading.

Rarefaction, Rank-abundance distributions
The OTU rank abundance curve makes it possible to show the species' richness and species' evenness visually 17 . According to figure 4, the species are ranked based on their richness . By obtaining the slope of fitted line, the species' evenness was counted. A steep gradient is an indication of low evenness 19 . The comparison of rank-abundance distributions reveals a shallow gradient highlighting high evenness since the abundances of various species are similar. When microbiota are gathered from a certain niche, it is encessary to assess how well a sample shows the true diversity of the specific niche. In other wods, it is identical with species' richness and relative abundance in time and space. Rarefaction curves measure the OTUs observed with a given depth of sequencing applied for comparing the observed richness among communities with an unequal samples 17 . As the table 1 shows the rarefaction result is in-line with Chao1, ACE, Shannon, and Simpson results. It can be concluded that the species richness and phylogenetic diversity have a similar alternation arrangement.
Regarding higher microbial diversity in the sludge cause enhancing microbial activity, so organic loading rate increasing unable to have considerable effect on Chemical Oxygen Demands (COD) efficiency 16,20  Here, HRTs and operation temperatures were similar; however, the organic loading rate had flactuation because of wastewater property.

Microbial Phylum and Class Identification
For bacteria, Firmicutes and Bacteroidetes were recurrently recognized as the chief phyla in many anaerobic digesters 23 . It is believed that Firmicutes and Bacteroidetes are significant in a hydrolytic/acidogenic digester fed with dried hay and straw 24 . It was also reported that these two phyla are vital during batch digestion. The phylum Bacteroidetes is a group of bacteria found in many different habitats, which live in soil or aquatic environments. It is also related to humans and animals that possess a high level of phenotypic and metabolic diversity 25 . Bacteroidetes are able to release more proteinaceous EPS that facilitate its colonization 16,26 . Their fimbriae helps them adhesion to other bacteria cells (colony forming) and to solid surfaces and improving their enrichment on the surfaces.
Several other studies also offered Bacteroidetes as one of the key microbial components in anaerobic reactors 27 . Bacteroidetes are fermentative bacteria that play a role in degradation of carbohydrates and proteins into acetate and NH 3 . The second chief bacterial phylum -Firmicutes represented 11% of the total bacterial clones in sludge sample. Firmicutes are hydrolytic bacteria that are capable of producing extracellular enzymes such as cellulases, lipases, and proteases 28,29 . Firmicutes is involved in the hydrolysis step of the anaerobic digestion (AD) sludge process by metabolizing a wide variety of compounds into acetate and butyrate (important VFA intermediates) 16 . According to Ding et al., Bacteroidetes, Chloroflexi, and Candidate with majority of 73.29% were the dominant phyla of the M-AnMBR sludge. As shown in figure 4, the abundance of Bacteroidetes in the current study is similar to the one studied by Ding et al. However, Chloroflexi in our study was lower (2%). As illustrated in figure 4, nine bacterial classes were detected in the sludge sample. Acidobacteria _ Gp 10, Anaerolineae, Bacteroidia, Betaproteobacteria, Clostridia, Deltaproteobacteria, Gammaproteobacteria, Synergistia, and Thermotogae were the main classes. The Clostridia, Thermotogae, Synergistia, and Gammaproteobacteria were the key four classes that accounted for 14%, 7.5%, 7%, and 5% of the total population in the sludge. Clostridia classes are bacilliform bacteria that respirate anaerobically. These class of bacterium was reported as one of the dominant bacterial class in the anaerobic sludge digesters 16 . The anaerobic sludge digestion can make a suitable environment for their growth 30 . The phylum Thermotogae includes anaerobic, thermophilic, as well as mesophilic bacteria. Some Thermotogae species are capable of surviving at high temperature. The bacteria of this class have metabolic ability to produce H 2 gas from complexcarbohydrates that makes them a likely candidate as a biotechnological source for CH 4 production 31 . Synergistia (phylum Synergistetes) that is one of the most abundant class of bacteria that are capable of decomposing complex organic matters efficiently. Spirochaetes is chemohetrotroph microorganism that have spiral shape. Regarding Spirochaetes can be found in AD reactor, so they able to recover methane through fermentation of carbohydrates or amino acids. Consequently, contribution of the above-mentioned three classes in the present study was possibly related to the organic matter removal within the digester. We also observed the Gammaproteobacteria is a class of the phylum Proteobacteria 32 . Some Gammaproteobacteria oxidize methane, photosynthetic, and hydrogen sulfide instead of water, producing sulfur as a waste product. Comparing our results with other systems with different in HRT and treating pattern showed that their microbial communities in the sludge were different with present study.

Microbial Genera Evolution
The literature reported high proportions of unidentified sequences at the genus level, which is confirmed in this study. Obviously, isolation and cultivation of uncultured microorganisms in the AD process 33-36 is a useful way to find out their role in the process, as recommended earlier. The genus Anaerobaculum, as the most overriding genus in genus level, comprises of Anaerobaculum thermoterrenum and Anaerobaculum mobile species and is a member of the phylum Synergistetes. The presence of Anaerobaculum strains in petroleum reservoirs, mining wastewater, solid waste digester, and methanogenic reactors were also reported. Anaerobaculum strains are fermentative organisms capable of producing acetate, CO 2 , and H 2 as their major metabolic products 37 . They are unable to reduce sulfate or sulfite to hydrogen sulfide. Sugihara et al. showed that Anaerobaculum-related species were dominant bacteria 38 . They concluded that Anaerobaculum-related species were able to degrade propionate. Consequently, one can say that higher abundance of Anaerobaculum in the studied digester possibly led to no accumulation of propionic acid and H 2 . According to the hydrolysis, acidogenesis, acetogenesis, and methanogenesis are the possible biogas production process stages 39 . At three first stages, the bacteria are principal organisms. CH 4 is formed through decarboxylation of acetate and methanisation of CO 2 and H 2 by acetogenotrophic and hydrogenotrophic archaea. Besides, synergistic acetate oxidation is involved in biogas formation. Moreover, during the hydrolysis of oils and fats, long-chain fatty acids (LCFA) are produced that are normally found in fattywastewaters. It is possible to convert these compounds effectively to methane in anaerobic bioreactors by applying the suitable condition. The LCFA converts to acetate and hydrogen by synchronizing the activity of syntrophic bacteria. In addition, the methanogenic archaea that utilizes these substrates is activated in this pathway. By prevailing this pathway, the complete LCFA degradation can be obtained. Here, the genus identification showed that Syntrophomonas sapovorans is the second most abundant species in the sludge sample. Meanwhile, studies show that 7 species of syntrophic bacteria are capable of growing on LCFA with more than 12 carbon atoms. Among these bacteria, only three mesophilic species are capable of utilizing the unsaturated LCFA, namely Syntrophomonas sapovorans (there was in studied sludge digester), Syntrophomonas curvata, and Syntrophomonas zehnderi. Generally, Syntrophomonas sapovorans is anaerobic, syntrophic and fatty acid-oxidizing, and obligately proton-reducing.
Based on our field data and microbial community data, a schematic diagram of the probable leading pathways can be made in the studied digester ( Figure  5). Acetate was metabolized via hydrogenotrophic methanogenesis by syntrophic acetate-oxidizing bacteria (SAOBs) in the studied digester. This is supported by the high abundance of Anaerobaculum in the digesters 40 , which was considered as syntrophic acetate oxidizers. Syntrophic acetate oxidation occurs under elevated temperatures and high ammonia concentrations, which was also represented in the studied digester (Table 2).   Issue (3), September 2020, 1091-102 Jehsd.ssu.ac.ir 1100 According to microbial data, it seems that acetoclastic methanogens is the key pathway for methane production. The gas analysis of produced biogas showed CH 4 63 ± 11, CO 2 35 ± 9%, and H 2 S 1800 ± 105 ppm. Presence of Syntrophomonas and Anaerobaculum in the digester help to have more effective SAO.

Conclusion
The sequencing study of this work was able to separate the microbial community and pathway of methane production in the full-scale AD of Tehran WWTP, successfully. According to our findings, 15 phyla were prevailing in the sludge of anaerobic digester, but Bacteroidetes and Firmicutes were the most abundance. However, the chief recognized genera were Anaerobaculum, Acinetobacter, Syntrophomonas, and Coprothermobacter, respectively. Exploring the leading bacteria and field data revealed that hydrogenotrophic methanogenesis acetate was metabolized by SAOBs. Generally, the findings of present study help the wastewater operators and process designers to use most effective method that able to treat waste sludge and methane production, simultaneously.