Microbiota alterations in proline metabolism impact on depression through GABA and ECM homeostasis


 The microbiota-gut-brain axis has emerged as a novel target in depression, a disorder with low treatment efficacy. However, the field is dominated by underpowered studies focusing on major depression not addressing microbiome functionality, compositional nature, or confounding factors. We applied a multi-omics approach combining pre-clinical models with three human cohorts including mild-depressed patients. Microbial functions and metabolites converging into glutamate/GABA metabolism, particularly proline, were linked to depression. Whole-brain dynamics revealed rich club network disruptions associated with depression and circulating proline. Proline supplementation in mice exacerbated depression along with microbial translocation. Human microbiota transplantation induced an emotional-impaired phenotype in mice and alterations in GABA-, proline-, and extracellular matrix-related pre-frontal cortex genes. Targeting the microbiome and dietary proline may open new windows for an efficient depression treatment.


Main
Depression affects more than 300 million people worldwide and is well known to constitute one of the main causes of disability. 1,2 Despite this, the underlying mechanisms of depression still remains a crucial unresolved research topic, which is evidenced by the lack of an appropriate treatment, with an overall efficacy below 50% and a relapse rate of 40% in responders. 3,4 In addition, benefits become clinically relevant only in the small minority of patient populations with severe major depression. Therefore, there is an urgent need for new insights into the pathophysiology of depression.
A promising novel area of investigation involves the microbiota-gut-brain axis, 5 which has recently shown to control cognitive function 6 and inhibitory behavior (Gut 2021, in press). The microbiome influences the gut-brain communication through neural, endocrine, immune and neuroactive pathways. The latter includes microbial-derived neurotransmitters (e.g., GABA, catecholamines) and metabolites (SCFA, bile acids), and brain-derived neurotrophic factors (e.g., BDNF). 7 However, a recent systemic review that identified 19 studies analyzing the gut microbiota of clinical groups with depression revealed large inconsistent findings. 8 Remarkably, this body of research is dominated by underpowered cross-sectional case-control studies comparing healthy controls to patients with major depressive disorders (MDD) focusing solely on taxonomic composition based on 16S rRNA sequencing. Therefore, it is mandatory to develop powerful longitudinal studies under controlled conditions to clarify the specific role of microbiota on depression. Subjects with mild depression should also be investigated in order to evaluate all the spectrum of depressive symptomatology and potentially provide preventive measures. In addition, different microbial signatures may result in the same phenotype due to functional redundancy. Going beyond taxonomic composition to include analyses of microbial functionality, such as metabolomics and shotgun metagenomics sequencing, is thus vital to accurately capture the host-microbiome interplay. Finally, with just two exceptions, all these studies did not consider all essential confounding variables known to influence depression or the gut microbiota composition, particularly diet 9 and psychotropic medication; and only one study applied appropriate statistical methods to take into account the compositional nature of the microbiome datasets, a key point to obtain reliable results.

Depression scores are associated with a specific microbial ecosystem
To overcome all these issues, we applied a unique integrative longitudinal, multi-cohort and multi-omics approach to reveal molecular mechanisms underlying the microbiome-gut-brain axis interplay in depression combining fecal shotgun metagenomics, plasma and fecal metabolomics, and whole-brain dynamic functional magnetic resonance imagining (fMRI) in three human cohorts with complementary mice experiments. We first assessed the relationships of bacterial composition and functionality with depression, diagnosed using the Patient Health Questionnaire 9 (PHQ-9), in a longitudinal discovery cohort (IRONMET, n = 116, Table S1) comprising nondepressed (n = 44, PHQ-9: 0-4), mild depressed (n = 47, PHQ-9: 5-9), and major depressed subjects (n =25, PHQ-9 >10). 10 With few exceptions, 11,12 previous studies analyzing the associations between the gut microbiome and depression have not used appropriate compositional data analysis, thereby rendering potential misleading conclusions. In order to take into account the compositional nature of the microbiome data 13 we modelled read counts using a Dirichlet distribution to deal with 0 count values and then applied a centered log-ratio (clr) transformation as implemented in the ALDEx2 R package. 14 Previous studies have reported inconsistent alpha diversity findings, with a limited number reporting lower indices in subject with MDD compared to controls. 8 In line with this, we did not observe significant differences in the alpha diversity measures such as species richness or the Shannon index among groups (Fig. 1a,b). However, non-depressed individuals had higher Fisher's alpha diversity indices than depressed subjects, while no differences were found between mild and major depression (Fig. 1c). The compositional alternative to the principal coordinate analysis plots of β-diversity is the principal component analysis (PCA). Therefore, we applied a PCA analysis to the clr-transformed data to reveal global variance patterns in the microbial profiles and identify outliers. This initial unsupervised exploratory analysis revealed significant differences in the microbiome composition between non-depressed and depressed subjects (Fig.   1b).
We next sought to identify those bacterial species contributing the alterations in the gut microbiome associated with depression. Importantly, the vast majority of previous studies were based on 16S rRNA gene sequencing and, consequently, they do not provide enough taxonomic resolution to report results at the species level. For each taxa, we fitted a robust linear regression model between the PHQ-9 scores and the clr-transformed data controlling for age, gender, BMI, education years, and antidepressant and anxiety medication. Using this approach, we were able to identify thirty bacterial species significantly associated with depression (padj<0.1) (Fig. 1d, Table S2). Patients with higher PHQ-9 scores had higher levels of Parabaceroides spp. and Acidaminococcus spp. but lower levels of Bifidobacterium pseudolongum and species from the butyrate-producing Lachnospiraceae family, including Roseburia spp. All previous case-control studies were also cross-sectional in nature. Longitudinal research is required to elucidate causality in the associations. When we evaluated baseline bacterial taxa predictive of the PHQ-9 score oneyear later, we found that lower levels of several Bifidobacterium spp. and Lachnospiracae species and higher levels of Prevotella and Enterobacter species were associated with higher PHQ-9 scores one year later (Fig. 1e, Table S3).

Microbial functions involved in glutamate family of amino acids metabolism/transport and neurotransmitter transporters are linked to depression
Due to functional redundancy, significantly different microbial compositions may produce the same functional outcome. 15 Different phylotypes can cover identical functions and produce the same metabolic outputs under the same environmental conditions. It is therefore imperative to include functional analyses to disentangle the complex host-microbiome interaction in health and disease. Unlike most previous studies based on 16S rRNA sequencing, the use of a shotgun metagenomics sequencing approach allowed us to study the functional role of the microbiome in depression. Reads originating from microbial genes were mapped to the Kyoto Encyclopedia of Genes and Genomes (KEGG) orthologs. The resulting KEGG counts were then used to identify differentially abundant microbial gene functions in the metagenome associated with depression.
Fitted generalized linear models to the clr-transformed KEGG counts partialling out the effects of age, gender, BMI, education years, and antidepressants and anxiolytics medication, revealed several pathways (KEGG Orthology level 3) associated with the PHQ-9 score (padj <0.1) (Fig.   1f). Of note, bacterial pathway involving arginine, proline and histidine metabolism were negatively associated with depression. Along with glutamine, these amino acids comprise the glutamate family of amino acids, as their catabolism converges into glutamate, which serves as a carbon source for fueling the TCA cycle and subsequent gluconeogenic reactions as well as the synthesis of γ-aminobutyric acid (GABA). Consistently, we also found significant associations of bacterial glutamate metabolism, glutamatergic synapse and GABAergic synapse with the host PHQ-9 scores. Deeper analysis of the functional terms (KEGG Orthology level 4) identified 327 out of 5714 bacterial functions significantly associated with depression scores (Table S4).
Notably, some of the bacterial functions with strongest associations were involved in TCA cycle; arginine, proline, histidine, and glutamate metabolism; as well as the transport of these amino acids (Fig. 1g). Remarkably, the solute carrier family 6 (SLC6) also had a negative association with PHQ-9 scores. SLC6 genes mainly encode transporters for neurotransmitters (e.g., GABA, monoamines, glycine and proline), but also proteinogenic amino acids, betaine, taurine and creatine, 16 which is consistent with alterations in the GABAergic synapse.

Metabolites from the histidine, arginine and proline degradation pathways converging into glutamate and GABA shunt are associated with depression
To further explore the microbiome functionally, we next performed a metabolic profiling ( 1 H-NMR and HPLC-ESI-MS/MS) of plasma and fecal samples. Unlike shotgun metagenomics, which only provides information about microbial genes and thus predicts the microbiome potential functionality, metabolic phenotyping provides a readout of the actual microbiota functional activity. 15 Applying a machine learning variable selection strategy based on multiple random-forest, 17 we were able to identify several metabolites linked to the PHQ-9 scores ( Fig.   2a-f). Consistent with our previous findings, several of these metabolites were involved in the TCA cycle (succinate, fumarate), histidine metabolism (urocanate) and proline and glutamate metabolism. We then aimed to validate these findings by performing an HPLC-MS/MS-based metabolic profiling of plasma samples from in an independent large validation cohort (IMAGEOMICS, n = 919, Table S5) consisting also of non-depressed, mild and major depressed subjects. The most consistent finding was the strong positive association of circulating proline with the depression scores (Fig. 2g). Other consistent findings implicated alterations in retinol metabolism (Fig. 2g) and TCA cycle intermediates (citric acid, Fig. 2h). To facilitate analysis and interpretation of these metabolomic results, we performed a pathway over-representation analysis mapping metabolites that were significantly associated with the PHQ-9 in the discovery cohort to the KEGG, Reactome, INOH and HumanCyc databases included in the ConsensusPathDB. 18 This enrichment analysis further highlighted a significant (qval<0.1) overrepresentation of pathways associated with the TCA cycle and oxidative phosphorylation, glutamate metabolism, and arginine, proline and histidine catabolism (Fig. 2i). In agreement with metagenomic findings, it also stressed a significant over-representation of SLC and amino acid transporters as well as the GABA synthesis/degradation pathway and its link with the TCA cycle through the GABA shunt, which is predominantly associated with neurotransmission in the mammalian brain. A summary of the main metabolites and bacterial functions associated with depression is shown in Fig. 2j.

High proline consumption is associated with increased depression scores and small intestine genes participating in glutamatergic and GABAergic synapse and extracellular matrix homeostasis
Diet has a strong impact in modulating the composition and metabolic activity of the gut microbiome. 19 Therefore, we next evaluated the associations of macronutrients, vitamins, minerals, amino acids and fatty acids derived from food frequency questionnaires with the PHQ-9 scores in the IRONMET cohort. Strikingly, partial spearman's rank correlation analysis revealed proline as the dietary factor with the strongest impact on depression (Fig. 3a). When subjects were categorized according the median levels of plasma proline and dietary proline consumption, those with both low circulating and dietary proline levels had the lowest PHQ-9 scores, while individuals with both high dietary proline consumption and circulating proline levels had the highest depression scores (Fig. 3b). As dietary proline is extensively metabolized by enterocytes in the small intestine (Fig. 3c), we performed an RNA sequencing of jejunum samples from a second independent cohort (INTESTINE cohort, n = 28, Table S6) to identified those transcripts associated with proline consumption. Differential gene expression analysis were performed following TMM normalization using the limma pipeline with the voom transformation, linear modelling and empirical Bayes moderation. 20 We identified 1,547 out of 15,144 significant gene transcripts associated with dietary proline (Fig. 3d and Table S7). To facilitate functional interpretation of differentially expressed genes, we performed over-representation analyses mapping those genes to both Reactome ( Fig. 3 Table S8) and KEGG (Fig. S1a,b and Table S9) pathways. As pathway information is inherently redundant, with genes often participating in multiple pathways, we collapsed redundant pathways into a single biological theme using EnrichmentMap 21 to overcome redundancy and further simplify interpretation ( Fig.   3f and Fig. S1b). Notably, Reactome-based analyses identified pathways involved in GABA receptor activation, synaptic interactions and axon guidance. It also revealed pathways participating in extracellular matrix (ECM), muscle contraction, MAPK signaling and GPRC signaling. In line with these results, KEGG-based analyses highlighted several pathway involved in neuron synapse (Fig. S1a,b), in particular GABAergic (Fig. 3g) and glutamatergic synapse ( Fig. S1c), which is consistent with metabolomics results. Our analyses in the jejunum also identified pathways participating in ECM and muscle contraction. We also mapped significant transcripts associated with dietary proline to the DisGeNET disease-based database, which contains one of the largest collections of genes associated with human diseases, using both expert curated and text mining data. Remarkably, enrichment analysis highlighted a significant overrepresentation of diseases associated with cognitive, neurodegenerative and CNS disorders, including schizophrenia and major depressive disorders (Fig. S2a-c and Table S10).

Proline supplementation exacerbates depressive-like behavior in mice in association with microbial translocation
To further investigate the direct effect of dietary proline on depression, we performed a supplementation study in mice (Fig. 3h). A total of n=40 C57BL/6J mice were fed either a standard diet (SD) or a crude fiber rich diet (CFD) and supplemented with either water or proline (36 g/L). Throughout the experiment, mice were chronically exposed to unpredictable mild stressors (UCMS model), an animal model with high face, constant and predictive validity, to develop and evaluate depressive-like behavior. After 6 weeks of UCMS exposure, mice supplemented with SD + proline had higher immobility times in the forced swim test (Fig. 3i), a well validated model of the despair behavior, and reduced sucrose intake (Fig. 3j), a behavioral model of anhedonia, a core symptom of depression, compared to stress-free control mice that only received water. Additionally, compared to control animals, SD + proline mice had higher levels of lipopolysaccharide binding protein (LBP) (Fig. 3k), which is indicative of microbial translocation and in agreement with our microbiome results. Interestingly, these behavioral and biochemical effects of proline did not appear in mice fed a CFD, suggesting a potential protective effect of this diet.

Brain iron deposition and whole-brain functional dynamics reveal rich-club network disruptions associated to both depression and circulating proline
Metagenomics, metabolomics and RNA-seq analyses were consistent in identifying alterations in glutamatergic and GABAergic systems. Since iron plays a crucial role in glutamate and GABA homeostasis, we assessed iron deposition in the brain of the IRONMET patients using magnetic resonance imaging (MRI) based on T2* and R2* relaxometry. In fact, iron overload in certain brain areas has been implicated in neurodegenerative disorders 22 and alterations in anxiety-like behaviour and mood, 23 but only few studies have characterized the impact of iron on depression.
Partial spearman's rank correlation analysis among relaxometry parameters in brain regions (AAL Atlas) and PHQ-9 scores highlighted low iron deposition (mean T2* values) of regions along the cingulum and frontal lobe positively associated with depression (Fig. 4a). These results were validated using a machine learning variable selection strategy (Fig. 4b). Moreover, the human brain is a complex network of structurally and functionally connected regions. Functional communication between these regions is thought to play a vital role in complex processes such as depression. However, most studies have focused on static functional connectivity, but functional connectivity among brain networks is not static over time. Studying the dynamics of resting-state brain activity across the whole-brain functional network might facilitate interpretation of brain functioning and provide better insights into the pathophysiology of disease. 24 Therefore, we analyzed resting-state fMRI data (Fig. 4c) in a subset of patients from the IMAGEOMICS cohort (n=591) and studied whole-brain functional dynamics applying a novel intrinsic-ignition framework across the whole-brain functional network (214 brain areas), 25 to assess the effect of spontaneous local activation events on local-global integration (Fig. 4d), followed by machine learning algorithms to identify those nodes predictive of the PHQ-9 scores. Using the Boruta algorithm we identified 72 intrinsic-ignition nodes associated with depression ( Fig. 4e). Similarly, it revealed 68 nodes linked to the circulating proline levels (Fig. 4f), 30 of which were also associated with the PHQ-9 scores ( Fig. 4g and Table S11). Notably, these shared intrinsicignition areas, including the superior frontal cortex, the precuneus, insula and subcortical areas such as the caudate, putamen, and hippocampus, mostly belong to the so-called "rich-club" (Fig.   4h), a set of high-degree nodes that tend to be more closely connected among themselves than with peripheral regions, i.e. lower degree nodes. 26,27 An emotional impairment is transferred to mice through the gut microbiota. Direct impact on the expression of mPFC genes participating in proline and GABA

transport and ECM and collagen metabolism in recipient mice
Finally, to evaluate a potential causal role of the microbiota in the development of emotional disorders, we transplanted microbiota from 20 human donors with different PHQ-9 scores into 20 antibiotic-treated mice (Fig. 5a). We then assessed whether an emotional impaired phenotype emerged in mice receiving microbiota from donors with higher depression scores using a fear conditioning induced freezing test, a well-recognized model of maladaptive response to stress.
Notably, donor's PHQ-9 scores were significantly correlated with the freezing time in recipient mice (Fig. 5b). One of the most consistent findings in MDD include decreased frontal lobe function, mainly involving the medial prefrontal cortex (mPFC). 28 In addition, results from a recent meta-analysis suggested that decreased levels of glutamatergic metabolites in the mPFC are linked with the pathophysiology of depression. 29 Considering our findings highlighting the impact of the microbiome on glutamatergic and GABAergic systems, we performed and RNAseq of the mPFC of recipient mice. We identified 59 out of 15,537 gene transcripts in recipient's mice mPFC significantly associated with donor's PHQ-9 scores (Fig 5.c and Table S12). To gain a better insight into the potential mechanisms underlying the microbial effects on depression, we built gene-gene interaction networks using the STRING database 30 (Fig. 5d,e). Furthermore, we mapped significant genes to KEGG (Fig. 5f) and Reactome (Fig. 5g) databases. Half of the transcripts negatively associated with the donor's PHQ-9 score clustered together and were involved in oxidative phosphorylation and neurodegenerative disease. Notably, we identified a cluster comprising gene transcripts encoding for transporters of GABA (slc6a12 and slc6a13) and Proline (slc6a20) (Fig. 5h), which is in agreement with our functional analyses. Consequently, we found an over-representation of Reactome pathways associated with GABA neurotransmission, which is also in consonance with the results from the jejunal RNA-seq analysis. These results are in line with metagenomics and metabolomics and highlight again the importance of the GABA shunt interconnecting the TCA cycle with GABA and glutamate metabolism. Not only that, but we also found another cluster of gene transcripts participating in ECM and collagen homeostasis and muscle contraction. Importantly, proline is essential for collagen biosynthesis, constituting 10% of its amino acid content. Therefore, these results validate and further highlight the importance of the proline-glutamate-GABA-microbiome interplay identified in our metabolomics and metagenomics analyses.

Discussion
The gut microbiota has emerged as a novel actor in the pathophysiology of depression. However, recent meta-analyses have revealed a strong inconsistency among studies in terms of gut microbiome signatures associated with depression. 8 These inconsistencies mostly arise from underpowered studies, methodological heterogeneity, varying depression diagnostic criteria, lack of consideration of confounding variables, inappropriate statistical analysis, and functional redundancy of the microbiota. Furthermore, a limitation of the vast majority previous studies includes the use of 16S-RNA sequencing, which does not have enough taxonomic resolution to report results at the species levels does not provide information about microbial functionality.
Here, we tackled these limitations by applying a multi-omics approach in three human cohorts and pre-clinical studies.
In concordance with some previous findings, 8 we found decreased levels of SCFA-producing bacteria such as species from the Lachnospiraceae family (including Roseburia spp.) and Bifidobacterium spp. 31 Notably, Bifidobacterium strains are amongst the most efficient GABAproducers. 32 The current functional metabolomics and metagenomics analyses highlighted several bacterial functions and metabolites involved in proline, histidine and arginine pathways that converge into glutamate and GABA metabolism linked to depression scores. In line with this, we found that patients with higher PHQ-9 scores had higher levels of Acidaminococcus spp., which use glutamate as the only carbon source and have shown to grow only in culture media containing arginine, glutamate or histidine. 33 Dysregulation of glutamate and GABA neurotransmission and increased circulating levels of glutamate and GABA have been reported in participants with MDD. [34][35][36] Therefore, the current findings point towards a potential role of the microbiome in depression through glutamate/GABA metabolism, compatible with the glutamate hypothesis of depression. 36 So far, only two studies have performed functional metagenomics analysis. 11,12 Although not significant, a recent study performing a targeted analysis pointed towards alterations in GABA shunt and glutamate degradation pathways in subjects with depression 11 . Similarly, in a small study (n=40) the GABA degradation pathway was prominent in the microbiome of individuals with MDD. 12 The most consistent finding was the positive association of circulating proline with the depression scores. In a recent meta-analysis of peripheral blood metabolites in major depressive disorders, a subgroup analysis revealed that antidepressant-free MDD patients had higher levels of Lproline. 37 It is worth noting that our analyses were controlled for antidepressant medication.
Hyperprolinemia has also been linked to epilepsy, schizophrenia, seizures and impaired cognitive function. Importantly, the proline degradative pathway can eventually generate glutamate and GABA. Hence, proline accumulation has shown to disrupt GABA production, glutamate release and impair and synaptic transmission, 38,39 while inactivation of proline transporter altered glutamatergic synapse and perturbed behaviours in mice. 40 In line with previous findings, we found that dietary proline was strongly associated with several pathways involving GABAergic and glutamatergic synapse. We also demonstrated a potential causal role of this amino acid in depression by supplementing mice with proline. Notably, circulating proline levels were consistently associated with nodes of intrinsic brain networks linked to depression. Alterations in the circulating proline levels were also linked to depressions cores through the so-called "richclub" of highly interconnected nodes. The rich-club of the brain network plays a critical role in global integration of neural information and is essential for efficient communication across multiple segregated and distant brain regions. 26,27 Therefore, network disruptions within the rich club have a direct impact on various behavioral and cognitive tasks. Connectivity-based neuroimaging studies have identified alterations in rich-club network in several neuropsychiatric disorders. 26 Recently, disruptions in the rich-club network organization has also been implicated in MDD (n=32). 41 Here, we extend these findings to mild-moderate depression using a much larger cohort (n=591).
Dietary proline was also strongly linked to the ECM homeostasis not only in the human intestine but also in the mouse brain. Collagen makes up 80% of ECM, with proline and hydroxyproline accounting for 25% of collagen amino acid content. Therefore, proline is stored in the ECM as collagen, which can serve both as a reservoir or dump for proline. 42 Of note, components of the ECM such as matrixmetalloproteinases, which degrade collagen and release proline, and perineuronal nets, have been recently suggested as new key players in the development of psychiatric disorders. 43  Hospital. All subjects were of Caucasian origin and reported a body weight stable for at least three months before the study. Subjects were studied in the post-absorptive state. The following exclusion criteria were considered: i) no systemic disease other than obesity; ii) free of any infections in the previous month before the study; iii) no liver diseases (specifically tumour disease and infections) and thyroid dysfunction, which will be specifically excluded by biochemical work-up. This protocol was revised, validated and approved by the Ethics committee of the Hospital Dr Josep Trueta. The purpose of the study was explained to participants and they signed written informed consent before being enrolled in the study.

Clinical and laboratory parameters
Completed medical history and anthropometric variables were collected from all participants. In fasting conditions, a blood sample was provided. Fasting plasma glucose (FPG) and lipid profiles were measured by standard laboratory methods using an analyzer (CobasR 8000 c702, Roche Diagnostics, Basel, Switzerland). High-sensitivity C-reactive protein (hsCRP) levels were determined by immunoturbidimetric method (CobasR 8000 c702, Roche Diagnostics, Basel, Switzerland). Glycated hemoglobin (HbA1c) was determined by high performance liquid chromatography (ADAMRA1c HA-8180V, ARKRAY, Inc., Kyoto, Japan).

Body composition
Fat mass (FM), fat free mass (FFM) and their distribution was measured by a dual energy xray absorptiometry (DEXA, GE lunar, Madison, Wisconsin).

Dietary pattern
The dietary characteristics of the subjects were collected in a personal interview using a validated food-frequency questionnaire. 44

Cognitive assessment
The Patient Health Questionnaire-9 (PHQ-9) is the depression module of the PRIME-MD diagnostic instrument for mental disorders. 45 It is a self-administered questionnaire and consists of 9 items of depression symptoms plus a question about functional impairment. The PHQ-9 can be scored either as a depression severity rating (range 0-27 points) or with an algorithm based on the DSM-IV criteria (major and minor episode). It can also be interpreted using a cut-off point applied the symptom severity score. Scores of 10 to 14 represent a moderate symptom severity level, 15 to 19 represent moderately severe symptoms, and 20 to 27 severe depressive symptoms.
Scores of 10 or more have an 88% sensitivity and specificity.

MRI acquisition and image pre-processing
All subjects were studied on a 1.5T Ingenia (Philips Healthcare, Best, The Netherlands) with eight channel head coils. As a part of a larger study protocol a multislice fluid attenuation inversion recovery (T2-FLAIR) with TR/TE/TI=6500/120/2200ms, flip angle 90º, in-plane resolution 0.78x0.78mm, slice thickness 5mm without gap and 20 axial slices was used to exclude preexisting brain lesions. MRI relaxometry was assessed by using a multi-echo gradient echo sequence with TR/1stTE/ΔTE=800/2.2/5ms, flip angle 80º, in-plane resolution2x2mm, slice thickness 5mm without gap and 20 axial slices. After acquisition, T2* and R2* maps were computed using Olea Sphere 3.0 (Olea Medical, La Ciotat, France) with Bayesian analysis algorithm. T2* maps were calculated by fitting the signal decay curve of the respective magnitude multiecho data and R2* maps were calculated as R2* = 1/T2*. In addition, a brain extraction tool was used to delete all non-brain tissues of calculated T2* and R2* maps. R2* were measured in s -1 .

Intrinsic-Ignition Framework
We applied the Intrinsic-Ignition Framework 25 to obtain the effect of naturally occurring activation events that reflect the capability of a given brain area to propagate activity to other brain areas. In brief, we transformed the BOLD time series to phase space by filtering the signals within the narrowband (0.04-0.07 Hz) and computed the Hilbert transform to obtain the phases of the signal between each pair of brain areas at each time point. Fig. 4c shows the representation of the Hilbert BOLD phase for a brain area over time in the complex plane. Fig.   4d shows the algorithm used to obtain the ignition value of each brain area evoked by an event within a set time window. The binary events were defined by transforming the time series into z-scores, z i (t), and fixing a threshold, θ. Then, a phase lock matrix P jk (t), was computed which describes the state of phase synchronization between brain areas j and k at time t as: where φ j (t) and φ j (t) correspond to the phases of the BOLD time series for brain areas and at time . Then, the integration was defined by measuring the length of largest connected component in the phase lock matrix ( ) and the value of integration was computed as the length of the connected component considered as an adjacent graph (i.e., the largest subcomponent). Finally, for each brain area, we averaged across the events the integration evoked at each time t within the set time window. A complete description of the method can be consulted in 25 .

Extraction of faecal genomic DNA and whole-genome shotgun sequencing
Total DNA was extracted from frozen human stools using the QIAamp DNA mini stool kit All 1 H-NMR spectra were recorded at 300 K on an Avance III 600 spectrometer (Bruker®, Germany) operating at a proton frequency of 600.20 MHz using a 5 mm PABBO gradient probe and automatic sample changer with a cooling rack at 4ºC.

Proline supplementation mice experiment
Animals. Male CD-1 mice (8 weeks old at the beginning of the experiment) were used. Animals were housed in reverse light-dark cycle (lights on from 20:00 to 8:00), standard temperature (21º +/-1º C) and with food and water available ad libitum. Mice exposed to the Unpredictable Chronic Mild Stress (UCMS) procedure were singled housed in one room. Control mice, not exposed to UCMS, were housed 2-4 per cage in a different room. Animals were divided into the following experimental groups based on the type of diet and water they were exposed to. Thus, 10 control mice were exposed to standard diet and water whilst another 10 control mice were exposed to an especial diet rich in crude fiber (RCFD; Altromin C1014, 20% crude fiber content) and to lproline (Sigma-Aldrich P0380) diluted in the drinking water (36 g/l) for the entire experimental procedure. Similar, mice exposed to the UCMS protocol were divided in the following 4 groups of 10 mice each: i) standard diet and water, ii) standard diet and l-proline (36 g/l) in drinking water, iii) RCFD and water, and iv) RCFD and l-proline (36 g/l) in the drinking water. L-proline was freshly prepared and replaced every 3-4 days.   were given twice weekly to throughout the study to reinforce donor microbiota phenotype. Animals were exposed to a fear conditioning test with nociception assessed by the hot plate test to ensure specificity. At the end of the study, animals were consecutively sacrificed and the brains were quickly removed and the medial prefrontal cortex was dissected according to the atlas of stereotaxic coordinates of mouse brain. 58 Brain tissues were then frozen by immersion in 2methylbutane surrounded by dry ice, and stored at -80ºC.

Unpredictable Chronic Mild
Emotional testing in mice. Fear conditioning was conducted as described previously with some modifications. 59  Metagenomics analysis. To take into account the compositional structure of the microbiome data and rule out possible spurious associations, raw counts were transformed using a centered logratio (clr) transformation as implemented in the "ALDEx2" R package. 14 It first uses a Dirichletmultinomial model to inter abundance from read counts and then applies a clr transformation to each instance. We used 128 Dirichlet Monte Carlo instances in the aldex.clr function. Bacterial species and functions associated with the PHQ-9 scores were identified using robust linear regression models as implemented in the Limma R package, 20 adjusting for age, body mass index, sex, education years, antidepressant and anxiety medication. Taxa and bacterial functions were previously filtered so that only those with more than 10 reads in at least five samples were selected. The p-values were adjusted for multiple comparisons using a Sequential Goodness of Fit 61 as implemented in the "SGoF" R package. Unlike FDR methods, which decrease their statistical power as the number of test increases, SGoF methods increase their power with increasing number of tests. SGoF has proven to behave particularly better than FDR methods with high number of tests and low sample size, which is the case of omics large datasets. Statistical significance was set at padj<0.1.
Metabolomics and whole-brain functional dynamics analysis. Metabolomics data were first normalized using a probabilistic quotient normalisation. Metabolomics and whole-brain functional dynamics data were analysed using machine learning (ML) methods. In particular, we adopted an all-relevant ML variable selection strategy applying a multiple random forest (RF)based method as implemented in the Boruta algorithm. 17 It has been recently proposed as one of the two best-performing variable selection methods making use of RF for high-dimensional omics datasets. 62 The Boruta algorithm is a wrapper algorithm that performs feature selection based on the learning performance of the model. 17 It performs variables selection in three steps: a) Randomization, which is based on creating a duplicate copy of the original features randomly permutate across the observations; b) Model building, based on RF with the extended data set to compute the normalized permutation variable importance (VIM) scores; c) Statistical testing, to find those relevant features with a VIM higher than the best randomly permutate variable using a Bonferroni corrected two-tailed binomial test; and d) Iteration, until the status of all features is decided. We run the Boruta algorithm with 500 iterations, a confidence level cut-off of 0.005 for the Bonferroni adjusted p-values, 5000 trees to grow the forest (ntree), and a number of features randomly sampled at each split given by the rounded down number of features/3 (the mtry recommended for regression). Pathway over-representation analysis was performed mapping metabolites that were significantly associated with the PHQ-9 in the discovery cohort to the KEGG, Reactome, INOH and HumanCyc databases included in the ConsensusPathDB. 18 RNA-seq analysis. Differential expression gene analyses were performed on gene counts using the "limma" R package. 20 First, low expressed genes were filtered, so that only gene with more than 10 reads in at least 2 samples were selected. After filtering, 15,144 genes out of 22,204 were retained for subsequent analyses. RNA-seq data were then normalized for RNA composition using the trimmed mean of M-value (TMM) as implemented in edgeR package. 63 Normalized counts were then converted to log2 count per million (logCPM) with associated precision weights to account for variations in precision between different observations using the "voom" function with donor's age, BMI, sex, education years, antidepressant and antxiety medication, and kcal intake as covariates. A robust linear regression model adjusted the previous covariates was then fitted to the data using the "lmFit" function with the option method = "robust", to limit the influence of outlying samples. Finally, an empirical Bayes method was applied to borrow information between genes with the "eBayes" function. Over-representation analyses were performed by mapping differentially expressed genes into the Reactome, KEGG, and DisGeNET databases. Pathway significance was assessed using a hypergeometric test and a Storey procedure (q-values) was applied for multiple testing correction. Statistical significance was set at qval<0.1.
For the FMT study, differentially expressed genes were also mapped to the Search Tool for Retrieval of Interacting Proteins/Genes (STRING) database (which integrates known and predicted protein/gene interactions) to predict functional gene-gene interaction networks. 30 Then, functional local clusters in the interaction network were determined using a Markov Cluster algorithm (MCL) with a inflation parameter = 3. Active interacting sources including text mining, experiments, databases, co-expression, and co-occurrence and an interaction score > 0.4 were used to construct the interaction networks. Pairwise differences between groups were assessed using the pairwise.adonis function adjusted for Bonferroni correction. *, P < 0.05; **, P < 0.01. c) Principal component analysis scores plot