Microbial Communities and Sulfate-Reducing Microorganisms Abundance and Diversity in Municipal Anaerobic Sewage Sludge Digesters from a Wastewater Treatment Plant (Marrakech, Morocco)

Both molecular analyses and culture-dependent isolation were combined to investigate the diversity of sulfate-reducing prokaryotes and explore their role in sulfides production in full-scale anaerobic digesters (Marrakech, Morocco). At global scale, using 16S rRNA gene sequencing, Proteobacteria, Bacteroidetes, Firmicutes, Actinobacteria, Synergistetes, and Euryarchaeota were the most dominant phyla. The abundance of Archaea (3.1–5.7%) was linked with temperature. The mcrA gene ranged from 2.18 × 105 to 1.47 × 107 gene copies.g−1 of sludge. The sulfate-reducing prokaryotes, representing 5% of total sequences, involved in sulfides production were Peptococcaceae, Syntrophaceae, Desulfobulbaceae, Desulfovibrionaceae, Syntrophobacteraceae, Desulfurellaceae, and Desulfobacteraceae. Furthermore, dsrB gene ranged from 2.18 × 105 to 1.92 × 107 gene copies.g−1 of sludge. The results revealed that exploration of diversity and function of sulfate-reducing bacteria may play a key role in decreasing sulfide production, an undesirable by-product, during anaerobic digestion.


Introduction
Anaerobic digestion is an effective way of energy production through the decomposition of waste. This technology is applied in the treatment of not only domestic wastewater [1] allowing an important production of energy, recovered as biogas, composed of 50-75% methane, 25-50% carbon dioxide, 0-10% nitrogen, 0-3% hydrogen sulfide, and traces of other gases [2]. Because of the wide variety of starting materials, a complex set of microbial populations, with high functional redundancy, is involved in the process of anaerobic digestion [3,4]. This is one of the reasons for the robustness of this process [5,6]. Nowadays, the application of next-generation sequencing technologies provides increased resolution for the study of microbial communities in large-scale anaerobic digesters [7,8].
The overlapping of both molecular-and culture-dependent approaches could be complementary to each other to give a clear vision of the microbial community present in an ecosystem [9].
The anaerobic digesters may have operational troubleshooting that results in increased operational costs, among them foam and scum production, dewatering difficulties, loss of efficiency, acidification, and production of toxic sulfides [10]. Indeed, sulfide production is essentially related to the dissimilatory reduction [11,12] of various sulfur compounds present in the sewage sludge by a heterogeneous metabolic group of sulfate-reducing prokaryotes (SRP). They have the ability to use sulfate, as well as other oxidized sulfur compounds (elemental sulfur, thiosulfate, sulfite), as a terminal electron acceptor [12] for the oxidation of simple organic compounds or hydrogen, in a process leading to a sulfide production [13,14]. Sulfides can be an undesirable by-product during biogas production that may reach very high concentrations ranging from 500 ppmv to 20,000 ppmv (2% v/v of the biogas; [15]). Sulfides induced the corrosion of metal and concrete [16,17], significantly shorten the life of gas engines in cogeneration units [18] and the emission of strong, smelly [19], and toxic emanations [20]. In addition, combustion of biogas with high concentrations of hydrogen sulfide causes pollution due to excess of SO x emissions. Furthermore, sulfide has an inhibitory effect on the anaerobic digestion process [21]. Indeed, for these reasons, the elimination of sulfides (purification of biogas) is a mandatory step before using biogas as a source of energy [15,21].
The wastewater treatment plant (WWTP) of Marrakech, Morocco brings an additional capacity of 33 million m 3 of water reuse, mainly for irrigation. The sludge is treated by anaerobic digestion and produces 18,000 m 3 .day −1 of biogas cogenerated into 30,000 KWh.day −1 supplying 45% of electrical energy requirements of the station. The produced biogas (40-65%) contains an average of 4000 ppm of sulfides originating from sulfate-reducing microorganism activities that must be treated (desulfurized) before use. In Marrakech, the elimination of sulfide from biogas is carried out through Thiopaq technology [22], a cost-effective process. The current remediation technologies thus treat symptoms rather than causes, as their mode of action is to reduce the H 2 S concentrations rather than decrease its production [15]. While major sulfidogenesis pathways have been well described, the microorganisms responsible for this metabolic activity in anaerobic digesters remain poorly characterized [23]. The critical step in elucidating the conditions that are favorable to H 2 S production in an environment would be to characterize its inhabiting population of sulfate-reducing bacteria (SRB) [24].
In this context, the present study describes an investigation of anaerobic digesters treating domestic wastewater in order to gain more insight on the factors that may affect H 2 S production and explore the potential involvement of SRB in this process. This study focuses on the diversity and abundance (through combination of molecular and culture dependent approaches) of the microbial communities found in the Marrakech WWTP. The high throughput sequencing analysis of the V4-V5 region of the 16S rRNA gene has been undertaken as well as quantitative PCR (qPCR) to reveal the distribution of Archaea and sulfate-reducing microorganisms during anaerobic digestion to evaluate their role in the production of methane and sulfides, respectively. In addition, some strains belonging to the sulfate-reducers have been isolated and eventually described [25] in order to gain a better view of the sulfate-reducing community, a group of economic interest for the WWTP from Marrakech.

Sampling Strategy and Operating Conditions
Sampling was carried out on 24 July 2015 from the municipal WWTP of Marrakech. Sludge samples were taken from the surface and the bottom of the four-anaerobic digesters designated DA, DB, D1 and D2, the primary sludge (PS, from primary decanters), and the flocculated sludge (FS, from thickeners). A total of 10 different samples were collected directly from the sampling valve into sterile-sampling vessels. To maintain anaerobic conditions of samples, the 120 mL sterile serum bottles were completely filled to displace air and then stored at 4 • C until further use. In addition, six replicates from each sample were stored in cryotubes at −80 • C for molecular analyses (MiSeq and qPCR). All anaerobic digesters were operated at variable mesophilic temperature conditions (Table 1) and had each a working volume of 6400 m 3 . The treatment capacity of WWTP of Marrakech are summarized in Table S1, whereas the data for biogas digesters are given in Table 1.

Extraction of Genomic DNA and High Throughput Sequencing
Total DNA extractions were performed in triplicate using PowerSoil ®® DNA Isolation kit (MO BIO) following the manufacturer instructions. The DNA was stored at −20 • C until further use. The V4-V5 hypervariable region of the 16S rRNA gene targeting Bacteria and Archaea was amplified using the primers set 515F and 928R [26][27][28]. The amplifications of the 16S rRNA gene were performed in a final volume of 50 µL of PCR mixture with 0.2 µM of each primer, 0.2 mM of deoxynucleoside triphosphates (dNTPs), 1.25 U Taq polymerase (Amp Taq Gold 360), and its buffer. The amplification reactions are carried out according to the following program: 1 cycle of 10 min at 94 • C (denaturation); 35 cycles of 3 steps: denaturation at 94 • C, hybridization at 65 • C, and elongation 40 s at 72 • C; final elongation cycle of 10 min at 72 • C. Amplification was conducted using a PCR thermal cycler Model PTC-100 (Bio-Rad, Richmond, CA). Forward and reverse primers contained the adapter sequences (5 GTGYCAGCMGCCGCGGTA) and (3 ACTYAAAKGAATTGRCGGGG), respectively, used in Illumina sequencing technology. The amplicons were sequenced using an Illumina MiSeq (250 bp paired-end reads) at the platform Genotoul Inc. (www.genotoul.fr, Toulouse, France).
Raw datasets of sequencing can be found into the NCBI Sequence Read Archive database (Bioproject accession number, PRJNA587578).

Processing of Sequencing Data and Statistical Analyses
The datasets from MiSeq sequencing were processed following the FROGS pipeline (Find, Rapidly, OTUs with Galaxy Solution) as described by Escudie et al., [29]. Operational taxonomic units (OTUs) were grouped at an identity threshold of 97%. From the taxonomic affiliation, a representative sequence of each OTU was extracted. Finally, a normalization of the abundances found for each OTU was carried out by the DESeq method [30]. Diversity and evenness indices [31] were calculated using the bacterial OTU data after normalization. The evolutionary history was inferred using the neighbor joining method [32]. An alignment was then carried out between OTUs; the alignment was used to build a phylogenetic tree with PhyML algorithm [33].
The statistical calculations and multivariate analyses were conducted using R software (https: //cran.r-project.org/) with the following libraries: phyloseq, vegan, ape, Venn diagram, heatmap2, gclus, ggplot2, GUniFrac, grid, and optparse. A one-way analysis of variance (ANOVA) was used to compare the mean values with respect to sampling sites. Duncan's multiple range tests at p = 0.05 were employed to compare the means and to group them accordingly. Pearson correlation coefficients were calculated from pairs of variables within and between process and bacterial community data. A principal component analysis (PCA) was performed to evaluate the microbial community composition of the samples.

Growth, Isolation of Sulfate-Reducing Microorganisms and Phylogenetic Identification
Enrichment cultures were inoculated with samples using basal medium [25] and supplemented with 10 mM sulfate as an electron acceptor and 10 mM lactate for the Desulfovibrio and Desulfomicrobium genera isolation or 10 mM propionate for the Desulfobulbus genus isolation. After three sub-culturing, the isolation of strains was performed using repeated agar serial dilutions [34]. The colonies were transferred from the highest dilutions to liquid media. The purity of the strains was checked routinely by phase-contrast microscopy and by performing growth tests in basal medium containing 1% (w/v) glucose, 1% (w/v) casamino acid, 1% (w/v) peptone, 1% (w/v) yeast extract and 1% (w/v) biotrypticase under oxic and anoxic conditions. All genomic DNA extractions from pure cultures were performed using a kit (UltraCleanTM Soil DNA Kit insulation MOBIO Laboratories, Carlsbad, CA, USA) according to manufacturer's recommendations. The amplifications of the 16S rRNA gene was carried out as described by Dias et al., [35]. Sequencing with the 8F and 1492R primers of the amplicons obtained was performed (GATC Biotech-Konstanz, Germany). An alignment was then carried out between OTUs sequences, sulfate-reducing bacteria reference sequences obtained from the SILVA database (version 128) and the isolated strain sequences. The alignment was used to build a phylogenetic tree with MEGA6 software (Molecular Evolutionary Genetics Analysis software, Version 6.0, 2013). The evolutionary history was inferred using the neighbor joining method [32] and conducted in MEGA6 software [36].

Quantitative PCR
Primer sets and reaction conditions were performed as described by Geets et al., [37] and Steinberg and Regan [38] for both genes dsrB and mcrA, respectively. The quantifications were performed in triplicate on a Mx3005P™ real-time PCR machine (Stratagene) in a final volume of 20 µL containing 1 µL of template DNA, 0.4 µM of each primer and the enzyme mixture of SYBR ®® Green DyNAmo™ Color Flash qPCR Kit (Finnzymes). The amplifications started with an initial denaturation step (10 min at 95 • C.) and then 40 cycles of 1/denaturation (1 min, 95 • C.), 2/hybridization of the primers to the hybridization temperature (54 • C for dsrB gene and 57 • C for mcrA gene, 30 s) and 3/elongation (40 s, 72 • C). A melting curve of the PCR products obtained was carried out at the end of each amplification.

Global Analysis of Microbial Communities
The total of 10 samples were analyzed in triplicates from the WWTP of Marrakech, including primary sludge (PS), flocculated sludge (FS), and the four anaerobic digesters (D1, D2, DA, DB) sampled at the surface (s) and the bottom (b). 1,243,519 sequences were founded in all the samples after MiSeq sequencing. 1,095,879 among them were retained, corresponding to 88.1% of the effective sequences obtained with a mean length of 250 bp. The number of OTUs found at 97% similarity was 984. The OTU number threshold for DB-s was significantly below all other anaerobic digester samples with values of 695 ± 13 OTU (Table 2), whereas, the number of OTUs for the two sludges PS and FS was 593 ± 43 and 566 ± 57 respectively. The diversity indexes listed in Table 2 were lower for the sample D1 (with the highest temperature), whatever the selected indexes. The low values found for exponential Shannon index, and reciprocal Simpson index, compared to Chao1 indicated a highly uneven distribution of the populations in the samples.
The comparison of the relative abundance of populations at surface and bottom of each digester showed that there was no significant difference between the percentages of the most abundant families. These findings were in accordance with by the Permanova test for analysis of variance using Bray-Curtis distance matrices showing no significant differences between surface and bottom of each digester (p-value: 0.001). These results were confirmed by the principal component analysis ( Figure S1) based on the weighted UniFrac β-diversity metric showing that the triplicates for all samples as well as the surface and bottom samples for each anaerobic digester were closely grouped together. As a result, from now on, the surface and bottom samples will be grouped and presented together (mean of six samples). Considering all the samples, the FS samples could be separated from the other samples (PS and the four anaerobic digesters) along axis 1 explaining 49.68% of the difference ( Figure S1a). The axis corresponds to a gradient of anaerobiosis (aerated FS separated from PS and anaerobic digesters). FS consists of limited representatives of the Euryarchaeota and Synergistes ( Figure 1). The flocculated sludge is obtained after the activated sludge processes undertaken under highly aerated conditions to promote the growth of the aerobic microorganisms and discriminate anaerobic microorganisms such as Euryarchaeota and Synergistes. The profile of the most abundant microbial communities in PS was similar to the anaerobic digesters. This was consistent with the results of principal component analysis ( Figure S1a). PS contain obligate anaerobes and Firmicutes (20.2%) compared to the FS (only 4.4% of Firmicutes).
Another principal Component analysis focused on the four anaerobic digesters ( Figure S1b) revealed a clear separation along the first axis, explaining 41.28% of the difference. The four digesters were placed along this axis corresponding to a gradient of operating temperature in each digester (Table 1). According to the physico-chemical parameters obtained from the WWTP (Table 1), D1 had a higher temperature, a bit higher than the mesophilic conditions (average temperature in July 2015 of 41.1 • C, max value of 44.2 • C), whereas DA has the lower temperature (average temperature of 29.7 • C). These results were consistent with the richness analysis ( Figure 2) shared between the four anaerobic digesters. The D1 presented 28 OTUs (representing 1.97% of cumulative relative abundance) that were unique for this digester as compares with less than 0.3% for all the others. Among these 28 OTUs, the most dominants belong to Clostridia. Moreover, D1 had higher relative contents of true thermophiles Firmicutes (OTU 1) represented by Coprothermobacter genus belonging to the Thermodesulfobacteriaceae family ( Figure 3).

Methanogens Communities Structure in Anaerobic Sewage Sludge
The domain Archaea was represented by only one phylum ( Figure 1). Among 984 OTUs, 31 OTUs were assigned the Euryarchaeota. Sequences from Euryarchaeota accounted for 0.3%, 3.1%, 4.9, 4.4%, 5.4%, and 5.7% in FS, PS, D1, D2, DA, and DB, respectively. The passage from the primary sludge with an anoxic gradient towards the anaerobic digesters having strict anaerobic conditions makes it possible to favor and enrich the anaerobic communities of which the Archaea is a part. The methanogens found out in the fours anaerobic digesters were classified into four well-established orders: Methanobacteriales, Methanomicrobiales, Methanosarcinales, and Thermoplasmatales. At the family level, the Archaeal OTUs could be classified into 7 different families that were detected present in the anaerobic digesters:  (Figure 4a). Based on qPCR analysis, it was also observed that the abundance of Archaea at the surface of all digesters (except DB) is greater than that at the bottom.

Diversity of Sulfate-Reducing Bacteria in Anaerobic Digesters
The OTUs that were identified and classified as sulfate-reducing are presented in a heatmap based on their abundance in the four digesters ( Figure 5). The analysis of the relative abundance of SRP into the four anaerobic digesters at the family level showed that Peptococcaceae, Desulfobacteraceae, Desulfobulbaceae, Desulfovibrionaceae, Desulfurellaceae, Syntrophaceae, and Syntrophobacteraceae were the different widespread sulfate-reducing families. The taxonomic affiliation at the genus level of the dominant Firmicutes (OTU 1 and OTU 25) found out in the anaerobic digesters were affiliated to their respective genera Coprothermobacter and Syntrophonas that are reported as non-sulfate-reducers ( Figure 3). However, the OTU 18 affiliated to Desulfosporosinus and genus member of the Peptococcaceae family was found out as the only of the sulfate-reducing candidate from Firmicutes.
The global profile of the relative abundances at OTUs level of SRP shown on the heatmap indicates that the D1 is different from the other digesters ( Figure 5). As for Archaea, the abundance of the SRP was lower in D1 (operating in higher temperature) compared to the other digesters. The Syntrophaceae family (represented by OTU 176, OTU 55, OTU 711, OTU 1051, OTU 30, OTU 244, OTU 103, and OTU 688) was the most abundant in digester 2, B, and A. The OTU 30 was abundant in DA and DB (0.53% and 0.50%; respectively) followed by the OTU 55 abundant in digester 2 compared to the DB, DA and D1 (0.50%, 0.23%, 0.39%, and 0,31%; respectively). However, the D1 was dominated by the Syntrophobacteraceae family represented by OTU 283, followed by OTU 55 (0.34%, 0.31%; respectively).
The family Desulfovibrionaceae (represented by OTU 349, OTU 1164, OTU 888) was not abundant in different digesters compared to the other families. However, The Desulfobulbaceae family (represented by OTU 264, OTU 226, OTU 882, OTU 330, OTU 587, OTU 107, OTU 897, and OTU 325) was the most abundant in DA compared to the other digesters. The OTU 107 (Desulfobulbaceae family) was the most abundant with following respective relative abundances: 0.28%, 0.28%, 0.34%, and 0.32%; D1, D2, DA, DB. As for the Archaea, the SRP abundance was more important at the surface of each digester compared to the bottom (Figure 4b). In addition, a large amount of dsrB gene copies (1.92 × 10 7 ) was detected in the primary sludge compared to the other samples. The abundance of sulfate-reducing microorganisms (dsrB gene quantification) was ranged from 9.57 × 10 5 to 9.25 × 10 6 copies.g −1 of sludge in anaerobic digesters.

Relationship of Sulfate-Reducing Communities in Anaerobic Digesters with Isolated Sulfate-Reducing Strains from Anaerobic Digesters Sludge
As an alternative approach to identifying active sulfate-reducing populations in anaerobic digesters, a phylogeny relationship was carried out between the OTU sequences founded in anaerobic sludge into all samples by high throughput sequencing, the sulfate-reducing bacteria reference sequences obtained from the SILVA database (version 128), and sequences of isolated sulfate-reducing strains from anaerobic digester sludge in this study ( Figure 6). The enrichment and isolation on lactate and propionate as electron donors in the presence of sulfate as an electron acceptor has made it possible to target three different genera of sulfate-reducing bacteria. Enrichment on lactate allowed the isolation of 12 strains affiliated with the genus Desulfovibrio and 5 strain with the genus Desulfomicrobium. However, in the enrichment on propionate, 5 Desulfobulbus strains were isolated. Desulfobulbus strains (Prop1, Prop4, Prop6 T , C3P4-2, and C3P4-3) were affiliated to the Desulfobulbus oligotrophicus strain Prop6 T isolated and previously described as new species within the genus Desulfobulbus.
The 13 Desulfovibrio isolated strains were affiliated to Desulfovibrionaceae family. One of these strains (DZ1) was affiliated to Desulfovibrio aminophilus strain ALA-3 T and was most related to the OTU 349 (99.3% of similarity); three isolated strains PBS3, PBS4, and PBS5 were affiliated to Desulfovibrio oxamicus strain Montecillo 2 T ; and the last 9 isolated strains were affiliated to Desulfovibrio vulgaris DP4 T . The isolated Desulfomicrobium strains (Halo1, Halo4, Halo6, DB4 and DB5 were affiliated to Desulfomicrobium escambiens esc1 T . However, no OTUs linked to the isolated Desulfomicrobium strains or Desulfovibrio vulgaris strains could be detected using molecular approaches.

Discussion
Primary sludge (PS), flocculated sludge (FS), and the four anaerobic digesters (D1, D2, DA, DB) were analyzed from WWTP Marrakech using molecular and culture-dependent approaches. Using MiSeq sequencing, 984 OTUs were found out present in all samples. These results are in the range of those found by other studies carried out with the same primers (V4-V5 region of the 16S rRNA gene) [3,39,40]. The study of dynamics of microbial association networks in mesophilic continuous anaerobic digesters revealed the presence of a total of 908 OTUs [40]. An average of 1036 ± 219 OTUs were identified when studying the microbial community structure and diversity in a municipal solid waste landfill [39]. The Hill numbers [31] were lower for the sample D1 (with the highest temperature) reflecting the lowest richness. These results were in accordance with previous findings revealing that thermophilic reactors have a lower number of species than mesophilic [41,42].
The microbial community dynamics were revealed to be heavily influenced by agitation [43] affecting the distance between syntrophic microbes and therefore increase the diffusion distance between synthrophs [44]. The analysis of distribution of microbial populations in anaerobic digesters based on principal component analysis ( Figure S1) and the Permanova test for analysis of variance using Bray-Curtis distance (p-value: 0.001) showed was no significant difference between the percentages of the most abundant families between surface and bottom samples for each digester. These results revealed that the mixing system applied in WWTP Marrakech digesters (agitation bay re-pumping biogas) works perfectly assuming a good distribution of microbial populations within bottom and surface of the digesters, the comparison of the relative abundance of the most abundant populations at surface and bottom of each digester based on the MiSeq sequencing showed that there was no significant difference between the percentages of the most abundant families between surface and bottom samples for each digester. The Permanova test for analysis of variance using Bray-Curtis distance matrices showed no significant differences between surface and bottom of each digester (p-value: 0.001). However, when considering all the samples, the FS samples were separated along axis 1 ( Figure S1a) corresponding to a gradient of anaerobiosis and contained limited representatives of the Euryarchaeota and Synergistes (Figure 1), phyla known to contain exclusively obligate anaerobes [45]. PS has a similar profile closed to the anaerobic digesters ( Figure S1a). The active Bacteroidetes and Firmicutes found out in PS and PS as well as anaerobic digesters intervene during the hydrolysis process of anaerobic digestion and were known to be involved in hydrolytic/acidogenic step in reactors operating under mesophilic and thermophilic conditions [46][47][48]. Previous studies reported that WWTPs was dominated by Bacteroidetes aerobic carbon degraders [49] dominant in aerated sludges [50]. Indeed, the FS samples contained a dominance of the Proteobacteria (38.9%), the Bacteroidetes (19.4%), The Actinobacteria (10.2%), the Chloroflexi (6.6%), and Planctomycetes (6.2%) as reported by [51].
The temperature and feeding substrate were the main driving forces of the anaerobic digestion microbiome diversification [52,53]. These results were in accordance with the richness analysis ( Figure 2) shared between the four anaerobic digesters showing the unique separation of the thermophilic D1 dominated by the Firmicutes (Thermodesulfobiaceae) and reported to dominate in methanogenic reactors alongside the Rikenellaceae and the Bacteroidia [54][55][56]. The four anaerobic digesters revealed a clear separation ( Figure S1b) depending on the gradient of operating temperature in each digester (Table 1) obtained from WWTP Marrakech.
The Proteobacteria, the Bacteroidetes, the Firmicutes, the Actinobacteria, the Synergistetes, the Euryarchaeota, and Cloacimonetes were the dominant phyla in WWTP Marrakech. These data are similar to those previously obtained using marker gene-based analysis [39,40,52,53,55,[57][58][59]. At the OTU level (Figure 3), OTU 1 (Firmicutes), OTU 2 (Cloacimonetes), and OTU 3 (Bacteroidetes) were the most abundant of the four digesters. The characterization of microbial community structure during continuous anaerobic digestion of straw and cow manure showing the Firmicutes, Bacteroidetes, and Cloacimonetes were identified as the three most dominant within the bacterial community [53,60]. Most of the members belonging to the Firmicutes phylum are syntrophic bacteria that can degrade various VFAs, which were often detected in both activated sludge systems and anaerobic digesters [61].
The Firmicutes were found as the most dominant phyla involved in polysaccharides degrading functions [53]. The Cloacimonetes, represented by OTU 2, has been suggested to be involved in the degradation of cellulose, in syntrophic propionate oxidation, and in the degradation of amino acid to acetate, carbon dioxide, and hydrogen, and thus depends on cooperation with hydrogen-consuming and/or acetate-utilizing partners such as Archaea [62,63].
The Euryarchaeota belonging to the Methanobacteriales, the Methanomicrobiales, the Methanosarcinales, and Thermoplasmatales were the main Archaea found out in the WWTP of Marrakech (Figure 1). The abundance of Archaea in anaerobic digesters was in concordance to those of previous studies [61,64,65] showed that Archaea was among the most abundant phyla in anaerobic digesters. The dissecting microbial community structure and methane-producing pathways, Methanobacteriales, Methanococcales, Methanomicrobiales, Methanosarcinales, and Methanopyrales were revealed to be present in a full-scale anaerobic reactor digesting activated sludge from wastewater treatment [64]. The study of the abundance of Archaea by quantification of mcrA gene (Figure 4) in sewage sludges of WWTP was similar to various types of sewage sludge, used in the continuous anaerobic digestion, and collected from eight WWTPs municipal domestic wastewater [66,67].
The analysis of the relative abundance of SRP at the family level showed that Peptococcaceae, Desulfobacteraceae, Desulfobulbaceae, Desulfovibrionaceae, Desulfurellaceae, Syntrophaceae, and Syntrophobacteraceae were the main sulfate-reducing microorganisms involved in sulfides production into the four anaerobic digesters ( Figure 5). These results were consistent with a study investigating the SRB populations between cattle manure and digestate in the elucidation of H 2 S production during anaerobic digestion of animal slurry [23]. The clone libraries analysis of dsrAB genes that revealed that the dominant genera among the biofilms from moving bed biofilm reactor wastewater treatment plants were Desulfomonile (Syntrophaceae), Desulfomicrobium (Desulfomicrobiaceae), Desulfococcus (Desulfobacteraceae), Desulfobulbus (Desulfobulbaceae), and Desulforhopalus (Desulfobulbaceae) [68].
The Firmicutes (OTU 1 and OTU 25) found out in the anaerobic digesters were affiliated to their respective genera Coprothermobacter and Syntrophonas that are reported as non-sulfate-reducers ( Figure 3). However, the OTU 18 was affiliated to Desulfosporosinus, a genus member of the Peptococcaceae family, and reported to be a candidate of sulfate-reducing Firmicutes [13]. Most SRB have been described from the Firmicutes phylum (represented in this study by OTU 18), including Desulfosporosinus acidiphilus, Desulfosporosinus acididurans, and Thermodesulfobium narugense [13]. It was reported that some candidates from Firmicutes (Desulfotomaculum subclusters, Moorella thermoacetica, and Ammonifex degensii) acquired dsrAB gene from a deltaproteobacterial donor, a gene involved in the sulfate-reducing pathway [69,70]. According to the microbial community resilience and regarding sulfides production, Syntrophobacter and Syntrophus were the most abundant genera in laundry wastewater and domestic sewage pilot-scale [71]. D1 is separated from the other digesters and dominated by Syntrophobacteraceae. However, D2, DA, and DB were dominated by Syntrophaceae. Sulfate-removal from wastewaters can be performed by means of anaerobic bioreactors, which contain Syntrophaceae and Desulfovibrionaceae [72,73].
The family Desulfovibrionaceae was not abundant in different digesters compared to the other families ( Figure 5). However, the Desulfobulbaceae family was abundant in DA compared to the other digesters. Desulfovibrio, Desulfobulbus, Desulfomicrobium were reported as the most widespread in the bioreactors of wastewater treatment plants [74]. Desulfobulbus sp. was the dominant member of the SRB community in all stages of biofilm growth [75,76], and were reported to grow on propionate and sulfate [25,77,78] and were able to survive under oxic conditions [79].
At a global aspect, the SRP as well as Archaea abundances were affected by high temperature (See Figures 3 and 5). Our results were consistent with another study showing that both methanogenesis and sulfidogenesis were largely suppressed under thermophilic relative to mesophilic conditions [80]. SRB were more abundant at the surface of each digester compared to the bottom (Figure 4b) and were most abundant in PS (1.92 × 10 7 dsrB gene copies). Desulfobulbus and Desulfovibrio species have been found to be the numerically important SRB in wastewater biofilms grown under oxic conditions [75,81] and activated sludge flocs [76,82]. The reason for the higher abundance and activity of these SRB species is due to their tolerance to O 2 [83,84]. This may explain the large amount of SRB detected in the primary sludge compared to the other samples. These results were similar to the average number of dsrA gene copies in biofilm samples from full-scale municipal wastewater treatment plants using a moving bed biofilm anaerobic reactor. The 16S rRNA gene pyrosequencing data of putative SRB-identified Desulforhopalus, Desulfomicrobium, and Desulfobacter as dominant members of the community. The family Syntrophaceae, (genera Desulfomonile, Desulfobacca, and other, non-SRB organisms) was absent from the 16S rRNA gene dataset, which reaffirms the importance of dsrAB gene-based analyses to characterize SRB communities [68]. SRB in anaerobic wastewater treatment facilities can contribute via continuous H 2 S formation to the production of undesired odorous volatile compounds [85]. Treatment of paper mill effluents involves an up-flow anaerobic sludge blanket reactor with active SRB participating in organic matter decomposition [21]. SRB of moving bed biofilm reactors possess a fixed biofilm and can be used in varying scales for the treatment of domestic and industrial waste [86].
The phylogeny relationship carried out between the OTU sequences, the sulfate-reducing bacteria reference sequences and sequences of isolated sulfate-reducing strains (obtained by culture-dependent approaches) were studied. From Desulfobulbus isolated strains ( Figure 6), Desulfobulbus oligotrophicus Prop6 T was isolated and described as new species within the genus Desulfobulbus [25] and were related to the OTU 882. In the previous study comparing the diversities of these genera at the beginning of sampling and after cultivation. Desulfovibrio, Desulfobulbus, and Desulfomicrobium genus as dominant among sulfate-reducers in the bioreactors were detected [74,[87][88][89]. Desulfobulbus genus can incompletely oxidize numerous organics including glucose, propionate, and pyruvate to acetate allowing the sulfidogenic anaerobic digestion of sulfate waste from anaerobic sludge [88]. It was reported that the Desulfobulbus are able to grow on propionate and sulfate [25,90,91]. However, in the absence of sulfate, they can ferment other organic acids and alcohols such as lactate, ethanol, and pyruvate [68].
The Desulfovibrionaceae family (represented by 13 Desulfovibrio isolated strains and affiliated to the OTU 349) was effective for sulfate removal and both hydrogen and sulfides production during anaerobic digestion [92]. Members of the completely oxidizing Desulfobacteraceae could be demonstrated by cultivation methods to dominate the SRB communities [93,94]. The Desulfomicrobium strains or Desulfovibrio vulgaris strains could not be detected using molecular approaches (no OTUs linked). However, the culture-dependent approaches allowed to isolate candidates from these both genera, showing the importance of the combination of molecular and culture-dependent approaches to gain more insight into the abundance of SRB. Desulfomicrobium and Desulfovibrio was widely discovered in natural environments and anaerobic reactors, which degraded lactic into acetate [95][96][97]. The Desulfobacteraceae and Desulfurellaceae OTUs found in anaerobic digesters were considered as key players in the carbon-and sulfur-cycling in organic carbon-and sulfate-rich reduction [98]. Their potential in treating industrial waste and drainage has been investigated under various settings [13]. The Syntrophobacteraceae OTUs were involved in the syntrophic degradation of propionate in an up-flow anaerobic sludge blanket (UASB) bioreactor granules [99]. Many other studies have shown that Syntrophobacter species are key organisms for the syntrophic oxidation of propionate in bio-granules [100][101][102].

Conclusions
Microbial sulfate-reduction is a process of great environmental and biogeochemical importance, which is associated with different environments because of their high levels of sulfate. Microbial sulfate-reduction initiates and supports the sulfur cycle and is one of the main biological processes. For these reasons, SRB play a key role in our understanding of the biogeochemical cycles of sulfur and carbon in different environments such as anaerobic digesters. The 16S rRNA high-throughput sequencing data, the qPCR quantification of dsrB and mcrA gene combined with the isolation approaches in this present study revealed a high diversity of microbial communities including SRB and methanogens in the anaerobic sludge from municipal wastewater treatment plant from Marrakech. The classification of SRB and communities based on 16S rRNA gene and functional gene was analyzed and compared according to the results of a previous study [57]. Overall, the present study is one of the few studies combined molecular to culture-dependent approaches to provide more comprehensive information regarding the SRB and Archaea community characteristics and suggested the existence of an active sulfur cycle in anaerobic sludge. In the last few years, our knowledge regarding archaeal microorganisms has increased from both an ecological standpoint [27,103]. The qPCR method proposed is a useful tool for studying active methanogenic populations and their influence on the biogas production rate. Our results showed an important abundance of methanogens compared to the previous study [61,64,65]. On a global scale, the Proteobacteria (that include the sulfate-reducing microorganisms responsible for sulfides production), the Bacteroidetes, the Firmicutes, the Actinobacteria, the Synergistetes, and the Euryarchaeota (methanogens microorganisms) are the most widespread microbial communities. The molecular and culture-dependent approaches allowed detected the Syntrophaceae, the Desulfobulbaceae, the Desulfovibrionaceae, the Syntrophobacteraceae, the Desulfurellaceae, and the Desulfobacteraceae as the main families responsible for sulfides production in anaerobic digestion of WWTP of Marrakech. Real-time quantitative PCR of dsrA and mcrA genes was used to better understand the abundance of SRB and active methanogenic communities compared to other microorganisms in the anaerobic digesters. A considerable abundance of SRB and Archaea were detected in anaerobic sludge. The results show that the temperature was among the key operating conditions that may affect the rate of methane and sulfide production. In addition, these results give more insight to understand the diversity and function of SRP in order to investigate the key role that they may play to affect sulfides production, an undesirable by-product, during anaerobic digestion.