Gut microbial carbohydrate metabolism contributes to insulin resistance

1We analysed 306 individuals (71% male) aged from 20 to 75 years (median age, 61 years), who were recruited during their annual health check-ups (Extended Data Fig. 1a). Individuals diagnosed with diabetes were excluded to avoid any long-lasting effects of hyperglycaemia,. Consequently, our study included relatively healthy individuals compared with most of the previous metagenomic studies of diabetes and obesity–,,; the median (interquartile range (IQR)) body mass index (BMI) and glycated haemoglobin (HbA1c) were 24.9 kg m−2 (22.2–27.1 kg m−2) and 5.8% (5.5–6.1%), respectively (Supplementary Table 1). The main clinical phenotype analysed in this study was insulin resistance (IR), which we defined as a homeostatic model assessment of IR (HOMA-IR) score of at least 2.5 (ref. ). We also analysed the associations between faecal metabolites and metabolic syndrome (MetS), an IR-related pathology. The clinical characteristics of IR and MetS largely overlapped except for blood pressure and sex ratio, for which there was no difference between individuals with IR versus normal insulin sensitivity (IS) (Supplementary Table 1). Untargeted metabolomics analysis using two mass spectrometry (MS)-based analytical platforms identified 195 and 100 annotated faecal and plasma hydrophilic metabolites, and 2,654 and 635 annotated faecal and plasma lipid metabolites, respectively (Extended Data Fig. 1a). To identify the overall difference in microbial functions, faecal metabolites and predicted genes were summarized into co-abundance groups (CAGs) and KEGG categories, respectively (Extended Data Fig. 1b). Transcriptomic information of peripheral blood mononuclear cells (PBMCs) was obtained using the cap analysis of gene expression (CAGE) method, which can measure gene expression at the transcription-start-site resolution.

2To examine how omics data of faecal samples can predict IR, we first compared the area under the curve (AUC) of receiver operating characteristic (ROC) curves on the basis of random-forest classifiers. Predictor variables for the models were selected using the minimum-redundancy maximum-relevance algorithm from the faecal 16S, metabolome, metagenome and their merged datasets (Supplementary Table 2). We found that the selected features of faecal metabolomic data generally outperformed those of 16S and metagenomics in predicting IR (Fig. 1a), suggesting that faecal metabolomics could be used to study IR pathogenesis.

3We next searched for the associations between clinical phenotypes and faecal metabolite CAGs (Fig. 1b and Supplementary Tables 3–8). Major confounding factors, namely sex and age, were adjusted throughout the correlation and regression analyses with clinical markers. Among the hydrophilic metabolites, most of the CAGs showing significant associations with IR were those of carbohydrate metabolites, mainly monosaccharides (hydrophilic CAGs 5, 12 and 15; Fig. 1b, top). Short-chain fatty acids (SCFAs), which are known as carbohydrate fermentation products, were also increased in IR (hydrophilic CAG 8). Hydrophilic CAG 18 remained unannotated as it included metabolites from different pathways (Supplementary Table 5). KEGG pathway enrichment analysis of the metabolites in these IR-related hydrophilic CAGs revealed that these metabolites were indeed involved in carbohydrate metabolism (Extended Data Fig. 2a). Specifically, we found that the major monosaccharides such as fructose, galactose, mannose and xylose significantly correlated with IR (Fig. 1c). Among the SCFAs, propionate was particularly increased in IR (Extended Data Fig. 2b), consistent with its role in gluconeogenesis. Faecal monosaccharides were similarly increased in MetS, obesity and prediabetes (Fig. 1d and Extended Data Fig. 2c,d). By contrast, disaccharides showed weak or no association (Extended Data Fig. 2b–d). These findings show that the end products of carbohydrate degradation—such as monosaccharides, which are readily absorbed and used by the host—are particularly increased in the faeces of individuals with IR and MetS. Supporting these findings, our analysis of previously published faecal metabolomics data from the TwinsUK cohort showed that faecal monosaccharides, notably glucose and arabinose, were positively associated with obesity and HOMA-IR, both of which relate to IR (Extended Data Fig. 3a–c and Supplementary Table 9). Similarly, the peak intensity of faecal fructose, glucose and galactose was associated with BMI in a small number of individuals without inflammatory bowel disease (IBD) from HMP2 data (Extended Data Fig. 3d). Together, these findings indicate that faecal carbohydrates are increased in IR and related pathologies and that this alteration is consistently observed across populations.

4In addition to hydrophilic metabolites, faecal lipid CAGs were also associated with IR (Fig. 1b). Lysophospholipids, bile acids and acylcarnitine were associated with IR and MetS as reported previously. Among them, a lipid CAG largely consisting of digalactosyl/glucosyldiacylglycerol (DGDG) (lipid CAG 11) came to our attention as DGDG is reportedly derived from bacteria,. These lipids contain glucose and/or galactose in their structures, although their biological functions in mammals are largely unclear. Most of the DGDGs in this cluster showed positive correlations with some of the precursor diacylglycerols and monosaccharides (that is, glucose and galactose) (Extended Data Fig. 4a). As diacylglycerols are deeply involved in IR pathogenesis, the biological functions of this metabolite class are of particular interest. Notably, DGDGs with different acyl chains in lipid CAG 41 showed no association with IR (Supplementary Table 7), implying that the differences in acyl chains of lipids may have a physiological importance as reported previously.

5We next investigated the alteration in gut microbiota and the functions of gut microbiota that are associated with IR. Gut microbiota diversity varied among individuals (Extended Data Fig. 5a–e). We then profiled the genus-level microbial composition of the study participants using 16S rRNA sequencing data and identified four bacterial groups (Extended Data Fig. 5f). Group 1 was dominated by the Lachnospiraceae family such as Blautia and Dorea, whereas group 2 was characterized by Bacteroidales (such as Bacteroides, Parabacteroides and Alistipes) and Faecalibacterium. Group 3 contained Actinobacteria genera. Group 4 did not form a distinct network. We could further classify the study participants into four clusters, A to D, on the basis of their taxonomic profiles (Fig. 2a). Individuals in cluster C distinctly harboured group 2 with Bacteroidales, whereas those in cluster D showed a higher abundance of group 1 and 3 bacteria (Extended Data Fig. 5g). Notably, the proportion of IR (Fig. 2a; P = 0.0071) was significantly lower in cluster C. Other metabolic parameters associated with IR and MetS such as HOMA-IR, BMI, triglycerides, HDL-cholesterol (HDL-C) and adiponectin were also different between cluster C (with the lowest proportion of IR) and the other three clusters (Fig. 2b and Supplementary Table 10). The proportion of IR among individuals with abundant group 1 and 3 bacteria was consistently higher than those with abundant group 2 bacteria, as identified on the basis of shotgun metagenomics data (Extended Data Fig. 5h). HOMA-IR showed negative associations with the genus Alistipes in the Rikenellaceae family and several species from Bacteroides, Bifidobacterium and Ruminococcus (Extended Data Fig. 5i and Supplementary Tables 11 and 12), partly recapitulating previous reports regarding individuals with obesity–. Notably, different genera and species correlated with other clinical markers, suggesting that the individual association between microbial taxa and clinical manifestation is not as robust as in the co-abundance analysis.

6We next constructed a microorganism–metabolite network on the basis of the significant positive or negative correlations (Supplementary Table 13). Although faecal SCFAs and lipids such as DGDG correlated with both IR- and IS-associated bacterial groups, IR-associated faecal carbohydrates predominantly correlated with genera in groups 1 and 4, the most prominent being Dorea in Lachnospiraceae (Fig. 2c,d). By contrast, the majority of these carbohydrates negatively correlated with IS-associated genera in group 2 bacteria such as Bacteroides, Alistipes and Flavonifractor (Fig. 2d and Extended Data Fig. 5j), with minimal correlations with bacteria in group 1. Accordingly, the faecal carbohydrate levels were distinctly different among the participant clusters (Extended Data Fig. 5k). Previous studies have suggested that several Lachnospiraceae species are involved in polysaccharide fermentation,, while Alistipes is increased on an animal-based diet rather than a polysaccharide-rich diet. These findings highlight a tight connection between carbohydrate-degradation products and IR- and IS-associated bacteria, suggesting that these bacteria may be involved in the aberrant faecal carbohydrate profile in IR.

7The IR-associated faecal carbohydrates were also correlated with KEGG pathways relating to carbohydrate metabolism and transportation, such as the phosphotransferase system (PTS), starch and sucrose metabolism, and galactose metabolism, while negatively associated with pathways relating to carbohydrate catabolism, such as glycolysis and pyruvate metabolism (Fig. 2e and Supplementary Tables 14 and 15). These pathways were also distinctly correlated with the participant clusters defined in Fig. 2a and the genera relating to carbohydrates defined in Fig. 2d. Amino acid metabolism was also different, particularly between clusters C and D, whereas lipid metabolism did not show distinct associations with microbiota (Extended Data Fig. 6a,b and Supplementary Table 16). Although carbohydrate pathways such as PTS and starch and sucrose metabolism showed strong positive associations with HbA1c and γ-GTP, the associations with other IR markers were generally sparse (Extended Data Fig. 6c and Supplementary Table 17), suggesting that metabolites are more sensitive to the clinical manifestations as shown in Fig. 1a. PTS is an essential component for bacteria to incorporate sugars into themselves as energy sources. Detailed analyses of KEGG orthologues revealed that faecal carbohydrates and participant clusters mainly correlated with PTSs relating to disaccharides and amino sugars (Extended Data Fig. 6d,e and Supplementary Table 18), suggesting that the preference of sugar use by microbiota through PTS may affect the metabolite levels. Glycosidases, which catalyse the breakdown of oligo- and disaccharides, were also associated with faecal monosaccharides (Extended Data Fig. 6f). Extracellular glucosidases such as β-fructofuranosidase (K01193, KEGG Orthology database), amylosucrase (K05341, KEGG Orthology database) and oligo-1,6-glucosidase (K01182, KEGG Orthology database), which were predicted to degrade sucrose and dextrin into glucose and fructose (Extended Data Fig. 6g,h), showed the highest positive correlations, especially with faecal glucose. By contrast, glucosidases relating to starch use such as α-amylases (K01176 and K07405, KEGG Orthology database) were negatively linked with faecal carbohydrates. Importantly, the abundance of these glycosidase genes was significantly different between participant cluster C and the other three clusters, suggesting that taxonomic profiles largely explain the variations of glucosidases (Fig. 2f, Extended Data Fig. 6h and Supplementary Table 18). Consistently, disaccharide-breakdown genes were predominantly conserved in the genomes of Blautia and Dorea abundant in cluster D, whereas they were almost lacking in Bacteroidales abundant in cluster C (Extended Data Fig. 6i). Together, our findings reveal four distinct populations with unique taxonomic profiles and carbohydrate metabolisms characterized by sugar use and degradation, which correlate with IR and its related markers.

8Consistent with previous reports,, the host cytokine, metabolomic and transcriptomic signatures were highly associated with IR (Supplementary Tables 19–21). Moreover, many of these PBMC genes were functionally involved in inflammation (Extended Data Fig. 7a) and possibly derived from monocytes (Supplementary Table 21). Several studies have suggested that microbial components such as lipopolysaccharides have a role in facilitating inflammation of metabolic diseases,. However, it remains unclear whether microbial metabolism is involved in low-grade inflammation. We therefore tried to infer possible associations between host inflammatory signatures of IR and faecal carbohydrates. First, the cross-omics correlation-based network with individual metabolites, bacteria, transcripts and cytokines associated with IR revealed that faecal carbohydrates were strongly tied with both bacteria and host IR-related signatures, especially cytokines, suggesting that these metabolites are the hubs of the host–microorganism network in IR (Fig. 3a, Extended Data Fig. 7b,c and Supplementary Table 22). Differential abundance, calculated as the ratio of their abundance in IR and IS, was most pronounced in the associations between faecal carbohydrates and cytokines. Notably, IL-10, a plasma cytokine, showed the most prominent associations with faecal carbohydrates and modestly with PBMC-derived transcripts, supporting recent studies showing its paradoxical effect to facilitate IR–. Faecal carbohydrates moderately explained the variance of IL-10 and, to a lesser extent, adiponectin, leptin and serpin E1, suggesting that faecal carbohydrates are particularly associated with these cytokines (Fig. 3b). Although the proportions of variance explained by faecal carbohydrates were lower than by plasma metabolites, they were much higher than those by genus-level abundance, highlighting the role of faecal metabolites linking gut microbiota and host inflammatory responses. We next sought to infer whether these cytokines mediated the effects of faecal carbohydrates on host metabolism using causal mediation analyses. We found that IL-10, serpin E1, adiponectin and leptin mediated most in silico causal relationships between faecal carbohydrates and host IR markers such as HOMA-IR (Fig. 3c, Extended Data Fig. 7d and Supplementary Table 23). Notably, there were unique correspondences between metabolites and cytokines; for example, IL-10 mediated the effects of fructose, mannose, xylose and rhamnose, but not other metabolites. Although the biological importance of these unique correspondences remains to be investigated, the combined analyses of faecal microbiota, metabolome and host inflammatory phenotypes in IR suggest a previously unrecognized interaction, whereby excessive monosaccharides may affect host cytokine expression.

9The above findings from human multi-omics analyses revealed an association between carbohydrate metabolites and IR pathology. To address the causal relationship between gut microbiota, faecal carbohydrates and metabolic diseases, we first analysed metabolites in the bacterial culture of 22 human faecal IS- and IR-associated bacteria. These bacteria were selected on the basis of the findings from the genus-level co-occurrence (Fig. 2a,b) and the species-level (Extended Data Fig. 5i) profiles. Principal component analysis plots of 198 metabolites indicated that Bacteroidales, a representative IS-associated bacterial order, showed a distinct metabolic profile along PC1 (Extended Data Fig. 8a,b and Supplementary Table 24). The top 10 metabolites contributing to the group separation included several amino acids and fermentation products such as succinate and fumarate, and the majority of these metabolites were preferentially produced by Bacteroidales (Extended Data Fig. 8b,c). We detected 13 out of 15 carbohydrates associated with IR (Fig. 1b) in the bacterial culture (Extended Data Fig. 8b). Most of these carbohydrates were plotted negatively along PC1, suggesting that these metabolites were negatively associated with Bacteroidales. Glucose, mannose and glucosamine were preferentially consumed by Bacteroidales compared with the other orders, whereas lactulose was mainly produced by Eubacteriales (Extended Data Fig. 8d). Alistipes indistinctus was the most potent in consuming a wide variety of carbohydrates (Extended Data Fig. 8e,f). These findings show that Bacteroidales species are potent consumers of several carbohydrates, driving the production of their fermentation products.

10We next tested the potential therapeutic effects of seven candidate bacteria shown to be associated with IS in human cohort findings. Postprandial blood glucose levels were particularly reduced in mice administered with A. indistinctus, Alistipes finegoldii and Bacteroides thetaiotaomicron that were fed a high-fat diet (Fig. 4a). Insulin tolerance tests also revealed that these strains ameliorated IR, most prominently by A. indistinctus administration (Fig. 4b,c). A. indistinctus administration ameliorated body mass gain, ectopic triglyceride accumulation in the liver and glucose intolerance (Extended Data Fig. 9a–d). Serum levels of HDL-C, adiponectin and, to a lesser extent, triglycerides, were also improved in mice that were treated with A. indistinctus (Extended Data Fig. 9e–g). The findings of the hyperinsulinaemic–euglycaemic clamp analysis indicated that A. indistinctus administration significantly improved IR and, particularly, whole-body glucose disposal (Extended Data Fig. 9h–j). Phosphorylation of AKT in the liver and epididymal fat was increased in mice treated with A. indistinctus and A. finegoldii mice (Extended Data Fig. 9k,l), suggesting that insulin signalling was improved in the liver and adipose tissue. These findings reveal a potency of A. indistinctus administration in ameliorating diet-induced obesity and IR.

11Mechanistically, metabolic measurement revealed that carbohydrate oxidation was significantly reduced in mice that were treated with A. indistinctus, implying that carbohydrate use is limited (Extended Data Fig. 9m and Supplementary Table 25). As dietary intake and locomotor activity remained unchanged (Extended Data Fig. 9n,o), we reasoned that host-accessible carbohydrates in the intestine were reduced by treatment with A. indistinctus. In this regard, A. indistinctus administration substantially altered caecal metabolites, characterized by a reduction in several carbohydrates including fructose, a lipogenic monosaccharide (Extended Data Fig. 10a–c and Supplementary Table 26). Fructose was similarly reduced in the serum (Extended Data Fig. 10d). Importantly, the AUC of insulin tolerance test was positively correlated with the caecal monosaccharides fructose, glucose and mannose (Fig. 4d). Collectively, these findings reveal that A. indistinctus ameliorates IR and affects intestinal carbohydrate metabolites in mice, supporting our observations in the human cohort.

12To deepen our understanding of the host–microorganism relationship in IR, we used multimodal techniques to conduct a comprehensive and extensive study investigating the interactions between the gut microbiome and metabolic diseases in humans. Although carbohydrate metabolism by the gut microorganisms has been suggested to influence the pathogenesis of obesity,, and prediabetes,, the actual mechanistic linkage has been elusive in humans owing to the lack of detailed metabolomic information. In this regard, the major strength of our approach is that we combine faecal metabolomics cataloguing more than 2,800 annotated metabolites with both microbiome and host pathology. This metabolome-based approach enabled us to identify the faecal metabolites related to IR, identify an association between faecal carbohydrates and low-grade inflammation of IR, and efficiently select candidate strains for functional validations in experimental settings (Extended Data Fig. 10e). Together, our study highlights the advantage of comprehensive omics strategy in exploring the involvement of microbial metabolism and their products in the pathogenesis of IR. Excessive monosaccharides have the potential to promote ectopic lipid accumulation while also activating immune cells, leading to low-grade inflammation and IR–. Fructose is a widely recognized risk factor for inflammation and IR due to its role in lipid accumulation, whereas galactose has been shown to participate in the energy metabolism of activated immune cells. Our in vivo studies confirm that A. indistinctus administration improves lipid accumulation and thereby IR, while simultaneously reducing intestinal monosaccharide levels (Fig. 4d). Nevertheless, we are aware that further mechanistic studies are needed to examine the kinetics of absorption and their effects on host metabolism. In particular, how Alistipes strains suppress carbohydrate metabolism is an intriguing question (for example, whether these bacteria per se inhibit carbohydrate metabolism, or whether they interact with other commensals), as it would directly open the possibility of a new therapeutic strategy. Given that A. indistinctus improved whole-body IS (Extended Data Fig. 9i), it would be important to investigate the involvement of insulin signalling not only in the liver but also in peripheral tissues, including skeletal muscle and adipose tissue, along with the accumulation of specific lipid molecules (such as ceramides and diacylglycerols) in these tissues. Such investigations hold the potential to shed light on the underlying mechanisms that contribute to A. indistinctus-mediated improvement of IR. Finally, two participants in the human study were unable to collect their faeces in the morning, which could potentially influence the outcomes due to the lack of stringent control over time-of-day and fasting conditions. We therefore believe that longitudinal studies incorporating a timely documentation of dietary habits are warranted to dissect the intricate impacts of microbial metabolism on the trajectory of diabetes and its complications while accounting for potential confounding factors.

13The study participants were recruited from 2014 to 2016 during their annual health check-ups at the University of Tokyo Hospital. The individuals included both male and female Japanese individuals aged from 20 to 75 years. The exclusion criteria were as follows: established diagnosis of diabetes, routine use of medications for diabetes and/or intestinal diseases, use of antibiotics within 2 weeks before sample collection and loss of 3 kg of body weight in the 3 months before sample collection. Written consent was obtained from the participants after a thorough explanation of the nature of the study at their health-checkups.

14To normalize the participants’ clinical characteristics, we planned to recruit around 100 healthy individuals, 100 individuals with obesity (BMI ≥25, based on the Japanese definition) and 100 individuals with a prediabetic condition (FBG ≥110 mg dl−1 and/or HbA1c ≥6.0%) on the basis of their clinical data, and stopped recruiting when the number of participants almost reached the goal. The sample size was determined on the basis of previous metagenomics studies showing microbial signatures of diabetic patients,. We enrolled 112, 100 and 101 individuals for the normal, obese and prediabetic groups, respectively. The participants were provided with instructions to fast overnight before their visits, and all clinical information and blood samples were collected in the morning during their hospital visit. Blood samples were immediately centrifuged to collect plasma and then stored at −80 °C until the sample preparation and analysis. The participants were also instructed to collect faecal samples in the morning and were provided with guidance on how to collect and preserve faecal samples, along with a kit comprising a sampling tube and an ice pack. The faecal samples were then transported to the hospital either by refrigerated shipping or by the participants themselves. In both scenarios, the samples were delivered in a chilled state within 24 h after collection and stored at −80 °C until sample preparation and analysis. Consequently, 256 participants collected their faeces in the morning on the day of their hospital visit. As for the remaining participants, they collected their faeces in the morning between 2 days before and 7 days after their hospital visit, with the exception of 5 individuals who collected their faeces in the morning more than 7 days after their hospital visit, 2 individuals who reported collecting their faeces in the evening 1 day before their hospital visit, and 5 individuals who did not provide faecal samples. Moreover, two individuals withdrew from the study after enrolment. Thus, 306 individuals who underwent physical examination, laboratory tests, faecal sampling for faecal 16S rRNA pyrosequencing and metabolomic analyses, and plasma sampling for plasma metabolomic analyses were included for the analysis. Owing to the limited samples, faecal metagenomics data were available for 290 individuals; CAGE analysis data for 298 individuals; and plasma cytokine and insulin data for 282 individuals. The number of samples included in each analysis is described in the figure legends. The clinical study was approved by the institutional review board of RIKEN and The University of Tokyo and performed in accordance with the institutes’ guidelines.

15Although we determined the criteria for enrolment, these criteria were not necessarily appropriate for subsequent analyses. For example, those in the prediabetes group were significantly leaner than those in the obese group (27.3 kg m−2 versus 25.2 kg m−2, P < 0.0001). Moreover, owing to the nature of the study participants (that is, those participated in regular health checkups), the blood glucose and HbA1c of the prediabetes group were significantly but only marginally higher than those of the obese group (FBG, 106 mg dl−1 versus 94 mg dl−1, P < 0.0001; and HbA1c, 6.2% versus 5.6%, P < 0.0001). We therefore reasoned that, in these subclinical conditions of diabetes, many metabolic traits may be overlapping between prediabetes and obesity groups and they do not necessarily capture their distinct features in metabolic and clinical continuums. This hinders us from distinguishing microbial and metabolomic characteristics directly related to human metabolic dysfunctions. We therefore considered that individual indices representing participants’ clinical conditions (that is, IR and MetS, as described below) may offer a better interpretation of the participants’ metabolic traits and data. Nevertheless, we observed consistent results even with the clinical criteria of obesity and prediabetes (Extended Data Fig. 2d).

16IR is defined as HOMA-IR ≥2.5, as has been set for the Japanese population. Similarly, normal IS was defined as HOMA-IR ≤1.6. HOMA-IR is calculated using the following formula: fasting plasma insulin (μU ml−1) × fasting plasma glucose (mg dl−1)/405. HOMA-IR values could be calculated for 282 individuals only, owing to the limited data of plasma insulin in some participants. MetS is diagnosed according to the Japanese criteria, which require an abdominal circumference of ≥85 cm for male and ≥90 cm for female individuals and at least two out of the following three clinical abnormalities: (1) dyslipidaemia, defined as triglyceride levels of ≥150 mg dl−1 and/or HDL-C levels of <40 mg dl−1; (2) elevated blood pressure, defined as systolic blood pressure of ≥130 mmHg and/or diastolic blood pressure of ≥85 mmHg; and (3) impaired fasting glucose, defined as FBG levels of ≥110 mg dl−1. Individuals who meet the criteria of abdominal circumference but only one clinical abnormality were defined as pre-MetS, as reported previously.

17Plasma cytokines were measured using Human Adipokine Magnetic Bead Panel 2 (Millipore, HADK2MAG-61K) and Human Obesity Premixed Magnetic Luminex Performance Assay Kit (R&D, FCSTM08) according to the manufacturers’ instructions. Measurements below the lower detection limits were considered to be zero, and those above the upper detection limits were considered to be the highest values of analysed cytokines.

18Aliquots (5 g) of faeces were blended with 30 ml methanol and filtrated with 100 µm of mesh filter to remove food residue after vigorous vortexing. The filtrate was centrifuged at 15,000g for 10 min at 4 °C and the supernatant (methanol extract) was used for metabolomics analysis. DNA of the faecal microbiome was extracted from the pellet.

19We followed the extraction process and gas chromatography-tandem MS (GC–MS/MS) measurement methods for water-soluble metabolites described previously with some modifications. In brief, a 10 μl aliquot of plasma was added to 150 μl methanol, 125 μl Milli-Q water, 15 μl internal standard solution (1 mM 2-isopropylmalic acid) and 60 μl CHCl3. For faecal samples, a 25 μl aliquot of methanol extract was added to 125 μl methanol, 150 μl Milli-Q water containing internal standard (100 μM 2-isopropylmalic acid) and 60 μl CHCl3. The solution was shaken at 1,200 rpm for 30 min at 37 °C. After centrifugation at 16,000g for 5 min at room temperature, 250 μl of the supernatant was transferred to a new tube and 200 μl of Milli-Q water was added. After mixing, the solution was centrifuged at 16,000g for 5 min at room temperature, and 250 μl of the supernatant was transferred to a new tube. The samples were evaporated dry using a vacuum evaporator for 20 min at 40 °C and lyophilized using a freeze dryer. Dried extracts were derivatized with 40 μl of 20 mg ml−1 methoxyamine hydrochloride (Sigma-Aldrich) dissolved in pyridine and shaken at 1,200 rpm for 90 min at 30 °C. The solution was then mixed with 20 μl of N-methyl-N-trimethylsilyl-trifluoroacetamide (MSTFA, GL Science) and incubated for 30 min at 37 °C with shaking at 1,200 rpm. After derivatization, the samples were centrifuged at 16,000g for 5 min at room temperature, and the supernatant was transferred to a glass vial. The analysis was performed using a GC–MS/MS platform on the Shimadzu GCMS-TQ8030 triple quadrupole mass spectrometer (Shimadzu) with a capillary column (BPX5, SGE Analytical Science). The GC oven was programmed as follows: 60 °C held for 2 min, increased to 330 °C (15 °C min−1), and finally 330 °C held for 3.45 min. GC was operated in constant linear velocity mode set to 39 cm s−1. The detector and injector temperatures were 200 °C and 250 °C, respectively. Injection volume was set at 1 μl with a split ratio of 1:30.

20We followed the SCFA extraction and GC–MS/MS measurement methods as previously described with some modifications. A 90 μl aliquot of plasma was added to 10 μl Milli-Q water containing internal standards (2 mM [1,2-13C2]acetate, 2 mM [2H7]butyrate and 2 mM crotonate). For faecal samples, a 25 μl aliquot of methanol extract was added to 10 μl Milli-Q water containing internal standards and then centrifugally concentrated at 40 °C and reconstituted with 100 μl of Milli-Q water. Then, 50 μl of hydrochloric acid (HCl) and 200 μl of diethyl ether were added to the solution and mixed well. After centrifugation at 3,000g for 10 min, 80 μl of the organic layer was transferred to a glass vial and 16 μl N-tert-butyldimethylsilyl-N-trifluoroacetamide (MTBSTFA, Sigma-Aldrich) was added to derivatize the samples. The vials were incubated at 80 °C for 20 min and allowed to stand for 48 h before injection. The analysis was performed using a Shimadzu GCMS-TQ8030 triple quadrupole mass spectrometer with a capillary column (BPX5). The GC oven was programmed as follows: 60 °C held for 3 min, increased to 130 °C (8 °C min−1), increased to 330 °C (30 °C min−1) and finally 330 °C held for 3 min. The detector and injector temperatures were 230 °C and 250 °C, respectively. GC was operated in constant linear velocity mode set to 40 cm s−1. Injection volume was set at 1 μl with a split ratio of 1:30. The data were processed and concentration was calculated by LabSolutions Insight (Shimadzu).

21Overall, 195 and 100 metabolites in the faecal and plasma samples, respectively, were detected by our GC–MS/MS platform and passed quality control. The values below the limit of detection were replaced with zero. Consequently, 110 faecal and 88 plasma metabolites that were detected (that is, above zero) in more than 75% of participants were included in subsequent analyses, for which they were combined into a common analysis pipeline and defined as hydrophilic metabolites.

22The lipidomics analysis was performed according to a previously reported study. Methanol, isopropanol, chloroform and acetonitrile of liquid chromatography (LC)–MS grade were purchased from Wako. Ammonium acetate and EDTA were purchased from Wako and Dojindo, respectively. Milli-Q water was purchased from Millipore (Merck). EquiSPLASH was purchased from Avanti Polar Lipids. Palmitic acid-d3 and stearic acid-d3 were purchased from Olbracht Serdary Research Laboratories.

23For plasma lipid extraction, an aliquot of 20 μl of human plasma sample was added to 200 μl of methanol containing 5 μl of EquiSPLASH, 10 μM palmitic acid-d3 and 10 μM stearic acid-d3, and vortexed for 10 s. Then, 100 μl of chloroform was added and vortexed for 10 s. After incubation for 2 h at room temperature, the solvent tube was centrifuged at 2,000g for 10 min at 20 °C. A total of 200 μl of supernatant was transferred to an LC–MS vial (Agilent Technologies). For faecal lipid extraction, 50 μl of the methanol extract was added to 145 μl of methanol containing 5 μl of EquiSPLASH, 10 μM palmitic acid-d3 and 10 μM stearic acid-d3 in a 2 ml glass tube, and vortexed for 10 s. Then, 100 μl of chloroform was added and vortexed for 10 s. After incubation for 1 h at room temperature, 20 μl of water was added and vortexed for 10 s. After 10 min incubation at room temperature, the solvent was centrifuged at 2,000g for 10 min at 4 °C, and the supernatant was transferred to the LC–MS vial. All of the samples were divided into four batches for plasma analyses and five batches for faecal analyses, with 70–80 and 55–60 samples per batch after randomization, respectively. For each batch, a series of samples was prepared, and subsequent LC–MS/MS measurements were performed. A quality control sample was prepared by mixing the same volume of plasma from the first batch subjects. A procedure blank was prepared by using the same volume of water instead of a biological sample. The blank sample was analysed at the beginning and the end of each analysis batch, and the quality-control sample was injected every ten study samples.

24The LC system consisted of a Waters Acquity UPLC system. Lipids were separated on an Acquity UPLC Peptide BEH C18 column (50 × 2.1 mm; 1.7 µm) (Waters). The column was maintained at 45 °C at a flow rate of 0.3 ml min−1. The mobile phases consisted of (A) 1:1:3 (v/v/v) acetonitrile:methanol:water with ammonium acetate (5 mM) and 10 nM EDTA; and (B) 100% isopropanol with ammonium acetate (5 mM) and 10 nM EDTA. A sample volume of 0.5−3 µl, depending biological samples, was used for the injection. The separation was conducted under the following gradient: 0 min, 0% B; 1 min, 0% B; 5 min, 40% B; 7.5 min, 64% B; 12 min, 64% B; 12.5 min, 82.5% B; 19 min, 85% B; 20 min, 95% B; 20.1 min, 0% B; and 25 min, 0% B. The sample temperature was maintained at 4 °C.

25MS detection of lipids was performed on a quadrupole/time-of-flight mass spectrometer TripleTOF 6600 (SCIEX). All analyses were performed in high-resolution mode in MS1 (~35,000 full width at half-maximum) and the high sensitivity mode (~20,000 full width at half-maximum) in MS2. Data-dependent MS/MS acquisition (DDA) was used. The parameters were MS1 and MS2 mass ranges, m/z 70–1,250; MS1 accumulation time, 250 ms; MS2 accumulation time, 100 ms; collision energy, +40/−42 eV; collision energy spread, 15 eV; cycle time, 1,300 ms; curtain gas, 30; ion source gas 1, 40(+)/50(−); ion source gas 2, 80(+)/50(−); temperature, 250 °C(+)/300 °C(−); ion spray voltage floating, +5.5/−4.5 kV; declustering potential, 80 V. The other DDA parameters were dependent product ion scan number, 16; intensity threshold, 100 cps; exclusion time of precursor ion, 0 s; mass tolerance, 20 ppm; ignore peaks, within m/z 200; and dynamic background subtraction, true. The mass calibration was automatically performed using an APCI positive/negative calibration solution through a calibration delivery system.

26MS-DIAL (v.4.48), was used with the following parameters: (data collection) retention time begin, 1.0 min; retention time end, 18 min; MS1 and MS2 mass range begin, 0 Da; MS1 and MS2 mass range end, 2,000 Da; MS1 tolerance, 0.01 Da; MS2 tolerance, 0.025 Da; (peak detection) minimum peak height, 3,000 amplitude; mass slice width, 0.1 Da; smoothing method, linear weighted moving average; smoothing level, 3 scan; minimum peak width, 5 scan; exclusion mass list, none; (identification) retention time tolerance, 1.5 min; MS1 accurate mass tolerance, 0.01 Da; MS2 accurate mass tolerance, 0.05 Da; identification score cut off, 70%; all lipid subclasses were used as the search space; (alignment) retention time tolerance 0.15 min; MS1 tolerance, 0.015 Da. The default values were used for other parameters. In faecal lipidomics, a total of 48,790 and 20,367 chromatographic peaks were detected in positive- and negative-ion mode data, respectively. Of these, 2,654 unique lipid molecules were annotated and semi-quantified in the MS-DIAL software program and used for further statistical analyses. Likewise, in plasma lipidomics, 1,469 and 2,167 chromatographic peaks were detected in positive- and negative-ion mode data, respectively, and 635 unique lipid molecules were annotated and semi-quantified. The semi-quantitative value of lipids was calculated by the internal standards according to the previous study. The abbreviations of lipids are listed in Supplementary Table 27. Details of lipid subclass characterization follow the previous study.

27To generate co-abundance clusters, 110 hydrophilic metabolites and 2,654 lipid metabolites detected in more than 75% of participants were included. These metabolites were clustered based on their co-abundance using the R package WGCNA (v.1.72-1). The following parameters were used for the analysis. For hydrophilic metabolites, soft thresholding β = 12, minimum cluster size = 3, deep split = 4, cut height = 0.9999, PAM clustering = F. For lipid metabolites, soft thresholding β = 14, minimum cluster size = 20, deep split = 4, cut height = 0.999, PAM clustering = F. As soft thresholding of WGCNA was not able to cluster all of the metabolites, the remaining metabolites that did not fit the criteria were subsequently clustered on the basis of biweight midcorrelation. The following parameters were used for the secondary clustering. For hydrophilic metabolites, minimum cluster size = 3, deep split = 4, cut height = 0.9999, PAM clustering = F. For lipid metabolites, minimum cluster size = 6, deep split = 4, cut height = 0.999, PAM clustering = F. The clusters with biweight midcorrelation above 0.8 were merged. The first principal component (PC1) of each cluster was calculated using the moduleEigengenes command of WGCNA and used as the representative value of the cluster for further analyses. The representative classes of the clusters were described in Supplementary Tables 2 and 3. KEGG pathway enrichment analysis of CAGs was performed on MetaboAnalyst (v.5.0) using 84 metabolite sets based on the KEGG pathway. Hypergeometric test and false-discovery rate (FDR)-adjusted P values were used to test significance. The enrichment ratio was calculated as the ratio of actual metabolite number to the expected value in each pathway.

28To validate the associations between clinical markers and faecal metabolites, we used the metabolomic data of TwinsUK and HMP2 (ref. ). The metabolome data of the TwinsUK cohort included 1,116 metabolites including 36 carbohydrates. The median (interquartile range) of age and BMI were 65 years (60–71 years) and 25.4 (22.8–28.8), and the proportion of males was 6.6%. As reported previously, the metabolite levels were scaled by run-day medians. The data were then log-transformed and scaled. For regression analyses, we filtered out the metabolites detected in less than 50% of participants; as a result, 759 metabolites including 29 carbohydrates were used for further analyses. The record of BMI and HOMA-IR were used for phenotypic outcomes. For BMI, we retrieved 786 samples measured on the same day of faecal collection. For HOMA-IR, plasma glucose and insulin obtained in the same year of the faecal collection were used for the following calculation: plasma glucose (mM) × insulin (pM)/6.945/22.5. We identified 550 individuals who underwent both faecal collection and glucose and insulin measurement in the same year and included them in the analysis. The HMP2 data were obtained from the Inflammatory Bowel Disease Multi’omics Database (https://ibdmdb.org/). Among the 26 out of 106 samples from non-IBD control, BMI data were available for 20 samples. We further excluded four individuals aged <10 years. As HMP2 is a longitudinal study, only the first faecal sampling for metabolomics was used for the current analysis to avoid redundancy. The intensity of fructose, glucose and/or galactose was log-transformed and scaled.

29DNA extraction was performed according to a protocol described previously with slight modifications. Before DNA extraction, the faecal pellet was washed once with PBS and suspended in a 10 mM Tris-HCl/20 mM EDTA buffer (pH 8.0). Lysozyme (Sigma-Aldrich), achromopeptidase (Wako) and proteinase K (Merck) were subsequently added to the samples for cell lysis. DNA was recovered by a phenol–chloroform extraction method. To purify the extracted DNA, RNA was digested with RNase (Nippon Gene). DNA was then precipitated in a solution containing polyethylene glycol 6000 (Hampton Research). The DNA concentration was quantified using Quant-iT PicoGreen (Thermo Fisher Scientific).

30The hypervariable V1–V2 region of the 16S rRNA gene was amplified by PCR using barcoded primers. PCR amplicons were purified using AMPure XP magnetic purification beads (Beckman Coulter), and quantified using the Quant-iT PicoGreen dsDNA Assay Kit (Life Technologies Japan). Equal amounts of each PCR amplicon were mixed and then sequenced using the MiSeq (Illumina) system.

31On the basis of sample-specific barcodes, reads were assigned to each sample using bcl2fastq. Next, the reads lacking both forward and reverse primer sequences were removed using BLAST and parasail followed by trimming of both primer sequences. Data were further denoised by removing reads with average quality values of <25 and possible chimeric sequences. Reads with BLAST match lengths of <90% with the representative sequence in the 16S databases (described below) were considered to be chimeras and were removed. The filter-passed reads were used for further analysis. The 16S database was constructed from three publicly available databases: the Ribosomal Database Project (RDP; v.10.27), CORE (http://microbiome.osu.edu/) and a reference genome sequence database obtained from the NCBI FTP site (ftp://ftp.ncbi.nih.gov/genbank/, December 2011).

32Operational taxonomic unit (OTU) clustering and UniFrac analysis from the filter-passed reads, 3,000 high-quality reads per sample were randomly chosen. All reads (the number of samples × 3,000) were then sorted according to their average quality value and grouped into OTUs using UCLUST (http://www.drive5.com/) with a sequence-identity threshold of 97%. The representative sequences of the generated OTUs were processed for homology search against the databases mentioned above using the GLSEARCH program for taxonomic assignments. For assignment at the phylum, genus and species levels, sequence similarity thresholds of 70%, 94% and 97% were applied, respectively.

33Metagenome shotgun libraries (insert size of 500 bp) were prepared using the TruSeq Nano DNA kit (Illumina) and sequenced on the Illumina NovaSeq platform. After quality filtering, reads mapped to the human genome (HG19) or the phiX bacteriophage genome were removed. For each individual, the filter-passed NovaSeq reads were assembled using MEGAHIT (v.1.2.4). Prodigal (v.2.6.3) was used to predict protein-coding genes (≥100 bp) in the contigs (≥500 bp) and singletons (≥300 bp). Finally, 6,458,217 non-redundant genes were identified in the 290 samples by clustering the predicted genes using CD-HIT with a 95% nucleotide identity and 90% length coverage cut-off. Functional assignment of the non-redundant genes was performed using DIAMOND (e-value ≤ 1.0 × 10−5) against the Kyoto Encyclopedia of Genes and Genomes (KEGG) database (release 2019-10-07) to obtain the KEGG orthologues. The genes with the best hit correlating to eukaryotic genes were excluded from further analysis.

34For taxonomic assignment of metagenomic reads, 1 million filter-passed reads were processed for mOTU analysis (v.2.5.1) to obtain the relative abundance at the species level. To functionally annotate the predicted genes, 1 million filter-passed metagenomic reads per individual were mapped to the combined reference gene set consisting of non-redundant genes identified in this study, JPGM and IGC using Bowtie2 with a 95% identity cut-off. Multi-mapped reads, that is, the reads that mapped to multiple genes with identical scores, were normalized to the proportion of the number of other reads that uniquely mapped to these genes, according to a strategy outlined in a previous report. The proportion of KEGG orthologues was calculated from the number of reads mapped to them. For the enrichment analysis of KEGG pathways, the significantly and positively (negatively) associated KEGG orthologue gave +1 (−1) for all of the upstream pathways linked to the KEGG orthologue, and the points were summarized as the ratio to the number of KEGG orthologues in the pathway. For the KEGG-orthologue-level analyses of PTS, those including ‘phosphotransferase system (PTS)’ in the KEGG pathway (02060) were selected for the following correlation analyses. In the analyses of glucosidases, ‘glycoside hydrolases’ defined in the CAZy database on the basis of EC numbers were selected. We further selected those included in ‘starch and sucrose metabolism’ in the KEGG pathway (00500). We defined intracellular glucosidase by their substrate described in the KEGG pathway map; those cleave phosphorylated carbohydrates were recognized as intracellular, and the rest of the genes were recognized to possess extracellular enzymatic activities. The pathways were further summarized into carbohydrate metabolism (09101), amino acid metabolism (09105), lipid metabolism (09103) and membrane transport (09131) on the basis of the KEGG Orthology database.

35The list of KEGG organisms used for this genome analysis is listed in Supplementary Table 28. All KEGG organisms from genera Alistipes, Bacteroides, Flavonifractor, Blautia, Dorea and Coprococcus, which showed the top three positive or negative correlations with faecal carbohydrates in Fig. 2d, were selected for this analysis. The lists of genes involving the ‘starch and sucrose metabolism’ pathway (00500) in these KEGG organisms were extracted using the R package KEGGREST (v.1.32.0). The representative protein sequences of Blautia hydrogenotrophica strain 2789STDY5608857 (taxonomy ID 53443), Dorea longicatena strain 2789STDY5608851 (taxonomy ID 88431) and Dorea formicigenerans strain ATCC 27755 (taxonomy ID 411461) were downloaded from the NCBI Datasets (https://www.ncbi.nlm.nih.gov/datasets/genomes/). KEGG annotation of these protein sequence files was performed using BlastKOALA (https://www.kegg.jp/blastkoala/) with ‘Bacteria’ used as the taxonomy group. The presence of KEGG orthologues relating to extracellular glycoside hydrolases in starch and sucrose metabolism pathways shown in grey in Extended Data Fig. 6f was summarized.

36Blood samples were collected in Vacutainer CPT tubes (Becton Dickinson) and mixed with the anticoagulant by gently inverting the tubes 8 to 10 times. After centrifugation of the blood for 30 min at 1,500g, PBMCs were isolated as a diffuse layer above the gel. The plasma was removed, and the PBMCs were collected in conical tubes with 500 μl RNAlater (Thermo Fisher Scientific). The conical tubes were centrifuged at 1,000g at room temperature for 3 min to pellet the cells and the supernatant was discarded. The RNA was then isolated using the Maxwell 16 LEV simplyRNA Blood Kit (Promega) according to the manufacturer’s instructions. The quality of the RNA was assessed using Bioanalyzer (Agilent), as recommended by the manufacturer. The RNAs were quantified using the GloMax plate reader (Promega) and Quant-iT RiboGreen RNA Assay Kit (Thermo Fisher Scientific).

37The CAGE libraries were constructed according to the dual-index nanoCAGE protocol, a template-switching-based variation of the standard CAGE protocol designed for low quantities of RNA,. cDNA libraries were prepared with RNA extracted from PBMC samples and sequenced on the Illumina HiSeq 2000 (50 bp paired-end). The sequenced reads were processed with the MOIRAI pipeline: low quality and rDNA reads were first removed, then the remaining reads were mapped to the human genome version hg38 patch 1 using BWA v.0.5.9 (r16). The mapped reads were overlapped with the FANTOM5 robust promoter set (http://fantom.gsc.riken.jp/5/datafiles/latest/extra/CAGE_peaks/) and mapped to the nearest GENCODE v.27 annotations within (500 bp),. The mapped reads falling under each FANTOM5 CAGE cluster were summed to produce the raw expression counts. Expression counts were converted to counts per million (CPM), and CAGE clusters expressed in less than 100 samples with at least 1 CPM and greater than 1 sample with at least 10 CPM were removed from further analysis. For each sample, the richness index was calculated using the R package vegan’s rarefy function with a subsample size of 100 on the filtered raw counts. Samples with a read library size of less than 1,000,000 and a number of unique CAGE clusters of <11,000 and richness less than 44 were removed as outliers, with the thresholds selected from visual inspection of the respective distributions. Cell type specificities of promoters of interest were determined using the FANTOM5 hg38 human promoterome view.11 in the ZENBU genome browser (https://fantom.gsc.riken.jp/zenbu/). Top-hit cells for analysed promoters were described. For cell-type gene set enrichment analysis of genes significantly associated with IR, annotated genes were analysed using Enrichr, and the Human Gene Atlas database, and the results of cell types with Padj < 0.05 were selected. The Enrichr combined score is defined as: c = log[p] z, where c is the combined score, p is the P value based on Fisher’s exact test and z is the z-score.

38The following strains were used for this culture analysis: A. indistinctus (JCM16068), A. finegoldii (JCM16770), Alistipes putredinis (JCM 16772), B. thetaiotaomicron (JCM 5827), Bacteroides xylanisolvens (JCM 15633), Bacteroides ovatus (JCM 5824), Bacteroides caccae (JCM 9498), Parabacteroides merdae (JCM 9497), Parabacteroides distasonis (JCM 5825), D. formicigenerans (JCM 31256), D. longicatena (JCM 11232), B. hydrogenotrophica (JCM 14656), Blautia producta (BP, JCM 1471), Coprococcus comes (JCM 31264), Faecalibacterium prausnitzii (JCM 31915), Flavonifractor plautii (JCM 32125), Clostridium spiroforme (JCM1432), Coriobacterium glomerans (JCM 10262), Roseburia hominis (JCM 17582), Adlercreutzia equolifaciens subsp. equolifaciens (JCM 14793), Eggerthella lenta (JCM 9979) and Collinsella aerofaciens (JCM 10188). All strains were obtained from RIKEN BioResource Research Center. All of the strains were cultivated in EG medium (JCM Medium No. 14) supplemented with 5% Fildes extract prepared by pepsin-digested horse blood instead of horse blood itself. For measurement of metabolites in bacterial culture, 60 μl of the bacterial culture grown in the EG medium was inoculated into 3 ml of the experiment medium (EG medium) and cultivated for 24 h. The samples were centrifuged, and the cell-free supernatant was collected for analysis. GC–MS was performed to measure hydrophilic metabolites as described above. We identified 261 metabolites by the analysis and used 198 metabolites observed in at least 30% of samples for the following analyses.

39C57BL6/N male mice aged 6 weeks were purchased from CLEA Japan. They were randomly assigned to either the control or treatment group and housed in a conventional animal facility of Yokohama City University. The mice were fed Quick Fat (CLEA Japan) for 3 weeks before bacterial administration and continued to be fed for 3 weeks during bacterial challenges. A. indistinctus (JCM16068), A. finegoldii (JCM16770), B. thetaiotaomicron (JCM 5827), B. xylanisolvens (JCM 15633), P. merdae (JCM 9497), F. prausnitzii (JCM 31915) and C. spiroforme (JCM1432) were used to broadly compare the efficacy of bacterial administration on the animal model. These strains were cultivated in EG medium overnight, and the concentration was adjusted to 2.5 × 108 CFU per ml by PBS. The bacteria and PBS, a negative control, were orally administered to the mice at a dose of 200 μl per mouse. The bacteria and PBS as the vehicle control were provided 3 times a week for 3 or 4 weeks. Body mass was measured before oral gavage. Postprandial blood glucose measurement and insulin tolerance test were performed 3 weeks after the initiation of bacterial challenges. After the insulin tolerance test, the mice were subjected to 5 h fasting before insulin injection, and 0.85 U kg−1 human regular insulin (Eli Lilly) was subsequently administered intraperitoneally. The intraperitoneal glucose tolerance test was performed 4 weeks after the initiation of bacterial challenges. The mice were subjected to 5 h fasting before glucose infusion, and 2.0 g per kg glucose (Nacalai Tesque) was administered intraperitoneally. In both experiments, the blood glucose was collected from the tail vein and serially measured using GLUCOCARD G Black (Arkray). For the necropsy, the mice were euthanized by isoflurane (MSD), and the fat mass of perigonadal and mesenteric fats was measured. Blood was drawn through cardiac puncture after the anaesthesia. HDL-C (Wako), triglycerides (Wako) and adiponectin (Otsuka) were measured in accordance with the manufacturers’ instructions. The Yokohama City University animal facility is maintained under a 12 h–12 h light–dark cycle at 24 ± 1.5 °C and 55 ± 10% humidity.

40To assess the metabolism, dietary intake and locomotor activity of mice, C57BL6/N male mice at the age of 6 weeks were purchased from CLEA Japan and were maintained in a vinyl isolator of SPF animal facility at RIKEN Yokohama branch. Using the same experimental protocol in the conventional condition, the mice were fed Quick Fat (CLEA Japan) for 3 weeks before bacterial administration and continued to be fed the diet during bacterial challenges and metabolic measurement. We gave three oral gavages of A. indistinctus or PBS (vehicle control) every other day and then placed the mice individually in acrylic cages. We further gave one oral gavage 2 days after the start of individual housing. Their metabolic activity, dietary intake and physical activity were subsequently monitored. There was no significant difference in body mass at the start of metabolic measurement (mean ± s.d. of body mass were 25.7 ± 2.6 g and 26.1 ± 1.4 g in the vehicle and A. indistinctus groups, respectively). Oxygen and carbon dioxide concentration was measured using the ARCO-2000 system, an open-circuit metabolic gas analysis system with a mass spectrometer (Arco Systems). VO2, VCO2, energy expenditure, fat oxidation rate, carbohydrate oxidation rate and respiratory quotient were calculated within the system. Dietary intake and physical activity were simultaneously monitored through ACTIMO-100M and MFD-100M (Shinfactory). The differences in diurnal variation were tested using two-way mixed ANOVA, and P values for interactions between time and group were reported. The RIKEN animal facility is maintained under a 12 h–12 h light–dark cycle at 23 ± 2 °C and 50 ± 10% humidity. The sample size was determined on the basis of our preliminary experiments. Bacterial administration and body mass measurements were performed by an independent researcher who was not involved in the grouping and outcome assessments. All experimental procedures were approved by the Institutional Animal Care and Use Committee of the RIKEN and Yokohama City University and performed in accordance with the institutes’ guidelines.

41To analyse phosphorylation of AKT (p-AKT) at Ser473, the mice administered with A. indistinctus, A. finegoldii and PBS (vehicle control) 3 times a week for 4 weeks were subjected to 6 h fasting before the insulin injection, and 0.85 U kg−1 human regular insulin (Eli Lilly) was subsequently administered from the inferior vena cava. The liver, epididymal fat (eWAT) and gastrocnemius muscle were subsequently collected 5 min after the insulin injection, weighed and snap-frozen by liquid nitrogen. To prepare the lysates for western blotting, the tissues were homogenized in buffer A (25 mM Tris-HCl, pH 7.4, 10 mM sodium orthovanadate, 10 mM sodium pyrophosphate, 100 mM sodium fluoride, 10 mM EDTA, 10 mM EGTA and 1 mM phenylmethylsulfonyl fluoride). Thereafter, the lysates were resolved on 10% SDS–PAGE. Phosphorylated or total protein of AKT was isolated by immunoblotting using specific antibodies after the tissue lysates were resolved by SDS–PAGE and transferred to a Hybond-P PVDF transfer membrane (Amersham Biosciences). Bound antibodies were detected with HRP-conjugated secondary antibodies using ECL detection reagents (Amersham Biosciences). Rabbit polyclonal antibodies directed against AKT and p-AKT (Ser473) were purchased from Cell Signaling Technology. Precision Plus Protein All Blue Standards (Bio-Rad) were used for the molecular mass markers.

42The protocol has been published elsewhere,. Mice administered with A. indistinctus or PBS (vehicle control) for 5 to 6 weeks were used for the experiment. Jugular vein catheterization was performed 1 day before the clamp test. In brief, a mouse was anaesthetized with isoflurane (MSD), and the right jugular vein was exposed. A double-channel catheter was subsequently inserted to the vein. The next day, the mice were subjected to 4 h fasting before the clamp test. Human regular insulin (Eli Lilly) was intravenously administered at 7.5 mU kg−1 min−1, and the blood glucose levels were monitored every 5 min for 120 min. 50% glucose solution containing 6,6-d2-glucose (Sigma-Aldrich) was simultaneously infused to keep blood glucose levels around 100 to 120 mg dl−1. To separate the plasma, approximately 25 μl of blood was also drawn from tail vein at 0, 90, 105 and 120 min, placed into a tube containing 2 μl of heparin (Mochida Pharmaceutical) and centrifuged at 12,000g at 4 °C for 5 min. The plasma levels of glucose and 6,6-d2-glucose were measured using GC–MS. In brief, a 5 µl aliquot of plasma was extracted and derivatized with methoxyamine hydrochloride (Sigma-Aldrich) and N-methyl-N-(trimethylsilyl)trifluoroacetamide (GL Sciences), as previously described. The analysis was performed using a GC–MS/MS platform on a Shimadzu GCMS-TQ8040 triple quadrupole mass spectrometer (Shimadzu) with a capillary column (BPX5) (SGE Analytical Science/Trajan Scientific and Medical). The programme of GC–MS/MS analysis was previously described with minor modifications. We integrated each derivative peak to obtain mass isotopomers of glucose for the following ions: m/z 319.1, 320.1 and 321.1. The glucose infusion rate was determined as the infusion rate at 90, 105 and 120 min. The rate of glucose disappearance was determined on the basis of the plasma levels of 6,6-d2-glucose and total glucose using a non-steady-state equation as described previously, and considered as the whole-body glucose disposal after insulin stimulation. Hepatic glucose production was determined as the subtraction of glucose disappearance rate and glucose infusion rate.

43For the necropsy, the mice were anesthetized using isoflurane (MSD), and the left half of liver was dissected, weighed and frozen in liquid nitrogen. The extraction of triglyceride contents from the liver tissue has been reported elsewhere,. In brief, the samples were homogenized in buffer A (25 mM Tris-HCl at pH 7.4, 10 mM sodium orthovanadate, 10 mM sodium pyrophosphate, 100 mM sodium fluoride, 10 mM EDTA, 10 mM EGTA and 1 mM phenylmethylsulfonyl fluoride) and mixed with chloroform/methanol (2:1, v/v). The mixture was shaken for 15 min, centrifuged and the organic layer was collected. The extraction step was repeated three times. The collected samples were evaporated and resuspended in 1% Triton X-100/ethanol. The triglyceride content was assessed using Triglyceride E-test Wako (Wako) according to the manufacturer’s instructions.

44For general statistical comparisons, two-sided Wilcoxon rank-sum tests were used for two-group comparisons, Kruskal–Wallis tests followed by Dunn’s post hoc analysis were used for comparisons of more than two groups, and Fisher’s exact tests were used for comparison of categorical variables. For general correlation analyses, Spearman’s rank correlation in the function corr.test of the R package psych v.2.1.6 was used. For partial correlation analyses, partial Spearman’s rank correlation in the function pcor.test of the R package ppcor v.1.1 was used. To predict the metabolite levels and their CAGs (Fig. 1b,d and Extended Data Figs. 2c,d and 3a), rank-based regression analyses were performed using the function rfit of the R package Rfit (v.0.24.2). For the ordinal independent variables (that is, IR, MetS, and original categories with obese and prediabetes), IS, no MetS, and healthy categories were considered as the references, respectively, and the coefficients and P values for other categories were calculated against these reference categories. For the analyses involving generalized linear models (GLM) such as Fig. 2b and Extended Data Figs. 5i and 6c, the dependent variables were assumed to follow a Gamma distribution and arcsine square root transformation was applied to the relative-abundance values of microbiota and KEGG orthologues. To enhance comparability, the standardized coefficient was also calculated by standard deviations of dependent and independent variables using the function lm.beta of the R package QuantPsyc v.1.5 in Extended Data Fig. 5i. In the reanalysis of TwinsUK data, we fitted generalized linear mixed-effects models with age, sex, zygosity and BMI as fixed effects and sample collection year as a random effect using the function glmer of R package lme4 v.1.1-27.1 to estimate the associations between HOMA-IR and faecal carbohydrate metabolites (Extended Data Fig. 3b,c). Similarly, in the reanalysis of HMP2 data, we fitted a generalized linear mixed-effects model with consent age and sex as fixed effects and sample collection site as a random effect to estimate the associations between BMI and faecal fructose, glucose and/or galactose (Extended Data Fig. 3d). To analyse the associations between the participants’ clusters and clinical markers in Fig. 2b, the clusters were reordered before regression analyses according to their proportion of individuals with IR, where cluster C showing the lowest proportion of IR was set as the reference. To calculate the KEGG pathway enrichment associated with the participant clusters (Fig. 2e and Extended Data Fig. 6a,b), the KEGG orthologues were compared between each cluster and the remaining three clusters using a two-sided Wilcoxon rank-sum test, and significant (Padj < 0.05) KEGG orthologues were summarized into the pathway level (Supplementary Table 16). For comparison of metabolites in bacterial cultures (Extended Data Fig. 8), one-way ANOVA followed by Tukey’s post hoc test was performed, followed by multiple testing corrections based on the Benjamini–Hochberg procedure. For comparisons of time-series data such as insulin tolerance test, two-way repeated-measures ANOVA was used and the between-group difference was analysed by estimated marginal means. P < 0.05 was considered to be significant. To analyse the body mass change in animal experiments, ANCOVA analysis was performed to adjust baseline body mass (that is, body mass change as a dependent variable and group and baseline body mass as independent variables). We also validated the assumption of this ANCOVA model, that is, homogeneity of regression slopes, homogeneity of variances and normality of residuals. For multiple-testing corrections, P values were corrected using the Benjamini–Hochberg procedure using the R function p.adjust. Padj < 0.05 was used as a cut-off unless otherwise specified. All data were collected using Microsoft Excel 2016. All statistical and graphical analyses were conducted using R v.4.1.1 using R studio v.1.4.1717, unless otherwise specified.

45To analyse ROC curves of omics datasets, the datasets of faecal metabolomics, including hydrophilic and lipid metabolites, faecal 16S rRNA gene sequencing at the genus level, faecal metagenome consisting of KEGG orthologues and clinical metadata, were included. We first selected feature variables in each dataset, that is, the best explaining variables in the given dataset, using the minimum redundancy maximum relevance (mRMR) algorithm. The function mRMR.classic of the R package mRMRe v.2.1.2.1 was used for the calculation. The datasets were square-root-transformed before mRMR calculation. We selected 5 to 50 variables in 5 increments as the maximum number of genera was 50. Using the selected variables, we next established random-forest models using the R package caret v.6.0-88 to classify the individuals into IR or not. Specifically, the results of mRMR were split into train and test datasets in a 3:1 ratio. The generated random-forest models were evaluated using a tenfold cross-validation method and applied to the test datasets to obtain probability scores. The accuracy of each classification model was described by the AUC of ROC curves using the R package pROC v.1.17.0.1.

46To construct the co-abundance networks of genus-level bacteria, we selected 28 genus-level microorganisms that were observed in more than 40% of the participants and calculated the correlations using the R package CCREPE (compositionality corrected by renormalization and permutation) v.1.28.0 with Spearman’s correlations and the default settings. Interactions with Padj < 0.05 were selected for further analysis. Bacteria that exhibited a positive correlation with one another were determined to be members of an independent co-abundance microbial group, except for the interaction between Bacteroides and Robinsoniella. We decided to categorize Robinsoniella into the Blautia and Dorea group owing to its stronger correlation with Blautia in comparison to Bacteroides, both of which showed the highest centrality within their respective networks. Those weakly associated with each other or negatively associated with the members of other CAGs were classified as miscellaneous (Extended Data Fig. 5f). To characterize the microbial profiles of the study participants, the individuals were clustered on the basis of the abundance of 28 genera, which includes 20 genera in co-abundance microbial groups identified with CCREPE and 8 unclustered genera, using the ward.D function of the R package pheatmap v.1.0.12. Four distinct clusters of participants were determined, and the proportion of IR was compared using Fisher’s exact tests. Microorganism–metabolite networks were constructed on the basis of the correlations between the 28 genera observed in at least 40% of samples and the faecal metabolites, including all hydrophilic metabolites (n = 110) and bacteria-related lipid metabolites (n = 259). Bacteria-related metabolites were defined according to previous reports,. The following classes were selected: DGDG, PE-Cer, MGDG O, FAHFA, Cer-AS, Cer-BDS, NAGly, NAGlySer, PI-Cer, SL, AcylCer, bile acids, DGDG O and AAHFA. Positive and negative Spearman’s correlations with Padj < 0.05 were separately depicted in the networks. The networks were visualized using Cytoscape (v.3.7.0).

47To construct and visualize a correlation-based network of omics data, we first analysed IR-associated host signatures using plasma cytokines, plasma metabolites and CAGE promoter expression data. We identified the significant host markers through the following models: (1) GLM with a gamma distribution: HOMA-IR as a dependent variable and host markers, age and sex as independent variables; (2) logistic regression model: IR (HOMA-IR ≥ 2.5 = 1, HOMA-IR ≤ 1.6 = 0) as a dependent variable and significant host markers in the model 1, age and sex as independent variables. In both models, host markers with Padj < 0.05 were considered to be significant. We finally identified 6, 21 and 36 significant associations from plasma cytokines, plasma metabolites and CAGE promoter expression data, respectively (Supplementary Tables 19–21). In terms of bacteria, 20 genera with significant interactions between each other, which were identified with CCREPE as shown in Extended Data Fig. 5f, were included. In terms of faecal metabolites, 15 carbohydrates associated with IR in the CAG analysis as shown in Fig. 1b were included. Pairwise partial Spearman’s rank correlations adjusted by age, sex, BMI and FBG between all given factors were calculated with the R package ppcor v.1.1. The correlations with Padj < 0.05 were selected for visualization. The size of nodes was determined as the ratio of median abundance in IR over IS. As the median values of genera Robinsoniella and Rothia were zero, these elements were removed from the visualization. The width of lines was determined as the absolute value of partial Spearman’s coefficient. The networks were visualized using Cytoscape v.3.7.0. as in the microorganism–metabolite networks described above.

48To assess the explained variance of ten plasma cytokines, we established random-forest models using the R package caret v.6.0-88 to predict the plasma cytokine levels using 15 IR-associated faecal carbohydrates identified in Fig. 1b; 20 genera with significant interactions with each other that were identified in Fig. 2a; 21 IR-associated plasma hydrophilic metabolites (Supplementary Table 20); or 36 IR-associated CAGE promoters (Supplementary Table 21). Plasma cytokines were log10-transformed and scaled before the regression analyses. The data were split into train and test datasets at a 4:1 ratio. The generated random-forest models were evaluated using a tenfold cross-validation method and applied to the test datasets to obtain predictions. The explained variance shown as R2 was calculated as its definition: 1 − sum(test − predict)2/sum(test − mean(test))2. The negative values were considered as zero.

49To infer the effects of plasma cytokines on in silico causal relationships between faecal carbohydrates and IR markers (HOMA-IR, BMI, triglycerides and HDL-C), we performed causal mediation analysis using the R package mediation (v.4.5.0). As previously reported, we first screened significant associations (Padj < 0.05) between 15 IR-associated faecal carbohydrates and four IR markers, and significant associations between ten plasma cytokines and four IR markers. Age and sex were included as independent variables in both models. We then performed causal mediation analyses with the following models: (1) Mediator models: cytokine ~ metabolite + age + sex; (2) outcome models: IR marker ~ metabolite + age + sex + cytokine. In both models, faecal carbohydrate and plasma cytokine values were scaled before the analyses, and GLM with Gaussian distribution was used. A nonparametric bootstrap procedure was used to calculate the significance, followed by multiple testing corrections using the R function p.adjust. Average causal mediation effects and average direct effects with Padj values from representative models are reported in Extended Data Fig. 7d, whereas all of the results including the total effects and proportion mediated are reported in Supplementary Table 23.

50Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

51Any methods, additional references, Nature Portfolio reporting summaries, source data, extended data, supplementary information, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at 10.1038/s41586-023-06466-x.