Human subjects
All study participants were part of the IDEO cohort, which has been previously described [18,25]. Briefly, IDEO consists of 25-65 year-old men and women of multiple ethnicities and across a wide BMI range (18.5–52 kg/m2) living in the San Francisco Bay Area; exclusion factors include smoking, unstable weight within the last 3 months (>3% weight gain or loss), a diagnosed inflammatory or infectious disease, liver failure, renal dysfunction, cancer, and reported alcohol consumption of >20 grams per day. Using IDEO, we recruited both lean and obese W and EA individuals into this study based on World Health Organization cut-offs: W/EA BMI≤24.9 kg/m2 (lean); W BMI≥30 kg/m2 (obese); and EA BMI≥27.5 kg/m2 (obese) [17,26,27].
Each participant consented to take part in the study, which was approved by the University of California San Francisco (UCSF) Committee on Human Research. We utilized demographic, medical, dietary, and lifestyle metadata on each participant that were part of their initial recruitment into IDEO, as previously reported [18,25,28]. Participants with Type 2 Diabetes (T2D) were classified in accordance with American Diabetes Association Standards of Medical Care guidelines [29], defined by having glycated hemoglobin (HbA1c) ≥6.5% or the combination of a prior diagnosis of T2D and the active use of an antidiabetic medication. For stool sample collection, participants took home or were mailed a stool sample collection kit and detailed instructions on how to collect the specimen. All samples were collected at home, stored at room temperature, and brought to the UCSF Clinical Research Center by the participants within 24 hours of defecation. Samples were aliquoted and stored at -80°C.
Anthropometric and body composition measurements
We leveraged host phenotypic and demographic data from IDEO, which was the focus of two previous studies [18,25]. For the convenience of the reader, we restate our methods here. Height and weight were measured using a standard stadiometer and scale, and BMI (kg/m2) was calculated from two averaged measurements. Waist and hip circumferences (to the nearest 0.5 cm) were measured using a plastic tape meter at the level of the umbilicus and of the greater trochanters, respectively, and waist-to-hip ratio (WHR) was calculated. Blood pressure was measured with a standard mercury sphygmomanometer on the left arm after at least 10 minutes of rest. Mean values were determined from two independent measurements. Blood samples were collected after an overnight fast and analyzed for plasma glucose, insulin, serum total cholesterol, high density lipoprotein (HDL) cholesterol, and triglycerides. Low density lipoprotein (LDL) cholesterol was estimated according to the Friedewald formula [30]. Insulin resistance was estimated by the homeostatic model assessment of insulin resistance (HOMA-IR) index calculated from fasting glucose and insulin values [31]. Two obese subjects on insulin were included in the HOMA-IR analysis (1 EA, 1 W). Body composition of the subjects was estimated by Dual-Energy X-ray Absorptiometry (DEXA) using a Hologic Horizon/A scanner (3-minute whole-body scan, <0.1 G milligray) per manufacturer protocol. A single technologist analyzed all DEXA measurements using Hologic Apex software (13.6.0.4:3) following the International Society for Clinical Densitometry guidelines. Visceral adipose tissue (VAT) was estimated from a 5 cm-wide region across the abdomen just above the iliac crest, coincident with the fourth lumbar vertebrae, to avoid interference from iliac crest bone pixels and matching the region commonly used to analyze VAT mass by CT scan [32–34] . The short version of the International Physical Activity Questionnaire (IPAQ) was used to assess the habitual physical activity levels of the participants. The IPAQ total score is expressed in metabolic equivalent (MET)-minutes/week [35].
Dietary assessment
IDEO participants completed two dietary questionnaires, as previously described [18,25], allowing for the assessment of usual total fiber intake and fiber from specific sources, as well as macronutrient, phytochemical, vitamin, and mineral uptake. The first instrument was an Automated Self-Administered 24-hour Dietary Assessment (ASA24) [36,37], which queries intake over a 24-hour period. The 24-hour recalls and supplement data were manually entered in the ASA24 Dietary Assessment Tool (v. 2016), an electronic data collection and dietary analysis program. ASA24 employs research-based strategies to enhance dietary recall using a respondent-driven approach allowing initial recall to be self-defined. The second instrument was the National Cancer Institute’s Diet History Questionnaire III (DHQIII) [38,39]. The DHQIII queries one’s usual diet over the past month [39]. Completing the DHQIII was associated with participant survey fatigue and completion rates were accordingly only 42% after 1 phone-based administration of the instrument, although they improved to 79% by the 2nd session and reached 100% within four sessions over a 5-month period. Due to the effort needed to achieve DHQIII completion, we modified our protocol to request completion of the simpler ASA24 at three separate times, at appointments where there were computers and personnel assistance for online completion, in addition to completion of the DHQIII questionnaire. By combining both instruments, we were able to reliably obtain complete dietary information on all participants.
DNA extraction
Human stool samples were homogenized with bead beating for 5 min (Mini-Beadbeater-96, BioSpec) using beads of mixed size and material (Lysing Matrix E 2mL Tube, MP Biomedicals) in the digestion solution and lysis buffer of a Wizard SV 96 Genomic DNA kit (Promega). The samples were centrifuged for 10 min at 16,000 g and the supernatant was transferred to the binding plate. The DNA was then purified according to the manufacturer’s instructions. Mouse fecal pellets were homogenized with bead beating for 5 min (Mini-Beadbeater-96, BioSpec) using the ZR BashingBead lysis matrix containing 0.1 and 0.5 mm beads (ZR-96 BashingBead Lysis Rack, Zymo Research) and the lysis solution provided in the ZymoBIOMICS 96 MagBead DNA Kit (Zymo Research). The samples were centrifuged for 5 min at 3,000 g and the supernatant was transferred to 1 mL deep-well plates. The DNA was then purified using the ZymoBIOMICS 96 MagBead DNA Kit (Zymo Research) according to the manufacturer's instructions.
16S rRNA gene sequencing and analysis
For human samples, 16S rRNA gene amplification was carried out using GoLay-barcoded 515F/806R primers [40] targeting the V4 region of the 16S rRNA gene according to the methods of the Earth Microbiome Project (earthmicrobiome.org) (Table S2). Briefly, 2 µL of DNA was combined with 25 µL of AmpliTaq Gold 360 Master Mix (Fisher Scientific), 5 µL of primers (2 µM each GoLay-barcoded 515/806R), and 18 µL H2O. Amplification was as follows: 10 min 95°C, 30x (30s 95°C, 30s 50°C, 30s 72°C), and 7 min 72°C. Amplicons were quantified with PicoGreen (Quant-It dsDNA; Life Technologies) and pooled at equimolar concentrations. Aliquots of the pool were then column (MinElute PCR Purification Kit; Qiagen) and gel purified (QIAquick Gel Extraction Kit; Qiagen). Libraries were then quantified (KAPA Library Quantification Kit; Illumina) and sequenced with a 600 cycle MiSeq Reagent Kit (250x150; Illumina) with ~15% PhiX spike-in. For mouse samples, 16S rRNA gene amplification was carried out as per reference protocol and primers [41]. In brief, the V4 region of the 16S rRNA gene was amplified with 515F/806R primers containing common adaptor sequences, and then the Illumina flow cell adaptors and dual indices were added in a secondary amplification step (see Table S3 for index sequences). Amplicons were pooled and normalized using the SequalPrep Normalization Plate Kit (Invitrogen). Aliquots of the pool were then column (MinElute PCR Purification Kit, Qiagen) and gel purified (QIAquick Gel Extraction Kit, Qiagen). Libraries were then quantified and sequenced with a 600 cycle MiSeq Reagent Kit (270x270; Illumina) with ~15% PhiX spike-in.
Demultiplexed sequencing reads were processed using QIIME2 v2020.2 (qiime2.org) with denoising by DADA2 [42]. Taxonomy was assigned using the DADA2 implementation of the RDP classifier [43] using the DADA2 formatted training sets for SILVA version 138 (benjjneb.github.io/dada2/assign.html). For Amplicon Sequence Variants (ASV) analysis, we utilized quality scores to set truncation and trim parameters. The reverse read of human 16S data suffered from low sequence quality and reduced the overall ASV counts, so we therefore analyzed only the forward reads, although a separate analysis using merged forward and reverse reads complemented the findings we report in this manuscript. For the manuscript, forward reads were truncated to 220 base pairs and underwent an additional 5 base pairs of trimming for 16S analysis of human stool. For gnotobiotic mice, forward and reverse reads were truncated to 200 and 150 base pairs respectively. ASVs were filtered such that they were present in more than one sample with at least a total of 10 reads across all samples. Alpha diversity metrics were calculated on subsampled reads using Vegan [44] and Picante [45] R packages. The PhILR Euclidean distance was calculated by first carrying out the phylogenetic isometric log ratio transformation (philr, PhILR [46]) followed by calculating the Euclidean distance (vegdist, Vegan [44]). Principal coordinates analysis was carried out using the pcoa function of APE [47]. ADONIS calculations were carried out (adonis, Vegan) with 999 replications on each distance metric. Centered log2-ratio (CLR) normalized abundances were calculated using the Make.CLR function in MicrobeR package [48] with count zero multiplicative replacement (zCompositions; [49]). ALDEx2 [50] was used to analyze differential abundances of count data, using features that represented at least 0.05% of total sequencing reads. Corrections for multiple hypotheses using the Benjamini-Hochberg method [51] were performed where applicable. Where described, a false discovery rate (FDR) indicates the Benjamini-Hochberg adjusted p-value for an FDR (0.1 unless otherwise specified). Analysis of distance matrices and alpha diversity mirror prior analyses developed in the Turnbaugh lab and were adapted to the current manuscript [9]. Calculations of associations between ASVs and ASA24 questionnaire data were completed by calculating a Spearman rank correlation and then adjusting the p-value for a Benjamini-Hochberg FDR using the cor_pmat function in the R package ggcorrplot [52]. The randomForest package [53] was employed to generate random forest classifiers. Given the total number of samples (n=46) we generated 46 classifiers trained on a subset of 45 samples and used each classifier to predict the sample left out. AUCs are visualized utilizing the pROC [54] and ROCR [55] packages.
Metagenomic sequencing and analysis
Whole-genome shotgun libraries were prepared using the Nextera XT DNA Library Prep Kit (Illumina). Paired ends of all libraries were sequenced on the NovaSeq 6000 platform in a single sequencing run (n=45 subjects; see Table S2 for relevant metadata and statistics). Illumina reads underwent quality trimming and adaptor removal using fastp [56] and host read removal using BMTagger v1.1.0 (ftp.ncbi.nlm.nih.gov/pub/agarwala/bmtagger/) in the metaWRAP pipeline (github.com/bxlab/metaWRAP) [57]. Metagenomic samples were taxonomically profiled using MetaPhlan2 v2.7.7 [58] and functionally profiled using HUMAnN2 v0.11.2 [59], both with default parameters. Principal coordinates analysis on MetaPhlan2 species-level abundances was carried out using Bray Curtis distances and the pcoa function of APE [47]. Metaphlan2 abundance outputs were converted to counts and subsampled to even sample depth. Differences between groups were determined utilizing the Aldex2 package as described above. Tables of gene family abundances from HUMAnN2 were regrouped to KEGG orthologous groups using humann2_regroup_table. Functional pathways relating to short-chain fatty acid production were manually curated from the pathway outputs from HUMANn2 and normalized by the estimated genome equivalents in each microbial community obtained from MicrobeCensus [60].
Quantification of bacterial load
Quantitative PCR (qPCR) was performed on DNA extracted from the human stool samples. DNA templates were diluted 1:10 into a 96-well plate. Samples were aliquoted in a 384-well plate, and PCR primers and iTaq Universal Probes Supermix were added utilizing an Opentrons OT-2 instrument then analyzed on a BioRad CFX384 thermocycler with an annealing temperature of 60˚C. The following primers including a FAM labeled PCR probe was used for quantification: 891F, TGGAGCATGTGGTTTAATTCGA; 1003R, TGCGGGACTTAACCCAACA; 1002P, [6FAM]CACGAGCTGACGACARCCATGCA[BHQ1]. Absolute quantifications were determined against a standard curve of purified 8F/1542R amplified Vibrio casei DNA. Reactions identified as inappropriately amplified by the instrument were rejected, and the mean values were used for downstream analysis. Absolute 16S rRNA gene copy number was derived by adjustments for dilutions during DNA extraction and template normalization dividing by the total fecal mass used for DNA extraction in grams.
Nuclear magnetic resonance (NMR) metabolomics
NMR spectroscopy was performed at 298K on a Bruker Avance III 600 MHz spectrometer configured with a 5 mm inverse cryogenic probe (Bruker Biospin, Germany) as previously described [61]. Lean and obese EA and W individuals (n=20 total individuals, five in each group) were selected and matched based on body composition and metabolic parameters. Stool samples from these subjects were subjected to NMR-based metabolomics. 50 mg of human feces were extracted with 1 mL of phosphate buffer (K2HPO4/NaH2PO4, 0.1 M, pH 7.4, 50% v/v D2O) containing 0.005% sodium 3-(trimethylsilyl) [2,2,3,3-2H4] propionate (TSP-d4) as a chemical shift reference (δ 0.00). Samples were freeze-thawed three times with liquid nitrogen and water bath for thorough extraction, then homogenized (6500 rpm, 1 cycle, 60 s) and centrifuged (11,180 g, 4 °C, 10 min). The supernatants were transferred to a new 2 mL tube. An additional 600 μL of PBS was added to the pellets, followed by the same extraction procedure described above. Combined fecal extracts were centrifuged (11,180 g, 4°C, 10 min), 600 μL of the supernatant was transferred to a 5 mm NMR tube (Norell, Morganton, NC) for NMR spectroscopy analysis. A standard one-dimensional NOESY pulse sequence noesypr1d (recycle delay-90°-t1-90°-tm-90°-acquisition) is used with a 90 pulse length of approximately 10s (-9.6 dbW) and 64 transients are recorded into 32k data points with a spectral width of 9.6 KHz. NMR spectra were processed as previously described [61]. First, spectra quality was improved with Topspin 3.0 (Bruker Biospin, Germany) for phase and baseline correction and chemical shift calibration. AMIX software (version: 3.9.14, Bruker Biospin, Germany) was used for bucketing (bucket width 0.004 ppm), removal of interfering signal, and scaling (total intensity). Relative concentrations of identified metabolites were obtained by normalized peak area.
Targeted gas chromatography mass spectrometry (GC-MS) assays
Targeted analysis of short-chain fatty acids (SCFAs) and branched-chain amino acids (BCAAs) was performed with an Agilent 7890A gas chromatograph coupled with an Agilent 5975 mass spectrometer (Agilent Technologies Santa Clara, CA) using a propyl esterification method as previously described [61]. 50 mg of human fecal samples were pre-weighed, mixed with 1 mL of 0.005 M NaOH containing 10 μg/mL caproic acid-6,6,6-d3 (internal standard) and 1.0 mm diameter zirconia/silica beads (BioSpec, Bartlesville, OK). The mixture was thoroughly homogenized and centrifuged (13,200 g, 4°C, 20 min). 500 μL of supernatant was transferred to a 20 mL glass scintillation vial. 500 μL of 1-propanol/pyridine (v/v=3/2) solvent was added into the vial, followed by a slow adding of an aliquot of 100 μL of esterification reagent propyl chloroformate. After a brief vortex of the mixture for 1 min, samples were derivatized at 60°C for 1 hour. After derivatization, samples were extracted with hexane in a two-step procedure (300 μL + 200 μL) as described [62]. First, 300 μL of hexane was added to the sample, briefly vortexed and centrifuged (2,000g, 4°C, 5 min), and 300 μL of the upper layer was transferred to a glass autosampler vial. Second, an additional 200 μL of hexane was added to the sample, vortexed, centrifuged, and the 200 μL upper layer was transferred to the glass autosampler vial. A combination of 500 μL of extracts were obtained for GC-MS analysis. A calibration curve of each SCFA and BCAA was generated with series dilution of the standard for absolute quantitation of the biological concentration of SCFAs and BCAAs in human fecal samples.
Targeted bile acid quantitation by UHPLC-MS/MS
Bile acid quantitation was performed with an ACQUITY ultra high pressure liquid chromatography (UHPLC) system using a Ethylene Bridged Hybrid C8 column (1,7 µm, 100 mm x 2.1 mm) coupled with a Xevo TQ-S mass spectrometer equipped with an electrospray ionization (ESI) source operating in negative mode (All Waters, Milford, MA) as previously described [63]. Selected ion monitoring (SIM) for non-conjugated bile acids and multiple reaction monitoring (MRM) for conjugated bile acids was used. 50 mg of human fecal sample was pre-weighed, mixed with 1 mL of pre-cooled methanol containing 0.5 μM of stable-isotope-labeled bile acids (internal standards) and 1.0 mm diameter zirconia/silica beads (BioSpec, Bartlesville, OK), followed by thorough homogenization and centrifugation. Supernatant was transferred to an autosampler vial for analysis. 100 µL of serum was extracted by adding 200 µL pre-cooled methanol containing 0.5 μM deuterated bile acids as internal standards. Following centrifugation, the supernatant of the extract was transferred to an autosampler vial for quantitation. Calibration curves of individual bile acids were drafted with bile acid standards for quantitation of the biological abundance of bile acids.
Gnotobiotic mouse experiments
All mouse experiments were approved by the UCSF Institutional Animal Care and Use Committee and performed accordingly. Germ-free mice were maintained within the UCSF Gnotobiotic Core Facility and fed ad libitum autoclaved standard chow diet (Lab Diet 5021). Germ-free adult male C57BL/6J mice between 6-10 weeks of age were used for all the experiments described in this paper. 10 lean subjects in our IDEO cohort were selected as donors for the microbiota transplantation experiments, including 5 EA and 5 W donors. The selected donors for gnotobiotic experiments were matched for phenotypic data to the degree possible (Table S4). Stool samples to be used for transplantation were resuspended in 10 volumes (by weight) of brain heart infusion media in an anaerobic Coy chamber. Each diluted sample was vortexed for 1 min and left to settle for 5 min, and a single 200 µL aliquot of the clarified supernatant was administered by oral gavage into each germ-free mouse recipient. In experiments LFPP1 and LFPP2, microbiome transplantations were performed for 2 donors per experiment (1 W, 1 EA) with gnotobiotic mice housed in sterile isolators (CBC flexible, softwall isolator) and maintained on ad libitum standard chow also known as low-fat, high-plant-polysaccharide (LFPP) diet. In LFPP1, 6 germ-free mice per colonization group received an aliquot of stool from a donor of either ethnicity and body composition (measured using EchoMRI) were recorded on the day of colonization and at 6 weeks post-transplantation (per group n=6 recipient mice, 1 isolator, 2 cages). In LFPP2, we shortened the colonization time to 3 weeks and used two new donor samples. For the third experiment (HFHS experiment), mice were weaned onto an irradiated high-fat, high-sugar diet (HFHS, TD.88137, Envigo) for four weeks prior to colonization and housed in pairs in Tecniplast IsoCages. The same 4 donors from LFPP1 and LFPP2 were included in the HFHS experiment, in addition to 6 new donors (per donor n=2 recipient mice, 1 IsoCage). Body weight and body composition were recorded on the day of colonization and again at 3 weeks post-transplantation. Mice were maintained on the HFHS diet throughout the experiment. All samples were sequenced in a single pool (Table S3). For comparisons between donors and recipient mice, donors and recipient mice were subsampled to even sequencing depth and paired between donor and recipient mice (range: 18,544-78,361 sequencing reads/sample).
Glucose tolerance tests
Food was removed from mice 10 hr (LFPP1 experiment) or 4 hr (HFHS experiment) prior to assessment of glucose tolerance. Mice received i.p. injections of D-glucose (2 mg/kg), followed by repeated collection of blood by tail nick and determination of glucose levels by handheld glucometer (Abbott Diabetes Care) over a 2-hour period.
Geographic analyses
Map tiles and distance data was obtained using GGMap [64], OpenStreet Maps [65], and the Imap R [66] packages. GGMap was employed using a Google Cloud API key and the final map tiles were obtained in July 2020 [64]. Spearman ranked correlation coefficients (rho) were calculated as embedded in the ggpubr [67] R package. 2018 US Census data for EA and W subjects was obtained (B02001 table for race, data.census.gov) for the ZIP codes available in our study and using the leaflet [68] package. The census data used is included as part of Table S2 to aid in reproduction. Each census region is plotted as a percentage of W individuals over a denominator of W and EA subjects. The leaflet package utilized ZIP Code Tabulation Areas (ZCTAs) from the 2010 census. We extracted all ZCTAs starting with 9, and the resulting 29 ZIP codes that overlap with IDEO subjects were analyzed (Table S2). Two ZCTAs (95687 and 95401) were primarily W when comparing W and EA subjects. There were two W subjects recruited from these ZTCAs. These ZIP codes are cut off based on the zoom magnification for that figure and as a result ZTCAs for 27 individuals are plotted. Distance to a central point in SF was calculated. The point of reference was latitude=37.7585102, longitude=-122.4539916.
Dietary questionnaire correlation analysis
DHQIII and ASA24 data were analyzed using a Euclidean distance matrix. These transformations were completed using the cluster package [69]. Subsequent analysis was completed using the vegan package [44]. Procrustes transformations were performed using 16S-seq data from human subjects, which was then subjected to a PhILR transformation. The resulting matrix was rotated against the distance matrix for ASA24 or DHQIII questionnaire data using the procrustes command in the vegan R package using 999 permutations. Mantel statistics were calculated utilizing the mantel command of the vegan package.
R packages used in this study
Picante [45], PhILR [46], MicrobeR [48], ALDEx2 [50], ggcorrplot [52], randomForest [53], GGMap [64], OpenStreetMap [65], IMap [66], ggpubr [67], leaflet [68], cluster [69], readxl [70], Rtsne [71], vegan [72], ape [73], tigris [74], lmerTest [75], qiime2R [76], gghighlight [77], Phyloseq [78], Janitor [79], table 1 [80], ggplot2 [81].
Statistical analyses
Statistical analysis of the human data was performed using the table1 package in R (STATCorp LLC. College Station, TX). Human data were presented as mean ± SD. Unpaired independent Student’s t tests were used to compare differences between the two groups in the case of continuous data and in the case of categorical data the χ2 test was utilized. These tests were adjusted for a Benjamini-Hochberg false discovery rate utilizing the command p.adjust in R, which is indicated as an adjusted p-value in the tables. In Tables S9 and S10 no values met an adjusted p-value cutoff of < 0.1. In Table S1, p-values indicated by numbers were pooled together for adjustments and those represented by symbols were separately pooled together for adjustment. All microbiome-related analyses were carried out in R version 3.5.3 or 4.0.2. Where indicated, Wilcoxon rank-sum tests were calculated. A Benjamini-Hochberg adjusted p-value (FDR) of 0.1 was used as the cutoff for statistical significance unless stated otherwise. Statistical analysis of glucose tolerance tests was carried out using linear mixed effects models with the lmerTest [75] R package and mouse as random effect. Graphical representation was carried out using ggplot2. Boxplots indicate the interquartile range (25th to 75th percentiles), with the center line indicating the median and whiskers representing 1.5x the interquartile range.