Introduction

Volatile fatty acids (VFAs) and alcohols, e.g., acetate, propionate, butyrate, and ethanol, have wide industrial applications but their current production mainly relies on fossil sources1. Acidogenic processes, e.g., homoacetogenesis, solventogenesis, fermentation, and acidification, can be alternative environmentally friendly approaches for producing VFAs and alcohols. In these processes, microorganisms catalyze CO2 or other organic wastes into VFAs and alcohols with complex interactions and metabolic pathways, simultaneously reducing CO2 emissions and achieving contaminant removal.

Understanding the distribution pattern of microbial communities and identification of functional microbes are key to improving the system performance of acidogenic processes. Currently, many studies have provided novel insights into the diversity, distribution, and functions of microbial populations in environmental biotechnologies through bioinformatic analysis. For example, Wu et al.2 examined the global diversity and biogeography of microbial populations in activated sludge collected from wastewater treatment plants and found that the bacterial composition was driven by stochastic processes. Among the high diversity of the microbial communities, 28 operational taxonomic units were identified as a global core microbial community which could be strongly linked to system performance2. Microbial community composition and diversity in anaerobic digestion systems were also investigated and could be linked to renewable biogas production3,4.

Despite the rapid advancement of biorefinery technologies, the microbial diversity and ecology of acidogenic systems are less understood. Some studies focused on the effects of one or two specific environmental variables on the performance of acidogenic systems. For example, Zhang et al.5 investigated the effects of fermentation substrate and zero-valent iron on the microbial community and metabolic function in a fermentation system. Park et al.6 investigated the metabolic flux and potential function of the microbial community in an acidogenic dynamic membrane bioreactor. Nevertheless, microbial profiles of acidogenic systems such as microbial diversity and the identification of the core microbial community are still far from being comprehensively characterized.

The microbial composition in acidogenic systems may be affected by many environmental factors such as pH, temperature, inoculum, and substrate type5,7,8. These environmental factors provide different microbial niches affecting the microbial community dynamics and metabolic pathways. For instance, solventogenesis usually occurs at a pH lower than 5.59, while propionate-producing bacteria prefer to grow at a neutral pH7. Besides, microbial electrosynthesis has been recently proposed as an efficient approach to control and stabilize fermentation processes10. Microbial electrosynthesis has been extensively used to produce different types of soluble molecules (e.g., alcohols or carboxylic acids) from the reduction of CO2 or the fermentation of organic matter, and to overcome some limitations of metabolic reactions11,12. Such bioelectrochemical systems exert significant advantages on microbial metabolism, interspecies interactions, and the selection of microbial populations10. Nonetheless, the detailed effects of electrode potential remain largely unrevealed. Accordingly, comprehensive profiling of microbial communities in acidogenic systems with different environmental variables may provide valuable information for predicting microbial ecosystems responding to environmental changes and benefit the optimization of biorefinery technologies.

In this study, 18 different inocula and 42 samples collected from lab-scale acidification reactors with different environmental variables were analyzed to explore the microbial community and diversity of acidogenic systems. The core microbial community and their metabolic pathways in acidogenic samples were identified. The acidogenic microbial network and subnetworks were established to highlight the distribution and topological structures of core features from the network aspect. Furthermore, the effects of key environmental variables, e.g., pH, substrate type, temperature, inoculum, and bioelectrochemical potential on the core microbial community and metabolic pathways were analyzed. Statistical analysis was conducted to investigate the potential correlation between the acetogenic performance and core microbial community or predicted functional pathways. The distribution pattern and metabolic pathways of core features revealed in this study can help guide the regulation of core microbial community to optimize acidogenic practices.

Results and discussion

Core microbial community

A total of 3,107,841 sequences were obtained from all 60 samples, with an average of 51,797 sequences per sample. The core microbial community was determined by the abundance and occurrence frequency of ASVs2. As shown in Fig. 1a, 23 out of 30732 ASVs were identified as the core features that accounted for 41% of the sequences in all acidogenic samples investigated. Most of the core features are acidogenic bacteria and are distinct from the core microbial community identified in domestic wastewater treatment plants or anaerobic digestion processes producing methane2,3, suggesting that the acidogenic system selects a unique core microbial community.

Fig. 1: Characterizations of core features identified.
figure 1

a Percentage and relative abundance of the core features and b the taxonomic composition of core features at the genus level and their frequency of abundance (green column), frequency of ubiquity (yellow column), and correlation with acetogenic performances.

At the phylum level, nine core features belonged to Firmicutes, followed by Bacteroidota (6), i.e., 65.2% of the core features belonged to these two phyla. At the genus level, Clostridium sensu stricto 12 was identified as the most abundant core feature, accounting for 13.9% of the sequence abundance. Clostridium is a well-known acetogen widely existing in acidogenic processes13,14,15. Surprisingly, the second most abundant core feature (4.2% of the sequence abundance) was a typical H2 utilizing methanogen Methanobacterium. Methanobacterium was ubiquitously detected in 71.4% and abundant in 59% of all acidogenic samples. One explanation is that many inoculum samples were collected from anaerobic digesters producing methane, and although bromoethane sulfonate (BES) was added into some of these samples, it seemed to fail to completely eliminate Methanobacterium9,16. Besides, it is reported that hydrogenotrophic methanogens have a higher H2 competitiveness than homoacetogens at a high H2 partial pressure17. Desulfovibrio was also identified as a core feature, reflecting the importance of incomplete oxidizing sulfate-reducing bacteria in acidogenic systems that develop out of inocula from full-scale anaerobic digesters.

The prevailing microorganisms are often suggested to be responsible for maintaining system performance3. To evaluate the correlation between the core features and acidogenic system performance, the highest acetate concentration and the maximum acetate-production rate were used to characterize acetogenic performance since acetate was the most common product in all samples. Three core features were positively correlated with both the acetate concentration and acetate-production rate (Fig. 1b), including Clostridium sensu stricto 12 (p ≤ 0.05), Bacteroidetes vadinHA17, and one ASV belonged to the phylum Clostridiaceae (p ≤ 0.05). The other 5 core features were positively correlated with at least one acetogenic performance while the remaining core features were negatively correlated with acetogenic performance. The core features that were negatively associated with acetate production might be critical for maintaining other acidogenic functions such as propionate and butyrate production. For example, butyrate can be produced via bioelectrochemical chain elongation from CO2 and acetate, or via acetyl-CoA reduction with H2 as electron donor18. All these processes consume acetate or precursors of acetate production. Consistent with this, the butyrate producers Oscillibacter19 and Caproiciproducens20 were negatively correlated with acetogenic performance. On the other hand, core features positively correlated with the acetogenic function might be the winners of substrate competition since acetate is the end product of homoacetogenesis, fermentation, and acidification. The correlation result shown in Supplementary Fig. 1 may support this hypothesis because many correlations among core features were negative. Particularly, Clostridium sensu stricto 12 and the ASV belonging to the phylum Clostridiaceae were negatively correlated with most (73–77%) of the other core features.

Microbial network

Microbial correlations were also presented by the network analysis (Fig. 2). As shown in Fig. 2a, a complex microbial network (855 genera involved) was observed with more than 15 thousand correlations (Pearson correlation coefficient >0.7, p ≤ 0.05) existing in the acidogenic community. The topological parameters of the acidogenic network are listed in Supplementary Table 1. The low network density (0.04) and high heterogeneity (0.9) suggested that the majority of genera tended to have few connections21. Besides, the low network centralization (0.1) implied that the centralities of genera are relatively similar22. Furthermore, the whole network consisted of many high-density circle-type subnetworks. There might be multiple relatively independent ecosystems within the acidogenic system investigated.

Fig. 2: Acidogenic network and its subnetworks.
figure 2

a Acidogenic network at the genus level (Pearson coefficient >0.7, p ≤ 0.05) and b Top 12 subnetworks of the acidogenic network. The orange nodes represent core features.

To investigate the in-depth structure of the acidogenic network, subnetworks were further extracted. A total of 28 subnetworks were formed in the acidogenic network. These subnetworks were marked from 1 to 28 according to the number of nodes. Figure 2b shows the structures of the largest 12 subnetworks, many of which were circle-type subnetworks, consisting of one or two circle-like structures. Based on the number of nodes, edges, and the average number of neighbors (Supplementary Fig. 2), the subnetworks can be divided into three scales. Subnetwork 1 had the largest scale, including 15.5% of nodes and 38% of edges, with a high average number of neighbors (81). Subnetworks 2–11 were medium-scale subnetworks, including more than half of the nodes and edges. Finally, subnetworks 12–28 were small-scale subnetworks, with much fewer nodes and edges. Supplementary Table 2 shows the average numbers of topological parameters of the subnetworks. The scales of subnetworks were positively correlated with the average clustering coefficient and were negatively correlated with the network heterogeneity. Another interesting result is that the subnetworks exhibited different characteristics from the whole acidogenic network. The acidogenic network was relatively discrete, while subnetworks were more concentrated. For example, the network densities of the subnetworks were much higher than the whole acidogenic network (Supplementary Table 2), which means that subnetworks were more populated with edges than the whole acidogenic network.

Supplementary Fig. 3 and Supplementary Table 3 show the distribution of core features in different subnetworks. Although subnetwork 1 had the highest number of nodes, only 3 core features were located in it. 15 core features were distributed in medium-sized subnetworks and the other 2 core features were located in small subnetworks. Specifically, the medium-sized subnetworks 5 and 8 had the largest number of core features (4), followed by the large subnetwork 1 (3). Subnetwork 5 mainly contained butyrate (Oscillibacter, the ASV belonged to the phylum Oscillospiraceae, and Caproiciproducens) and propionate (Paludibacter) oxidizers19,20, while subnetwork 8 included bacteria able to produce acetate, butyrate, and propionate (Petrimonas, Bacteroides, the ASV belonging to the family Lachnospiraceae, and NK4A214 group). On the other hand, core features in subnetwork 1 include fermenters and syntrophic bacteria that are able to produce acetate (Mesotoga, the ASV belonging to the order Syntrophales, and the ASV belonging to the family Spirochaetaceae). The core features that co-occurred in the same subnetwork might have similar metabolic preferences, which however should be confirmed by more solid evidence. In addition, Clostridium sensu stricto 12 and the ASV belonging to the phylum Clostridiaceae did not have a strong (Pearson correlation coefficient >0.7) correlation with the other core features and were located in two small subnetworks.

Table 1 summarizes the topological properties of the core features. Core features presented obvious variations in some of these topological properties, reflecting that they might have different interactions with other microbes and might play different roles in the microbial network. This is in line with the fact that core features were distributed within different subnetworks rather than clustered together. For example, the ASV belonging to the family Spirochaetaceae had 103 neighbors, connecting most of the genera in subnetwork 1, therefore was likely to be a hub in the subnetwork. In contrast, Clostridium sensu stricto 12 and Christensenellaceae R-7 group had only 3 and 2 neighbors, respectively. Thus, they did not interact with most of the other nodes in the network23.

Table 1 Topological properties of core features in the microbial network at the genus level.

Previous studies reported that inoculum samples had higher diversity than acclimated samples, but none of the studies had clarified this through the aspect of the microbial network. To assess this, the inoculum network was also established and was compared with the acidogenic network (Supplementary Fig. 4 and Supplementary Table 1). Generally, the inoculum network presented a distinct structure. However, the difference between the acidogenic network and the inoculum network was not in the number of nodes (i.e., number of genera), but in the interacting relationships. Specifically, the inoculum network had 109% more edges and the average number of neighbors was also 132% higher than that of the acidogenic network. Therefore, genera in the inoculum samples had much more complex interactions than the acidogenic samples. The higher clustering coefficient and network density of the inoculum network also supported this result (Supplementary Table 1). These results indicate that acidogenic acclimation processes decreased the microbial interactions. One possible explanation is that compared with the inoculum samples, acidogenic samples have fewer microbial functions which were interlinked with interactions such as syntrophic oxidation of substrates.

Rare community

In contrast with the dominant species, rare species greatly contribute to the richness but account for a low percentage of abundance in the microbial ecosystem4. To understand whether there was a rare community in the acidogenic systems investigated, rare features were defined and identified. Surprisingly, only one ASV belonging to the phylum Firmicutes was found to be the rare feature when the criterion of ubiquitous features was defined as the occurrence frequency in more than 70% of all samples. Rare features in acidogenic systems are more likely to exist stochastically, unlike the existing pattern of the core microbial community.

Decreasing the criterion of the occurrence frequency from 70 to 50% led to the discovery of 5 rare features, including Candidatus Caldatribacterium, Phaselicystis, Thermovirga, one ASV belonging to the family Hungateiclostridiaceae, and one belonging to Firmicutes (Fig. 3). Candidatus Caldatribacterium, Phaselicystis, Thermovirga, and family Hungateiclostridiaceae are fermenting bacteria, able to function in VFAs production24,25,26,27.

Fig. 3: Characterizations of rare features.
figure 3

a The taxonomic composition of rare features at the genus level and their frequency of abundance (green column), frequency of ubiquity (yellow column), and correlation with acetogenic performance and b Pearson correlation among rare features.

Rare species play different roles in ecosystems. Generally, rare species serve as a seed bank with huge ecological potential which can respond to environmental changes and maintain system stability28. Besides, rare species are responsible for functions such as the removal of micropollutants29. Although the relative abundance of rare features was low, the fact that 4 of 5 rare features were able to function in VFAs production implied a possibility that they might to some extent contribute to acidogenesis which however needs further investigation.

Top metabolic pathways

PICRUSt2 was applied to predict metabolic pathways of acidogenic systems based on 16 S rRNA marker genes. The NSTI values of acidogenic samples ranged from 0.03 to 0.20 with a mean value of 0.12. Low NSTI values suggested that these samples can be tracked for PICRUSt2 prediction30. A total of 410 pathways from the MetaCyc database were predicted in all acidogenic samples. Figure 4a shows the top 40 pathways, including 33 biosynthesis pathways, 2 pathways associated with degradation/utilization/assimilation, and 5 pathways associated with the generation of precursor metabolites and energy. Biosynthesis pathways constitute the cells’ spectrum of biosynthetic capacities, including the pathways of synthesis of small molecules, macromolecules, and cell structure components. Most (85%) of the top 33 biosynthesis pathways belonged to amino acid biosynthesis, fatty acid and lipid biosynthesis, and nucleoside and nucleotide biosynthesis. Similar pathways were reported to be dominant in fermentation systems5,31. It should be noted that only 2 fermentation pathways were found in the top 40 pathways, including one pathway associated with “pyruvate fermentation to acetate and lactate II” and another one associated with “pyruvate fermentation to isobutanol”.

Fig. 4: Top metabolic pathways and acid-producing pathways.
figure 4

a Relative abundance of top 40 predicted pathways and b relative abundance of fermentation pathways.

Acid-producing pathways

Apart from the 2 fermentation pathways mentioned above, another 13 pathways associated with acidogenesis were further identified in the acidogenic samples (Fig. 4b), which in total accounted for 3.6% of the total number of all metabolic pathways identified. The acidogenic samples were clustered into two groups according to the abundance of these 15 pathways (Fig. 4b), with one group containing the high abundance of pathways associated with “pyruvate fermentation to isobutanol (engineered)”, “pyruvate fermentation to acetate”, “pyruvate fermentation to acetate and lactate II”, “homolactic fermentation”, and “acetyl-CoA fermentation to butanoate II”, and another group only containing high abundance of pathways associated with “pyruvate fermentation to isobutanol (engineered)”. The acidogenic samples in these two groups, however, could not be clearly distinguished by the category of substrates or products. Pyruvate, which was correlated with 5 of 15 pathways, seemed to be a key intermediate in the acidogenic process. According to the MetaCyc database, all CO2, lactate, and components of brewery spent grain (e.g., glucose) could be converted to pyruvate via multiple reactions. Therefore, one strategy to promote the targeted acid or alcohol yield is to increase the pyruvate pool by knocking out genes that contributed to by-product formation using pyruvate as the precursor. For example, to increase the isobutanol production from pyruvate in engineered Escherichia coli, the genes that regulated the by-product formation from pyruvate were deleted and levels of pyruvate available for isobutanol production were increased32.

As shown in Supplementary Fig. 5, four CO2 fixation pathways were identified in the acidogenic samples investigated. Supplementary Table 4 summarizes the reactants, products, and subsequent pathways of CO2 fixation. The “Calvin Benson Bassham cycle” turned out to be the most abundant CO2 fixation pathway, being responsible for generating 3-phospho-D-glycerate, glycerone phosphate, or D-glyceraldehyde 3-phosphate which could be used in the biosynthetic metabolism of autotrophic bacteria. The “reductive tricarboxylic acid (TCA) cycles I and II” from autotrophic eubacteria and archaea were also found in acidogenic samples. The “reductive TCA cycles I and II” are commonly found in autotrophic eubacteria and archaea. The basic cycle results in the fixation of 2 molecules of CO2 and the generation of 1 molecule of acetyl-CoA33. Another molecule of CO2 can be fixed by the carboxylation of acetyl-CoA to produce pyruvate34. The “reductive acetyl coenzyme A pathway”, known as the Wood-Ljungdahl pathway, is the most effective nonphotosynthetic CO2 fixation pathway by acetogens35,36. Two molecules of CO2 can be converted to acetyl-CoA in the Wood-Ljungdahl pathway. Acetate, ethanol, and butyrate could be produced from CO2 through the detected Wood-Ljungdahl pathway coupled with the acetyl-CoA reduction in acidogenic samples9,16,20.

Roles of core features in acidogenic metabolism

To examine the contribution of core features to metabolic pathways, genes involved in the top 40 pathways as well as the 15 fermentation pathways were further assigned to ASVs. As shown in Supplementary Fig. 6a, 4.6–74.3% of the top 40 pathways belonged to core features, with an average of 26.2% per pathway. The highest average percentage of core feature assignment was found in samples producing ethanol and acetate, accounting for an average of 45.2% per pathway, while the lowest one was found in samples producing butyrate and acetate (20%). Similar results were found in fermentation pathways where core features accounted for 4.9–66.1% of the fermentation-pathway genes (Supplementary Fig. 6b), with an average of 27.0% per pathway. This result implies that the core features participated not only in the top 40 pathways, but also in the 15 fermentation pathways in the acidogenic samples.

Microbial diversity

Alpha diversity indexes associated with evenness (Pielou’s evenness and Shannon’s index) and richness (Observed ASVs and Faith’s phylogenetic diversity) were calculated to investigate the microbial diversity of the acidogenic samples (Fig. 5). When grouped by substrates, samples fed with NaHCO3 showed the highest level of both evenness and richness, while samples fed with H2–CO2 demonstrated the lowest diversity. Particularly, diversity of samples fed with NaHCO3 was significantly higher than that of samples fed with H2–CO2, CO2, and H2–CO2–NaHCO3 (p ≤ 0.05) (Supplementary Table 5). On the other hand, when grouped by products, sludge samples producing butyrate and acetate showed the highest level of diversity, while those with ethanol and acetate as main products presented the lowest diversity (Supplementary Table 6). Compared with acidogenic samples, the diversity indexes of inoculum samples were relatively higher, of which the observed ASVs index was significantly different (p ≤ 0.01). This might be due to the washout or decay of some microorganisms that were not tolerant to the low pH in the acidogenic systems or the immigrated microorganisms from the substrates lead to the high microbial diversity of the inoculum samples3.

Fig. 5: Alpha diversity indexes (Pielou evenness, Shannon, observed features, and Faith’s phylogenetic diversity) of acidogenic samples grouped by substrate, product, inoculum, and bioelectrochemical characteristic.
figure 5

The line and square in the box mark the mid-point and average values of the data, respectively. The upper and lower edges of the box mark that the seventy-five and twenty-five of the data fall below the upper quartile and the lower quartile, respectively. The upper and lower whiskers represent the maximum and minimum values respectively, excluding outliers. The black diamond represents the outlier.

It is generally accepted that a high diversity increases functional resilience under stress conditions and results in better system performance37,38. However, as shown in Fig. 6, all alpha diversity indexes were negatively correlated with the acetogenic performance. This may be because that acetogenesis is less dependent on the diversity. For example, when CO2 and H2 are fed, homoacetogenesis occurs spontaneously without interacting with other microbes, thus increasing the number of other microorganisms involved. Tao et al.3 also reported that only a weak correlation between diversity indices and performance parameters was found in anaerobic digestion samples.

Fig. 6: Pearson correlation between alpha diversity indexes and acetogenic performances.
figure 6

The size of the circle represents the size of the absolute value of the index and the color of the circle represents the value of the index. Red is positive and blue is negative. The star represent the significant value.

To understand community assembly patterns, PCoA was applied to investigate the variations in community composition (beta diversity) based on taxonomic (Jaccard distance and Bray–Curtis distance) and phylogenetic (unweighted and weighted UniFrac distance) diversity (Fig. 7). Considerable variations between acidogenic samples were observed at both the taxonomic and phylogenetic levels. Among them, the beta diversity results based on the weighted Unifrac distance explained the highest variation (38.36%). Microbial communities of acidogenic samples grouped by substrate or products exhibited significant differences (p = 0.001) (Fig. 7), suggesting that these two parameters also affected the microbial community composition.

Fig. 7: PCoA of beta diversity.
figure 7

a Bray–Curtis distance, b Jaccard distance, c unweighted UniFrac distance, and d weighted UniFrac distance grouped by product (color) and substrate (shape).

Effects of environmental variables on the core microbial community

As demonstrated by CCA (Fig. 8a), environmental variables (e.g., electrode potential, substrate, and pH) and acetogenic performances had a clear impact on the core features. Together, environmental variables and acetogenic performance explained 73.5% of the core feature variation (p ≤ 0.05) (Supplementary Table 7). Surprisingly, the electrode potential was the key contributing factor (p ≤ 0.01), accounting for 48.2% of the explanation. Core features responded to these factors variously. Many core features were located near the arrow of the bioelectrochemical factor. This may be because extra electrons provided by the electrodes can motivate energetically unfavorable reactions11, and may thus favor the growth of some acidogenic bacteria. Specifically, the distribution of 16 core features was significantly different (p ≤ 0.005) in samples with or without the electrode potential (Supplementary Table 8). For example, Clostridium sensu stricto 12, Bacteroides, Bacteroidetes vadinHA17, and the ASV belonging to the family Lachnospiraceae were more abundant in samples without the electrode potential (Supplementary Table 8). In contrast, Methanobacterium, Desulfovibrio, Oscillibacter, and Caproiciproducens were dominant in the samples with the electrode potential (Supplementary Table 8). Methanobacterium species could produce H2 instead of CH4 when confronted with an excess supply of reducing equivalents on a cathode in an acetate-producing bioelectrochemical system39. Furthermore, Desulfovibrio and Pseudomonas (acetogen) are commonly recognized as electroactive bacteria, able to accept electrons from cathodes in microbial electrolysis cells40,41. In addition, from the subnetwork aspect (Fig. 2 and Supplementary Table 3), all core features located in subnetwork 8 were favored by the electrode potential. Environmental variables may also play a role in the establishment of subnetworks.

Fig. 8: Effects of environmental variables on core features and predicted pathways.
figure 8

a Canonical correlation analysis (CCA) of core features, environmental variables, and acetogenic performance as well as b redundancy analysis (RDA) of predicted pathways, environmental variables, and acetogenic performance. The full name of each pathway can be found on MetaCyc database (https://metacyc.org/).

On the other hand, Clostridium sensu stricto 12 and the ASV belonging to the phylum Clostridiaceae were located far away from other core features and near the arrow of the acetate-production rate (Fig. 8a), confirming their dissimilarity of distribution and their important roles in maintaining acetogenic performance. Besides, although Methanobacterium was identified as one of the core features in acidogenic samples, the CCA result shows that it was not closely correlated with most of the other core features and had a distinct pH preference.

Effects of environmental variables on the metabolic pathways

The RDA further revealed that the top 40 metabolic pathways, 15 fermentation pathways, and 4 CO2 fixation pathways responded differently to the environmental variables and gave a different acetogenic performance (Fig. 8b) which explained a total of 65.3% of the pathway variation (p ≤ 0.05, contribution of acetate concentration excluded). The bioelectrochemical factor was the highest contributing variable (p ≤ 0.01), accounting for 49.3% of the explanation, followed by pH (25.2%) (Supplementary Table 9). Furthermore, more than half of the fermentation pathways were positively correlated with the bioelectrochemical factor even though the bioelectrochemical factor was negatively correlated with acetogenic performance, implying that the electrode potential might benefit the production of other VFAs rather than acetate, which needs further confirmation. For example, propionate generation from pyruvate (P108-PWY), which can be first produced during lactate fermentation42, was positively correlated with the bioelectrochemical factor. Another intriguing result is that most of the top 40 pathways were clustered together and had a high positive correlation with the acetate-production rate, but were negatively correlated with the bioelectrochemical factor. Since 33 of the top 40 pathways are related with biosynthesis, it is also necessary to confirm in the future studies whether the electrode potential application negatively impacts the growth of microbes.

Effects of environmental variables on the diversity

Alpha diversity indexes of acidogenic samples were distinguished by their pH (p ≤ 0.01). As shown in Supplementary Fig. 7, samples with a low pH tended to harbor low microbial diversity and vice versa. Some samples with pH lower than 4.5 were those producing ethanol via solventogenesis. This is in line with the fact that solventogenesis usually occurs at a pH ranging from 4.5 to 5.59. Furthermore, these samples were fed with H2–CO2 but not NaHCO3, consistent with the result from Dessi et al.16 that changing carbon source from bicarbonate to CO2 reduced the alpha diversity of acidogenic systems. This may be also because bicarbonate was able to act as a buffer, avoiding the pH drop and eventually retaining more microbes that are sensitive to pH. On the other hand, the temperature which is also a key factor affecting the microbial population8 did not show a significant effect on diversity (p > 0.05) (Supplementary Fig. 8). Acidogenic samples cultivated at 25 and 37 °C might share relatively similar diversity.

Another intriguing result is that the electrode potential exerted significant effects on the diversity (p ≤ 0.05). Sludge samples collected from the bioelectrochemical systems showed a higher diversity than those without the electrode potential. It is recognized that electrochemical control affects the extracellular and intracellular metabolism of fermentative microorganisms10. Figure 5 further shows that the electrode potential could benefit the diversity of acidogenic samples, which might be partially due to the enrichment of electroactive bacteria. Moreover, when considering beta diversity, samples classified as with or without electrode potential showed a clear distinction that could be observed (p = 0.001) (Supplementary Fig. 9). This may be because the electrode potential on the cathodes can provide a unique microbial niche by establishing interactions either between microorganisms or between microorganisms and the electrode surface10, both of which can result in a shift of the community composition.

Implications

The above results provide insights into the optimization of acidogenic processes. A core microbial community was identified in the studied acidogenic samples. Since most of these core features are acidogenic bacteria and harbored a large percentage of predicted genes involved in acid-producing pathways, they are likely to be the main acid producers in many acidogenic systems for anaerobic biorefinery of valuable bioproducts. From the network aspect, the decentralized distribution of core features in different subnetworks suggested that many of them did not have co-occurring relationship and tended to play different roles in the acidogenic network. Although the exact function of each subnetwork remains unknown, the existence of many high-density subnetworks might imply a potential subnetwork-based regulation strategy of the acidogenic system.

Although diversity is suggested to play important roles in system performance37,38, in this study, all alpha diversity indexes were negatively correlated with the acetogenic performance. This may be because acidogenic processes are usually operated on relatively low pH conditions which are not suitable for many microbes. Therefore, the importance of the high microbial diversity to acidogenic performance should be reevaluated. The distinct differences between the inoculum network and the acidogenic network also supported the hypothesis that microbial diversity plays a minor role in acidogenic systems. The acclimation did not obviously affect the number of microbes (diversity) that have a strong positive correlation but significantly decreased their interacting relationships. Instead, core features and their response to environmental changes should be focused on when optimizing the acidogenic performance.

The bioelectrochemical factor has the greatest influence on not only core features but also the fermentation pathways. More than half of the core features and fermentation pathways were positively correlated with the electrode potential. Therefore, when aiming to enrich the associated core features or produce the associated products, microbial electrosynthesis may be a promising approach. Since some of the core features that positively affected by the electrode potential were located in the same subnetwork, a hypothesis is that the bioelectrochemical factor is able to regulate the microbial community at the subnetwork level but this needs further verification. Nevertheless, it should also be noted that most of the top 40 metabolic pathways, and some other core features were negatively correlated with the bioelectrochemical factor. Thus, in future studies, the accurate effects of the bioelectrochemical factor on shaping the microbial community should be further examined as well.

This study investigated the core microbial community, microbial network, metabolic pathways, and microbial diversity of acidogenic samples. A total of 23 microbes were identified as the core features, with Clostridium as the most abundant one. The non-centralized distribution of core features in different subnetworks and their topological parameters suggested that they play different roles in the microbial network. On average, 27% of the predicted fermentation-pathway genes belonged to the core microbial community, implying their important roles in acidogenesis. Electrode potential contributed the most to the distribution of both the core microbial community (48.2%) and metabolic pathways (49.3%). More than half of the core features and fermentation pathways were positively correlated with the electrode potential. Posing an electrode potential, as e.g., in microbial electrosynthesis, may thus be a promising tool to shape an anaerobic digestion microbiome towards a VFAs and alcohols one.

Methods

Sample description

A total of 60 biomass samples were collected from four different experiments with acidifying bioreactors conducted at National University of Ireland Galway, Ireland. The reason to choose samples from these experiments was that they were conducted under different environmental variables including different inoculum, temperature, pH, substrate, and electrode potential. The wide variety of environmental variables allowed us to have a more general profile regarding how they affect the microbial community. Detailed sample information is summarized in Table 2. Specifically, twelve samples were collected from group I9. The inoculum of group I was collected from an up-flow-anaerobic sludge bed (UASB) reactor producing methane from dairy wastewater. Acidogenic samples were collected from batch-test reactors operated at 25 °C or 37 °C, fed with H2–CO2, H2, CO2, or NaHCO3 as the substrate. Samples of group I were collected at the end of the batch-test incubation. The main products of these samples were ethanol and/or acetate. Fifteen samples were collected from group II20. The inoculum of group II was collected from an anaerobic digester treating cheese processing wastewater (Dairygold, Cork, Ireland). Acidogenic samples were collected from microbial electrosynthesis reactors producing propionate and acetate from lactic acid at 25°C and an applied cathodic potential of −1.0 V vs. Ag/AgCl. Samples of group II were collected at the end of the electrosynthesis experiment. Another fifteen samples were collected from group III16. The inoculum was digested sludge from a dairy processing industry (Dairygold, Ireland). Acidogenic samples were collected from microbial electrosynthesis reactors fed with CO2 or NaHCO3 at 25 °C and the applied cathodic potential of −1.0 V vs. Ag/AgCl. Samples of group III were collected on days 47, 183, 206, and 225 of reactor operation (total operation of 225 days). The main products of these samples were acetate and/or butyrate. Eighteen samples were collected from group IV43. The inoculum samples included three from rumen fluid, three from rhino manure, and three from tiger manure. Acidogenic samples were collected on day 21 from batch-test reactors operated at 37°C degrading brewery spent grain (total operation of 30 days)43.

Table 2 Details of the 60 samples collected from the four experiments.

DNA was extracted using DNA extraction kits following the manufacturer’s protocol9,16,20,43. 16S rRNA genes were amplified using the universal primer pair 515F/806R. Library preparation and high throughput sequencing were conducted on the Illumina Miseq platform.

Microbial community analysis

Amplicon sequence variants (ASVs) were constructed using the Qiime2 workflow44. Briefly, sequence merging, quality filtering, ASVs table construction, and taxonomic assignment were performed44. Sequence merging was conducted using PEAR with a minimum overlap of 10 bps and a minimum assembled length of 200 bps. Data2 was applied for sequence quality filtering and ASVs table construction45. Taxonomic analysis was performed by comparing with the SILVA database (silva-138-99-nb-classifier). The definition of the core microbial community was based on Wu et al.2 with a slight modification. The criterion of “overall abundant features” was defined as ASVs with a total relative abundance higher than 0.1% of all acidogenic samples. The criterion of “ubiquitous features” was defined as ASVs with an occurrence frequency in more than 70% of all acidogenic samples. The criterion of “frequently abundant features” was defined as ASVs abundant in at least 50% of the acidogenic samples. ASVs fulfilling all the above three criteria were defined as core features in acidogenic systems. The rare community was defined similarly except that the criteria of “overall rare features” was defined as ASVs with the total relative abundance lower than 0.1% and the criteria of “frequently rare features” was defined as ASVs rare in at least 50% of the samples.

Microbial metabolic pathways prediction using PICRUSt2

PICRUSt2 was used to predict the metabolic pathways of the microbial communities based on 16S rRNA marker gene profiles46. The functional annotation of PICRUSt2 predictions was obtained based on the MetaCyc database47. The nearest sequenced taxon index (NSTI) was calculated to evaluate the reference genome coverage for accurate PICRUSt prediction. Estimated MetaCyc pathways related to the functional potential involved in acidogenesis were manually classified, based on their MetaCyc identification. Predicted pathways were assigned to ASVs for the clarification of microbes involved in metabolic pathways.

Statistical analysis

Alpha and beta diversities were calculated by Qiime2 diversity plugin using the extracted ASVs44. The Kruskal–Wallis test was used to examine the significance of the alpha diversity48. The significance of the beta diversity was evaluated by the PERMANOVA test (pairwise comparison and 999 permuations). Principle coordinates analysis (PCoA) was applied to investigate the beta diversity. Heatmaps of metabolic pathways were generated using the STAMP software49. Detrended correspondence analysis (DCA), canonical correlation analysis (CCA), and redundant analysis (RDA) were conducted by the Canoco5 software. To normalize the types of substrates for statistical analysis, carbon substrates were calculated based on the molecular mass (brewery spent grain excluded). Since the length of axis 1 of DCA is lower than 2, RDA was selected to evaluate the correlation between the metabolic pathways and environmental variables50. Pearson correlation analysis was performed using the Origin software. Network analysis was performed in Cytoscape and Origin software based on the Pearson correlation coefficients (Pearson correlation coefficients >0.7, p value ≤ 0.05). Subnetworks were extracted using the CytoCluster plugin51. A p value ≤ 0.05 was considered statistically significant.