Methanogenic community during the anaerobic digestion of different substrates and organic loading rates

Abstract Three anaerobic reactors using pig manure (PM), maize straw (MS), and a mixture of the two as substrates were compared for archaeal community structure and diversity, and for methanogens response to increased organic loading rate (OLR, expressed in the mass of volatile solid (VS)). Methanogenic archaeal richness during codigestion of pig manure with maize straw (ACE: 2412) was greater than that during the others (ACE: 1225, 1467) at an OLR of 4 g L−1 day−1, accompanied by high specific methane yield. Euryarchaeota and Crenarchaeota predominated during overall digestion of different substrates; with relative abundances of 63.5%–99.0% and 1.0%–36.3%, respectively. Methanosarcina was the predominant genus that accounted for 33.7%–79.8% of the archaeal community. The diversity in the PM digester decreased with increase in OLR, but increased in the MS digester. The diversity was stable during the codigestion with increased OLR. The relative abundances of hydrogenotrophic methanogens increased by 2.6 and 2.1 folds; the methanogenic community shifted from acetoclastic to hydrogenotrophic methanogens during digestion of MS, and of the mixture of MS and PM. Canonical correspondence analysis revealed a strong relationship between reactor parameters and methanogenic community.

over other forms of waste treatment, such as less biomass sludge, minimal odor emissions, and so on (Smet, Van Langenhove, & De Bo, 1999;Ward, Hobbs, Holliman, & Jonesw, 2008). The characteristics of agricultural waste, such as imbalance in nutrition and high proportion of proteins or lignocellulosic biomass pose a challenge for process engineering. Furthermore, they affect the microbial community involved in AD, since the high concentration of ammonia likely inhibits methanogens activity, and high fiber content can cause blockage of pumps or sinking or floating layers (Munk, Guebitz, & Lebuhn, 2017).
AD is a complicated microbial process which involves four sequential steps: hydrolysis, acidogenesis, acetogenesis, and methanogenesis. Its performance and stability is strongly related with the microbial community structure. Archaea, especially methanogens, are key players during methanogenesis, thus attracting much attention of researchers (Li, Rui, et al., 2015). Methane is produced through hydrogen oxidation and acetate cleavage by hydrogenotrophic and acetoclastic methanogens, respectively . Munk et al. (2017) evaluated the potential utility of grass silage, which is rich in nitrogen. The result showed that the reactors could be operated stably as sole substrate at low OLR, although the ammonia concentration was high and hydrogenotrophic methanogens were dominant in the thermophilic and mesophilic reactors. Despite some novel findings, these results were very limited.
Due to the accumulation of inhibitory intermediates such as volatile fatty acids (VFAs) and ammonia, AD of agricultural waste is more prone to failure at high OLR compared to that at a low OLR. Investigation of the structure and dynamics of microbial community during the process of elevating the OLR should enable elucidation of the syntrophic interactions in AD, which may be used to optimize operational conditions. In this study, the typical agricultural wastes like MS, PM, and mixture of MS and PM were digested in laboratory-scale completely mixed anaerobic reactors, at different OLR of 2 and 4 g L −1 day −1 for 219 days.
Methanogenic community was investigated using high-throughput sequencing technology. The objective of this study was to compare the structure of methanogenic community in AD with different substrates (PM, MS, and PM+MS), analyze the changes in methanogens due to elevated OLR, and elucidate the link between methanogenic community and reactor performance. It is expected that the results presented herein would enable the optimization of operational conditions in order to achieve high efficiency AD for agricultural applications.

| Substrates and inoculum
Both PM and MS were collected from Yi Lilai Breeding Co. Ltd.
(Xiqing district, Tianjin, China). MS was smashed to fractions of nearly 1 mm in size by laboratory blender (Waring Commercial, USA). PM and MS were stored at 4 ± 1°C. The inoculum sludge was taken from a lab-scale completely mixed reactor that ran AD of PM, at 35°C. The characteristics of the substrates and inoculum sludge are shown in Table 1.

| Anaerobic reactor operation
The experiment was carried out in three completely mixed anaerobic digesters with 7 L working volume and mesophilic condition (35 ± 0.5°C). Three reactors were seeded with inoculum sludge (4 L) and adjusted to 7 L with tap water. The substrates were PM (R1), mixture of PM and MS (R2), and MS (R3), respectively. The VS ratio of PM/MS in R2 was 2:1, which maintained the ratio of C/N at around 26. The substrates were manually added to the reactors after discharging the digestate from the outlet each day. The reactors were run with a start-up OLR of 2 g L −1 day −1 , this was increased to 3 and 4 g L −1 day −1 on the 69th and 150th day at corresponding HRTs (60, 40, 30 days).

| Physicochemical analysis
Biogas was collected in an aluminum bag (20 L), the volume of which was measured with a wet gas meter every day. Gas and slurry samples were taken from the reactors at intervals of 3 day. Biogas composition was determined by GC (gas chromatography, Thermal Trace-1300) equipped with TCD (thermal conductivity detector).
SMY (specific methane yield, ml g −1 day −1 ) was the daily methane yield of per gram VS, the value equals daily biogas volume times methane content and divide by VS mass.

| DNA extraction and sequencing
Biomass was sampled at the end of each stage, running with OLR of 2 and 4 g L −1 day −1 , respectively, and frozen at −80°C

| Diversity analysis and statistical analysis
The raw sequences were sorted based on the unique sample barcodes, quality control for sample sequence using Prinseq, and removed nontarget region sequences and chimeric sequences using the Mothur program (Schloss et al., 2009) and Chimeras.uchime.
The sequences were clustered into operational taxonomic units (OTUs) defined by 97% similarity (3% max distance) using Uclust (Edgar, 2010). Venn diagram was charted according to the distribution of OTUs in the samples by Mothur and R language package "VennDiagram". Alpha and beta diversity were calculated using the QIIME software package (Caporaso et al., 2010). For alpha diversity analysis, Chao1, Shannon, ACE, Coverage, and Simpson indices were calculated by Mothur program, and rarefaction curves were generated. Species were classified by RDP classifier (Wang, Garrity, Tiedje, & Cole, 2007) according to the silva database, and statistics the relative abundance at phylum, class, order, and genus level. For beta diversity analysis, cluster analysis was preceded by 2D Principal Coordinates Analysis (PCoA) using the QIIME software package. QIIME calculates both unweighted and weighted UniFrac distance (Kuroda et al., 2015;Zhang, Sun, Zeng, Chen, & Sun, 2015).
Canonical correspondence analysis (CCA) was performed to discern the correlations between the archaeal community and operational parameters, including the substrate, OLR, NH 4 + -N, and VFA.

| Bioreactor performance
The physical and chemical parameters of the three reactors at different OLRs were shown in Table 2. These values were the average of every stage. The average specific methane yield (SMY) of R2 (PM and MS) was 218 and 254 ml g −1 day −1 at OLR of 2 and 4 g L −1 day −1 , respectively, which were 7.4%, 10.4%, and 33.7%, 192.0% higher, respectively, compared with that in R1 and R3. This difference could be explained by the optimal C/N ratio of 26 in R2. In R1 (PM), the average SMY at an OLR of 4 g L −1 day −1 reached 230 ml g −1 day −1 , which was 8.0% higher than that at OLR 2 g L −1 day −1 . This was consistent with the report by Bolzonella, Pavan, Battistoni, and Cecchi (2005) which demonstrated that a higher OLR (up to 4 g L −1 day −1 ) for shorter SRTs (15 days) generated higher methane levels under both mesophilic and thermophilic conditions. However, the value decreased by 46.6% in R3 (MS), this was supported by the research of Zhou et al. (2017), who reported that biogas yield improved when OLR increased at levels below 2 g L −1 day −1 during AD of rice straw.
Ammonia is an important nitrogen source for the growth of methanogens as well as a key pH-stabilizing agent for the neutralization of VFAs. However, high concentration of ammonia inhibit the activity of methanogens (Li, Liu, et al., 2015). Ammonia concentration The concentration of VFAs is about 50-250 mg L −1 in a stable anaerobic digester, and the inhibitory concentration is around 1,500 mg L −1 (Khanal, 2008). Other studies report that AD didn't fail until up to 10,000 mg/L of acetic acid or butyric acid (Khanal, 2008), this large difference may be attributed to the buffer capacity, which maintains the pH at the desired level, as well as to the methanogens' tolerance to toxic agents such as ammonia. In this study, the VFAs in three reactors changed when OLR was elevated: unlike the observed accumulation of VFAs in R1 and R3, acetate, propionate, and VFAs concentrations in R2 decreased by 18.0%, 83.3%, and 36.7%, respectively, as OLR increased from 2 to 4 g L −1 day −1 , and propionic acid concentration was below the inhibitory threshold value of 900 mg L −1 (Wang et al., 2016). The ratio of propionate/acetate is also an indicator of the stability of AD process. A ratio over 1.25 was reported to induce biomethanation process failure (Wang, Wang, Cai, & Sun, 2012). In this study, the ratio was under 1.25 over the experimental period.

TA B L E 2 Process parameters of the three reactors used in this study
Although the concentration of VFAs in R1 (965.5 mg L −1 ) was more than ten times that in R2 (94.4 mg L −1 ), the system was capable of sustaining its performance since the higher alkaline capacity of manure allows higher OLR without accumulation of VFAs. Cheng and Zhong (2014) studied the effect of codigestion on biogas production from AD of cotton stalk at a feed-to-inoculum ratio of 4:1 and 5% TS concentration of substrates, and determined that the concentration of VFAs was 320 and 3,000 mg/L during the monodigestion of cotton stalk and codigestion with PM, respectively.  Figure 1).

| Richness and diversity analysis of OTUs
The diversity estimators for each sample based on a species level of 97% similarity were showed in Table 3. Rarefaction curves of OTUs profiled the change in archaeal species richness in different feedstock and OLRs (Figure 2). R3 showed the highest community diversity, followed successively by R2 and R1, as confirmed by the Simpson and Shannon index (Table 3). Although none of the rarefaction curves approached a plateau, the Shannon's diversity index rarefaction curves approached asymptotes, indicating that the sampling depths were sufficient to capture the overall microbial diversities in all six samples (Zhang, Sun, et al., 2015). The mixing of different substrates improves the nutrients in the feedstock, thereby increase the growth rate of microbial organisms, and additionally enhances digestion efficiency of the system (Li, Cheng, et al., 2016). Anaerobic codigestion of PM with MS (R2) enhances the nutritional balance and reduces the possibility of inhibition induced by lipids and ammonia.
Interestingly, when the OLR was elevated, archaeal species richness decreased in R1, increased in R3, whereas that of R2 was relatively stable during the experiment.
To compare the methanogens community in the three samples at OLR of 4 g L −1 day −1 , the principal coordinates analysis (PCoA) were performed to cluster community. Based on unweighted UniFrac PCoA, the cluster community had the maximum variation in 63.1% (PC1) and 36.9% (PC2) (Figure 3a). It demonstrated a clear regional reparation, and samples from R1 and R2 tended to cluster together, whereas R3 was clearly different from them. These clustering results suggested that R1 and R2 shared similar methanogens community, which were clearly different from those in R3. Weighted UniFrac PCoA (Figure 3b) represented that R1 and R3 were grouped together along PC2 which only accounted for 9.8% of the total variations; however, the three samples were separated from each other according to PC1, which accounted for 90.2% of the total variations.

| The methanogenic archaeal community for different substrates
Archaeal community compositions of different substrates at OLR of 4 g L −1 day −1 were compared (Figure 4) Methanobacteriales accounted for 31.3% of total sequences in R2, which was higher than that for R1 (17.0%) and R3 (9.6%). The relative abundance of Methanomicrobiales was 1.1%, 5.8%, and 16.0% in R1, R2, and R3, respectively (Figure 4c). Although Methanobacteriales and Methanomicrobiales were both hydrogenotrophic methanogens, they showed different correlation pattern with regard to AD performance; Methanobacteriales appeared to play an important role in high loading AD, in contrast to the negative correlation of the Methanomicrobiales to biogas and OLR (Vrieze et al., 2015).
Among methanogens, Methanosarcina and Methanosaeta are well-known for utilizing acetate for methanogenesis. Methanosaeta is specialized in producing methane by acetate cleavage, whereas Methanosarcina is a relative generalist whose metabolic features are diverse and include both acetoclastic and hydrogenotrophic pathways (Li et al., 2013;Liu & William, 2008). Compared with  (Liu & William, 2008). In this study, the concentration of acetic acid in the three reactors was not in the optimal concentration range for Methanosaeta (though not so for Methanosarcina) and induced selective proliferation of Methanosarcina (Guo, Wang, Sun, Zhu, & Wu, 2014). In addition, the intermittent stirring in completely mixed reactor may be responsible for conferring an advantage to Methanosarcina since Methanosarcina was reported to be frequently predominant in fixed and stirred tank digesters (Liu & William, 2008).
Other methanogens, such as Methanosphaera, Methanospirillum, Methanoculleus, Methanobrevibacter, and Methanobacterium detected in this study are hydrogenotrophic methanogens. In R1, the relative abundance of acetoclastic methanogens was 80.2%, which was much higher than that of hydrogenotrophic methanogens (total relative abundance was 18.1%); therefore, it can be concluded that acetoclastic methanogens represent a major pathway for the digestion of PM. However, for R2 and R3, the relative abundances of acetoclastic methanogens were 46.9% and 34.8%, and that of hydrogenotrophic methanogens 37.2% and 23.5%, respectively. Therefore, acetoclastic and hydrogenotrophic pathways for methanogens occur at approximately equal extents during the digestion of the PM and MS mixture, or that of MS solely. This result was consistent with the finding that the methane-producing microbial community is involved in the anaerobic codigestion of pretreated wheat straw with cattle manure and solid state codigestion of kitchen waste, pig manure, and excess sludge (Li et al., 2013;Song & Zhang, 2015).
In addition, a large proportion of archaeal sequences belonging to class Thermoprotei, phylum Crenarchaeota were unclassified, especially in R3 (34.4%); this proportion was approximately equal to that of Methanosarcina (33.7%). Unfortunately, no further data on the predominance of this microorganism were available, owing to the lack of information regarding archaea in recent reports; however, F I G U R E 5 Relative abundance of methanogens 16S rDNA gene sequences of R1 (a), R2 (b), and R3 (c) at OLR of 2 g L −1 day −1 (left) and 4 g L −1 day −1 (right). The sequences showing a percentage of reads below 1.0% in all samples were grouped into 'Others' this microorganism may be significant for methanogens, especially for digestion with MS as mono or multiple substrates (Francisci, Kougias, Treu, Campanaro, & Angelidaki, 2015).

| The methanogenic archaeal community at different OLRs
The genus Methanosarcina represented the predominant phylotype under different OLRs ( Figure 5) Methanoculleus has been reported to show high identity with methanogenic archaea in stable anaerobic cellulose-degrading reactors (Chin, Lukow, Stubner, & Conrad, 1999). The total abundance of hydrogenotrophic methanogens increased to 37.2% from 17.7% when OLR was elevated. In R3, the relative abundance of Methanosarcina decreased from 46.9% to 33.7%, that of Methanoculleus decreased from 7.0% to 1.3%, and that of Methanobacterium increased from 0.8% to 9.4%. Methanospirillum, which was also present, showed a relative abundance of 12.6%. The total abundance of hydrogenotrophic methanogens increased from 8.9% to 23.5%.
Methanosarcina was predominant in the three reactors, and its advantage in R1 was enhanced with the increase in ammonia concentration from 2,482 to 4,384 mg L −1 . In previous studies (Li, Liu, et al., 2016), the stable CH 4 production that accompanies an increase in ammonia level may be explained by the increasing activity of hydrogen-utilizing methanogens. This is because hydrogenotrophic methanogens are capable of tolerating ammonia concentrations of 6,000 mg L −1 , which is two folds higher than the threshold ammonia concentration for Methanosarcina.
In this study, although the ammonia concentration increased strongly, the relative abundance of hydrogenotrophic methanogens decreased from 31.3% to 18.1%, and the SMY increased from 303 to 330 ml g −1 day −1 . Since the inoculum was cultured with high ammonia during AD feeding with PM, the ammonia in R1 (2,482-4,384 mg/L) had no inhibition for Methanosarcina. In addition, acetate has a crucial impact on the presence and relative abundance of acetoclastic methanogens (Yu et al., 2014); the acetate concentration in the three reactors was higher than the threshold value that favors the growth of Methanosaeta. With increased in OLR, the relative abundance of Methanosarcina decreased, whereas that of hydrogenotrophic methanogens in R2 and R3 increased by 2.1 and 2.6 folds, respectively.
Furthermore, the ratio of acetoclastic/hydrogenotrophic methanogens declined to 1.26 and 1.48 from 3.46 to 5.30, which indicated that the methanogenic pathway apparently shifted from mainly acetoclastic to the coexistence of acetoclastic and hydrogenotrophic methanogens pathways in the system (Goux et al., 2015). The coexistence of methanogenic pathways seems necessary for response to perturbation and maintenance of stable process performance (Lerm et al., 2012); this was evidenced by the performance of R2, in which the SMY increased by 13.3%, whereas it decreased by 46.6% in R3.

| The correlation between archaeal population and reactor parameters
Canonical correspondence analysis (CCA) was used to highlight the influence of altered process parameters on the archaeal community (Goux et al., 2015). As shown in Figure 6, the first and second  (Jang, Kim, Ha, & Park, 2014). These results agree with those for the taxonomic distribution of archaeal community ( Figure 5).
CCA revealed that substrate, ammonia, and VFA play key roles in determining the community structure; these factors are positively F I G U R E 6 Canonical correspondence analysis (CCA) of the archaeal community and the operational parameters in the three reactors. Blue vectors represent the influence of the process parameters such as OLR, SMY, substrate type, VFA, and NH 4 + -N; red vectors represent methanogenic archaea identified by highthroughput sequencing at the genus level. Black points represent the substrate type at different OLRs correlated with the abundance of Methanosarcina. On the contrary, the presence of Methanosaeta was negatively correlated with these factors, and with the abundance of Methanosarcina. These findings may be attributed to the greater tolerance of Methanosarcina for ammonia and VFA relative to that of Methanosaeta (Hao et al., 2016).
The analysis also showed a clear positive relationship between OLR and hydrogenotrophic methanogens; with an increase in OLR, the relative abundance of hydrogenotrophic methanogens increased, whereas that of Methanosarcina decreased in R2 and R3, and the concentrations of ammonia and VFA decreased (Table 2, Figure 5).

CO N FLI C T O F I NTE R E S T
The authors declare no conflict of interest.

DATA ACCE SS I B I LIT Y
All data generated or analyzed during this study are included in this published article.