Proteobacteria and Erysipelotrichaceae species are negatively associated with insulin sensitivity

Few studies have investigated the relationship between gut microbiota and insulin sensitivity using the hyperinsulinaemic-euglycaemic clamp.6–11 Importantly, all these studies assessed microbial composition using 16S rRNA sequencing, which did not provide information on microbial composition at the species levels or microbial function and did not take into account the compositional nature of the metagenomics datasets. To overcome these limitations, we analysed the faecal microbiota of three independent studies (ADIPOINST, n=15; IRONMET, n=121 and FLORINASH, n=67) by shotgun metagenomics. We applied linear regression models in the log-scale taking into account the underlying compositional structure of the metagenomics data using the analysis of compositions of microbiome with bias correction (ANCOM-BC)12 adjusting for age, sex and body mass index (BMI). We found a consistent negative association of insulin sensitivity with several genera and species (figure 1A-E and online supplemental tables 2–7) from the Erysipelotrichaceae family and the Proteobacteria phylum, specifically from the Enterobacteriaceae family, across three independent studies. We further validated these results on a subset of patients from the FLOROMIDIA cohort with available hyperinsulinaemic-euglycaemic clamps (n=24) using faecal 16S rRNA sequencing. Once again, the genus Desulfovibrio from the Proteobacteria phylum was strongly negatively associated with insulin sensitivity (online supplemental figure 1 and online supplemental table 8).

• Download figure
• Open in new tab
• Download powerpoint

Figure 1

Associations of the faecal microbiota composition and functionality with the hyperinsulinaemic-euglycaemic clamp across studies. (A–F) Volcano plots of differential microbial (a-c) genera and (d-f) species associated with insulin sensitivity (hyperinsulinaemic-euglycaemic clamp) in the ADIPOINST (n=15), IRONMET (n=121) and FLORINASH (n=67) studies, respectively; identified using the analysis of microbiomes with bias correction compared. The log₂ (Fold Change) and the −log₁₀ (p values) adjusted for multiple testing are plotted for each taxon. Significantly different taxa are coloured according to phylum. (g–i) Dot plots of the KEGG pathway over-representation analyses (q value <0.1) mapping the KEGG orthologues significantly associated with insulin sensitivity in the ADIPOINST, IRONMET and FLORINASH studies, respectively. Dots are coloured according to the q value. KEGG, Kyoto Encyclopaedia of Genes and Genomes.

Conversely, insulin sensitivity was positively associated with microbial species from the Bifidobacterium and Prevotella genera. In the ADIPOINST and FLORINASH cohort, we also found a positive association between short-chain fatty acid (SCFA)-producing species from the Blautia and Faecalibacterium genus and insulin sensitivity (figure 1D-E). Recent studies also demonstrate the causal role of gut microbiota in insulin sensitivity regulation. Faecal microbiota transplantation (FMT) from lean donors has been shown to enhance peripheral insulin sensitivity,7 9 whereas FMT from metabolic syndrome donors has had the opposite effect.6 Consistent with our results, FMT from lean male donors to male subjects with obesity increased SCFA-producing bacteria like Bifidobacterium pseudolongum, while reducing Escherichia coli levels compared with the autologous group.9 Additionally, increases in Desulfovibrio spp. predicted a ≥10% decrease in insulin sensitivity, whereas Prevotellaceae spp predicted non-deterioration.6 Similarly, treatment with vancomycin in male patients with obesity decreased peripheral insulin sensitivity and increased gram-negative bacteria, mainly Proteobacteria.8 Elevated Proteobacteria levels have also been observed in individuals with T2D,13 14 although these studies often did not account for medication use. Conversely, other studies found no clear link between increased Proteobacteria and insulin sensitivity after vancomycin treatment in obese individuals,15 and under certain conditions, a decrease in Proteobacteria has been associated with insulin sensitivity or with improved metabolic outcomes, such as after gastric bypass16 17 or in toll-like receptor 2 knockout mice.18

Medication is known to impact the gut microbiota.19 20 Therefore, to rule out the confounding effects of medication on our current associations, we conducted additional analyses controlling for hypertension and dyslipidaemia medication, as well as the consumption of proton pump inhibitors, which were the primary medications in our studies. Additionally, although the ADIPOINST and IRONMET studies did not include patients with T2D (HbA1c <6.5% and fasting glucose <126 mg/dL), nine patients (13%) in the FLORINASH cohort had T2D. To further eliminate the confounding effect of diabetes, we controlled the analyses for the presence of T2D in the FLORINASH cohort. Remarkably, after controlling for these additional covariates, not only did most of the associations remain significant, but genera from the Proteobacteria phylum also exhibited the strongest negative fold changes and most significant (p_adj) associations with insulin sensitivity (online supplemental figure 2 and online supplemental tables 9–11).

Previous studies have overlooked microbial functionality analyses, which are essential to accurately capture the host-microbiome interactions due to microbial functional redundancy.21 Hence, we next performed functional analyses based on Kyoto Encyclopaedia of Genes and Genomes (KEGG) orthologs in the three studies (online supplemental tables 12–14) and a pathway enrichment analysis using those microbial genes associated with insulin sensitivity (p_adj<0.05). In the ADIPOINST, IRONMET and FLORINASH studies, we found a consistent over-representation of pathways involved in energy metabolism such as the tricarboxylic acid cycle, oxidative phosphorylation, glycolysis/gluconeogenesis, pyruvate metabolism and the pentose phosphate pathway (figure 1G-I and online supplemental tables 15–17). Notably, both the IRONMET and FLORINASH studies showed a strong association with lipoic acid and histidine metabolism. Alpha-lipoic acid has been demonstrated to increase insulin sensitivity measured with hyperinsulinaemic-euglycaemic clamp in patients with T2D.22–24 A functional analysis based on KEGG modules identified histidine degradation as one of the most significant pathways (online supplemental figure 3 and online supplemental tables 18–20). Notably, imidazole propionate, a microbially produced histidine metabolite, is elevated in subjects with prediabetes and diabetes and has been shown to impair insulin signalling.25 26

A T cell-related jejunal transcriptomic signature is positively associated with systemic insulin sensitivity

Recent evidence suggests the gut is a key tissue in obesity-related insulin resistance through the gut microbiota-induced inflammation27 by disrupting intestinal permeability and enhancing absorption of microbial products such as LPS.3 4 Most of the digestion of lipids and sugars takes place in the jejunum and is regulated by insulin.28 However, little is known about the jejunal transcriptome in relation to insulin sensitivity. Considering this, we performed an RNAseq of jejunal samples from the ADIPOINST cohort (n=26) and applied robust linear regression models adjusted for age, BMI and sex, to identify jejunal transcriptome signatures associated with the hyperinsulinaemic-euglycaemic clamp (figure 2A and online supplemental table 21). A pathway over-representation analysis highlighted that most of the upregulated genes played a key role in both the CD4 and CD8 helper T cell receptor signalling (TCR) and differentiation (Th1 and Th2), innate immunity, modulation of the inflammatory response, the RHO GTPase cycle and maintenance of the intestinal epithelium integrity (figure 2B,C and online supplemental table 22). CR2, involved in the complement cascade, and CCL11, encoding for an eosinophil-specific chemokine, were the two jejunal genes that had the highest effect size. These findings hint at the preservation of intestinal immunity and integrity with increased insulin sensitivity.

• Download figure
• Open in new tab
• Download powerpoint

Figure 2

Cross-cohort and cross-omics associations in the small intestine with insulin sensitivity and resistance. (a) Volcano plot of differentially expressed genes associated with insulin sensitivity (hyperinsulinaemic-euglycaemic clamp) in the jejunum of patients from the ADIPOINST cohort (n=26) identified by limma-voom analysis controlling for age, sex and BMI. The log₂ fold change associated with a unit change in the clamp and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (b) Dot plot of pathways significantly associated (q value <0.1) with insulin sensitivity in the jejunum identified from a pathway over-representation analysis mapping significantly upregulated genes to the Reactome, Kyoto Encyclopaedia of Genes and Genomes, Wikipathways, PID and NetPpath databases. Dots are coloured by the q value. (c) Over-representation analysis results were mapped as a function network of pathways using an enrichment map. Edges connect overlapping gene sets, while node size reflects the total number of genes in each pathway. Overlapping gene sets tend to cluster together, making it easy to identify functional modules. Functionally related pathways are clustered based on the Markov Cluster Algorithm and coloured with the same colour. (d–g) Scatter plot of the partial Spearman’s rank correlations (adjusted for age, sex and BMI) between the fasting glucose or HbA1c levels and the number of T cells or the percentage of T cells in the ileum of patients from the SIMMUNIDIA cohort (n=43). The ranked residuals are plotted. (h) Volcano plots of differential microbial genera associated with HbA1c and (i) the number (#) of T cells in the ileum of the SIMMUNDIA cohort (n=42) identified using ANCOM-BC controlling for age, sex and BMI. The log₂ (Fold Change) and the −log₁₀ (p values) adjusted for multiple testing are plotted for each taxon. Significantly different taxa are coloured according to phylum. (j) Boxplots of the normalised variable importance measure for the metabolites/lipids associated with insulin sensitivity in the jejunum of the ADIPOINST cohort. The red dot represents the mean and the colour bar above each plot indicates the sign of the association between the metabolites/lipids and insulin sensitivity, with red indicating negative correlation and green positive correlation. Significant metabolites were identified using the Boruta algorithm with 5000 trees and 500 iterations. ANCOM-BC, analysis of microbiomes with bias correction; BMI, body mass index; HbA1c, glycated haemoglobin.

As the TCR signalling pathway had the strongest association with insulin sensitivity in the jejunum, we sought to measure the content of T cells in the ileum of another independent cohort (SIMMUNIDIA, n=42). We found that both the numbers and proportion of ileal T cells were negatively associated with fasting glucose and HbA1c (figure 2D-G). While interpreting these results, it is important to consider that the ileum has a slightly distinct immune microenvironment than the jejunum.29 Overall, these findings indicate that the presence of T cells seems to exercise a protective function in the small intestine. Similarly, the frequency of Th2 cells has also been negatively correlated with insulin resistance in mice and humans in other tissues like VAT.30 31

Ileal Proteobacteria and deoxycholic acid (DCA) exhibit a negative association with insulin sensitivity

The gut microbiota and derived metabolites are key in intestinal barrier disruption and local immune responses that contribute to inflammation and lead to metabolic disease.32 However, only a limited number of studies have described the composition of the small intestine microbiota, and there is even little evidence about its relationship with obesity and insulin resistance.33 Thus, we analysed the microbiota of the ileum mucosa in the SIMMUNIDIA cohort (n=42) using 16S rRNA gene sequencing and identified bacterial genera associated with HbA1c as an indirect marker of insulin resistance (figure 2H). In agreement with our findings in faeces, bacterial genera from the Proteobacteria phylum (Escherichia, Shigella and Enterobacter) were associated with insulin resistance (figure 2H and online supplemental table 23) and negatively with the total number of T cells (figure 2I and online supplemental table 24), while the Bifidobacterium genus and SCFA-producing genera such as Faecalibacterium were negatively associated. In line with our results, a high-fat diet (HFD) resulted in alterations in the ileum microbiota that impaired the immune system of the small intestine by decreasing the number of interleukin (IL)-17-producing CD4 T cells (Th17) and induced T2D in mice.34

We next profiled the lipidome and metabolome of the jejunum in the ADIPOINST cohort (n=26) to elucidate the small intestine biochemical pathways linked to insulin sensitivity. We applied a random forest-based variable selection machine learning algorithm to identify lipids and metabolites predictive of the hyperinsulinaemic-euglycaemic clamp. DCA, a microbial-derived secondary bile acid, was identified as the most important metabolite negatively associated with insulin sensitivity (figure 2J and online supplemental table 25). Valine and isoleucine were also strongly negatively associated with insulin sensitivity, whereas taurine had a positive association.

Consistently, jejunal levels of DCA were positively associated with the presence of proinflammatory bacterial species in faeces such as Gammaproteobacteria (E. cloacae, E. coli and Citrobacter) and negatively associated with species from the Bifidobacterium and Faecalibacterium genera (online supplemental figure 4 and online supplemental table 26). In line with our results, FMT from donors with metabolic syndrome resulted in an increase in faecal levels of DCA and other secondary bile acids while decreasing the insulin sensitivity of recipient patients.6 DCA is known to induce intestinal inflammation, disrupting the epithelial barrier and increasing gut permeability.35–37 Elevated levels of intestinal DCA have been linked to HFDs and western diets38 and increased abundance of the genus Clostridium, the main producer of DCA, promoting colonic inflammation.38 Plasma levels of DCA were associated with insulin resistance and significantly elevated in patients with T2D.39

Proteobacteria and a T cell-related transcriptomic signature in the colon are also associated with insulin sensitivity

We also performed a colon transcriptomics analysis in the FLOROMIDIA cohort to identify genes associated with insulin sensitivity (figure 3A and online supplemental table 27). Consistent with our findings in the jejunum, pathways involved in the CD4/CD8 helper T cell lineage commitment and the Th17 type immune response were enriched (figure 3B, online supplemental table 28). We found that insulin sensitivity was positively associated with interleukin 17 receptor E (IL17RE) and interleukin 6 receptor (IL6R), but negatively with IL-23 subunit alpha (IL23A) (figure 3A,C). Similarly, the fasting glucose levels in the SIMMUNIDIA cohort were associated with the expression of IL-17 receptor B (IL17RB) and tumour necrosis factor receptor (CD40) (figure 3D and online supplemental table 29). We then analysed the microbiota of the colon in the SIMMUNIDIA cohort (n=55) using 16S rRNA gene sequencing. In agreement with our results at the ileum, we found that genera from the Proteobacteria phylum (Escherichia, Shigella and Gammaproteobacteria) were positively associated with insulin resistance assessed with the HbA1c and negatively associated with the number of T cells in the colon, whereas Faecalibacterium were again negatively associated (figure 3E,F and online supplemental tables 30,31).

• Download figure
• Open in new tab
• Download powerpoint

Figure 3

Cross-cohort and cross-omics associations in the colon with insulin sensitivity and resistance. (a) Volcano plot of differentially expressed genes associated with insulin sensitivity (hyperinsulinaemic-euglycaemic clamp) in the colon of patients from the FLOROMIDIA cohort (n=22), identified by limma-voom analysis controlling for age, sex and BMI. The log₂ fold change associated with a unit change in the clamp and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (b) Dot plot of gene ontology-biological processes significantly associated (q value <0.1) with insulin sensitivity in the colon identified from a gene ontology over-representation analysis using significant genes associated with insulin sensitivity. Dots are coloured by the q value. (c) Gene-concept network depicting the linkage of significant genes associated with insulin sensitivity participating in Th17 immune response and CD4+ or CD8+, alpha-beta T cell lineage commitment. (d) Volcano plot of differentially expressed genes associated with the fasting glucose levels in the colon of patients from the SIMMUNIDIA cohort (n=22), identified by limma-voom analysis controlling for age, sex and BMI. The log₂ fold change and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (e) Volcano plots of differential microbial genera associated with HbA1c and (F) the number (#) of T cells in the colon of the SIMMUNDIA cohort (n=55) identified using ANCOM-BC controlling for age, sex and BMI. The log₂ (Fold Change) and the –log₁₀ (p values) adjusted for multiple testing are plotted for each taxon. Significantly different taxa are coloured according to phylum. ANCOM-BC, analysis of microbiomes with bias correction; BMI. body mass index; HbA1c, glycated haemoglobin.

An adipose tissue signature involved in T cell and Rho GTPase signalling is associated with systemic insulin sensitivity and resistance

Chronic inflammation in the adipose tissue has become a hallmark of obesity-related insulin resistance. Hence, we performed an RNA-seq analysis of the subcutaneous adipose tissue (SAT) of patients with and without obesity from the ADIPOMIT cohort (n=740). After fitting robust linear regression models controlling for age, BMI and sex, we identified several genes associated with HbA1c (figure 4A and online supplemental table 32). An enrichment analysis highlighted that genes negatively associated with HbA1c were involved in pathways that included insulin signalling, insulin resistance, IRS activation or leptin signalling, whereas genes positively associated with HbA1c were mainly involved in the TCR signalling, the immune system and the signalling by Rho GTPases (figure 4B,C and online supplemental tables 33,34), consistent with our findings in the jejunum and colon. However, while higher expression of the genes associated with these pathways in the intestine appears to be linked to greater insulin sensitivity, in adipose tissue, it seems to be associated with increased insulin resistance. This disparity might indicate different stages of development and/or asynchronous responses to environmental stimuli. In the intestine, higher insulin sensitivity could correspond to a stronger immune response, suggesting that an adequate intestinal immune response might prevent gut dysbiosis-induced insulin resistance. Conversely, in adipose tissue, an increased immune response could indicate insulin resistance-associated adipose tissue dysfunction. Chronic activation of this response in adipose tissue can become detrimental. However, it should be noted that these interpretations remain speculative.

• Download figure
• Open in new tab
• Download powerpoint

Figure 4

SAT transcriptomic signatures associated with HbA1c. (a) Volcano plot of differentially expressed genes associated with HbA1c in the SAT of patients from the ADIPOMIT cohort (n=740), identified by limma-voom analysis controlling for age, sex and BMI. The log₂ fold change and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (b) Manhattan-like plot of pathways significantly associated (q value <0.1) with insulin sensitivity in the jejunum identified from a pathway over-representation analysis mapping significantly downregulated genes and (c) significantly upregulated genes to the Reactome, Kyoto Encyclopaedia of Genes and Genomes, Wikipathways, PID and NetPath databases. HbA1C, glycosylated haemoglobin; SAT, subcutaneous adipose tissue.

We next performed an RNA-seq analysis of the VAT of patients from the ADIPOINST cohort (n=29). After applying robust linear regression models adjusted for age, BMI and sex, we identified a large number of VAT genes associated with insulin sensitivity (figure 5A and online supplemental table 35). Remarkably, a pathway enrichment analysis based on the genes negatively associated with insulin sensitivity revealed associations with innate immunity, specifically an over-representation of the TCR signalling pathway; the signalling by Rho GTPases; the IL4-, IL8-, IL3-, IL-18-, IL2-mediated signalling events and anticipated metabolic pathways such as leptin, adiponectin, glucagon and phosphatidylinositol 3-kinase (PI3K)/protein kinase B signalling (figure 5B and online supplemental table 36).

• Download figure
• Open in new tab
• Download powerpoint

Figure 5

VAT transcriptomic signatures associated with insulin sensitivity and cross-omics and cross-tissue integration in the ADIPOINST cohort. (a) Volcano plot of differentially expressed genes associated with insulin sensitivity (hyperinsulinaemic-euglycaemic clamp) in the VAT of patients from the ADIPOINST cohort (n=29) identified by limma-voom analysis controlling for age, sex and BMI. The log₂ fold change associated with a unit change in the clamp and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (b) Manhattan-like pathways significantly associated (q value <0.1) with insulin sensitivity in the VAT identified from a pathway over-representation analysis mapping significantly downregulated genes to the Reactome, Kyoto Encyclopaedia of Genes and Genomes, Wikipathways, PID and NetPath databases. (c) Correlation circle plot for the integration of the jejunum and VAT genes, jejunum metabolites and faecal microbial species associated with insulin sensitivity in the ADIPOINST cohort using a multiblock sparse PLS model. Strongly positively associated variables or groups of variables are projected close to one another on the correlation circle (~0° angle). The variables or groups of variables strongly negatively associated are projected diametrically opposite (~180° angle) on the correlation circle. Variables not correlated are situated ~90° from one another. (d) Dot plot of significantly enriched (q value <0.1) gene ontology-biological processes and (e) pathways from jejunal genes strongly positively associated with insulin sensitivity included in cluster 2 with a heatmap displaying the gene participating in each biological term. Dots are coloured by the q value and genes in the heatmap are coloured by the log₂ Fold Change of the association with insulin sensitivity. (f) Dot plot of significantly enriched (q value <0.1) gene ontology-biological processes and (g) pathways from VAT and (h) jejunal genes strongly negatively associated with insulin sensitivity included in cluster 1 with a heatmap displaying the gene participating in each biological term. Dots are coloured by the q value and genes in the heatmap are coloured by the log₂ Fold Change of the association with insulin sensitivity. (i) Volcano plots of differential microbial genera associated with insulin sensitivity in the VAT of the IRONMET (n=12) identified using ANCOM-BC controlling for age, sex and BMI. The log₂ (Fold Change) and the –log₁₀ (p values) adjusted for multiple testing are plotted for each taxon. Significantly different taxa are coloured according to phylum. ANCOM-BC, analysis of microbiomes with bias correction; BMI, body mass index; VAT, visceral adipose tissue

A multiomics and multitissue integration links Proteobacteria with jejunal DCA, jejunal and VAT genes involved in actin cytoskeleton, insulin and T cell signalling

Single-omics analyses only provide information about one layer of biological regulation. Therefore, we integrated the metagenomic, metabolomic and transcriptomic (jejunal and VAT) profiles linked to insulin sensitivity in the ADIPOINST cohort to obtain a holistic picture of the mechanisms underlying the gut microbiome-host interactions related to insulin sensitivity. Due to the large number of VAT genes associated with insulin sensitivity, we applied a multivariate method based on the multiblock-sparse projection to latent structures (MB-sPLS) using an L1 regularisation to select the 50 most relevant VAT genes and improve interpretation.

Cluster 1 depicted a strong negative association with insulin sensitivity and jejunal metabolites and lipids (DCA, valine, isoleucine, dehydroepiandrosterone sulfate (DHEAS)), microbial genera from the Proteobacteria phylum (Escherichia, Shigella and Legionella) and the Erysipelotrichaceae family (Allobaculum and Solobacterium) and transcriptomic signatures in the jejunum and VAT (figure 5C). Notably, this cluster was strongly negatively associated with cluster 2, with a strong correlation among metabolites (taurine), microbial genera (Faecalibacterium), jejunal and VAT genes strongly positively associated with insulin sensitivity.

DHEAS), microbial genera from the Proteobacteria phylum (Escherichia, Shigella and Legionella) and the Erysipelotrichaceae family (Allobaculum and Solobacterium) and transcriptomic signatures in the jejunum and VAT (figure 5C). Notably, this cluster was strongly negatively associated with cluster 2, with a strong correlation among metabolites (taurine), microbial genera (Faecalibacterium), jejunal and VAT genes strongly positively associated with insulin sensitivity.

To explore the role of the genes of both clusters in relation to insulin sensitivity, we conducted gene ontology-biological process (GO-BP) and pathway enrichment analyses. The most significant BP from jejunal genes in cluster 2 (positively associated with insulin sensitivity) involved lymphocyte migration and the regulation of the immune system (figure 5D). Consistently, the most significant jejunal pathways from this cluster included modulators of TCR signalling and T cell activation, the complement cascade, the innate immune system and the Rho GTPase cycle (figure 5E), whereas the most significant VAT pathway was precisely involved in the regulation of insulin-like growth factor (IGF) transport and uptake by IGFBPs (q value=0.0008).

Additionally, supramolecular fibre organisation (CCL1, CHRFAM7A, CYFIP2 and RHOH) and actin cytoskeleton organisation (CCL11 and RHOH) were among the GO-BP over-represented in cluster 2. Actin disassembly in both adipose and muscle cells has been shown to inhibit insulin-induced events such as glucose transporter recruitment to the cell surface and enhanced glucose transport40; actin is indispensable for the insulin-stimulated translocation of the glucose transporter GLUT4.41 Plasma CCL11, has been positively associated with steatosis severity42 and CCL11 knockout mice fed a high-fat, high-calorie diet improved insulin sensitivity.43 Additionally, RHOH is a key factor for the development and proper function of T cells.44 45

Of note, the GO-BP and pathways from the VAT most strongly associated with cluster 1 (negatively associated with insulin sensitivity) comprised again the regulation of the actin cytoskeleton organisation and regulation but also the insulin signalling (figure 5F,G) Similarly, the over-represented BP from these cluster in the jejunum also included the musculoskeletal movement and the supramolecular fibre organisation (figure 5H).

Finally, to assess whether the negative associations between faecal Proteobacteria and insulin sensitivity translated to the adipose tissue, we performed 16S rRNA sequencing of the VAT in a subset of patients from the IRONMET cohort (n=12). In line with our results, we found that genera from the Proteobacteria phylum (Bradyrhizobium and Moraxellaceae families) in the VAT were negatively associated with insulin sensitivity (figure 5I and online supplemental table 37).

A multiomics and multitissue integration links proteobacteria, liver DCA and liver genes involved in actin cytoskeleton, insulin secretion and the immune system

The liver plays a central role in modulating glucose and insulin. Therefore, we next profiled the liver transcriptome of obese patients from the FLORINASH cohort (n=80). After fitting robust linear regression models adjusted for age, BMI and sex, we identified 1356 genes associated with the hyperinsulinaemic-euglycaemic clamp (figure 6A and online supplemental table 38). Consistent with our findings in the jejunum and the VAT, these genes were involved in the TCR signalling, the complement system, the Ras and Rap1 GTPase signalling and the cytokine–cytokine receptor interaction (figure 6B and online supplemental table 39). We also profiled the lipidome of the liver of a subset of these patients. We found several lysophosphatidylcholines and lysophosphatidylethanolamines associated with insulin sensitivity (figure 6C, online supplemental figure 5A and online supplemental table 40). Importantly, we found again a negative association of DCA with insulin sensitivity (figure 6C and online supplemental figure 5A). We also analysed the microbiota of the liver of a subset of participants from Spain (n=17) and Italy (n=20) in the FLORINASH study. Consistent with our previous results in faeces, ileum, jejunum and colon, we found that genera from the Proteobacteria phylum in the liver exhibited the strongest and/or most significant negative associations with insulin sensitivity, after controlling for age, BMI, sex, dyslipidaemia and hypertension medication and diabetes both in Spain (figure 6D and online supplemental table 41) and Italy (online supplemental figure 5B and online supplemental table 42). A MB-sPLS-based multiomics integration of the faecal microbiota and liver lipids and genes (only the 50 most relevant using an L1 regularisation) associated with insulin sensitivity revealed the presence of a cluster of highly interconnected features including genera from the Proteobacteria phylum (Deltaproteobacteria_uc, Desulfovibrionaceae_uc, Bilophila and Desulfovibrionales_uc) and the Erysipelotrichaceae family, primary (taurocholic acid) and secondary (DCA) bile acids and a transcriptomic signature negatively associated with insulin sensitivity (figure 6E).

• Download figure
• Open in new tab
• Download powerpoint

Figure 6

Liver transcriptomic signatures associated with insulin sensitivity and cross-omics and cross-tissue integration in the FLORINASH cohort. (a) Volcano plot of differentially expressed genes associated with insulin sensitivity (hyperinsulinaemic-euglycaemic clamp) in the liver of patients from the FLORINASH cohort (n=80) identified by limma-voom analysis controlling for age, sex and BMI. The log₂ fold change associated with a unit change in the clamp and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (b) Manhattan-like plot of pathways significantly associated (q value <0.1) with insulin sensitivity in the liver identified from a pathway over-representation analysis mapping significant genes to the Reactome, KEGG, Wikipathways, PID and NetPpath databases. (c) Boxplots of the normalised variable importance measure for the lipids associated with insulin sensitivity in the liver of the FLORINASH cohort (n=41). The bar above each plot indicates the sign of the association between the lipids and insulin sensitivity, with red indicating a negative correlation and green a positive correlation. Significant metabolites were identified using the Boruta algorithm with 5000 trees and 500 iterations. (d) Volcano plots of differential microbial genera associated with insulin sensitivity (clamp) in the liver of a subset of patients from Spain of the FLORINASH cohort (n=17) identified using ANCOM-BC controlling for age, sex and BMI. The log2 (Fold Change) and the −log10 (p values) adjusted for multiple testing are plotted for each taxon. Significantly different taxa are coloured according to phylum. (e) Correlation circle plot for the integration of the liver lipids, liver genes and faecal microbial species associated with insulin sensitivity in the FLORINASH cohort using a multiblock sparse partial least squares model. Strongly positively associated variables or groups of variables are projected close to one another on the correlation circle (~0° angle). The variables or groups of variables strongly negatively associated are projected diametrically opposite (~180° angle) on the correlation circle. Variables not correlated are situated ~90° from one another. (f) Dot plot of significantly enriched (q value <0.1) pathways (based on Reactome, KEGG, Wikipathways, PID and NetPath) and (g) heatmap of genes participating in these pathways, identified from liver genes strongly negatively associated with insulin sensitivity from the transcriptome signature clustering with Proteobacteria and deoxycholic acid and negatively associated with insulin sensitivity. Dots are coloured by the q value and genes in the heatmap are coloured by the log₂ Fold Change of the association with insulin sensitivity. KEGG, Kyoto Encyclopaedia of Genes and Genomes.

Consistent with our previous results in the jejunum and VAT, a pathway enrichment analysis on this transcriptome identified an over-representation of the Rap1 signalling pathway, the leucocyte transendothelial migration, the regulation of the actin cytoskeleton, the chemokine signalling pathway, the immune system, the cytokine signalling in the immune system and the secretion of insulin, among others (figure 6F,G and online supplemental table 43). Again, a GO-BP analysis highlighted an enrichment of pathways mainly involved in lymphocyte activation, leucocyte cell–cell adhesion, regulation of leucocyte migration and T cell activation (online supplemental table 44).

HFD and mono associations with Enterobacter cloacae induce transcriptional changes associated with insulin resistance in D. melanogaster intestine and fat body

To further explore the multiorgan transcriptional changes associated with insulin sensitivity, we used the model organism D. melanogaster. To model insulin resistance in Drosophila, we administer a HFD that causes an obese-like phenotype in flies, including increased circulating glucose levels and insulin-like resistance.46 To further assess the role of Proteobacteria in insulin resistance, E. cloacae was selected for monocolonisation experiments in D. melanogaster due to its dual relevance as both a commensal of the Proteobacteria phylum in Drosophila47 48 and its established association with obesity-related phenotypes in humans and animal models (figure 7A).49–51 We depicted the Drosophila orthologues of the jejunal genes identified in cluster 2 involved in the regulation of the immune system CR2 (hig), DOCK8 (Zir), GRAP2 (drk), CYFIP2 (Cyflip), SPN (Deaf) and RHOH (Rac1) (figure 5D) and assessed relative gene expression by qRT-PCR in the fly’s intestine fed with standard diet (SD) or HFD either raised in axenic conditions or mono-associated with E. cloacae. HFD-fed flies exhibited a significant decrease in the expression of hig, Rac1, Cyflip and Deaf, which aligns with our findings in humans, as these genes were positively associated with insulin sensitivity in the jejunum (figure 7B-G). Similarly, SD-fed flies mono-associated with E. cloacae showed a significant reduction in the intestinal expression of hig, Rac1 and a tendency for Deaf (figure 7B-G). These findings collectively suggest that both HFD and mono-associations with E. cloacae lead to alterations in intestinal gene expression patterns reminiscent of those observed in insulin-resistant humans.

• Download figure
• Open in new tab
• Download powerpoint

Figure 7

Enterobacter cloacae mono-association and HFD supplementation produce transcriptional changes in orthologues of insulin sensitivity-associated genes in the Drosophila intestine and fat body. (a) An experimental scheme was followed to generate Drosophila wild-type flies under sterile (germ-free) or monocolonisation conditions. Sterile flies can be easily generated by egg sterilisation. Subsequently, mono-associations were established by fly food supplementation with E. cloacae or the vehicle (sterile flies). On the 10th day of adulthood, dissections were conducted. 25–30 fly intestines or fat bodies were collected per sample. For complete flies, eight adults were collected per sample. Quantitative reverse transcription-PCR results; bars represent relative gene expression of b) hig, c) Rac1, (d) Cyflip, (e) Deaf, (f) drk, (g) Zir, (h) Col4a1, (i) Moe, (j) siz, (k) Past1, (l) cta, (m) ssh, (n) vib, (o) Ac76E, (p) Dsor1, (q) mys, (r) dlip2 and (s) Hpd in flies fed with SD and HFD non-colonised with E. cloacae or left sterile for the. P values were determined using the one-way analysis of variance combined with Fisher’s (least significant difference) multiple comparisons test when unequal variances Kruskal-Wallis non-parametric test with multiple comparisons was conducted. (t) Larval haemolymph glucose clearance after glucose feeding at 0, 10, 20, 30 and 90 min. Statistical significance at each time point was tested using an unpaired, two-sided t-test. n=3 independent replicates pulling eight larvae each. Error bars represent SE of the mean (#p<0.1, *p<0.05, *p<0.01 and *p<0.001). (u) Schematic representation of FOXO activation through insulin signalling in Drosophila. HFD, high-fat diet; SD, standard diet

Next, we conducted gene expression analysis in the Drosophila fat body, which serves as the equivalent of human adipose tissue. We assessed the expression of Cluster 1 genes, Past1, ssh, siz, cta, Dsor1 and Col4a1, which correspond to the human homologues of EHD2, SSH1, IQSEC1, GNA12, MAP2K1 and COL4A2, respectively (figure 5F,G). HFD-fed flies displayed a significant increase in the expression of Col4a1, siz and Past1 compared with SD-fed flies (figure 7H-M). In humans, these genes were negatively associated with insulin sensitivity in the VAT (figure 5F,G). Mono-associations with E. cloacae did not affect FB expression of the selected genes.

Fruit flies do not contain an organ equivalent to the liver; its functions are assumed by the fat body.52 53 Therefore, we assessed fat body gene expression of vib, Moe, Ac76E and mys (figure 7N-Q), the Drosophila orthologues of PITPNA, MSN, ADCY7, and ITGB2, all of them associated with the hyperinsulinaemic-euglycaemic clamp in liver transcriptomics (figure 6G). HFD-fed flies mono-associated with E. cloacae exhibited a significant increase of vib expression in the FB (figure 7N). Remarkably, we observed that PITPNA, its corresponding human orthologue, was depicted as one of the genes associated with insulin resistance together with species from the Proteobacteria phylum (figure 6E). Furthermore, Moe, the Drosophila orthologue of MSN, which displayed an inverse association with insulin sensitivity in the human liver, exhibited a significant upregulation in HFD-fed flies (figure 7O).

Finally, in order to investigate the role of Proteobacteria in insulin resistance and validate our model, we assessed the expression of dilp2 (Drosophila insulin-like peptide 2), one of the Drosophila orthologues of human insulin. We observed increased expression of dilp2 in HFD-fed flies compared with SD-fed flies, resembling insulin resistance as previously described.46 Mono-association of flies with E. cloacae further elevated the levels of dilp2 expression of SD-fed flies, indicating that this bacterium exacerbates insulin resistance also in D. melanogaster (figure 7R). Additionally, flies monocolonised with E. cloacae showed a significant decrease in the expression of hydroxyphenylpyruvate dehydrogenase (Hpd) compared with axenic flies (figure 7S). Hpd expression is regulated by the transcription factor dFOXO, which is activated downstream of PI3K in response to insulin.54 In conditions of insulin resistance, there is a decrease in dFOXO signalling and the expression of its downstream genes55 (figure 7U). Finally, when challenged for glucose tolerance,56 larvae monocolonised with E. cloacae did not clear glucose from haemolymph as efficiently as sterile larvae. Ninety min after ingesting a glucose solution, the E. cloacae-colonised larvae had significantly higher levels of glucose in the haemolymph (figure 7T). Given that Drosophila larvae primarily use trehalose as a circulating sugar,56 we also measured the levels of circulating trehalose in addition to glucose. Our findings indicated that larvae monocolonised with E. cloacae consistently exhibited higher levels of trehalose and glucose in the haemolymph compared with sterile larvae (online supplemental figure 6).

Overall, these experiments align with the transcriptional changes observed in the human intestine, VAT and liver, correlating with insulin sensitivity and further supporting the role of the Proteobacteria in the exacerbation of insulin resistance pathophysiology when present in the intestine.

Fasting glucose levels are linked to antiviral response genes

Our results identified consistent associations between insulin sensitivity and glycaemic state with specific microbiota and T cell numbers and signalling in different tissues. Given the importance of T cells in responding to viral infections, and our observation of intriguing associations between insulin sensitivity and viral infections such as Epstein-Barr virus infection, viral myocarditis and human cytomegalovirus infection (figure 2B and figure 6F), we also analysed the associations of fasting glucose with the transcriptome in the jejunum, ileum, VAT and liver (figure 8A-C,G and online supplemental tables 45–48). Strikingly, a GP-BP analysis revealed the defence response to the virus (figure 8D,H and online supplemental table 49) and viral transcription (figure 8E,H and online supplemental table 50) as the most over-represented processes in the jejunum and ileum. It also highlighted the type I interferon (IFN) production and the cytoplasmic pattern recognition receptor pathway in response to the virus (figure 8D). A KEGG over-representation analysis in the colon and VAT identified a significant enrichment of the viral protein interaction with cytokine and cytokine receptor (figure 8F,H and online supplemental table 51), and the viral infection, including the hepatitis C, influenza A and cytomegalovirus infection, which have been associated with an increased risk for the development of prediabetes and diabetes.57 A pathway over-representation analysis in the jejunum also highlighted an enrichment of the IFN signalling, the antiviral mechanism by IFN-stimulated genes and the oligoadenylate synthetase antiviral response (online supplemental table 52), while the type II IFN signalling was also over-represented in the colon (figure 8H and online supplemental table 53). A recent study also identified that downregulation of glucose metabolism promoted retinoic acid-inducible gene I-like receptor (RLR)-induced type I IFN production along with a reduction in virus replication and identified elevated lactate dehydrogenase and/or glycolysis as a potential mechanism for viruses to evade host defence by inhibiting RLR-triggered type I IFN production.58 Type I IFNs play a crucial role in protecting cells from viral pathogens. A screening of over >380 human IFN-stimulated genes revealed eight broadly acting effectors, including IRF7 and IFIH1,59 which we also found associated with glucose levels in the ileum (figure 8H). Genes identified in the jejunum included OASL, a gene involved in targeted antiviral specificity and well-known antiviral effectors such as EIF2AK2, or the IFN-induced genes, OAS1, OAS2 and OAS3, which play a key role in innate cellular antiviral response and have recently been identified as associated with critical illness caused by COVID-19.60

• Download figure
• Open in new tab
• Download powerpoint

Figure 8

Cross-cohort and cross-tissue transcriptomic associations with fasting glucose levels. Volcano plot of differentially expressed genes associated with the fasting glucose levels in the (a) jejunum, (b) ileum, (c) VAT and (g) colon from the ADIPOINST (n=26), SIMMUNIDIA (n=16), ADIPOINST (n=29) and FLOROMIDIA (n=37) studies, respectively. The log₂ fold change associated with a unit change in the fasting glucose levels and the log₁₀ p values adjusted for multiple testing are plotted for each gene. (d–f) Dot plot of gene ontology biological processes (d,e) and Kyoto Encyclopaedia of Genes and Genomes pathways (F) significantly associated (q value <0.1) with the fasting glucose levels in the jejunum, ileum and VAT, respectively. (h) Gene-concept network depicting the linkage of significant genes associated with fasting glucose levels participating in pathways involved in the antiviral response in the jejunum, ileum, VAT and colon. VAT, visceral adipose tissue.