GutMIND: Using machine learning to link gut bacteria and brain function in mental health disorders across multiple groups

We conducted a search for the keywords “gut microbiota”, “shotgun sequencing”, “gut-brain axis”, “mental disease”, among others, in PubMed and Google Scholar (as of June 2024). A primary study exclusion criterion: (1) corresponding metagenomics data available for download from the National Center for Biotechnology Information (NCBI) and China National GeneBank (CNGB) databases. (2) relevant disease information. Metadata were extracted through manual curation of publications, database records, and author communications. The GutMIND database comprised 34 datasets derived from 31 independent studies, stratified by disease type and sequencing approach (details shown in Table S1).

Raw sequence files from NCBI/CNGB were processed using a bioinformatic pipeline. First, fastp42^,43 (v0.23.2) was used for quality control to remove low-quality reads and bases. Subsequently, Bowtie244 (v2.5.3) was employed for host sequence removal against the human reference genome (GRCh38). The resulting clean reads were then analyzed for taxonomic profiling using MetaPhlAn445 (v4.0.6) with the mpa_vOct22_CHOCOPhlAnSGB_202212 database, and for functional pathway analysis using HUMAnN46 (v3.9).

We applied a series of stringent filters to ensure the quality and relevance of our dataset. Specifically, we excluded any sample with a read count below 10 million to ensure sufficient sequencing depth and minimize technical variation.47 Additionally, samples with more than 30% unclassified sequences at the species level were removed to effectively exclude a small number of potential low-quality outlier samples, while retaining over 95% of the dataset for the final meta-analysis. We also excluded samples from individuals with a BMI greater than 30 kg/m² to reduce the confounding influence of obesity on the microbiota-gut-brain axis.48 Furthermore, if a study had fewer than 10 samples remaining after these filters, the entire study was discarded. After implementing these comprehensive filters, a total of 3,492 samples remained in our final dataset. To minimize technical and phenotypic heterogeneity across studies, we categorized the database into three groups: the high-quality dataset comprising Illumina paired-end sequenced samples with matched case-control designs (8 disease types), the platform-expanded dataset containing non-Illumina samples covering 3/8 overlapping diseases for platform comparison, and the disease-expanded dataset including samples representing additional disease phenotypes beyond the core set, enabling robust cross-platform validation while maximizing disease coverage (details shown in Table S1). Species profiling was refined using a relative abundance cutoff of 1 × 10⁻⁴ and an occurrence rate threshold of 1%. Specifically, samples with missing metadata were coded as “NA” and removed for statistical analyzes. The applied abundance and occurrence thresholds, along with sample inclusion criteria, reduce stochastic noise, minimize inter-sample heterogeneity, and improve the robustness and reproducibility of downstream analyzes, including differential abundance and predictive modeling.

Batch effect correction was applied to each disease cohort using three distinct methods: MMUPHin49 (adjust_batch function from R package MMUPHin v1.18.1), DEBIAS-M50 (DebiasMClassifier from Python module debiasm v0.0.1), and ComBat51 (ComBat function from R package sva (v3.52.0)), all run with default parameters while explicitly including disease status as a covariate to preserve the biological signal of interest. These methods were selected as they represent complementary strategies that respectively address global linear shifts, compositional and group-dependent effects, and signal–batch confounding. ComBat is a classical empirical Bayes method that models batch effects as linear additive shifts, providing a conventional baseline correction. MMUPHin extends this framework to microbiome compositional data, incorporating group structure and covariate effects to better handle compositional constraints. DEBIAS-M, in contrast, employs a debiased representation learning framework that removes batch effects confounded with biological signals. Principal Coordinates Analysis (PCoA) plots based on Bray-Curtis distances were generated to visualize microbial community distributions before and after batch correction. Comparative analyzes included pre-correction and post-correction results for each method across all disease cohorts.

The effectiveness of correction was quantitatively evaluated by calculating R² values using PERMANOVA (Permutational Multivariate Analysis of Variance) (999 permutations). PERMANOVA was conducted using the adonis2 function in the vegan52 (v2.7.1) R package, with Bray-Curtis distance as input. The R² values represent the proportion of variance explained by the tested factor (e.g., disease status or batch effect, depending on the analysis).

Alpha diversity was calculated using the Shannon index by vegan. Beta diversity was evaluated by Bray-Curtis dissimilarity and Jaccard index using vegan.52

To comprehensively delineate microbiota-neuropsychiatric associations with cross-regional validity, we conducted the Gut Microbiome in Multinational Integrated Neuropsychiatric Disorders (GutMIND) database. We acquired shotgun metagenomic data from NCBI repositories and manually curated clinical metadata from original publications. Uniform bioinformatics processing was performed on all sequencing data, followed by systematic organization and integration of the metadata. All datasets underwent standardized processing with stringent quality control (methods), MetaPhlAn4 and HUMAnN3 were used for taxonomic profiling and functional profiling (Figure 1a, Figure S1a).

The final GutMIND resource integrates 3,492 human gut microbiome samples (21.1 Terabases (Tb) of sequencing data) derived from 31 studies across 12 countries, representing a major integrated resource for gut microbiome studies in neuropsychiatric disorders. Geographically, China contributed the largest subset (n = 1,194; 34.2%), followed by the United Kingdom (n = 721; 20.6%), with remaining samples distributed across North America, Europe, Asia, and Oceania. The GutMIND cohort demonstrates broad demographic diversity, encompassing participants aged 1–95 years with balanced gender representation (1,083 males [47.6%] and 1,190 females [52.4%]) based on available phenotypic metadata (Figure 1a,b). This integrated resource covers 14 major neuropsychiatric conditions, with Parkinson's disease (PD; n = 1,350; 38.7%) and schizophrenia (SCZ; n = 536; 15.3%) representing the two largest diagnostic groups (Figure 1c, Table S1). Taxonomic profiling revealed Firmicutes as the dominant phylum (53.0 ± 19%), with Bacteroides (11.1 ± 11.1%) and Phocaeicola (9.1 ± 11.0%) representing the most abundant genera, and Phocaeicola vulgatus (4.9 ± 8.0%) together with Faecalibacterium prausnitzii (4.8 ± 5.0%) constituting the prevalent species (Figure 1d). For a specific disease, the dominant phylum is still Firmicutes (57.8 ± 9.7%) (Figure 1e). HUMAnN3-based functional analysis identified 118 conserved metabolic pathways as the core functional repertoire of the gut microbiome (Figure 1a), collectively providing multi-dimensional characterization of gut microbiome features.

Our analytical framework leverages this resource through four integrated modules: (1) Identification of cohort-stable microbial signatures using multi-modal approaches (diversity metrics, differential abundance testing, and variance partitioning); (2) Disease condition-specific differential analysis; (3) Development of the MetaClassifier workflow for constructing diagnostic models and deriving the Microbiota-Gut-Brain Axis Health Index (MGBA-HI); and (4) Systems-level analysis revealing a neuroprotective taxon consortium consistently associated with mental health outcomes, functionally enriched for neurotransmitter biosynthesis and immunomodulation pathways (Figure S1b).

To systematically investigate the gut microbiota-disease interactions across neuropsychiatric disorders, we established a comprehensive analytical framework operating at both cohort and disease levels. Our analysis utilized a high-quality dataset (n = 2,734) generated through Illumina paired-end sequencing (Materials and methods, Figure 2a). To enhance analytical robustness, we used both relative abundance and prevalence metrics. Taxonomic profiling conducted using MetaPhlAn4 demonstrated consistent microbial community structures across all cohorts (Figure S2).

Ecological community structure was quantified via α-diversity indices and β-diversity dissimilarity matrices. Across cohorts, PD patients exhibited a significant elevation or upward trend in α-diversity compared to healthy controls, whereas other neuropsychiatric disorders demonstrated inconsistent or non-significant alterations. Integration of multiple PD cohorts revealed a significant increase in microbial alpha diversity (p = 2.47 × 10⁻⁶) (Figure 2b, Table S2). These findings challenge the “dysbiosis-as-depletion” paradigm.58 Furthermore, consistent β-diversity perturbations emerged in disease groups, with statistically significant increases in Jaccard and Bray-Curtis dissimilarities observed in 15 out of 22 cohorts (68% of studied populations) (Figure 2c, Table S2). Such patterns corroborate prior observations of disrupted microbial community architecture and reinforce the “dysbiosis” paradigm in neuropsychiatric disorders.59^,60 PERMANOVA was performed with disease status as the explanatory factor to quantify its contribution (R²) to overall community variation. Five neuropsychiatric disorders - PD, SCZ, depressive disorder (DD), major depressive disorder (MDD), and bipolar disorder (BD) - demonstrated significant gut microbiome associations (PERMANOVA p < 0.05 for each disorder). Notably, autism spectrum disorder (ASD) showed significance only in the Wang_2019_ASD study, with observed significance specifically when using Jaccard distance-based R² calculations. Alzheimer’s disease (AD, including prodromal stages) and social anxiety disorder (SAD) had no significant microbiota-associated R² relationships across cohorts (Figure 2d, Table S2).

To elucidate shared microbial signatures across neuropsychiatric disorders, we identified taxonomic and functional variations across cohorts and diseases using MaAsLin357 on paired abundance and prevalence profiles. For disorder-specific differential analysis, we systematically evaluated three batch correction methods for taxonomic relative abundance data: ComBat, DEBIAS-M, and MMUPHin (Figure S3, 4) and PERMANOVA was performed with batch effect as the explanatory factor to quantify its contribution (R²) to overall community variation. Based on this comparative assessment, we ultimately implemented both DEBIAS-M and MMUPHin for integrated multi-cohort analysis following batch effect correction. These methods demonstrated effective batch effect reduction of 74.4% and 57.5%, respectively.

Integrating these analytical approaches, we identified significant microbial alterations in four distinct neuropsychiatric disorders: ASD, MDD, PD and SCZ. Meta-analysis revealed 294 taxon-specific recurrence events with disease-associated enrichment patterns (FDR-adjusted p < 0.1), representing cumulative microbial associations with disease status across participating cohorts (Table S3). Within this set, 23 species exhibited health-associated enrichment across various disease conditions, while 6 species showed enrichment in disease cohorts, and 16 context-dependent species exhibited disorder-specific enrichment patterns (Figure 3, Table S3). Notably, PD and SCZ demonstrated the highest concordance in microbial alterations, sharing 15 health-enriched and 4 disease-enriched species (Figure S5). Most control-enriched species were short-chain fatty acid (SCFA)-producing species. For example, Faecalibacterium prausnitzii, a SCFA-producing bacterium,61 exhibited significant enrichment in healthy individuals across three independent diseases encompassing MDD, PD and SCZ patients. In contrast, Ruthenibacterium lactatiformans demonstrated consistent enrichment in patient cohorts across three independent diseases (ASD, PD and SCZ). These cross-disease conservation patterns highlight core microbial regulators of brain health, suggesting broad-spectrum therapeutic targets for restoring gut-brain axis homeostasis. We also found context-dependent microbial patterns, particularly for Lacrimispora amygdalina and the unclassified GGB3746_SGB5089, which exhibited diametrically opposite enrichment profiles across disorders. These species were significantly depleted in healthy controls relative to MDD patients, yet showed the reverse pattern in PDs. This bidirectional dysregulation suggests these microbes may play disease-specific roles - potentially as protective commensals in PD that become pathogenic in MDD, or alternatively, may reflect disorder-specific alterations in gut microenvironment favoring their proliferation.

Our meta-analysis identified significant alterations in 108 metabolic pathways across four neuropsychiatric disorders (ASD, MDD, PD, SCZ; FDR-adjusted p < 0.1), including 20 consistently health-enriched pathways, 27 disease-enriched pathways, and 16 context-dependent pathways showing disorder-specific patterns (Figure S6, Table S4). PD and SCZ exhibited the strongest pathway overlap, sharing 29 disease-enriched and 18 health-enriched pathways (Figure S5). Key findings included depletion of adenosylcobalamin salvage (COBALSYN.PWY) in MDD/PD and enrichment of GABA shunt (GLUDEG.I.PWY) in ASD/PD patients.

To address the challenge of inconsistent diagnostic accuracy and poor generalizability of microbiome-based diagnostic models, we developed a standardized framework termed MetaClassifier to optimize model construction and evaluation (see Materials and Methods). For our disease-specific analyzes, we employed the comprehensive workflow illustrated in Figure 4a. This workflow begins with an initial data harmonization step, involving the integration and batch correction of multi-study data where applicable. Following this, the core modeling and validation stages are executed using the MetaClassifier framework, which encompasses nested CV for feature selection and the final external validation of the model's generalization and specificity (see Materials and Methods).

A comparison of five models across individual studies indicated that abundance-based CatBoost achieved the best overall predictive performance mean nested CV AUROC = 0.71 (Figure S7a,b). However, these single-study models exhibited limited generalizability and dropped to near-random levels when applied to independent cohorts of the same disease, a finding that held for models based on both abundance and prevalence data (Figure S7c-f, Table S5). To overcome this, we integrated studies of same disease and applied batch correction, subsequently identifying a DEBIAS-M-corrected CatBoost model as optimal for abundance data and a CatBoost model for prevalence data through rigorous re-benchmarking (Figure S8, Table S6).

Our framework successfully generated robust predictive models for a majority of the diseases studied. For abundance data, six of the eight diseases achieved nested cross-validation AUROCs at or exceeding 0.7. Across all eight diseases, the tuning cross-validation AUROCs were consistently high at approximately 0.8, and the final fit AUROCs approached 1.0, reflecting the model's fit to the full high-quality dataset (Figure 4b). Therefore, to assess robust generalization performance, we prioritized the nested cross-validation results within the high-quality dataset (Figure 4b) and independent external validation using the platform-extended dataset (Figure 4f,g). Prevalence-based models demonstrated a comparable level of performance (Figure 4c). Critically, all predictive models proved to be highly specific, showing only random-chance performance when applied to other diseases (mean cross-disease AUC ≈ 0.5) (Figure 4b,c).

The integration of disease-specific studies was crucial for improving model generalizability over single-project approaches. Our single-study model analysis of diseases ASD, PD, and SCZ, each with at least three cohorts, quantitatively demonstrated this necessity (Figure 4d–e). While models trained on a single study performed poorly when tested across studies (mean AUC dropping from 0.74 to 0.59), our integrated models maintained consistently high and robust performance across all individual studies (mean AUROC = 0.73). This proves that our integration strategy is essential to overcome cross-study heterogeneity within the same disease. Notably, the disease-specific classifiers also demonstrated considerable transferability to a platform-expended dataset comprising independent external cohorts, achieving AUROCs of approximately 0.70 for PD and SCZ (Figure 4f,g). Besides, shared microbial features emerged across diagnostic classifiers. For instance, Bifidobacterium longum enrichment emerged as a discriminative feature in SCZ discrimination models and was also identified as a key differentiator in abundance-based classifiers for AD, ASD, MDD, PD, and SCZ, as well as prevalence-based classifiers for ASD, AD, PD and SCZ (Figure S9, 10).

To inform preventive interventions for neuropsychiatric disorders and promote long-term brain health, early detection and risk stratification of neuropsychiatric status are essential. We developed the Microbial Gut-Brain Axis Health Index (MGBA-HI), a data-driven indicator to provide a single, quantitative score reflecting the overall microbial signature that distinguishes the neuropsychiatric and disease state from controls (Materials and methods).

The MGBA-HI demonstrated exceptional discriminative power in the high-quality dataset, with patients scoring significantly lower than healthy controls (Wilcoxon rank-sum p < 2.22 × 10⁻¹⁶, Cliff’s delta = –0.45), achieving significant separation for every individual disorder (Figure 5a,b, Table S7). This diagnostic capability extended to the more heterogeneous, platform-extended dataset (p = 4.25 × 10⁻⁴, Cliff’s delta = –0.21), where it maintained significance for PD (p = 2.06 × 10⁻⁴, Cliff’s delta = –0.34) and SCZ (p = 3.65 × 10⁻³, Cliff’s delta = –0.20) (Figure 5c,d). Remarkably, the MGBA-HI also showed powerful generalizability to disorders entirely excluded from model training, successfully differentiating patients with anorexia nervosa (AN; p = 8.12 × 10⁻⁴, Cliff’s delta = –0.28), obsessive-compulsive disorder (OCD; p = 3.19 × 10⁻³, Cliff’s delta = –0.35) in a disease-extended dataset (see Materials and Methods, Figure S11). In stark contrast, the conventional Shannon alpha-diversity index failed to distinguish patients from controls across these same scenarios. It showed no significant differences between the overall patient cohort and healthy controls, nor within the DD and MDD subgroups (Figure 5e,f). This limitation persisted in the platform-extended dataset for PD and SCZ (Figure 5g,h) and for all tested disorders in the disease-extended dataset (Figure S11). These consistent, opposing outcomes underscore a key advantage of the MGBA-HI: its superior sensitivity in detecting disease-associated microbial signatures that are missed by traditional diversity measures.

Beyond its diagnostic accuracy, MGBA-HI scores showed robust correlations with key clinical biomarkers across multiple conditions. In MDD, MGBA-HI scores exhibited a strong inverse correlation with depression severity (Kovtun_2022_MDD, CES-D: rho = –0.61, p = 7.62 × 10⁻⁷) and motor deficits in PD (Mao_2021_PD, UPDRS-III: rho = –0.34, p = 3.76 × 10⁻²). Furthermore, in SCZ (Zhu_2020A_SCZ), higher MGBA-HI scores were associated with better overall neurocognitive performance, as measured by the NAB (rho = 0.20, p = 2.16 × 10⁻²), and lower KYNA (rho = –0.17, p = 4.11 × 10⁻²), suggesting a potential association with neurotransmitter regulation, a hypothesis that warrants further mechanistic investigation (Figure 5i-l, Table S8).

The cross-disease conservation patterns highlighted core microbial species as regulators of brain health, critical for maintaining neurological health. We integrated health-enriched or high-prevalence species in healthy individuals compared with patients of neuropsychiatric disorders, and MGBA-HI-related species, to systematically delineate the core microbiota signature as the protective indicator of neuropsychiatric health. Differential abundance analysis identified 23 health-associated species (Figure 3). The MGBA-HI subsequently prioritized 107 species with robust discriminative capacity for healthy individuals, while prevalence screening revealed 73 microbial features consistently colonizing ≥50% of healthy individuals (Table S9). The intersection of these species ultimately defined 9 core neuropsycho-protective species (Figure 6a, Table S9).

To explore the neuroactive potential of the core protective species, we conducted genome-scale metabolic reconstruction of 21 genomes within the core consortium. Functional interrogation revealed that most genomes harbor complete gut–brain modules (GBMs)56 including SCFAs synthesis, amino acid metabolism, and the synthesis of several neurotransmitters. Notably, extensive functional overlap has been revealed among the 9-core gut bacterial species. All 9 species maintained genetic capacity for both glutamate and acetate biosynthesis, while 8 species (88.9%) encoded S-adenosylmethionine (SAM) synthesis pathways. Firmicutes bacterium AF16-15 exhibits genetic potential for butyrate synthesis among the core species (Figure 6b, Table S10).

MetaClassifier is publicly available at https://github.com/juyanmei/MetaClassifier↗. These related public data were derived from the following resources available in the public domain: NCBI (https://www.ncbi.nlm.nih.gov/sra↗) and CNSA (https://db.cngb.org/cnsa/↗). The related tables are available at https://figshare.com/s/416eae4ed99b8787304b↗.