Clostridial Butyrate Biosynthesis Enzymes Are Significantly Depleted in the Gut Microbiota of Nonobese Diabetic Mice.

Increasing evidence suggests that the intestinal microbiota is involved in the pathogenesis of type 1 diabetes (T1D). Here we sought to determine which gut microbial taxa and functions vary between nonobese diabetic (NOD) mice and genetically modified NOD mice protected from T1D (Eα16/NOD) at 10 weeks of age in the time window between insulitis development and T1D onset. The gut microbiota of NOD mice were investigated by analyzing stool samples with a metaproteogenomic approach, comprising both 16S rRNA gene sequencing and microbial proteome profiling through high-resolution mass spectrometry. A depletion of Firmicutes (particularly, several members of Lachnospiraceae) in the NOD gut microbiota was observed compared to the level in the Eα16/NOD mice microbiota. Moreover, the analysis of proteins actively produced by the gut microbiota revealed different profiles between NOD and Eα16/NOD mice, with the production of butyrate biosynthesis enzymes being significantly reduced in diabetic mice. Our results support a model for gut microbiota influence on T1D development involving bacterium-produced metabolites as butyrate.IMPORTANCE Alterations of the gut microbiota early in age have been hypothesized to impact T1D autoimmune pathogenesis. In the NOD mouse model, protection from T1D has been found to operate via modulation of the composition of the intestinal microbiota during a critical early window of ontogeny, although little is known about microbiota functions related to T1D development. Here, we show which gut microbial functions are specifically associated with protection from T1D in the time window between insulitis development and T1D onset. In particular, we describe that production of butyrate biosynthesis enzymes is significantly reduced in NOD mice, supporting the hypothesis that modulating the gut microbiota butyrate production may influence T1D development.

OBSERVATION

Type 1 diabetes (T1D) is an autoimmune disease characterized by the specific destruction of pancreatic insulin-producing β-cells (1). Although there is evidence of a strong genetic basis for T1D (2), a remarkable increase in T1D incidence has been measured in the last decades, suggesting a significant contribution from the environment (3). Among possible environmental factors influencing T1D development, many may be significantly related to the gut microbiota, including hygiene, antibiotic use, and diet (4). Gut microbial colonization starts at birth, playing a key role in priming the immune system (5), as well as in complementing the host's metabolism (6). Alterations of the gut microbiota during host development have been hypothesized to have an impact on T1D autoimmune pathogenesis (7), although a full comprehension of the host-microbiota interactions related to T1D is far from being achieved.

The most widely studied animal model of autoimmune diabetes is the nonobese diabetic (NOD) mouse. Diabetes in NOD mice is a T cell-dependent disorder, with the major histocompatibility complex (MHC) locus being the dominant genetic determinant, as in human T1D (8). NOD mice do not express the MHC class II E complex due to deletion of the Eα promoter (9), and genetically modified NOD mice expressing the Eα molecule are completely protected from T1D (10). We recently demonstrated that this protection operates via modulation of the composition of the intestinal microbiota during a critical early window of ontogeny (11). In particular, we observed that the microbiota of Eα16/NOD mice (protected from T1D) starts to differentiate from that of NOD mice, both in terms of alpha- and beta-diversity, before initiation of insulitis.

Here, we carried out a deep metaproteogenomic characterization of fecal microbiota collected from 10-week-old Eα16/NOD and NOD mice (6 Eα16/NOD mice versus 6 NOD mice) in order to identify which microbial functions specifically discriminate the two microbiotas before T1D onset and therefore possibly correlate to T1D protection/development. Materials and methods are detailed in Text S1 in the supplemental material.

Supplemental methods. Download Text S1, DOCX file, 0.02 MB.

Metaproteomic (MP) analyses led to the detection of 153,044 microbial peptide-spectrum matches, mapped to 12,790 microbial sequences, and assigned to 1,205 functions and 448 taxa/243 genera. In parallel, 2,139,608 reads were obtained upon 16S rRNA gene sequencing (16S), mapped to 14,022 operational taxonomic units (OTUs), and assigned to 639 taxa/297 genera.

First, we sought for taxonomic differences in the gut metaproteome profiles between Eα16/NOD and NOD mice. The application of a paired sample test (taking into account the cage effect, with each cage housing two littermates with different genotypes) enabled the identification of 19 differential taxa (Data Set S1). A global depletion of Firmicutes proteins was found in the NOD mice microbiota compared to the level in Eα16/NOD mice, mainly attributable to the order Clostridiales. At lower levels, 7 clostridial genera, namely, Blautia, Coprococcus, Dorea, Johnsonella, Lachnoclostridium, Oribacterium from Lachnospiraceae, and Oscillibacter from Oscillospiraceae, were detected with significantly lower abundances in the NOD mouse fecal microbiota, as displayed in Fig. 1. On the other hand, Lactobacillus proteins were found to be significantly higher in abundance in the NOD metaproteome than in the Eα16/NOD counterpart. MP results were globally consistent with 16S data (Data Set S2), with the large majority of differential OTUs and genera that were more abundant in Eα16/NOD mice being assigned to the Clostridiales (Fig. S1 and S2). Interestingly, all OTUs higher in abundance in NOD mice were classified as belonging to the order Bacteroidales. Our taxonomic results are consistent with those of previous investigations carried out with human subjects at the onset of T1D. Specifically, several 16S rRNA studies of the gut microbiota from children matched for age, sex, early feeding history, and human leukocyte antigen risk genotype and differing in T1D-associated autoantibodies demonstrated a striking depletion of butyrate-producing clostridial species and an increase in Bacteroidetes spp. in children with β-cell autoimmunity (12–15).

FIG 1

Microbial genera showing significantly differential abundances between Eα16/NOD and NOD mouse fecal metaproteomes. Genera detected in all samples of at least one group are shown. Each dot represents a different mouse, with Eα16/NOD and NOD littermates reared in the same cage having the same shape. A paired sample test based on an inverted beta binomial model was used, with each pair of samples comprising two littermates with different genotypes. *, q value < 0.05; **, q value < 0.01.

Bar graph illustrating the distribution of differential OTUs based on their taxonomic annotation (red, taxa higher in NOD mice microbiotas; blue, taxa higher in Eα16/NOD mice microbiotas). A paired sample (inverted beta binomial) test was applied to find the differential OTUs. Download FIG S1, TIF file, 1.4 MB.

Clostridial genera showing significantly differential abundances between Eα16/NOD and NOD mouse fecal microbiotas upon 16S rRNA gene analysis. Each dot represents a different mouse, with Eα16/NOD and NOD littermates reared in the same cage having the same shape. A paired sample test based on an inverted beta binomial model was used, with each pair of samples comprising two littermates with different genotype. *, q value < 0.05; **, q value < 0.01. Download FIG S2, TIF file, 0.6 MB.

Data set containing identification, annotation, and differential analysis information from metaproteomic data. Download Data Set S1, XLSX file, 2.1 MB.

Data set containing identification, annotation, and differential analysis information from 16S rRNA gene sequencing data. Download Data Set S2, XLSX file, 1.4 MB.

We then focused on the protein functions expressed by the members of the microbiota, as enabled by the MP approach. As a result, we found 59 functions with significantly differential abundances between the microbiotas of NOD and Eα16/NOD mice (Data Set S1). When associating phylum and function information, we found 105 phylum-specific functions showing differential abundances between the two genotypes (Data Set S1); interestingly, all 87 of those organisms more abundant in the microbiota of Eα16/NOD mice were associated with Firmicutes, while most of those more abundant in the microbiota of NOD mice belonged to Bacteroidetes. Protein functions that were more abundant in the NOD mouse microbiota and assigned to Bacteroidetes were implicated primarily in polysaccharide transport and degradation. Intriguingly, we observed the increase of a few Firmicutes proteins in the NOD mouse microbiota (a trend opposite to that of all other proteins from the same phylum), namely, a small acid-soluble spore protein, enolase, and a cell surface glycoprotein (S-layer protein), all possibly related to bacterial survival and host immune evasion (16), and phosphoketolase. On the other hand, proteins specifically enriched in the microbiotas of T1D-protected mice (mostly from the Clostridia) were related to many different activities, including enzymes involved in acetogenesis (carbon monoxide dehydrogenase), butyrogenesis (see below), pyruvate metabolism (including pyruvate synthase), polysaccharide metabolism (e.g., glycogen synthase), and amino acid biosynthesis (such as cysteine and tryptophan synthases), as well as proteins associated with carbohydrate transport, response to stress, and chemotaxis.

Remarkably, we found that all key enzymes involved in butyrate biosynthesis and expressed by members of the Clostridia were significantly depleted in the gut microbiotas of NOD mice, compared to those of T1D-protected Eα16/NOD mice. Figure 2 shows the expression profiles of six butyrogenic enzymes (mapping almost entirely to the butyrate biosynthesis pathway, from pyruvate to butyrate via butyrate kinase), along with the taxonomic assignments of the related peptides. As a result, the increase of butyrogenic activity in the microbiotas of T1D-protected mice was revealed to be due mainly to members of the Clostridiales, including Clostridium, Eubacterium, Roseburia, and Oscillibacter. Throughout the metabolic pathway, some catalytic functions were found to be family or even genus specific, while others were shared by several members of the microbiota or could not be mapped at the lowest taxonomic levels, being highly conserved sequences. Butyrate represents the main source of energy for enterocytes and is able to enhance mucus production in the gut, thus having a profound impact on the integrity of the gut barrier (3). Impaired butyrate production provides a possible explanation for the loss of tight barrier function which has been related to T1D pathogenesis according to the so-called "leaky gut" hypothesis (17). Additionally, butyrate is known to contribute significantly to the induction of colonic regulatory T cells (18, 19) and has been demonstrated to exert anti-inflammatory effects via several mechanisms, including histone deacetylase (20) and NF-κB (21) inhibition.

FIG 2

Expression of butyrate biosynthesis enzymes in the fecal microbiotas of Eα16/NOD and NOD mice according to metaproteomic results. (Top) Schematic of butyrate biosynthetic pathway (from pyruvate to butyrate), with Enzyme Commission numbers corresponding to each enzyme. EC 2.3.1.54, formate acetyltransferase; EC 2.3.1.9, acetyl-CoA acetyltransferase; EC 1.1.1.157, 3-hydroxybutyryl-CoA dehydrogenase; EC 4.2.1.150, short-chain-enoyl-CoA hydratase; EC 1.3.8.1, butyryl-CoA dehydrogenase; EC 2.3.1.19, phosphate butyryltransferase; EC 2.7.2.7, butyrate kinase. (Middle) Scatterplots illustrating the relative abundance of each butyrogenic enzyme in the fecal microbiotas of Eα16/NOD and NOD mice; each dot represents a different mouse, with Eα16/NOD and NOD littermates reared in the same cage having the same shape. A paired sample test based on an inverted beta binomial model was used, with each pair of samples comprising two littermates with different genotypes. *, q value < 0.05; **, q value < 0.01; ***, q value < 0.001. (Bottom) Heatmaps showing the relative distribution of the taxonomic assignments (according to the lowest-common-ancestor principle) associated with each enzyme. Taxa associated with at least one enzyme and found in at least 4 mice are reported. Boxes: 1, Firmicutes; 2, Firmicutes, Clostridia; 3, Firmicutes, Clostridia, Clostridiales; 4, Firmicutes, Clostridia, Clostridiales, Clostridiaceae, Clostridium; 5, Firmicutes, Clostridia, Clostridiales, Eubacteriaceae, Eubacterium; 6, Firmicutes, Clostridia, Clostridiales, Lachnospiraceae; 7, Firmicutes, Clostridia, Clostridiales, Lachnospiraceae, Roseburia; 8, Firmicutes, Clostridia, Clostridiales, Oscillospiraceae, Oscillibacter.

Indirect evidence that microbial production of butyrate also exerts an important role on preventing T1D was previously reported with different experimental approaches. The proportion of butyryl coenzyme A (butyryl-CoA) dehydrogenase genes was shown to be lower in the metagenomes of T1D cases than in those of healthy controls (12). Another study reanalyzed 16S rRNA data from 6-month-old children of a prospective cohort by applying microbial cooccurrence networks, suggesting a protective role of butyrate in T1D pathogenesis (22). Moreover, the direct administration of short-chain fatty acids (SCFAs) or SCFA-yielding diets to rodent models of T1D was proven to provide a high degree of protection from T1D-related autoimmune responses (23, 24). Specifically, butyrate was demonstrated to promote the abundance and the function of T regulatory cells.

In conclusion, we report here for the first time that butyrate biosynthesis enzymes are depleted in a T1D model. Hence, our data further support the hypothesis of involvement of SCFA production in T1D pathogenesis and specifically underline a possible key role of butyrate in the complex series of events leading to autoimmune insulitis and preceding T1D onset.

Data availability.

16S rRNA gene sequencing data have been deposited into the European Nucleotide Archive with the study accession number PRJEB25325.

The mass spectrometry proteomics data (including the above-mentioned sequence database) have been deposited to the ProteomeXchange Consortium via the PRIDE (25) partner repository with the data set identifier PXD003616.