Associations of the Fecal Microbial Proteome Composition and Proneness to Diet-induced Obesity

Individuals exhibit marked heterogeneity in the degree of adiposity that results from consumption of an obesogenic diet. We hypothesized that markers of inflammation, fecal microbiome, and/or proteome composition might associate with heterogeneity in diet-induced obesity. We found the fecal proteome, which was comprised of microbial proteins that changed markedly upon exposure to the obesogenic diet associated with extent of obesity and associated dysglycemia. These results suggest the fecal microbial proteome plays a role in diet-induced obesity.


INTRODUCTION
Obesity is an emerging 21 st century epidemic. Obesity, and the disease states it drives, including type 2 diabetes, cardiovascular disease, and liver disease threaten to overwhelm healthcare systems 1 . Thus, obesity is a contemporary medical concern that poses a grave public health crisis in dire need of a solution. The increased incidence in obesity is thought to have been driven by broad societal changes that have resulted in reduced physical activity and increased availability of palatable low-cost energy-rich foods 2 . Yet the extent to which individuals develop obesity in such an environment is highly heterogeneous. Variations in individual genetics contribute to, but are insufficient to fully explain, such heterogeneity. For example, studies characterizing weightdiscordant monozygotic twins has shown that individuals with shared environmental, physical activity, and genetic factors display heterogeneity in adiposity 3 . Similarly, ratbased studies show marked heterogeneity in weight gain and adiposity in response to obesogenic diets even when using highly inbred animals in a well-controlled environment 4,5 . Better understanding non-genetic factors that influence proneness to obesity might aid the identification of individuals at-risk for development of obesity and can yield modifiable factors to ameliorate this disease state.
A number of factors that are at least partially independent of genetics are proposed to influence proneness to diet-induced obesity (DIO). One potential central nexus of such factors is inflammation, impacting metabolic signaling pathways including insulin and leptin 6 , which have well-established roles in feeding behavior. Inflammation is also suggested to promote behavioral patterns such as anxiety-like and anti-social behaviors that can impact food consumption 7 . While numerous elements impact by guest on January 10, 2021 https://www.mcponline.org Downloaded from inflammation, one increasingly appreciated factor is the gut microbiota [8][9][10][11][12] , which is the collective term for the large diverse community of microorganisms that inhabit the gastrointestinal tract. Indeed, in humans, gut microbiota composition is associated with obesity. One way microbiota composition influences metabolic signaling is via lipopolysaccharide (LPS), which activates pro-inflammatory signaling via Toll-like receptor 4 (TLR4) resulting in production of molecules including tumor necrosis factor alpha (TNF-α), and interleukin-6 (IL-6). These molecules interfere with leptin and insulin signaling, wherein LPS derived from gamma-proteobacteria is a particularly potent activator of TLR4 13 . Another host-microbiota interaction implicated in inflammation and obesity is the sensing of flagella through TLR5, which keeps motile bacteria in-check by a range of mechanisms including production of antimicrobial peptides and promoting production of anti-flagella immunoglobulins that help regulate the microbiota in the healthy gut 14 . In addition to its impacts on inflammation, microbiota composition is also reported to influence energy harvest from ingested food 12,15 . Hence, in light of its ability to impact inflammation, metabolism, and behavior, gut microbiota composition might provide a means of identifying host proneness to obesity when presented with an obesogenic diet.
Here, we sought to identify microbiota-based markers that might predict proneness to diet-induced obesity, specifically exposing mice to a western-style, lowfiber high fat diet (HFD). Both targeted and untargeted approaches were utilized including 16S rRNA gene amplicon sequencing for microbial community profiling and a Tandem Mass Tag (TMT) based multiplexed mass spectrometry (MS) approach for analysis of the fecal metaproteome. Additionally, we measured behavior, inflammatory markers, and metabolic parameters. Notably, we show that the fecal metaproteome appears to be a promising candidate for distinguishing mice with differential responses to obesogenic diets. Collectively, this study provides insight into potential mechanisms regarding the host-microbiota interactions mediating response to HFD exposure, and highlights putative biomarkers for predicting DIO.

Mice and high-fat diet administration
Female, [3][4][5] week old C57BL/6 mice were purchased from Jackson Laboratory (Bar Harbor, ME) and maintained at Georgia State University, Atlanta, Georgia, USA under institutionally approved protocol under approved protocols (IACUC # A14033), housed 5 mice per cage, were subjected to metabolic monitoring, including behavior analysis and sample collection over a 3-week period. During this time, the mice were fed standard grain-based chow (GBC), which is comprised of relatively unrefined ingredients. The cohort of mice was then switched to a diet composed of 60% kcal from fat (Research Diet, D12492) for 8 weeks. Mice were then euthanized, and colon length, colon weight, spleen weight and adipose weight were measured. Serum, feces, and organs were collected for downstream analysis.

Fecal metaproteome data acquisition
Fecal samples were measured out to ~0.2 g and suspended in 10 mL of ice-cold, sterilized TBS. A 20 µM vacuum, steriflip (Milipore) filter was used to remove particulate from the samples. Cells were pelleted through centrifugation at 4000 rpm for 10 min.
Next, cells were lysed in 2 mL of buffer containing 75 mM NaCl (Sigma), 3% sodium dodecyl sulfate (SDS, Fisher), 1 mM NaF (Sigma), 1 mM beta-glycerophosphate (Sigma), 1 mM sodium orhtovanadate (Sigma), 10 mM sodium pyrophosphate (Sigma), 7 through two 10-second intervals of probe sonication at 25% amplitude. Proteins were then reduced with dithiothreitol (DTT, Sigma), alkylated through iodoacetamide (Sigma), and quenched as previously described 17 . Proteins were then precipitated via chloroform-methanol precipitation and protein pellets were dried 18 . Protein pellets were re-suspended in 1M urea in 50 mM HEPES, pH 8.5 and digested overnight at room temperature with LysC (Wako) 19 . A second, 6-hour digestion using trypsin at 37 ºC was performed and the reaction was stopped through addition of 10% trifluoroacetic acid (TFA, Pierce). Samples were then desalted through C18 Sep-Paks (Waters) and eluted with a 40% and 80% Acetonitrile solution containing 0.5 % Acetic Acid 20 . Concentration of desalted peptides was determined with a BCA assay (Thermo Scientific). 50 µg aliquots of each sample were dried in a speed-vac, additional bridge channels consisting of 25 µg from each sample were created and 50 µg aliquots of this solution were used in duplicate per TMT 10-plex (Thermo Scientific) as previously described 21 .
These bridge channels were used to control for labeling efficiency, inter-run variation, mixing errors and the heterogeneity present in each sample 22 . Each sample or bridge channel was resuspended in 30% dry acetonitrile in 200 mM HEPES, pH 8.5 for TMT labeling with 7 µL of the appropriate TMT reagent 23 . Reagents 126 and 131 (Thermo Scientific) were used to bridge between mass spec runs. Remaining reagents were used to label samples in random order. Labeling was carried out for 1 hour at room temperature, and quenched by adding 8 µL of 5% hydroxylamine (Sigma). Labeled samples were acidified by adding 50 µL of 1% TFA. After TMT labeling, each 10-plex experiment was combined and desalted through C18 Sep-Paks and dried in a speedvac. Each 10-plex experiment was fractionated using a High pH Reversed-Phase by guest on January 10,2021 Peptide Fractionation Kit (Pierce) per manufacturer instructions. All LC-MS 2 /MS 3 experiments were carried out on an Orbitrap Fusion (Thermo Fisher Scientific) with an in-line Easy-nLC 1000 (Thermo Fisher Scientific) and chilled autosampler. Separation and acquisition settings were as previously defined 24 .

Metaproteome data processing
Data was processed using Proteome Discoverer 2.1 (Thermo Fisher Scientific). MS 2 data was searched against a catalog of mouse gut genes 25 25 . A peptide and protein false discovery rate of 1% was enforced using a reverse database search strategy 30-32 . by guest on January 10,2021 TMT reporter ion intensities were extracted from MS 3 spectra for quantitative analysis and signal-to-noise ratios were used for quantitation. Additional stringent filtering was used removing any moderate confidence peptide spectral matches (PSMs), or ambiguous PSM assignments. Additionally, any peptides with a spectral interference above 25% were removed, as well as any peptides with an average signal to noise ratio less than 10. Normalization occurred as previously described 24 . Briefly, relative abundances are normalized first to the pooled standards for each protein and then to the median signal across the pooled standard. An average of these normalizations was used for the next step. To account for slight differences in amounts of protein labeled, these values were then normalized to the median of the entire dataset and reported as final normalized summed signal-to-noise ratios per protein per sample.

Behavioral analysis
Three weeks after arrival in the animal facilities, behavior in the open field and in the home cage was assessed in a counter-balanced fashion over the course of two days.
Behavioral testing occurred within the last 4 h of the light and quiescent phase and was conducted under illumination of overhead white lighting (between 300-400 lux). Open field arenas were cleaned with 70% ethanol between trials, and home cage bedding was changed after each trial. Behavioral tests were videotape using a Sony camcorder for later analysis by The Observer version XT11 (Noldus Information Technology Inc., Wageningen, The Netherlands). An experiment blinded as to the treatment conditions scored behavioral tests in the Observer. by guest on January 10, 2021

Open Field Test
Locomotor behavior was assessed in a 43.2 X 43.2 X 30.5cm (WxLxH) Plexiglas arena (Med Associates, Inc., St. Albans, VT) containing 2 arrays of infrared transmitters strips (16 beams each) located on the bottom of the arena (in the X and Y plane). The center zone of the arena was defined as square containing the center 8 beams (e.g., beams 4-12) in the X and Y plane. Each mouse was placed into the arena with its nose facing the wall and allowed to freely investigate for 10 min. The total distance traveled, the total time spent in the center of the arena, and circling behaviors, which are defined as movements below a preset ambulatory threshold, were calculated by Activity Monitor (Med Associates, Inc.) on a computer connected to the open field arenas.

Home Cage Behavior
Mice were placed into a clean housing cage containing 2 cm deep Alpha-dri bedding (Shepherd Specialty Paper, Fibercore, Cleveland, OH, USA) and video recorded for 10 min. An experimenter blinded as to condition scored the occurrence and duration of i) the time spent walking, as defined by locomotion along the bottom of the enclosure, around the arena, ii) grooming, as defined by stroking or scratching the face of body iii) digging, as defined as using the fore-or hind paws to displace the bedding, and iv) the rears, as defined by standing on the hind legs with either the forepaws unsupported or when the forepaws were supported by the walls of the enclosure, were quantified using the Observer.

Fasting blood glucose measurement and body composition measurement.
by guest on January 10, 2021 For fasting blood glucose tolerance test, mice were fasted for 5 hours, and baseline blood glucose were measured by using a Nova Max plus Glucose meter and expressed in mg/dL. Measurement of percent fat mass and lean mass was performed via MRI (Bruker MiniSpec) at day 0, prior to diet treatment, and day 28 and 56, after diet treatment.

Fecal sample preparation for immunoglobulin quantification.
Fecal sample preparation of enzyme-linked immunosorbent assay (ELISA) has been previously described 33
Briefly, 96-well microtiter plates (Costar, Corning, New York) were coated with 100 ng/well of laboratory-made flagellin in 9.6 pH bicarbonate buffer overnight at 4° C.
Serum samples from mice were then applied either pure or at a 1:100 dilution for 1 hour at 37° C. After incubation and washing, the wells were incubated with either horseradish peroxidase-linked anti-mouse IgG (GE Healthcare Life Sciences, Pittsburgh, Pennsylvania) or horseradish peroxidase-linked anti-IgA (Southern Biotech, by guest on January 10,2021 Birmingham, Alabama). Quantification of immunoglobulin was then developed by the addition of 3,3',5,5'-Tetramethylbenzidine and the optical density was calculated by the difference between readings at 450nm and 540nm.

Fecal flagellin and lipopolysaccharide load quantification
We quantified flagellin and lipopolysaccharide (LPS) as previously described using frozen extruded feces using a PowerSoil-htp kit from MoBio Laboratories (Carlsbad, California, USA) with mechanical disruption (bead-beating). The 16S rRNA genes, region V4, were PCR amplified from each sample using a composite forward primer and a reverse primer containing a unique 12-base barcode, designed using the Golay errorcorrecting scheme, which was used to tag PCR products from respective samples 40 . We

16S rRNA gene sequence analysis
Forward and reverse Illumina reads were joined using the fastq-join method 41,42 , sequences were demultiplexed, quality filtered using Quantitative Insights Into Microbial Ecology (QIIME, version 1.8.0) software package 43 . QIIME default parameters were used for quality filtering (reads truncated at first low-quality base and excluded if: (1) there were more than three consecutive low quality base calls (2), less than 75% of read length was consecutive high quality base calls (3), at least one uncalled base was present (4), more than 1.5 errors were present in the bar code (5), any Phred qualities were below 20, or (6) the length was less than 75 bases). Sequences were clustered to operational taxonomic units (OTUs) using UCLUST algorithm 44  Size) was used to investigate bacterial members that drive differences between groups 49 . Unprocessed sequencing data will be deposited in the European Nucleotide Archive.

Study design
The overall study included 50 mice which were followed longitudinally for monitoring of weight gain and other measures. The sample size was determined for statistical power based on previous publications 4, 5 and experience. The metaproteome analysis included four mice with the highest weight gain and the lowest weight gain. These samples sizes were determined to be sufficient based on previous reports of strong differences in related animal models with similar sample sizes 50 .
Extra files associated with the analysis within the notebooks are deposited as supplementary files in the MassIVE (https://massive.ucsd.edu) repository for this study (Study ID: MSV000083891). All analysis was performed on the proteins identified in both TMT 10-plex experiment. Qiime2, version 2019.1 (https://qiime2.org/), was used for principle coordinates analysis through the command "qiime diversity core-metrics" as well as for determining significance of beta-diversity clustering through the command "qiime diversity beta-group-significance". K-means clustering was performed through by guest on January 10, 2021 Morpheus (https://software.broadinstitute.org/morpheus). Enriched and depleted proteins were determined by π-score, which accounts for both fold change and pvalue 51 . A statistical cutoff for highly ranked associations was set to |π| > 1 (α ~ 0.05), which provided an adequate number of proteins for functional and taxonomic assessment 52 while maintaining a moderate stringincy. Volcano plots were visualized using GraphPad Prism (version 7.0b). Mouse protein gene functional enrichment analysis was performed using DAVID 53 , with all mouse proteins identified as a background list. The python package, Seaborn (version 0.9.0) was used for boxplots, swarmplots, and catplots. Statistical analysis between groups within the boxplots was performed using ANOVA with Dunnett corrected p-values through GraphPad Prism (version 7.0b).

Statistical analysis
Significance was determined using unpaired two-tailed t-test or linear regression analysis (GraphPad Prism software, version 6.01). Differences were noted as significant *p≤0.05 for t-test or linear regression analysis. For clustering analyzing on principal coordinate plots, categories were compared and statistical significance of clustering was determined via Permanova 54 .
by guest on January 10, 2021

RESULTS
Stratification and characterization of mice prone, or resistant, to HFD-induced metabolic syndrome.
The primary goal of this study was to elucidate factors that predict and possibly govern, susceptibility to developing obesity in response to administration of an obesogenic diet. Hence, we designed a prospective study wherein 50, 3-week old female C57BL/6 mice, housed 5 mice per cage, were subjected to metabolic monitoring, including behavior analysis and sample collection over a 3-week period. During this time, the mice were fed standard grain-based chow (GBC), which is comprised of relatively unrefined ingredients. The cohort of mice was then switched to a diet compositionally low in fiber (5%) and high in fat (35% by mass, 60% by calories), herein referred to as an obesogenic diet or high-fat diet (HFD) for an 8-week period. Prior to, during and after administration of the high-fat diet, sample collection and monitoring was performed as outlined in Fig. 1A. In accord with our previous rodent-based studies 5 , the extent of obesity following the obesogenic diet was quite heterogeneous with many mice weighing between 20-25 grams, which is the approximate weight of age-matched GBC-fed mice of this strain/gender. In contrast, some mice appeared to dramatically gain weight over the course of the experiment with final weights over 30 grams.
Therefore, based on their final body weight, mice were stratified into tertiles as being prone, intermediate, or relatively resistant to being obese following exposure to an obesogenic diet (Fig. 1B). First, we examined if mice prone or resistant to DIO clustered within cages but did not observe a distribution pattern to support this possibility (Fig.   1C). Nor were these groupings significantly related to the initial weight of the mice (Fig. by guest on January 10, 2021 1D). The total weight gain over the period of exposure to the obesogenic diet for resistant mice was about 40%. This observation is approximately the expected agerelated weight gain of GBC-fed mice during this period, while prone mice increased in weight by about 70% during this period (Fig. 1E). Fat mass, as determined by magnetic resonance imaging (MRI), prior to, during, or at the end of exposure to HFD, was highly correlated with body weight within the cohort (Fig. 1F). Accordingly, post-euthanasia weight of the periovarian fat pad, which has classically been used to assess adiposity in mice correlated closely (r 2 = 0.8229) with final body weights confirming our stratifications reflected degree of adiposity (Fig. 1G). Final body weights were also correlated with fasting glucose concentration (r 2 = 0.2218), suggesting mice that were prone to diet-induced obesity were also prone to its downstream consequences (Fig.   1H). Lastly, in light of the appreciation that low-grade intestinal inflammation can promote adiposity and its consequences, we measured weight/length ratio of the colon 55,56 . This measurement was also correlated with final body weight (r 2 = 0.2781; Fig. 1I), supporting the notion that the obese mice were in a state of low-grade gut inflammation.

Associations of inflammatory markers/mediators and proneness to obesity.
Low-grade inflammation is reported to associate with, and promote obesity 56,57 .
Accordingly, we investigated levels of pro-inflammatory mediators to determine if they might mark mice that would be prone to becoming obese following exposure to an obesogenic diet. Hence, we measured levels of fecal lipocalin-2 (Lcn-2), which is a broadly dynamic marker of gut inflammation 37 . Levels of fecal Lcn-2 did not correlate by guest on January 10, 2021 with final body weights when measured 14 days prior (r 2 = 0.0156) to exposure or 4 weeks after the initiation of the diet (r 2 = 0.0074; Fig. 2A, B). Additionally, the levels of serum pro-inflammatory cytokines CXCL1 and IL-6 when measured 7 days prior to administration of the obesogenic diet were also not correlated to final body weight (r 2 = 0.0177, 0.022 respectively, Fig. 2C, D). However, 4 of the 5 mice that displayed detectable serum IL-6 at this time point were in the top tertile of obesity following the diet suggesting the subset of mice displaying this parameter might be more prone to DIO. To further investigate this subset of mice, we tested for differences within various parameters associated with diet-induced obesity between the subsets of mice with or without detectable IL-6. At day -7, several of these parameters were consistent with the possibility that detection of IL-6 can discriminate high or low responders, but ultimately, none reached statistical significance (Fig. S1). Nevertheless, such elevations in IL-6 were not maintained when assayed after 8-weeks of diet (Fig. 2E). Other findings supporting the notion that obesity is associated with low-grade inflammation included a modest correlation after 8-weeks of diet between body weights and CXCL1, which is a chemokine expressed by many cell types and often used as a general serum marker of low-grade inflammation (r 2 = 0.0352, Fig. 2F). In contrast, there was no correlation between final body weights and levels of intestinal myeloperoxidase (MPO), which is a widely used marker of classic inflammation in the intestine 58 (Fig. 2G).
Gut bacterial components, flagellin and lipopolysaccharide (LPS), are well known for their inflammatory properties 59,60 . Fecal levels of each were measured 3 weeks prior to administration of the obesogenic diet with neither correlating with final body weight ( Fig. 2H, I). However, flagellin, but not LPS, were modestly correlated with final body by guest on January 10, 2021 weight when measured 4 weeks following initiation of the obesogenic diet (r 2 = 0.0984, 0.0151 respectively, Fig. 2J, K). Moreover, there was a correlation in levels of antiflagellin antibodies at the time of diet administration (for IgG but not IgA) and 8 weeks following exposure to obesogenic diet (Fig. 2L-O). Levels of anti-flagellin antibodies likely reflect exposure of the immune system to this molecule, which can be influenced by both levels of flagellin in the gut, bacterial-epithelial distance, and intestinal permeability 35,61 . Together, these studies did not reveal a reliable predictive marker of proneness to diet-induced obesity but suggest exposure to bacterial products, such as flagellin, might have some predictive power.

Quantitative measures of behavior did not predict obesity proneness
The gut-brain axis is increasingly appreciated to play a role in the pathogenesis of many neurological and metabolic diseases 62 . Hence, we investigated the extent to which certain behavioral parameters are able to predict proneness to weight gain.
Compulsive behavior and activity level were measured in a home cage behavior test, and time spent digging (Fig. S2A), time spent grooming (Fig. S2B), and total distance travelled ( Fig. S2C) were quantified. Additionally, anxiety-like behavior was assessed using the open field test, represented by time spent in the center zone (Fig. S2D) and distance travelled in the center (Fig. S2E) of the open field arena. Ultimately, none of these measures had a significant ability to predict extent of obesity in response to the obesogenic diet.

Impact of DIO on fecal metaproteome.
by guest on January 10, 2021 We next turned to a contemporary metaproteomics approach to study the fecal protein composition of our cohort of mice. While administration of an obesogenic diet is well known to rapidly alter gut microbiota species composition 63 , whether it might also impact the fecal metaproteome, let alone whether the fecal metaproteome might predict responsiveness to such a diet, has not been described. While metaproteomic analysis presents the challenges of discriminating host and bacterial proteins from potentially millions of proteins, the field is an area of rapid growth currently developing standard methodology 64 .
To this end, we applied our recently developed TMT-based metaproteomic methods 52 , in combination with a two-step database search method 29 on feces from mice that developed the highest and lowest degree of obesity (n= 4 mice per condition).
Our analysis included specimens collected before (day 0) and after 56 days of exposure to the obesogenic diet. The final data included quantitation of 13,975 proteins of which 1,108 were derived from mice (Table S1).
For a broad scale perspective of the data, an unsupervised Principle Coordinates analysis of the metaproteome data using the Bray-Curtis distance metric exhibited clear separation of samples before and after diet administration, reflecting a dramatic impact of the obesogenic diet on the overall fecal metaproteome (Fig. 3A). This analysis also exhibited clustering at the 56 day time point discriminating high and low response to diet (Permanova pseudo-F = 1.99, p = 0.058). Using K-means clustering, we identified 6 protein clusters, some of which were associated with increased representation of particular taxa and functions (Fig. 3B). These groupings include Group 4, which appeared to show an increased presence of Clostridiales and lipid transport and by guest on January 10, 2021 metabolism proteins in high responder mice after exposure to HFD (Fig. 3B). These groupings provide putative taxonomic associations to the functional differences observed before and after administration of HFD.
Comparing all samples before and after HFD exposure made evident that there were widespread changes in the fecal metaproteome. By using a statistical ranking method accounting for both fold change and t-test p-values, π-score 51 , we observed that 58% (3670/6311) of proteins displayed a high level of association to diet exposure (|π| > 1, Fig. 3C). The proteins associated with the dietary intervention contained large differences in their taxonomic and functional annotations. Taxonomic differences included a larger portion of proteins from Clostridiales and Bacteroidales before HFD exposure while a large portion (~40%) of proteins enriched after 8-weeks exposure were derived from Lactobacillales (Fig. 3D).
Functional categorization of the proteins associated to the dietary intervention was performed through the Evolutionary Genealogy of Genes: Non-supervised Orthologous Groups (eggNOG) database. These studies revealed a very strong association between proteins related to motility and HFD exposure as expression levels of 141 proteins were reduced following exposure to HFD with only 3 proteins increased after HFD, resulting in a 32-fold difference (Fig. 3E). In accord, we note that, on average, levels of fecal flagellin decreased by about 5-fold when measured 21 days preceding or 4 weeks following administration of the obesogenic diet (Fig. 2H, J).
Further, when subsetting all 680 flagellin proteins from the metaproteome dataset, we observed statistically significant decreases in abundance for both high and low responding mice (p < 0.0001; Fig. S3A). Functional assessment of the proteins by guest on January 10, 2021 enriched after HFD exposure resulted in weaker associations, the strongest of which was a 1.5-fold increased representation of Transcription proteins (Fig. 3E).

Functional and taxonomic characterization of fecal metaproteome in low-and
high-responder mice fed the obesogenic diet.
We next focused the analysis on discovering patterns in the fecal metaproteome that might have preceded or accompanied degree of responsiveness to the obesogenic diet. Toward this end, we examined the broad-scale functional composition of each sample's metaproteome through the eggNOG database. This revealed only modest variance amongst the samples (Fig. 4A). In contrast, viewing the composition of taxonomic orders in this manner revealed differences, both preceding and following diet exposure, that associated with a high-and low-response to the diet (Fig. 4B). There were 424 highly ranked proteins (|π| > 1) differentiating high and low responders at the initial time point (Fig. 4C). These proteins had large differences in their taxonomic origins with all proteins distinguishing the low responders belonging to Clostridiales while high responders had over 50% of proteins derived from Bacteroidales and Lactobacillales (Fig. 4D). Functionally, the proteins distinguishing high responders had a 14-fold enrichment in "Posttranslational modification, protein turnover, and chaperones" proteins, and a 5.6-fold enrichment in "Cell motility" proteins (Fig. 4E).
Many of the posttranslational modification, protein turnover, and chaperone proteins with the largest differences between high and low responders were chaperone proteins (Fig. S3B), a potential indication of a microbial stress response. The increased representation of cell motility proteins was a result of a subset of flagella, mostly derived by guest on January 10, 2021 from the order Clostridiales, that were significantly increased in high responders (p < 0.0001; Fig. S3C).
We also noted unique sets of carbohydrate metabolism and transport proteins differentially expressed at the initial time point. High responders had increased expression of Bacteroidale metabolism proteins including isomerases, kinases and aldolases, while Clostridiales uniquely had an increased expression of sugar transporters within low responders (Fig. S3D). This could be an indication of unique energy harvesting capacities in the microbiome present before the onset of HFD treatment 12 .
In the samples collected following administration of the obesogenic diet, there were 970 proteins distinguishing high and low responders (Fig. 4F). In contrast to proteins discriminating responses at the onset of HFD exposure, the proteins corresponding to high response were derived entirely from Clostridiales, while nearly 60% of proteins in low responders were derived from Lactobacillales (Fig. 4G).
Changes following the obesogenic diet associated with a high response to the diet included an 11-fold increased representation of lipid transport and metabolism proteins ( Fig. 4H). It is possible that the increase in lipid metabolism proteins from Clostridiales mediates more efficient harvesting of energy from lipids in high responding mice.

Analysis of fecal mouse proteins.
In addition to analyzing microbial proteins from the metaproteome, we next subset the data to determine associations within the fecal mouse proteome. In total, 699 host proteins were quantified within all samples and therefore included in the statistical by guest on January 10, 2021 analyses. A large portion (77%) of mouse proteins were strongly influenced by HFD exposure, and 92% of those proteins were increased after HFD (Fig. 5a). Using DAVID functional enrichment 53 , we identified significant (Bonferroni adjusted p-values < 0.05) enrichment of digestion and myosin proteins before HFD treatment and mitochondrial proteins after HFD (Table S2). Notably, myosin proteins occupied 6 of the top 10 proteins associated with the initial samples (Fig. 5b). Phosphorylated myosin light chains have previously been linked to intestinal permeability after HFD exposure 5 . Thus, our observed decrease in heavy chain myosin proteins may be related to changes in intestinal permeability. In regards to the increase of mitochondrial proteins, it was shown that HFD results in mitochondrial dysfunction 65 , and our data likely reflects this phenomena.
We next looked for mouse fecal proteins that might discriminate high and low responders prior to HFD administration. Here, 109 mouse proteins strongly differed between the groups, all of which were enriched in the DIO prone mice (Fig. 5c). Of these proteins, 40 (37%) were related to immunoglobulin. This strong enrichment for immunoglobulin genes was confirmed through DAVID, which showed a significant, 3fold enrichment (Bonferroni p-value = 7.0E-11) for Immunoglobulin V-set proteins (Table S2). This enrichment of immunoglobulin variable domains was also illustrated in the top 20 proteins associated with a heightened response to HFD (Fig. 5d). These findings further illustrate the link between low-grade inflammation and DIO, as this is a potential indication of increased immune activity in high responder mice, before administration of HFD.
by guest on January 10, 2021 Applying the same analysis to samples collected after 8-weeks exposure to HFD also revealed interesting insight into proneness to HFD exposure. After the dietary intervention, 63 mouse proteins were highly ranked in their ability to discriminate between high and low responders. All but two of these proteins were enriched within the low responders (Fig. 5e). Functional analysis showed a significant 6-fold enrichment (Bonferroni p-value = 4.0E-8) for keratin within the proteins enriched within low responders (Table S2). This increase of keratin could be an indication of greater colonic stress 66 in low responders at the final time point. High responders at this time had several immunoglobulin proteins within the top discriminatory proteins, though most were only modestly associated (Fig. 5f). Of note, many of the immunoglobulin proteins were among the strongest discriminatory proteins in high responders at both the initial and final day.

Analysis of microbiota composition vs. proneness and severity to DIO.
Lastly, we examined the potential of fecal microbiota composition, as analyzed by 16S rRNA gene sequencing to identify and/or reflect proneness to DIO. Visualization of fecal microbiota composition of all 50 mice at all time points by unweighted UniFrac revealed the expected dramatic difference in microbiota composition before and following administration of the obesogenic diet (p = 0.001; Fig. 6A). This analysis also showed clear, but much more modest differences between the 5 and 8-week postdietary change time points (Fig. 6A). In contrast, using this approach to examine differences in beta diversity did not identify differences in microbiota composition in high or low-responders either prior to (p = 0.977; Fig. 6B), or following administration of the by guest on January 10, 2021 obesogenic diet (p = 0.323; Fig. 6C). Rather, in accord with other diets, 8-week administration of the diet, which provided mice an additional 8 weeks to share their microbiota with their cage-mates, we observed strong cage clustering of microbiota compositions (p = 0.001; Fig. 6D). Nonetheless, levels of alpha-diversity, prior to administration of the obesogenic diet were moderately but significantly correlated (r 2 = 0.0873, p = 0.0394) with final body weights (Fig. 6E) suggesting that microbial community structure had some ability to predict proneness to DIO. An analogous but not significant trend was observed 8-weeks post-dietary change (Fig. 6F).
That overall assessment of microbiota composition lacked ability to identify highand low-responder mice does not preclude the possibility that select OTUs might provide such power. Hence, we selected specific OTUs whose abundance was enriched or depleted at time 0 in the mice that developed the greatest degree of obesity in response to HFD. This yielded an array of bacterial groups, those of which had the ten lowest p values represented here, Fig. S4. However, determining whether these differences are reproducible and/or biologically significant will require further experimentation. Additionally, we examined the ability of bacterial candidates generated by the metaproteomic analysis to predict (day 0), or reflect (day 56), proneness to the obesogenic diet. Of the 12 taxa analyzed in the day 0, 3 groups showed correlations predicted by the proteomic analysis with p values lower than 0.1 (Fig. 6G-I) while 9 did not ( Fig. S5 A-I). Regarding the taxa proteomic analysis identified as correlating with extent of obesity in the day 56 samples, 2 taxa correlated as determined with p-values lower than 0.1 (Fig. 6J-K) while the 10 others analyzed did not meet this criteria (Fig.   S6 A-J). Thus, overall, while development of approaches to predict proneness to by guest on January 10, 2021 obesity via analysis of the fecal proteome and/or microbiome remains a work in progress, these findings support its potential to contribute to such prognostications.

DISCUSSION
The goal of this study was to improve understanding of non-genetic determinants of DIO, focusing on parameters that might be impacted by gut microbiota, which is known to play a role in dictating severity of DIO. As obesity is promoted by low-grade, systemic inflammation, which can be driven by exposure to microbiota products 5  To find microbial factors that may correlate to DIO, we next turned to an untargeted metaproteomic approach. Our results confirm prior research showing large shifts in the overall structure of the metaproteome after administering HFD 50 . These broad shifts seem to be driven by proteins derived from Clostridiales and Bacteroidales, which decreased upon exposure to HFD, while the composition of Lactobacillales proteins expands. Interestingly, this difference was not observable in our analysis of the microbiota by 16S sequencing. One possible explanation is the known discrepency between genomic and proteomic technologies 52 which is supported by the notion that differences in protein abundance are not directly associated with species composition due to complex regulatory processes. However, other DNA sequence-based studies have also shown significant alterations in Clostridiales and Bacteroidales upon exposure to HFD 5, 69 , further suggesting a role for these taxonomies in HFD response.       After sub-setting the mouse derived proteins from the metaproteome data, differentially expressed proteins were determined using a statistical cut-off of |π-score| > 1. Volcano plots are shown demonstrating the log2 fold change (x-axis) and log10 p-value (y-axis) for (A) differences between final and initial samples, (C) differences between high and low responders and the initial time point, and (E) differences between high and low responders at the final time point. The π-score of the most significant proteins from each analysis are shown below each volcano plot in bar plots (B,D,F).