Multi-omic analysis reveals the anti-aging impact of sulforaphane on the microbiome and metabolome

Dietary factors may modulate many complex interactions between the microbiome, metabolome, and immune system and can have an impact on the functional status of older adults. Sulforaphane (SFN), a natural compound and Nrf2-related activator of cytoprotective genes, provides a wide range of biological effects from cancer prevention to reducing insulin resistance. We have shown that SFN increased survival and improved cardiac and skeletal muscle function in a mouse model of aging. This study aims to investigate the anti-aging effects of SFN on the gut microbiome and metabolome. Young (6-8 weeks of age) and old (21-22 months of age) male C57BL/6J mice were provided regular rodent chow or chow containing SFN for 2 months. Fecal samples were collected right before and at the completion of SFN administration. We profiled the gut microbiome and applied global metabolomic profiling to fecal samples. Multi-omics datasets were analyzed individually and integrated to investigate the relationship between SFN diet, the microbiome, and metabolome. Microbial diversity, composition and functional capacity varied substantially across different age groups. On a global level, in old mice we observed that the SFN diet restored the gut microbiome to mimic that in young mice. In old mice, the SFN diet enriched bacteria associated with an improved intestinal barrier function and the production of anti-inflammatory compounds. In addition, the tricarboxylic acid cycle, central in cellular respiration, was decreased and amino acid metabolism-related pathways were increased. SFN diet induced metabolite biomarkers in old mice that are associated majorly with the genera, Oscillospira , Ruminococcus , and Allobaculum . In mice, directed the metabolic potential to that of animals. Integrated microbiome and metabolome analyses revealed metabolite biomarkers that could be modulated by bacteria and contribute to the anti-aging effects of SFN. Collectively, our results provide evidence in support of a novel hypothesis that SFN exerts anti-aging effects by influencing the gut microbiome and metabolome. Although further investigations to identify precise mechanisms, modulating the gut by SFN may have the potential to promote healthier aging.


Background
The gut microbiome is an extremely diverse and complex ecosystem of bacteria, viruses, and fungi inhabiting the intestinal tract, which interact with each other and their host [1]. The human gut microbiome is unique to each individual, but is dominated by anaerobic bacteria belonging to two phyla, Firmicutes and Bacteroidetes [2]. It is becoming apparent that there are many complex interactions between the gut microbiome, immune system, and metabolism in health and disease [3,4]. Accumulating data support the contention that microbial metabolites play major roles in the regulation of the immune system [5]. Unsurprisingly, diet is a major modifiable factor shaping the gut microbiome structure and metabolic activity both in the short and long terms as reported in observational and intervention studies [6,7]. Given these associations, dietary factors may have significant therapeutic utility in modulating many interactions between the gut microbiome, metabolism, and immune system.
The human gut microbiome fluctuates over the individual's life span, undergoing the most prominent deviations during infancy and old age [8]. Interestingly, our immune health shows similar fluctuation patterns like the gut microbiome and the most unstable state during infancy and old age [8]. Although the causal relationship between the gut microbiome and aging is unclear, the gut microbiome may serve as a target for anti-aging intervention. We hypothesize that diet can revert the gut microbiome to a younger composition. In this study, we tested this hypothesis with a sulforaphane (SFN) containing diet in a mouse model of aging. SFN is a compound in cruciferous vegetables that has been investigated for its anti-cancer, anti-aging, antioxidant, antimicrobial, anti-inflammatory, anti-diabetic, and neuroprotective properties [9,10]. Many studies have investigated SFN for its role in aging. SFN mechanisms of action have been linked to Nrf2 and NF-κB, key cellular transcription factors. SFN-induced activation of Nrf2 and inhibition of NF-κB result in the induction of redox-modulating genes and the inhibition of inflammation, respectively [11]. However, the effect of SFN on the gut microbiome has not yet been investigated in aging.
We hypothesize that administration of SFN reshapes the old gut microbiome into a younger composition. This may lead to an enrichment of beneficial bacteria that produce short chain fatty acids (SCFAs), known for their anti-inflammatory actions that may help improve aging-related pathologies. In this study, we focus on the anti-aging impact of SFN on the gut microbiota in a mouse model. To disentangle the contribution of aging and non-aging related variables to the gut microbiome, we divided mice into four groups: an old control group receiving regular rodent chow, an old group receiving the same chow but supplemented with SFN, a young control group, and a young SFN diet group ( Table 1). Fecal samples were collected immediately before (indicated in text and figures as 0 months) and 2 months after the start of SFN administration. We conducted multi-omic profiling to investigate the relationship between SFN diet, the gut microbiome, and metabolome. Our results identified global relationships and highlighted novel associations between SFN diet and microbiome structure and function, and metabolite biomarkers that may be modulated by gut bacteria.
measures. SFN-induced beta diversity separation of old group was not obvious with unweighted-UniFrac distance ( Figure 2D). Nonetheless, differences between OS2 and YC2 with all diversity measures considered, and the difference between OS2 and YS2 with weighted UniFrac were statistically significant (Table S2). Beta diversity analysis implies again that SFN diet alters and may restore the young gut microbiome in old mice.
Taxonomic profiling. To further uncover microbial composition characteristics in different groups, we analyzed ASVs assigned for phylum, class and genus in Figure 3. At the phylum level ( Figure  3A), Firmicutes, Bacteroides and Verrucomicrobia were the most predominant phylum groups. At the class level ( Figure 3B), Clostridia (belonging to Firmicutes at the phylum level), Erysipelotrichia (Firmicutes), Verrucomicrobiae (Verrucomicrobia), Bacteroidia (Bacteroidetes) and Bacilli (Firmicutes) were the most predominant classes. At the genus level ( Figure 3C), Akkermansia (belonging to Verrucomicrobiae at the class level), Allobaculum (Erysipelotrichi), Bacteroides (Bacteroidia), Lactobacillus (Bacilli), Odoribacter (Bacteroidia), Oscillospira (Clostridia), Parabacteroides (Bacteroidetes), Prevotella (Bacteroidetes), Ruminococcus (Clostridia) and Turicibacter (Erysipelotrichi) were predominant genera. The abundances of these microbial taxa changed significantly across groups. At the phylum level, SFN diet increased Firmicutes and decreased Bacteroidetes in old mice. This effect was not noticeable in young mice as compared to old mice. Within the phylum Firmicutes, two classes, Erysipelotrichi and Clostrida showed the opposite trend in the SFN diet group. Erysipelotrichi occupied a large portion of the gut microbiota in stool sample equally in OS2, YC2 and YS2, whereas Clostrida was decreased in these groups. Bacteroidia belonging to the phylum Bacteroidetes accounted for a smaller proportion of gut microbiota in OS2 compared to all other groups, which might be due to SFN administration in old mice. We observed that the YC2 and YS2 groups contained many more sequences annotated at the genus which might be a characteristic of the gut microbiome in young mice. Interestingly, OS2 showed a similar number of sequences annotated at the genus to young mice, unlike OC2. Within the class Erysipelotrichi, Allobaculum was greatly increased in OS2 compared to OC2 ( Figure 3C). But, Allobaculum was also increased equally in young control and SFN groups at 2 months. Therefore, the genus Allobaculum might be a characteristic of the gut microbiome in young mice, and the SFN diet restored Allobaculum in the old gut microbiome. Akkermansia decreased and Oscillospira was increased in abundance in old and young mice due to the SFN diet. It is noted that young groups at the baseline (YC0 and YS0) somehow showed significant differences in abundance of Prevotella and Bacteroides (genera of the phylum Bacteroidetes). Prevotella are beneficial bacteria that degrade dietary fiber into short chain fatty acids (SCFAs). A high Prevotella-to-Bacteroides ratio is associated with a loss of body weight and body fat [14]. This evidence might support the contention that there is no unique optimal young gut microbiota composition. Based on all the evidence shown above, we conclude that SFN diet restored the young gut microbiome in old mice by increasing the abundance of specific taxa. This recovery of the microbiome may contribute to reversing functional and physiological effects of aging.
Linear discriminant effect size analysis (LEfSe). To identify distinctive features between groups, a linear discriminant analysis (LDA) effect size (LEfSe) analysis [15] was performed (Figure 4). LEfSe analysis showed that the genus Allobaculum, the order Erysipelotrichales, the class Erysipelotrichi and the phylum Firmicutes were significantly enriched in OS2 compared to OC2. In addition, LEfSe analysis also showed that the genus Adlercreutzia, the order Coriobacteriales, the class Coriobacteriia, and the phylum Actinobacteria were significantly higher in OS2 compared to OC2. However, the genus Adlercreutzia is from the class Actinobacteria and its abundance is negligible compared to other classes and genera in Figure 3B and Figure 4B. None of the significantly enriched features in OS2 compared to OC2 were identified as enriched features in YC2 and YS2 compared to OS2 (Figure S1). On the other hand, LEfSe analysis identified the family S24-7, the order Bacteroidales, the class Bacteroidia, and the phylum Bacteroidetes as enriched in OC2 compared to OS2.
Potential functional annotations of gut microbiota in two age groups with and without SFN diet. We used the PICRUSt2 software [16] for stringent prediction of metagenomes and functional metabolic pathways from the 16S survey data, yielding output as MetaCyc [17] pathway abundances. PCoA analysis was executed with predicted functional profiles by KEGG Orthology [18] (Figure 5A) and metabolic pathways based on the MetaCyc database ( Figure 5B). The analyses included only the microbiome data collected at 2 months after initiation of diet treatment. The first two principal components of PCoA explain 43% of variance with the KEGG Orthology profile and 51% of variance with the metabolic pathway profile. Interestingly, both profiles showed an orthogonal relationship between OS2 and OC2 ( Figure 5). Especially, the MetaCyc pathway profile showed that three groups (OS2, YC2 and YS2) were horizontally positioned along the first principal component whereas OC2 was positioned along the second principal component ( Figure  5B). Figure 5C represents the LEfSe result which displays the effect size of each differentially abundant MetaCyc metabolic pathway between OC2 and OS2 with LDA score cutoff of 3. Comparing OC2 to OS2, we identified 43 significant metabolic pathways among a total of 323. Prominent observations are as follows: amino acid metabolism related pathways, L-lysine biosynthesis II, superpathway of L-phenylalanine biosynthesis, superpathway of L-tyrosine biosynthesis, superpathway of L-lysine L-threonine and L-methionine biosynthesis I, superpathway of L-methionine biosynthesis (transsulfuration), superpathway of S-adenosyl-Lmethionine biosynthesis were significantly enriched in OS2 compared to OC2. Two central metabolism related pathways, incomplete reductive tricarboxylic acid (TCA) cycle and pyruvate fermentation to propanoate I were enriched in OC2 compared to OS2.
Microbiome and metabolome data integration analysis reveals microbiome dependent metabolic changes. We explored two orthogonal approaches, MIMOSA [19] and MelonnPan [20] which are fundamentally different. MIMOSA uses a metabolic model framework that integrates metabolic potential from bacterial genomes and metabolome composition into a unified analysis. Community-wide metabolic potential (CMP) scores are predicted based on the microbial metabolic genes for each metabolite and sample, and then are compared to the actual metabolome data obtained experimentally. MelonnPan calculates the correlation between the microbiome and metabolome composition to identify bacterial populations that might be responsible for metabolite patterns using an elastic regularization technique. Our focus was to identify well-predicted metabolites by changes in the microbial community due to the SFN diet in old mice, which resulted in 565 metabolites among 4158 putative metabolites in negative ion mode and 192 metabolites among 4197 putative metabolites in positive ion mode. First, to perform MIMOSA, we mapped the corresponding metabolite names to KEGG identifiers by mapping the compound IDs to the Human Metabolome Database (HMDB), which includes crossreferences to KEGG compound identifiers. Only 10 metabolites out of 565 significant metabolites in negative ion mode and 6 metabolites out of 192 significant metabolites in positive ion mode were identified with unique KEGG compounds. With log-transformed metabolite values, MIMOSA identified only one metabolite, Phosphatidylcholine. Figure S2A represents the relationship between the CMP score estimated by a rank-based method and the experimentally measured metabolite value for each sample. Of the metabolite values, 23% (model R-squared of 0.233) was explained by taxa that can contribute to phosphatidylcholines. For phosphatidylcholine, MIMOSA identified 29 contributing taxa, among which 6 taxa were selected based on an abundance cutoff of 10, sample frequency of 20% and contribution of 5% to model R-squared ( Figure S2B). MIMOSA identified the genus Parabacteroides as the major contributing taxon for phosphatidylcholine, followed by the family S24-7.
To run MelonnPan, we limited our analysis to OTUs which are prevalent in at least 20% of the samples, with relative abundance larger than the mean of abundances of OTUs (=4.209726) and with a variance larger than the mean of variances of OTUs (=63.0079). This resulted in the analysis of 186 OTUs out of a total of 929 OTUs. We normalized the abundance of the selected OTUs and log-transformed metabolite values into relative abundance. MelonnPan fitted a permetabolite elastic net model to the normalized data and determined the tuning parameters (the elastic net mixing parameter and sparsity parameter) in the elastic net model using 10-fold cross validation, which led to an optimal subset of OTUs whose abundances predict a given metabolite. Metabolites with a Spearman correlation coefficient between experimentally measured metabolite values and predicted metabolite values across samples > 0.3 were defined as well-predicted by MelonnPan. For each well-predicted metabolite, we kept only OTUs whose coefficients in the fitted elastic net model were ≥ 0.0002. Phosphatidylcholine by MIMOSA, was also predicted by MelonnPan (data not shown because coefficients in the fitted elastic net model < 0.0002). Among well-predicted metabolites, 61 metabolites had coefficients of organisms larger than 0.0002 in the fitted elastic net models. We kept only organisms which are annotated at the genus level where the same genus occurred more than once for a given metabolite, then we kept the one with the largest absolute coefficient. Figure 6 represents coefficients of organisms for well-predicted metabolites which have no known involvement in drug and plant metabolism. The genus Oscillospira was the most associated with well-predicted metabolites. The genera which were also majorly associated with well-predicted metabolites were Ruminococcus gnavus, Allobaculum, and Akkermansia muciniphila. The genera Anaeroplasma, Lactococcus, Mucispirillum, Schaedleri,Pasabacteroides, Rc4-4, Smb53, and Turicibacter were not associated with any of these endogenous metabolites (Figure 6), nor were they found in any abundance in OS2 ( Figure  3C).

Discussion
To our knowledge, this study is the first attempt to investigate the impact of SFN on structural and functional changes in the gut microbiota and metabolome in an animal model of aging.
Age-dependent microbial signatures of the mouse gut microbiome. First, in comparing OC2 to YC2, we determined that the gut microbiome of old mice was composed of more bacteria which were not annotated at the genus level, showing much lower taxonomic resolution compared to the microbiota in young mice. Gut microbial alpha diversity, a holistic estimator generally decreases when people age, likely due to changes in physiology, diet, medication, and lifestyles. Decreased diversity is considered an indicator of an unhealthy microbiome and has been linked to different chronic conditions such as obesity and type 2 diabetes [21]. A decline in gut microbiota diversity was also observed in our dataset ( Figure 1). As in the human gut microbiome, Firmicutes and Bacteroidetes were predominating bacterial phyla in our data. We observed a decrease of the ratio of Firmicutes to Bacteroides with age ( Figure 3A), which is also seen in aging humans [22] and is known to be associated with several conditions such as weight gain, obesity, insulin resistance, gut permeability, inflammatory bowel disease, and depression [23,24]. At the class level, gram-negative Bacteroidia (Bacteroidetes phylum) and gram-positive Clostridia, Bacilli, and Erysipelotrichi (all from the Firmicutes phylum) were dominant classes and altered in old mice. Especially, a decreased abundance of Erysipelotrichi in the old gut microbiome was the most visible. At the genus level, beneficial bacteria Allobaculum, Lactobacillus, Odoribacter, and Ruminococcus were decreased in the old mouse gut microbiome. A study of calorie restriction [25] found that Allobaculum and Lactobacillus made a significant contribution to the decrease of Firmicutes during aging. At the genus level, we found SCFAs-producing taxa, Allobaculum, Lactobacillus, Odoribacter, and Ruminococcus enriched in YC2 compared to OC2. Among them, Allobaculum (from the Erysipelotrichi class) is a SCFA-producing genus in the gut which has recently been reported to be an important functional phylotype [26,27] and especially, protects intestinal barrier function by producing SCFAs. Ruminococcus is related to polysaccharide fermentation into SCFAs and bile acid dihydroxylation. The members of the genus Allobaculum are already known to be inversely correlated with dietary-induced inflammation markers [28].
Odoribacter is a butyrate producer that belongs to the phylum Bacteroidetes [29][30][31][32]. The decrease of these SCFA producers in old-control mice might suggest an increase in the possibility of inflammaging, i.e. increased chronic, low-grade inflammation. The beneficial microbe, Akkermansia (Phylum Verrucomicrobia) produces both propionate and acetate [33,34] and is inversely correlated with several disease states. On the other hand, certain members of this genus (Akkermansia muciniphila) may also exacerbate infection via mucin degradation. Contrary to our expectations, Akkermansia was enriched in old mice ( Figure 3C). SFN diet-dependent microbial signatures in the mouse gut microbiome. Gut microbiota have emerged as an attractive therapeutic target to promote healthy aging and anti-aging effects which may result in improved quality of life [35]. The gut microbiota respond to dietary interventions very quickly, and short-term changes in the diet can alter the overall structure of the gut microbiota. We have recently shown that an SFN-containing diet reduces the effects of aging on cardiac and skeletal muscle function (unpublished data). In this study, we investigated the impact of SFN diet on the mouse gut microbiota, which may contribute to these anti-aging effects. SFN diet increased the diversity of the gut microbiota in old mice (Figure 1). Further, SFN diet encouragingly reshaped the gut microbial community structure of old mice to where they approach the one of young mice for all four beta diversity measures (Figure 2). Together, these results demonstrate the possibility of anti-aging effects of SFN on the gut microbiota.
At the phylum level, the SFN-diet enriched Firmicutes depleting Bacteroides in old mice which are known to be associated with a healthy gut microbiome. The ratio of Firmicutes to Bacteroidetes in the human gut microbiota changes with age, undergoing a significant increase from birth to adults and a significant decrease from adults to elderly, to where there is no significant difference between infants and elderly [36]. In our study, the ratio of Firmicutes to Bacteroidetes in young mice fed SFN was similar to the one in young control mice. In old mice, on the other hand, SFN diet increased the ratio of Firmicutes to Bacteroidetes. Our study pointed to the class Erysipelotrichi and the genus Allobaculum as major contributors to the SFN-induced enrichment of the phylum Firmicutes in the old gut microbiota. Increased numbers of Erysipelotrichi are associated with a phenotype of impenetrable mucus layer in the mouse [37]. Moreover, exercise increases butyrate-producing Erysipelotrichaceae along with several taxa in animals and adults regardless of diet [38]. Allobaculum, one of the most prevalent genera in young gut microbiota (Figure 3C), showed a significant change in abundance with SFN administration in the old mice (Figure 4). Allobaculum is a genus of gram-positive, non-spore-forming bacteria, strictly anaerobic and non-motile with tryptophan-catabolizing functions and suggested to be the most active glucose utilizer [39]. Allobaculum has been suggested to be beneficial for host physiology, and its increase was associated with low-fat feeding compared with high-fat diet feeding in a mouse model [40]. Allobaculum was depleted in mice with age-related mitochondrial dysfunction [41]. The end products of Allobaculum fermentation are SCFAs, specifically butyrate which has been shown to have epigenetic consequences as a histone deacetylase inhibitor, possess antiinflammatory and anti-carcinogenic properties [42]. Therefore, increases in Allobaculum due to the SFN diet may contribute to an increased intestinal barrier function and reduced inflammation in old mice. Oscillospira, one of the most abundant genera, is an under-studied anaerobic bacterial genus from Clostridial cluster IV in the Firmicutes phylum. Oscillospira showed an increased abundance in OS2 compared to OC0 mice. A study of Oscillospira metabolism [43] inferred that Oscillospira species are butyrate producers and found that their abundance was decreased in inflammation-related diseases. Given the strong evidence linking Oscillospira to human leanness or body mass index in both children and adults [44][45][46][47], further studies are necessary to ascertain the clinical significance of Oscillospira in human health and aging.
SFN diet-dependent microbial functional signatures in the mouse gut microbiome. PCoA analysis based on Euclidian distance demonstrated that the old group fed the SFN diet showed metabolic pathway profiles that were more similar to the ones of young animals than the mice in the old control group (Figure 5). Therefore, the SFN diet could be reasonably presumed to shape the metabolic capability of the aging gut microbiome toward a younger gut microbiome. Noticeably, the TCA cycle, one of the essential functions of mitochondria [48] and a central hub for energy metabolism and macromolecule synthesis, and pyruvate fermentation to propanoate I were significantly enriched in the old control group. A comparison study of the human gut microbiota between centenarian, elderly and young individuals, based on the premise that centenarians are a model of healthy aging, identified these two metabolic pathways enriched in the centenarian group compared to the young and elderly groups [49]. Since pyruvate is the major precursor for the synthesis of three major SCFAs, acetate, propionate, and butyrate that are related to a strong gut epithelium, this study interpreted the enrichment of the TCA cycle and pyruvate fermentation to propanoate I as a signature of longevity and healthy microbiome. On the other hand, increased metabolism of SCFAs such as pyruvate, butanoate and propanoate in obese rodent models may provide an extra energy source and induce insulin resistance [50]. In the current study, metabolic pathways of amino acids such as L-lysine, L-tyrosine, L-methionine, L-threonine, and Lphenylalanine were significantly enriched in the old group on SFN diet compared to the old control group. Among those amino acids, lysine, tyrosine, methionine, and phenylalanine are suggested to play an important role in the regulation of aging [51]. Interestingly, L-lysine is an essential ketogenic amino acid and critical building block for proteins. An increase in L-lysine in the human gut was also observed in young compared to centenarian and elderly groups [49]. Our metabolic pathway analysis was done with inferred functional profiles. There are questions that still remain unanswered and may need further study. For example, does SFN diet activate all genes in TCA cycles or L-lysine biosynthesis II or is there a selectivity for particular ones? SFN diet-dependent microbiome-dependent metabolites. To identify microbe-associated metabolites, we applied metabolic model-based (MIMOSA) and multiple regression with penalty term-based (MelonnPan) methods for our paired data. Both methods identified phosphatidylcholine (C00157) as a metabolite associated with microbiome structure where phosphatidylcholine showed an increased fecal concentration in the old group fed the SFN diet compared to the old control group. Phosphatidylcholine is one of the major phospholipids comprising the cellular membrane and is known to have several health-promoting properties. Phosphatidylcholine showed a lifespan-extending effect under oxidative stress (one of the major causal factors of aging) and delayed age-related decline of motility in a worm model [52]. MIMOSA identified Parabacteroides and several members from the family S24-7 as major contributors to phosphatidylcholine. Note that Parabacteroides did not show a significant differential abundance between the SFN-treated and control old mice in the current study. Among 61 differentially abundant metabolites between old SFN-treated and old control mice as identified by MelonnPan, several are not related to drug or plant metabolism. LysoPE(18:4(6z,9z,12z,15z/0:0) and phaseolic acid showed an increased concentration with SFN treatment in old mice and are in positive association with Allobaculum, Oscillospira, Ruminococcus/Ruminoccus gnavus, where Allobaculum is one of the SFN diet-induced anti-aging microbial signatures. Phaseolic acid is a hydroxycinnamic acid. Hydroxycinnamic acids and their derivatives have been reported to possess properties of antioxidant, anti-inflammatory, antimicrobial, and anti-tyrosinase activities [53] which might suggest that phaseolic acid can be exploited as an anti-aging agent. Glycinamide ribonucleotide is an intermediate in de novo purine biosynthesis [54]. It showed an increased concentration in SFN-treated old mice and positive associations with Adlercreutzia, Allobaculum, and Oscillospira, where Adlercreutzia and Allobaculum were both SFN-diet induced anti-aging microbial signatures. Lastly, Glaucarubin concentration was increased with SFN treatment in old mice, and this metabolite was positively associated with Allobaculum. Glaucarubin is found in fats and oils. Zarse et al reported that Glaucarubin may promote metabolic health and lifespan in mammals and possibly humans [55]. Many studies reported protective effects of SFN against brain diseases [56]. Among the SFN-diet induced metabolite signatures, taurocyamine is known as an endogenous alkaline shifter which effectively reduces the extent of brain intracellular lactic acidosis brought about by anoxic insult [57], and is an inhibitor of taurine transport and a glycine receptor antagonist in the brain [58]. Taurocyamine is positively associated with Allobaculum and Ruminococcus. MelonnPan does not utilize prior information of genomic metabolic capability and solely depends on abundance of species/functions and metabolite concentration. In the analysis of integrated metabolome/microbiome data, it might be a good strategy to use MelonnPan first to identify a subset of metabolites within untargeted metabolomics data, followed by the use of targeted metabolomics data of those metabolites identified by MelonnPan to investigate the relationship between the microbiome and metabolome using MIMOSA.

Conclusions
The immune system is likely influenced by the gut microbiota, and their interaction plausibly contributes to the process of inflammation. Ongoing studies suggest that diet has an effect on both the gut microbiota and systemic inflammation, with an impact on functional status of older adults. Manipulating the intestinal microbiota may be beneficial for maintaining health and treating disease and potentially mitigating aging related effects on systemic metabolism. SFN is known to have a wide range of biological effects including anticancer, anti-inflammatory, antioxidant, and anti-aging [10]. Here, we have used a mouse model of aging that allowed us to control for confounding factors from human data and investigated the impact of diet supplemented with SFN on the structure and function of the gut microbiota and metabolome. We have observed structural and functional changes in the microbiome that correlate with age. We demonstrated that SFN diet restored the structure and function of young microbiota in old mice, restoring beneficial microorganisms associated with health. In particular, SFN diet enriched bacteria associated with improved intestinal barrier function and anti-inflammatory pathways. Moreover, inferred metagenome-based data analyses revealed that the SFN diet decreased abundance of metabolites of the TCA cycle and increased amino acid metabolism related pathways in old mice. Integrated microbiome and metabolome analysis revealed putative metabolite biomarkers of SFN-induced in old mice that could be modulated by the microbiome. Probiotics have exerted very modest effects on the microbiome in the mice [59,60] relative to that in humans, and even in humans the effectiveness of probiotics is controversial [61]. Consequently, the beneficial modification of the microbiome reported by SFN here should be investigated in humans as an alternative to probiotic use. Our observations support a novel finding that SFN diet exerts its antiaging effect by influencing the composition and function of the gut microbiota. Many observations we made in our study are consistent with earlier studies of the human gut microbiome and aging, suggesting there may be some parallel shifts that occur in aging human and mouse populations. However, translational benefits of SFN in the human gut microbiome remain to be demonstrated. In addition, further studies are needed to investigate potential gender differences in the effects of SFN diet on shaping the gut microbiota. In conclusion, according to our knowledge, our study represents the first effort to investigate the impact of SFN on the gut microbiome and metabolome during aging, providing a new prospective for potential targets for microbiota-targeted intervention.

Methods and Materials
Study design and sample collection. This study conformed to the Guide for the Care and Use of Laboratory Animals of the National Institutes of Health and the work was performed in accordance with a protocol (IACUC#646767-5) approved by the Central Arkansas Veterans Healthcare System Institutional Animal Care and Use Committee. Animals were housed in the Veterinary Medical Unit at the Central Arkansas Veterans Healthcare System in Little Rock, Arkansas. Young and old male C57BL/6 mice were obtained from aged rodent colonies of the National Institutes of Health (Bethesda, MD 20892). Young (6-8 weeks of age) and old (21-22 months of age) mice were fed TD 96163 diet (Teklad, Madison, WI) (control groups) or TD 96163 diet supplemented with sulforaphane (SFN) (442.5 mg per kg diet; treated groups), for 2 months. Immediately before and at the completion of the 2 months-administration of SFN and control diet administration, fecal pellets were collected from 5 mice per age and treatment. For this purpose, mice were individually placed in a Plexiglas box to obtain fresh fecal pellets. Fecal pellets were collected and stored at −80°C until processing. Each pellet was divided into two parts under liquid nitrogen. one half was shipped to the University of California Los Angeles for 16S rRNA amplicon sequencing and the other half to Georgetown University for metabolomics as described in Casero et al [62] and also described below.
Sample preparation for microbiome analysis. Sample preparation and analysis were performed as described by us previously [62]. Briefly, genomic DNA was extracted using the PowerSoil DNA Isolation Kit (MO BIO Laboratories, Carlsbad, CA, USA) as per manufacturer's instructions. 16S rRNA amplicon bacterial gene sequencing was performed using extracted genomic DNA as the template using a standard protocol. The PCR primers used in the study (F515/R806) targeted the V4 hypervariable region of the 16S rRNA gene, with the reverse primers including a 12-bp Golay barcode. Thermal cycling was performed in an MJ Research PTC-200 (Bio-Rad Inc., Hercules, CA, USA) with the following parameters: 94°C for 5 min; 35 cycles of 94°C for 20 s, 50°C for 20 s, and 72°C for 30 s; 72°C for 5 min. PCR products were purified using the MinElute 96 UF PCR Purification Kit (Qiagen, Valencia, CA, USA). DNA sequencing was performed using an Illumina HiSeq 2500 (Illumina, Inc., San Diego, CA, USA), in paired-ended mode. Clusters were created using template concentrations of 4 pM and PhiX at 65 K/mm 2 . Sequencing primers targeted 101 base pair reads of the 5′ end of the amplicons and 7 base pair barcode reads. Reads were filtered using the following parameters: minimum Q-score: 30, maximum number of consecutive low-quality base calls allowed before truncating: 3, and maximum number of N characters allowed: 0. All filtered V4 reads had a length of 150 bp.
Sample preparation for metabolomics analysis. Fecal samples were processed by initially homogenizing in extraction solvent cocktail containing 50% methanol, 30% isopropanol, 10% water and 10% chloroform and internal standards to allow a broad based extraction, followed by protein crash using 1:1 acetonitrile [63]. The samples were centrifuged and the supernatant was dried and resuspended in water containing 50% methanol for MS analysis. The samples were resolved using reverse phase chromatography using an Acquity UPLC (Waters Corporation, USA) system online Xevo-G2-QTOF-MS ( Waters Corporation USA) operating in positive and negative ion mode, the details of tune page parameters have been described before [64][65][66]. Several measures including randomizing the sample queue, use of pooled quality controls and standards were used to monitor data quality, retention time drifts and signal intensity.
Bioinformatic processing with 16S rRNA amplicon sequencing data. We first removed adapters (barcodes, forward primer:GTGYCAGCMGCCGCGGTAA, and reverse primer: GGACTACNVGGGTWTCTAAT) in demultiplexed paired-end sequencing data of 40 sequencing libraries from one batch using Cutadapt [67], which resulted in about 3% of reduction in sequencing depth. For microbiome bioinformatic processing, raw sequence data were then imported into QIIME2 [68]. Sequencing depth ranged from 67,885 to 283,004 with a mean of 199,196 and a median of 219,054. After inspecting the quality profiles, it was clear that the reverse read quality dropped off more severely than in the forward read. Accordingly, we trimmed the reverse reads at position 136 and removed chimeric reads by consensus method followed by denoising and merging with DADA2 [69] (via q2-dada2) which resulted in a table of 1,746 amplicon sequence variants (ASVs) for 40 samples with the total frequency of 3,341,551.
For taxonomic assignment, we used a classifier that has been pretrained on Greengenes database (v13_8) with 99% operational taxonomic units (OTUs) and primers used for amplification and the length of sequence reads. Taxonomy was assigned to ASVs using the q2feature-classifier [70] classify-sklearn naïve Bayes taxonomy classifier against the Greengenes (v13_8) 99% OTUs reference sequences [71]. ASVs were filtered out through the two steps: first, ASVs which belong to chloroplast at the class level or mitochondria at the family level were filtered out, and then only those ASVs with at least phylum-level assignment were only kept.; second, low abundance ASVs with frequency less than 0.0005% of reads in the total dataset were removed as recommended for Illumina amplicon data [72] which allowed us to perform differential abundance analysis with increased detection sensitivity. The two filtering steps resulted in a feature-table of 1,377 ASVs with the total frequency of 3,336,130 which was used for downstream analysis. All ASVs were aligned with Mafft [73] (via q2-alignment) and a phylogeny was constructed using Fasttree2 [74] (via q2-phylogeny). For alpha and beta diversity analysis, alphadiversity metrics (observed OTUs and Faith's Phylogenetic Diversity [75]) and beta diversity metrics (weighted UniFrac [76], unweighted UniFrac [77], Jaccard distance, and Bray-Curtis dissimilarity) were used. Sample ordination based on principle coordinate analysis (PCoA) were estimated using q2-diversity after samples were rarefied (subsampled without replacement) based on rarefaction analysis by phylogenetic diversity index and species richness to the number of sequences of the sample with the least number of sequences (24,788).
For alpha diversity significance test, the Kruskal-Wallis test was used to evaluate the effect of the experimental factors on the relative abundance at ASV level. Similarly, significant differences in beta diversity between groups were tested using permutational multivariate analysis of variance (PERMANOVA) with 999 Monte Carlo permutations. We generated an OTU-table by clustering all ASVs into OTUs using Vsearch [78] (via q2-Vsearch) with a similarity threshold of 97% against the Greengenes (v13_8). The number of OTUs detected in 40 samples was 1,047, and the mean frequency was 70,387. After removing rare OTUs with frequency less than 0.0005% of reads in the total dataset and then rarefying the OTU table at a depth of 21,221 sequences per sample based on rarefaction analysis by phylogenetic diversity index and species richness, the resulting OTU-table provides a highly replicated, deeply sequenced dataset with 928 OTUs and mean frequency of 70,310. For functional analysis of the microbiome, we predicted functional profiles from the rarefied OTU-table using PICRUSt2 v.2.1.4-b software [16]. With the maximum NSTI cut-off of 2, PICRUSt2 predicts abundance of Kyoto Encyclopedia of Genes and Genomes (KEGG) Orthology (KO) gene family and Enzyme Commission (EC) number based on the KEGG database [18]. Afterwards, we inferred the abundances of metabolic pathways based on MetaCyc database [17] from the predicted abundances of EC numbers. The linear discriminant analysis (LDA) effect size (LEfSe) method [15] was used to identify statistically significant differences between groups of experimental design in taxonomic and metabolic pathway features. For taxonomic feature, the ASV table was collapsed at the species level. LEfSe algorithm uses Kruskal-Wallis rank sum test to detect features with significantly different abundances between the two groups, followed by LDA to estimate the effect size of each feature. A strength of the LEfSe method compared with standard statistical approaches is that in addition to providing p values, it estimates the magnitude of the association between each feature and the grouping categories. A significance alpha level of 0.05 for the Kruskal-Wallis test and an effect size (LDA score) threshold of 3 were used for all biomarkers. The cladogram from the LEfSe method indicates the phylogenetic distribution representing differentially abundant taxonomic groups. The size of each node represents its relative abundance.
Computation frameworks to integrate microbiome and metabolome for the identification of potential mechanistic links. We explored two different approaches to integrate microbiome and metabolome for the identification of potential mechanistic links: MIMOSA [19] using the updated software (https://borenstein-lab.github.io/MIMOSA2shiny/) and MelonnPan [20]. MIMOSA integrates metabolic potential from bacterial genomes and metabolome into a unified analysis. PICRUSt2-predicted metagenome from a reference-based OTU table rarefied (described above) was used for the MIMOSA analysis. MIMOSA first performs metabolic network modelling using the Predicted Relative Metabolic Turnover framework [79] derived from KEGG enzymatic reactions [18]. For each metabolite in each sample, community-wide metabolic potential (CMP) scores were computed as the matrix multiplication of a stoichiometric enzyme reaction matrix (M) and PICRUSt2-predicted metagenome matrix represented by KO-relative abundances (G) so that CMP scores represent the relative capacity of the community in a given sample to generate or deplete each metabolite. The integration with metabolomics data was performed by comparing CMP scores to actual LC-MS normalized metabolite abundances, by matching metabolite putative ids with KEGG compound ids. Due to this algorithmic step, only putative metabolites which are annotated by KEGG compound id, can be examined for integration of microbiome and metabolome by MIMOSA. Finally, MIMOSA identifies contributions of taxa based on the amount of variation in a metabolite explained by CMP scores analyzing the regression model fit. For metabolites that are significantly associated with predicted metabolic potential with a regression model p-value less than 0.1, the taxa with the largest contributions are hypothesized to be the main drivers of change of metabolite pattern across samples. However, it is difficult to apply or validate MIMOSA results in a data-driven manner because of its algorithmic limitations of applicability to metabolites not annotated with KEGG compound. In addition, MIMOSA depends on accurate characterization and annotation of species-and even strainspecific metabolites so that it cannot scale well to complex communities with partially referenced taxa or metabolites. We applied another computational framework named MelonnPan which does not rely on a limited number of well-characterized taxa, enzymes, metabolites, and functional annotation unlike MIMOSA. MelonnPan uses elastic net regularization [80] to identify which OTUs are predictive for a given metabolite based on only their relative abundance profiles. For each metabolite, MelonnPan fits the elastic net model and optimizes the tuning parameters (i.e., both the elastic net mixing parameter and sparsity parameters) based on cross-validation. MelonnPan evaluates the predictability of each metabolite based on Spearman correlation coefficient (r) between the experimentally measured and predicted metabolite concentrations across samples. Metabolites with r > 0.3 are defined as well predicted. For significance testing on well-predicted metabolites, MelonnPan repeatedly shuffles the sample labels in both a metabolite table and an OTU table, and applies the MelonnPan to the randomized data by shuffling to link OTUs to metabolites and finally compares well-predicted metabolites between the original data and the randomized data. For multiomics integration, we focused on statistically significant metabolites between OC2 and OS2 since we are interested in bacterial populations that might be responsible for changes in metabolite pattern due to SFN diet in old mice.

Author Contributions
S.P.S. conceptualized, obtained funding, designed the studies for the project and contributed to manuscript preparation. A.C. and M.B. contributed to the study design and obtained funding for the microbiome and metabolomic study. S.J. and S.P.S. wrote the manuscript. S.J., A.C., C.B., performed experiments, analyzed data, performed bioinformatics analyses and contributed to manuscript preparation. P.T.P., M.B., and S.A. contributed to manuscript preparation and analyses of literature. All authors read and approved the final manuscript. Figure 1. Alpha diversity in the different experimental groups. (A) phylogenetic diversity index, and (B) species richness. The boxes denote interquartile ranges (IQR) with the median as a black line and whiskers extending up to the most extreme points within 1.5 times IQR. Outliers are noted as points. Note that O means Old and Y stands for Young; C means control diet and S indicates SFN diet; 0 means 0 months (immediately before diet administration) and 2 stands for 2 months after initiation of diet administration. Alpha diversity significance test results are summarized in Table S1.      Figure S1. Identification of bacterial biomarkers using linear discriminant effect size analysis. (A) Comparison of old SFN-treated mice (OS2) with young SFN-treated mice (YS2). (B) Comparison of OS2 and young mice that received control diet for 2 months (YC2). The distribution bar chart of LDA values shows the species with LDA scores greater than 3 and the species with significantly different abundances in different groups. The length of the histogram represents the size of the impact of significantly different species. The circle radiating from inside to outside represents the classification from the phylum to the genus level. Each small circle represents a classification at that level at different classification levels. Those taxa in each level are colored by farm for which it is more abundant. The diameter of the small circle is proportional to the relative abundance. Note that no differentially abundant species between YC2 and YS2 were identified.   Supplemental Tables   Table S1. Alpha diversity significance test. The Kruskal-Wallis test was used to calculate pvalue between experimental groups for phylogenetic diversity index and species richness. A pvalue <0.05 is flagged with an asterisk. P-value for Phylogenetic Diversity Index P-value for Species Richness