Berberine and its structural analogs have differing effects on functional profiles of individual gut microbiomes.

The understanding of the effects of compounds on the gut microbiome is limited. In particular, it is unclear whether structurally similar compounds would have similar or distinct effects on the gut microbiome. Here, we selected berberine (BBR), an isoquinoline quaternary alkaloid, and 16 structural analogs and evaluated their effects on seven individual gut microbiomes cultured in vitro. The responses of the individual microbiomes were evaluated by metaproteomic profiles and by assessing butyrate production. We show that both interindividual differences and compound treatments significantly contributed to the variance of metaproteomic profiles. BBR and eight analogs led to changes in proteins involved in microbial defense and stress responses and enrichment of proteins from Verrucomicrobia, Proteobacteria, and Bacteroidetes phyla. It also led to a decrease in proteins from the Firmicutes phylum and its Clostridiales order which correlated to decrease proteins involved in the butyrate production pathway and butyrate concentration. Three of the compounds, sanguinarine, chelerythrine, and ethoxysanguinarine, activated bacterial protective mechanisms, enriched Proteobacteria, increased opacity proteins, and markedly reduced butyrate production. Dihydroberberine had a similar function to BBR in enriching the Akkermansia genus. In addition, it showed less overall adverse impacts on the functionality of the gut microbiome, including a better maintenance of the butyrate level. Our study shows that ex vivo microbiome assay can assess differential regulating effects of compounds with subtle differences and reveals that compound analogs can have distinct effects on the microbiome.

Introduction

Therapeutic drugs can interact with the gut microbiome, leading to changes in drug efficacy and changes in the microbiome which in turn can affect the host.[1] Although there is a growing interest in studying drug–microbiome interactions, our understanding of these complex interactions remains limited. It has been proposed that structurally similar compounds would interact with the same enzymes in microbiomes.[2] Maier et al. have shown that drugs that are structurally similar had more similar anticommensal activity compared with drugs that were structurally different.[3] Similarly, Dutta et al. found that L-captopril and its derivatives are all potential inhibitors of microbial enzyme DapE.[4] However, chemically similar compounds can also have markedly different biological actions and activities.[5] For example, Wiggers et al. demonstrated that while sulfasalazine inhibited bacterial diguanylate cyclase inhibitor, its two structurally related molecules sulfadiazine and sulfathiazole did not.[6] Notably, current structure–activity studies are based on single strains of bacteria. The human gut microbiome is composed of different bacteria, and the composition of the gut microbiome is different between individuals. Therefore, it remains unclear whether structurally similar compounds will affect the gut microbiome in a similar way.

Berberine (BBR) is an isoquinoline quaternary alkaloid used extensively in Asia as a nonprescription drug to treat diarrhea, dysentery, stomatitis, and hepatitis.[7,8] Numerous studies[9-18] have reported the mechanism of actions of BBR on the host glucose and lipid metabolism, cardiovascular functions, gastrointestinal tract, inflammation, etc. Recent studies have shown that BBR can affect the human gut microbiota, including increase of Akkermansia spp.,[19] fewer members of Firmicutes (Lactobacillus spp. and Clostridium spp.),[20,21] and Bacteoidetes-to-Firmicutes ratio and reduction of the gut microbiota diversity.[22]

Several studies have suggested that BBR’s structural analogs could have similar or improved functions.[23] For example, studies have shown that BBR analogs also increase the activity of the low-density-lipoprotein receptor.[24] Two BBR analogs were found to be good acetylcholinesterase inhibitors and more potent than BBR as radical scavengers.[25] Three synthesized BBR derivatives were found to induce a stronger effect of cell cycle arrest and cell death through apoptosis.[26] Moreover, a study reported the hypoglycemic activity of modified BBR.[27] BBR is known to have mild antibacterial activity,[28] and similar antibacterial activity of its analogs was reported.[29] However, how the BBR analogs modulate the human gut microbiome is unexplored.

In this study, we assessed whether structurally similar BBR analogs have similar or different effects on the human gut microbiota. We have previously optimized an in vitro culture model which maintains the composition and functions of individual gut microbiome and recapitulate in vivo microbiome responses to compounds.[30] Microbiome functionality is difficult to assess by 16S rDNA sequencing or metagenomic technologies as they only predict the potential function of a gut microbiome.[31] Instead, here, we used metaproteomics to accurately quantify proteins that are actually expressed,[32] a more accurate representation of the function of the microbiome. Briefly, we cultured ex vivo individual gut microbiomes in the presence of BBR and 16 analogs and used metaproteomics to analyze alteration of protein expressions in the gut microbiomes in response to treatment. We evaluated major functional alterations and identified major bacterial contributors to the functional shifts. Moreover, we observed changes in enzymes involved in butyrate production pathways, which was further validated through direct measurement of butyrate in the cultures.

Results

Metaproteomic response of individual gut microbiome to BBR analogs

We tested the effect of BBR and 16 BBR analogs (Figure 1a,b) on seven individuals’ ex vivo gut microbiome. The 16 compounds were structurally similar to BBR. We used a previously optimized workflow for culturing gut microbiome[30] combined with drug stimulation (Figure 1c). We have previously validated that our in vitro model maintains microbiome composition and functional expressions over 24 hours of culture.[30] Briefly, fresh human stool samples were collected in pre-reduced phosphate-buffered saline (PBS) and were homogenized and gauze-filtered immediately before inoculation. The inoculants were cultured in each well of the 96-well DeepWell plates containing 1 ml of optimized culture medium and one of the following: (1) BBR or one of the BBR analogs at 250 µM of concentration, pre-dissolved in 5 µl dimethyl sulfoxide (DMSO) before being added to 1 ml medium; or (2) 5 µl of DMSO (as the negative control). The DMSO control samples were cultured in technical triplicates for each individual gut microbiota. Following 24 hours of culture, the bacterial cells were pelleted, and proteins were extracted and digested for metaproteomic analysis, as previously described.[33] Samples were analyzed by HPLC-ESI-MS/MS using an Agilent 1100 Capillary LC system (Agilent Technologies, San Jose, CA) coupled to an LTQ-Orbitrap XL mass spectrometer (Thermo Electron, Waltham, MA). In total, 180 MS raw files were obtained with a total of 4,127,910 MS/MS spectra. We quantified 70,319 peptides corresponding to 19,123 protein groups with a false discovery rate (FDR) threshold of 1%. The ratio of MS/MS identification was 24.6 ± 8.0%.

Figure 1.

Screening berberine and its analogs against the gut microbiome. (a) Structures, chemical names, and abbreviations of berberine and its analogs involved in this study; (b) Analysis of compounds by structural and property similarity. Multidimensional scaling (MDS) was performed using ChemMine, http://chemmine.ucr.edu/; (c) In vitro culturing and metaproteomics-based approach to study microbiome response to berberine analogs.

Principal component analysis (PCA) using all quantified protein groups (Figure 2a) revealed different effects of the BBR analogs on the gut microbiota. Although the effects of some BBR analogs (dots) on individual microbiota by PCA clustered close to the DMSO control (in squares), other BBR analogs were separated from the control. Permutational multivariate analysis of variance (PerMANOVA) suggested that both individual features and compound treatments had significant effects on the variations of protein abundances (p = .001; Supplementary information S1), and the variance explained by individual feature (R2 = 0.37) is greater than that of the drug effects (R2 = 0.28). Considering the possibility that the individual differences in the microbiome may overshadow drug responses, we applied an empirical Bayesian-based approach (ComBat) to reduce the individual variance (Figure 2a,b). After data processing, the PerMANOVA result showed that the R2 of individual feature decreased to 0.10 and R2 of the drug effects increased to 0.47 (Supplementary information S1). Non-parametric ANOVA (Kruskal–Wallis tests) across all drug treatments and the control showed that the processed dataset doubled the observation of statistically significant differences (Supplementary information S1, Supplementary Table S2, and Supplementary Figure S1). Using the processed dataset, PCA based on individual compounds revealed that the responses to BBR derivatives were microbiome dependent. Nine compounds generally affected the microbiomes compared to the DMSO controls (Figure 2c and Supplementary Figure S2). Despite the existence of nonresponders in four of the compounds (COBA, PMTB, SANGR, and CLTR), the nine compounds showed valid partial least-squares discriminant analysis (PLS-DA) models (Figure 2d), while the other 11 compounds did not, suggesting strong impacts of the compounds on individual gut microbiomes. In agreement with PCA and PLS-DA results, Bray–Curtis distance between DMSO control and drug-treated microbiome shows the extent of the effect of drugs (Figure 2e). Significantly higher impacts on the microbiomes were observed for EOSANGR, SANGR, CTS, and BBR. The remaining five compounds (lighter blue) also suggested strong impacts on subsets of microbiomes.

Figure 2.

Berberine and its analogs showed marked effects on individual gut microbiomes’ metaproteomic profile. (a) PCA plots of the dataset before and (b) after ComBat transformation. Different colors indicate different individual microbiomes (V20, V22, V24, … are numbers of volunteers). (c) PCA plots of individual compounds based on ComBat-corrected data. Nine of the compounds with better separation are shown; PCA of the other compounds is shown in Supplementary Figure S1. (d) PLS-DA cross-validation results based on individual compounds. (e) Bray–Curtis distance between DMSO control and drug-treated individual gut microbiomes (n = 7). Different letters indicate statistically significant differences at p < .05 level by Tukey’s b test. Box spans interquartile range (25th to 75th percentile), and line within box denotes median.

Functional annotation of the protein groups using Clusters of Orthologous Groups (COGs) revealed clustering of eight compounds with DMSO, while the remaining nine compounds showed changes in abundances for specific COG categories (Figure 3a,b). Notably, functions such as defense mechanisms, cell wall/membrane/envelope biogenesis, signal transduction mechanisms, replication, recombination and repair, and transcription were significantly increased by compounds SANGR and EOSANGR. As well, BBR and CTS also induced a significant increase of defense mechanisms and replication, recombination, and repair functions.

Figure 3.

Microbiome functional alterations in response to berberine analogs. (a) Heatmap of COG categories. Sixteen significantly differently abundant COG categories are shown (nonparametric ANOVA, heat colors are based on averages of all tested microbiomes, n = 7). (b) Significantly increased functions found in subgroups of compounds (Mann–Whitney test).

Identifying microbial contributors to functional alterations

The PLS-DA models (Figure 2c) revealed 490 protein groups, with variable importance in projection (VIP) scores greater than one, increasing in abundance in response to these nine compounds. 332 of these protein groups were assigned at the phylum level, among which, 270 (81%) were matched to Bacteroidetes (105 protein groups), Proteobacteria (102 protein groups), and Verrucomicrobia (63 protein groups). Heatmap based on VIP scores of protein groups corresponding to these three phyla (Figure 4a) revealed a division of the nine compounds into three clusters: cluster I consists of compounds 13MBBR, DHBBR, COBA, and PMTB; cluster II includes SANGR, CLTR, and EOSANGR; and cluster III includes BBR and CTS. Moreover, the row clustering indicated that protein groups assigned to the same phylum tend to cluster together, indicating a similar pattern of functional response of the phyla to the stimuli. Cluster I compounds showed a major number of increased proteins from the Bacteroidetes phylum, cluster II compounds showed an increase of proteins from the Proteobacteria phylum, while cluster III had a major increase in proteins from the Verrucomicrobia phylum. Taxonomic enrichment analysis (Figure 4b) on these differential proteins showed that the Bacteroides genus was increased by compounds 13MBBR, DHBBR, COBA, and PMTB, whereas the Enterobacteriaceae family and Escherichia coli were increased by SANGR, CLTR, and EOSANGR. While there is a high number of significantly enriched (p < .01) proteins in the Akkermansia genus by CTS and BBR, they both also exerted increases of Bacteroides and Enterobacteriaceae. Interestingly, although DHBBR belonged to cluster I due to enrichment of Bacteroides and less abundance of Enterobacteriaceae, it also showed increased proteins from Verrucomicrobia.

Figure 4.

Taxonomic contributors to functional alterations. (a) Heatmap based on VIP scores of protein groups corresponding to these three phyla in each PLS-DA model. (b) Taxonomic enrichment analysis of increased protein groups identified by PLS-DA models. Nonsignificant results (p > .01) were marked with a “×.” (c) Taxon-specific functional enrichment analysis of the increased protein groups in response to 13MBBR, SANGR, and BBR. See also Supplementary Figure 2 for COBA, PMTB, EOSANGR, CLTR, DHBBR, and CTS. (d) Taxonomic enrichment analysis of decreased protein groups identified by PLS-DA models.

We analyzed the taxon-specific function enrichment using iMetaLab.[34] The enriched taxa and functions for the proteins with PLS-DA VIP >1 are increased following drug stimulation. Different patterns of taxon-specific functional responses were found (Figure 4c and Supplementary Figure S3). 13MBBR, COBA, and PMTB (“Bacteroides pattern”) had high number of enriched functions correlated to Bacteroides genus; SANGR, CLTR, and EOSANGR (“Enterobacteriaceae pattern”) had more than half of the enriched functions correlated to Enterobacteriaceae; while DHBBR, CTS, and BBR (“Akkermansia pattern”) had relatively high numbers of enriched functions of Akkermansia genus. We identified commonly increased proteins in each of the three patterns and found that in “Bacteroides pattern” and “Akkermansia pattern,” the differential proteins are enriched in basic metabolism pathways such as oxidative phosphorylation, glycolysis/gluconeogenesis, galactose metabolism, and fructose and mannose metabolism. Notably, in the “Enterobacteriaceae pattern,” COG3637, i.e., opacity protein and related surface antigens, is the most highly enriched protein among the three compounds. In agreement with the overall functional analysis (Figure 2), proteins involved in posttranslational modification, protein turnover, chaperones, signal transduction mechanisms, and cell wall/membrane/envelope biogenesis were found to be enriched in the “Enterobacteriaceae pattern.” In addition, we also examined decreased protein groups with VIP scores greater than 1. Taxonomic enrichment analysis (Figure 4d) shows that these protein groups are enriched in Firmicutes phylum and Clostridiales order.

Altered butyrate synthesis pathways in response to BBR analogs

Using metaproteomics, we observed enzymes from three major butyrate synthesis pathways[35] (Figure 5a,b), i.e., 4-amonobutyrate/succinate pathway, acetyl-CoA pathway, and lysine pathway. Correlation based on all cultured individual microbiomes in our dataset (Figure 5a) showed that enzymes of the acetyl-CoA and lysine pathway are clustered separately, with enzymes of each pathway significantly correlated (p < .05). Abfh/Isom of the 4-amonobutyrate/succinate pathway is significantly correlated to the acetyl-CoA pathway. Heatmap of these enzymes (Figure 5c) suggests that different butyrate production pathways responded differently to the compounds: in cluster III containing the control (DMSO) as well as BBR analogs that did not show overall functional impact (in agreement with Figure 2a), enzymes in all the butyrate synthesis pathways were increased; cluster II compounds showed a generally mild decrease in these enzymes, among which, BBR, 13MBBR, DHBBR, and COBA selectively inhibited the acetyl-CoA pathway. In contrast, cluster I showed a strong decrease in these enzymes. Taxon–function correlation analysis showed that these enzymes were mainly from Firmicutes phylum, especially the Clostridiales order. This was in agreement with Figure 4d showing an overall decrease of Firmicutes and Clostridiales in response to the nine compounds. We further analyzed the level of butyrate, iso-butyrate, and 2-methyl-butyrate in an individual microbiome cultured with these compounds. Results suggested that cluster I compounds labeled in Figure 4c led to a decrease in butyrate concentration. In addition, compound CLTR showed a mild decrease in both acetyl-CoA and lysine pathways also led to a decrease in butyrate level.

Figure 5.

Effect of BBR analogs on bacterial butyrate synthesis pathways and butyrate concentration in the culture. (a) Correlation of enzymes involved in butyrate synthesis pathways. All samples are used to calculate the Pearson’s correlation coefficient r; nonsignificant results (p > .05) were marked with a “×.” (b) Enzymes involved in three butyrate synthesis pathways observed in our dataset. (c) Heatmap showing alteration of enzymes involved in butyrate synthesis by the tested BBR analogs, and mean values from all individual microbiomes are shown (n = 7). (d) Top 30 links between taxon and function among butyrate synthesis-related protein groups in our dataset. (e) Concentration of butyrate, iso-butyrate, and 2-methyl butyrate in an individual’s gut microbiome cultured in the presence of the BBR analogs (mean ± SD, n = 3).

Correlating microbiome functional responses to compound features

Next, we attempted to correlate compound features with their biological activities through analysis of their structural or property similarity. From the multidimensional scaling (MDS) plot based on compound features (Figure 1b), the compound that had an effect on the gut microbiota, i.e., BBR, EOSANGR, CTS, DHBBR, 13MBBR, COBA, SANGR, CLTR, and PMTB, tends to appear on the left side of axis V1. Compounds CTS, SANGR, and EOSANGR are closely clustered, suggesting that their structural similarities could result in similar functional responses. Interestingly, the physicochemical properties of DMBBR and JATZ are closely clustered with COBA and PMTB (Supplementary Figure S4 and Figure 1b). However, DMBBR and JATZ did not show overall effects on the microbiomes. In order to explore the specific compound features that are correlated to the functions, we extracted 293 molecular descriptors from each compound using the R package “Rcpi.” The top 10 features that were different between functionally clustered compounds (Figure 3a) are shown in a heatmap (Figure 6). CTS, SANGR, and EOSANGR were clearly separated from other compounds. Interestingly, these compounds have a true value of the feature “Lipinski Failures.” This feature suggests that they satisfy Lipinski’s rule of five,[36] which predicts that poor absorption or permeation is likely when there are more than 5 H-bond donors, 10 H-bond acceptors, the molecular weight is greater than 500, and the calculated Log P (CLogP) is greater than five.[36] In addition, a higher XLogP indicates lower water solubility and higher lipid solubility of these compounds. BBR, CTS, 13MBBR, DHBBR, COBA, and PMTB had an overall difference compared in higher C2SP2, XLogP, nAromBond, and naAromAtom features. The remaining compounds (gray) were not separated based on their features.

Figure 6.

Top 10 different features between clustered compounds.

Discussions

In this study, we examined the effects of 16 BBR analogs on seven ex vivo gut microbiomes using in vitro culture, metaproteomics, and butyrate analysis. Our results showed that 9 out of the 16 compounds, including BBR, showed marked effects on the seven individual gut microbiomes. These compounds resulted in valid PLS-DA models in comparison with the DMSO control, suggesting globally modulating effects of these compounds on all tested individual’s gut microbiome.

Interestingly, BBR, CTS, DJBBR, CLTR, 13MBBR, and COBA showed enriched functions corresponding to the Akkermansia genus (Figure 4b). Studies have shown that Akkermansia is decreased in several conditions, including obesity, diabetes, intestinal inflammation, liver diseases, or chronic alcohol consumption.[37-39] Higher level of Akkermansia was found in modernized populations compared to ancestral populations.[40] It is associated with restored gut barrier function, increased mucus layer thickness, and improved  metabolic disorders.[41] Moreover, the capability of BBR to increase Akkermansia has been discussed to be a potential contributor to the antiatherosclerotic and metabolic protective effects of BBR.[19] Nevertheless, the roles of Akkermansia are still not fully clear. In our study, we found that analogs of BBR could also increase Akkermansia. The compound CTS, which is closely clustered with BBR in all functional analyses (Figures 2a and 3a) and in molecular features (Figure 1b), showed the closest similarity with BBR in upregulating Akkermansia functions. Notably, BBR also enriched functions in the Enterobacteriaceae family. Many species from this family are members of the normal intestinal flora; however, this family also includes overt and opportunistic pathogens responsible for a wide range of infections. The COG3637, opacity (Opa) protein and related surface antigens, was highly enriched by compounds SANGR, CLTR, and EOSANGR. This protein is associated with the pathogenesis of members in the Proteobacteria phylum.[42] Opa-protein interactions with host receptor can lead to bacterial attachment and invasion.[43] Therefore, the increase in Opa-protein highlights a potential risk for these three compounds. In addition, these three compounds showed a significant increase in bacterial defense and envelope-related pathways, suggesting stress responses of the microbiome.[44] Molecular structural analysis showed that these distinct functional impacts from SANGR, CLTR, and EOSANGR are correlated to their poor absorbability, low water solubility, and high lipid solubility. In addition, we found that proteins that are downregulated by BBR and its analogs are mainly enriched in the Firmicutes phylum and Clostridiales order. Interestingly, we found that proteins of the butyrate production pathways are also enriched in these taxa. Butyrate is produced by bacterial fermentation of polysaccharides[45] in the colon, and the consumption of polysaccharides, such as resistant starch, can increase butyrate production in the colon.[46-48] Once produced by bacteria, butyrate is released and provides energy for the host colonocytes. Butyrate is a an important mediator of human health[49] owing to its anti-inflammatory effects.[50] Decreased butyrate production has been reported in disease.[51] BBR and several analogs significantly inhibited pathways leading to butyrate productions and decreased the butyrate concentration. As BBR preserved the lysine pathway, it showed a slight decrease in butyrate concentration; however, several analogs showed a much stronger impact on butyrate production.

Dihydroberberine (DHBBR) also showed an increase in Akkermansia. Moreover, functional analysis showed that in contrast with BBR, DHBBR did not lead to significant shifts of functional categories related to defense mechanisms and stress responses. DHBBR also showed less increase in Enterobacteriaceae and weaker impact on butyrate concentration. DHBBR is a more biologically available derivative of BBR.[52] In vivo studies have shown that compared with BBR, DHBBR improved efficacy of counteracting increased adiposity, tissue triglyceride accumulation, and insulin resistance in high-fat diet-fed rodents;[52] DHBBR has higher activity in inhibiting pancreatic lipase, while BBR may have an adverse influence on ligand–pancreatic lipase affinity. Our results also suggest that DHBBR might be a better substitute to BBR in regulating the gut microbiome. Further animal study would be necessary to comprehensively evaluate its effects in vivo.

Conclusion

We examined the effect of BBR and 16 analogs on the ex vivo human gut microbiome. Eight of BBR analogs showed a global shift in the metaproteomic profiles of the tested microbiomes, leading to alterations of functional compositions mainly associated with microbial defense and stress responses. According to the major bacterial responders, we divided the compounds into three subgroups, enriching proteins from Verrucomicrobia, Proteobacteria, and Bacteroidetes, respectively. BBR and these eight analogs inhibited the proteins of Firmicutes phylum and its Clostridiales order, which was also correlated to different extents of decrease in butyrate pathway and concentration. Compounds SANGR, CLTR, and EOSANGR showed generally negative activities by activating bacterial protective mechanisms, increasing opacity proteins, and significantly decreasing butyrate production. Molecular feature comparison suggested that these changes may correlate to poor absorbability of these three compounds. Both BBR and DHBBR showed enrichment in the Akkermansia genus, and DHBBR showed less overall adverse impacts on the functionality of the gut microbiome. Our study provided a step toward the discovery of new BBR substitutions targeting the gut microbiome. Moreover, our study demonstrated that structurally similar compounds can have different effects on the microbiome.

Materials and methods

Compound information

Berberine chloride was obtained from Sigma-Aldrich (catalog no. PHR1502). The 16 BBR analogs were obtained from Chengdu DeSiTe Biological Technology Co., Ltd., and Chengdu Biopurify Phytochemicals Ltd., and all compounds were with purity ≥98% (Supplementary Table S1). Each compound was dissolved using DMSO to obtain a 50 mM stock solution. 5 µl of the stock solution will be added to each 1 ml of culturing medium, resulting in a 250 µM incubation concentration for each compound.

Stool sample collection

The study was performed in compliance with the Research Ethics Board protocol (# 20160585-01 H) for stool sample collection, approved by the Ottawa Health Science Network Research Ethics Board at the Ottawa Hospital. Stool samples were collected from seven health volunteers (age range 22–39 y; males and females). Exclusion criteria were irritable bowel syndrome (IBS), inflammatory bowel disease (IBD), or diabetes diagnosis; antibiotic use or gastroenteritis episode in the last 3 months; use of pro-/pre-biotic, laxative, or antidiarrheal drugs in the last month; or pregnancy. A fresh stool sample was collected in sterile PBS pre-reduced with 0.1% (w/v) L-cysteine hydrochloride. The sample was immediately weighed and transferred into an anaerobic workstation (5% H2, 5% CO2, and 90% N2 at 37°C). Then, samples were homogenized with a vortex mixer and filtered using sterile gauzes to obtain the microbiome inoculum.

In vitro culturing and drug treatments

We used a previously optimized in vitro model named MiPro to examine the responses of the individual microbiome to BBR and its analogs.[30] Briefly, 96-well DeepWell plates were used to culture the individual gut microbiomes. The composition of the optimized culture medium was as reported previously.[30] pH of the culture medium was adjusted to 7 using 1 M HCl. 5 µl stock solution of each compound was added into one well containing 1 ml of sterile and pre-reduced culture medium, and 5 µl of DMSO was used as the negative control. The microbiome inoculums were inoculated into the wells at a concentration of 2% (w/v). Then, each 96-well DeepWell plate was shaken at 500 rpm in an anaerobic incubator for 24 hours before the analyses.

Metaproteomic sample analysis

Metaproteomic sample processing was based on a 96-well plate-based workflow. Briefly, the bacterial cells were pelleted and washed three times with cold PBS in the 96-well DeepWell plates. Then, samples were stored overnight at −80°C before bacterial cell lysis in 150 µl 100 mM Tris-HCl buffer (pH = 8.0) containing 8 M urea, Roche PhosSTOP™, and Roche cOmplete™ Mini tablets. Sonication was performed using a sonicator (Q125 Qsonica, USA) with an 8-tip-horn probe. The lysed samples were reduced and alkylated with 10 mM dithiothreitol (DTT) and 20 mM iodoacetamide (IAA) and were digested using trypsin (Worthington Biochemical Corp., Lakewood, NJ) at 37°C for 18 hours. Digested samples were desalted using house-made 96-channel filter tips packed with 10-μm C18 resin (Dr. Maisch GmbH, Ammerbuch, Germany).

Tryptic peptides were dissolved in 0.1% formic acid, and 8 μg of protein was loaded for LC-MS/MS analysis with an Agilent 1100 Capillary LC system (Agilent Technologies, San Jose, CA) and an LTQ-Orbitrap XL mass spectrometer (Thermo Electron, Waltham, MA). Peptides were separated with a tip column (75 μm i.d. × 15 cm) packed with 1.9 μm/120 Å ReproSil-Pur C18 resin (Dr. Maisch GmbH, Ammerbuch, Germany) with a 240-min gradient of 5–25% acetonitrile (v/v) at a flow rate of 300 nL/min, with 0.1% formic acid (FA) in water as solvent A and 0.1% FA in acetonitrile as solvent B. Other analysis settings were as previously described.[53] All samples were run on LC-MS/MS in a randomized order.

Butyrate analysis

A series of butanoic acid standard solutions was prepared at a concentration range of 1–3000 µM. Samples and calibration standards were mixed with 2-nitrophenylhydrazine (20 mM, in ethanol, 30 µL) and 1-ethyl-3-(3-dimethylaminopropyl)-carbodiimide (250 mM solution in ethanol, with an equal volume pyridine (3% (v/v), in ethanol), 30 µL). The mixture was heated at 60°C for 45 min, and sodium hydroxide buffer (5% (w/v) sodium hydroxide aqueous solution:methanol (80:20, v/v)) was added. After heating at 60°C for 25 min, the reaction mixture was mixed with hydrochloric acid (1 N, in water, 300 µL), and the butanoic acid derivative was extracted with 400 µL of methyl tert-butyl ether. The methyl tert-butyl ether layer was evaporated to dryness under a stream of nitrogen at room temperature. The residue was dissolved in methanol (100 µL) and an aliquot of 1 µL for analysis. Samples were analyzed by triple quadrupole 6500 mass spectrometer (Sciex, USA). The ACQUITY UPLC BEH C8 (1.7 μm, 2.1 mm × 100 mm, Waters, USA) was used for the analysis. Gradient elution was used with a mobile phase composed of solvent A (water containing 5 mM ammonium acetate) and solvent B (acetonitrile:isopropanol (1:1, v/v)).

Data processing and analysis

For metaproteomic data processing, protein/peptide identification and quantification, taxonomic assignment, and functional annotations were performed using the MetaLab software (version 1.1.1). Spectra clustering strategy was applied to generate a sample-specific database from all raw files. The human gut microbiome gene catalog database comprising 9,878,647 sequences was obtained from the http://meta.genomics.cn/.[54] The identified protein lists were generated with a target-decoy strategy at an FDR cutoff of 0.01, and quantitative information of proteins was obtained with the maxLFQ algorithm[55] on MaxQuant (version 1.5.3.30). Carbamidomethyl (C) was set as a fixed modification and oxidation (M) and N-terminal acetylation (Protein N-term) were set as variable modifications. Razor and unique peptides were included for protein quantification with the minimum ratio count of 1. Then, LFQ protein group intensities were processed by a ComBat process[56,57] using iMetalab.ca[58] to remove possible batch effects between individual microbiomes. Using the ComBat-corrected data, PCA was performed using R function prcomp(), and PLS-DA was performed on MetaboAnalyst.ca.[59] PerMANOVA tests were performed using R packages “vegan” and “BiodiversityR.” PLS-DA model was evaluated by cross-validation of R2 and Q2. Bray–Curtis distance was calculated based on original LFQ intensities using the R package “vegan.” For statistical analysis, we applied nonparametric statistical hypothesis tests, Wilcoxon rank-sum tests were performed in R, and nonparametric ANOVA was carried out on MetaboAnalyst.ca. Taxon–function enrichment analysis was carried out on iMetaLab.ca. Correlation analysis was performed using R function cor() and visualized with package “ggcorrplot.”

Data accessibility

Metaproteomic dataset in this study has been deposited to the ProteomeXchange Consortium via the PRIDE (https://www.ebi.ac.uk/pride/archive/) partner repository with the dataset identifier PXD015934.

Click here for additional data file.