Prediction of sporulating Firmicutes from uncultured gut microbiota using SpoMAG, an ensemble learning tool

Sample collection and bioinformatics processing

Samples were collected from the GUARANI One Health Brazilian Group Network (Lemos et al., 2022). Between February and April 2020, rectal swabs were obtained in triplicates from cattle (n = 30), swine (n = 15), poultry (n = 30), and human fecal samples (n = 32) (Lemos et al., 2022). The sampling encompassed five Brazilian geographic regions: Northern (Castanhal-Pará/PA), Southern (Blumenau-Santa Catarina/SC), Southeastern (Bragança Paulista-São Paulo/SP), Midwestern (Dourados-Mato Grosso do Sul/MS), and Northeastern (Fortaleza-Ceará/CE), as reported previously (Lemos et al., 2022; Machado et al., 2024).

DNA extraction, sequencing, and bioinformatics analyses followed established protocols described by Lemos et al. (2022). MAGs were evaluated for quality using the Minimum Information about a Metagenome-Assembled Genome (MIMAG) criteria, and only high-quality MAGs (>90% completeness, <5% contamination) (Bowers et al., 2017) were retained for downstream analyses. To assess species-level similarity between MAGs from different hosts, we applied the FastANI method to calculate Average Nucleotide Identity (ANI), by considering a fragment length of 1,020 bp (–fragLen 1020) (Hernández-Salmerón, Irani & Moreno-Hagelsieb, 2023).

Selection of sporulation-associated genes

A total of 160 genes associated with bacterial sporulation (Table S1) were selected based on their conserved presence across Firmicutes species and established functional roles in sporulation processes (Galperin et al., 2022). Genes primarily involved in housekeeping functions during vegetative growth were excluded, as they are ubiquitously present in vegetative cells and do not specifically contribute to sporulation (Galperin et al., 2022). The selected genes represent key steps in sporulation, including sporulation onset and checkpoints, Spo0A regulon, engulfment, sigma factor regulons (SigF, SigE, SigG, SigK), spore cortex formation, spore coat assembly, and germination.

Selection of publicly available genomes of sporulating and non-sporulating bacteria

For supervised machine learning model development, we compiled a dataset of 136 high-quality complete bacterial genomes from a published study (Galperin et al., 2022). The dataset was balanced, comprising 68 sporulating and 68 non-sporulating species (Table S2), to minimize bias during model training and evaluation. Spore-forming genomes included 34 Bacilli and 34 Clostridia species, representing the two major lineages known for sporulation within the Firmicutes phylum, thus ensuring phylogenetic diversity. The non-sporulating group included 30 Bacilli and 38 Clostridia genomes, selected to closely match phylogenetic diversity while maintaining a clear functional distinction (i.e., absence of sporulation capability).

Although certain sporulation-associated genes are partially conserved across Firmicutes, others are class-specific (e.g., Bacilli-exclusive regulators). Thus, including both Bacilli and Clostridia allowed the modeling approach to capture both universal and lineage-specific genetic signatures, enhancing generalizability across diverse Firmicutes species, even non-described-yet species from metagenomics studies.

All genomes were uniformly processed using the same annotation pipeline applied to the MAGs to ensure consistency. Open reading frames (ORFs) prediction was performed using Prodigal v.2.6.3 in single-genome mode (-p single) with a bacterial translation table (-g 1) (Hyatt et al., 2010). Functional annotations for each genome were performed using eggNOG-mapper v.2.0.8 (e-value ≤ 1e−5, identity > 60%, query/subject coverage > 60%) (Huerta-Cepas et al., 2017) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) database (Kanehisa et al., 2025).

Data preparation for supervised machine learning models to predict sporulation probability

The target variable was encoded as a binary factor with the levels “Sporulating” and “Non-Sporulating”. To guarantee reproducibility, we implemented random stratified partitioning of the dataset into training (70%) and testing (30%) subsets using the “createDataPartition” function from the caret package in R (Kuhn, 2008), version 7.0.1. A binary presence/absence matrix of the sporulation genes was constructed prior to model training, representing fixed biological characteristics for each genome. Genomes included in the test set were completely withheld from the training process to prevent any form of data leakage. Given that only high-quality complete isolated genomes were used, missing values were treated as absences and encoded as zero, reflecting the assumption that unobserved genes were absent in the genomes.

Model training and tuning

We evaluated four supervised machine learning algorithms: Random Forest (RF), support vector machine (SVM), eXtreme Gradient Boosting (XGBoost), and neural network (NN). All models were implemented using caret (Kuhn, 2008) within a 10-fold cross-validation framework. Model performance was primarily assessed using the area under the receiver operative characteristic (ROC) curve (AUC).

The RF model was tuned by varying the number of variables randomly sampled at each split (mtry = 2, 3, or 4), while the number of trees (ntree) was fixed at 500. SVM was implemented with a radial basis function kernel. Hyperparameters were tuned over a grid of cost values (C = 0.25, 0.5, 1, 2, 4, 8, 16, 32, 64, and 128), while the kernel coefficient (sigma) was held constant at 0.0378. The XGBoost model was trained with learning rate (eta = 0.1), maximum tree depth (max_depth = 6), subsample ratio (subsample = 0.8), and column subsample by tree (colsample_bytree = 0.8). Early stopping was applied with a patience of 10 rounds, and the model was trained for a maximum of 100 rounds. Lastly, the NN was implemented as a feedforward multilayer perceptron using the nnet package, version 7.3-20 (Venables & Ripley, 2002). The architecture consisted of a single hidden layer, with the number of neurons tuned via repeated 10-fold cross-validation. Training was performed using backpropagation with a maximum of 200 iterations.

Development of the SpoMAG meta-classifier using an ensemble strategy

Based on their superior classification performance, RF and SVM algorithms were selected as base learners for the development of the SpoMAG, an ensemble classifier designed to predict probability of sporulation potential in MAGs. SpoMAG employs a stacked generalization strategy, where the probabilistic outputs of the RF and SVM models are used as input features for a meta-classifier. This meta-classifier was implemented as an RF model trained with 10-fold cross-validation (Ghasemieh et al., 2023), enhancing predictive robustness.

The meta-classifier’s hyperparameters, including the number of features considered at each split (mtry parameter), were optimized to improve classification performance. Both the base models and the ensemble were evaluated using standard metrics, including accuracy, recall, specificity, precision, F1-score, and AUC. To further assess predictive performance, we applied a non-parametric bootstrap procedure with 300 resampling iterations of the training set (Singh et al., 2021). Medians and 95% confidence intervals were calculated for all evaluation metrics, providing estimates while accounting for variability in the training data (Mokhtar, Yusof & Sapiri, 2023; Huang & Huang, 2023).

Before prediction, SpoMAG processes the functional annotation file of each MAG to identify genes based on gene names and KEGG Orthology (KO) identifiers. It then represents each MAG as a binary matrix encoding the presence or absence of 160 sporulation-associated genes and outputs a probabilistic estimate of sporulation potential for each genome.

Model interpretation and feature importance

To determine the most predictive genes for sporulation capability, we computed SHapley Additive exPlanations (SHAP) values (Lundberg & Lee, 2017) for both RF and SVM models using the iml package (Molnar, 2018), version 0.11.4, in R. The SHAP analysis measures the marginal contribution of each gene to individual predictions, offering an interpretable and consistent metric of feature importance.

Only genes identified as positive predictors in both RF and SVM models were selected for final biological interpretation. For visualization of feature importance, a bar plot was generated showing the mean absolute SHAP values of the contributing genes, arranged in descending order to highlight the most significant predictors of sporulation. The visualization was generated using ggplot2 (Wickham, 2016), version 3.5.1, and ggpubr (Kassambara, 2023), version 0.6.0, packages in R.

Ordination and multivariate dispersion analyses of sporulation profiles in Clostridia and Bacilli

To investigate variation in sporulation-associated gene profiles, we performed principal component analysis (PCA) using the prcomp function from the stats package (R Core Team, 2024), version 4.3.3. The analysis was based on binary matrices encoding the presence or absence of sporulation-associated genes across both reference genomes with known sporulation phenotypes and MAGs predicted as spore-forming by SpoMAG.

To assess potential differences in gene composition between the predicted sporulating and non-sporulating groups within the Clostridia and Bacilli classes, we performed a multivariate dispersion analysis using the betadisper function from the vegan package (Oksanen et al., 2024) (version 2.6.8). First, binary distance matrices were computed using the dist function (stats package), based on presence/absence profiles. The dispersion, i.e., within-group variability, was then calculated as the distance of each MAG to the centroid of its respective phenotypic group. Statistical significance of group dispersion differences was evaluated using 1,000 permutations and an analysis of variance (ANOVA) test at a 5% significance level.

The complexity of sporulation gene profiles challenges linear separation between sporulating and non-sporulating genomes

To explore whether the presence and absence of sporulation-associated genes alone could distinguish spore-forming from non-spore-forming organisms, we performed a PCA across the 136 Clostridia and Bacilli genomes, including both phenotypes (Fig. 2). The first principal component (PC1) accounted for 41.99% of the total variance, followed by PC2 with 14.48%. While PCA captured discernible clustering patterns, we observed significant phenotypic overlap between spore-forming and non-spore-forming genomes in both phylogenetic classes. This substantial overlap demonstrates that dimensionality reduction of binary gene presence-absence data alone provides insufficient resolution for reliable sporulation phenotype prediction.

Comparison of supervised machine learning models for classifying genomes as sporulating or non-sporulating

Among the four supervised classification models evaluated, RF and SVM exhibited the highest median performance in distinguishing sporulating from non-sporulating bacterial genomes (Table 1). The RF model showed the highest median accuracy (87.5%) and specificity (85.7%), reflecting robust generalization and a low false-positive rate, with confidence intervals indicating consistent predictions. Meanwhile, SVM exhibited the highest median recall (95.2%), indicating a high sensitivity for correctly identifying spore-forming genomes, albeit with lower specificity (75%). Both models reached comparable median AUC values (RF 89.7%, SVM 93%), reflecting their discriminative capabilities.

To leverage the complementary strengths of RF and SVM, we developed SpoMAG, a stacking ensemble classifier that integrates these two base learners. SpoMAG preserved the high median accuracy of RF (87.5%) and recall (90.5%), also achieving the highest median F1-score (88.2%) among all models. Importantly, it showed AUC of 92.2%, confirming that the ensemble approach enhanced predictive balance and stability without compromising classification performance.

SpoMAG demonstrates high specificity in predicting the absence of sporulation in non-Firmicutes genomes

To evaluate the specificity of SpoMAG, we applied the model to a dataset of 496 high-quality MAGs from bacterial phyla outside Firmicutes, including genomes derived from cattle (n = 136), poultry (n = 126), human (n = 105), and swine (n = 129). As expected, SpoMAG assigned a predicted sporulation probability of zero to all MAGs (Fig. 3, Table S3), demonstrating SpoMAG’s high specificity and its ability to correctly distinguish non-sporulating genomes, even in complex metagenomic datasets. This performance underscores SpoMAG’s robustness and affirms its applicability for large-scale genome annotation across heterogeneous microbiomes.

Sporulation potential across Bacilli-class MAGs reveals order-specific patterns

A total of 51 high-quality MAGs from the Bacilli class were recovered across cattle (n = 10 MAGs), poultry (n = 15 MAGs), human (n = 20 MAGs), and swine (n = 6 MAGs). SpoMAG predictions revealed distinct order-level sporulation patterns consistent with established phenotypic characteristics of Bacilli lineages (Fig. 4, Table S4).

In cattle microbiota, one MAG from the order Paenibacillales (family Paenibacillaceae, genus Cohnella) was predicted to be spore-forming with 99% probability (Fig. 4), consistent with the established sporulation capability of this order. MAGs from other orders, including RF39, ML615J-28, Izemoplasmatales, Erysipelotrichales, and Acholeplasmatales, were all assigned 0% sporulation probability, corroborating their non-sporulating status (Table S4).

Poultry-derived MAGs exhibited a similar trend, with four MAGs from Paenibacillales (family Paenibacillaceae with three from Paenibacillus and one from Cohnella) predicted as spore-formers with 99% probability (Fig. 4). Two Bacillales MAGs (family Planococcaceae, genus Rummeliibacillus) showed an average sporulation probability of 81.1%. All remaining MAGs from Erysipelotrichales and Lactobacillales scored 0%.

Among human-derived MAGs, two genomes were predicted as sporulating: one from Bacillales (Planococcaceae, genus Ureibacillus) and another from DSM-1321 family (genus Pradoshia), both with average sporulation probability of 93%. All other human MAGs, including those from Erysipelotrichales and Lactobacillales, were predicted as non-spore-forming.

Notably, swine samples lacked Bacillales and Paenibacillales representatives. All six swine MAGs belonged to orders traditionally considered non-sporulating (e.g., Acholeplasmatales, Erysipelotrichales, Lactobacillales, and RFN20), and all were assigned a 0% probability by SpoMAG.

Analysis of gene content further supported SpoMAG’s predictions, with MAGs from Paenibacillales and Bacillales containing a higher number of sporulation-associated genes compared to other Bacilli orders (Fig. 5, Table S5). All MAGs from other Bacilli orders scored 0% of sporulation probability, reinforcing the tool’s biological consistency and specificity.

To explore within-group variability in gene profiles, we performed a multivariate dispersion analysis (Betadisper) in the predicted spore-formers and non-spore-formers from Bacilli. The analysis revealed significant differences in gene profile dispersions (F = 12.62, p =0.0007), indicating distinct within-group variability in gene composition.

Sporulation potential across Clostridia-class MAGs reveals order-specific patterns

We recovered 262 high-quality MAGs from the Clostridia class, distributed across cattle (n = 53 MAGs), poultry (n = 58 MAGs), humans (n = 95 MAGs), and swine (n = 56 MAGs). SpoMAG predictions revealed consistent and order-specific patterns in predicted sporulation potential (Fig. 6, Table S6).

In human samples, seven MAGs from the order Acetivibrionales were classified as spore-forming, with an average predicted probability of 98.5%. All belonged to the family Acetivibrionaceae, including two from the genus Ruminiclostridium and five unclassified at the genus level. One MAG from the same order, classified under (DSM-8532 family) was predicted to be non-sporulating (0% probability). Within Christensenellales, only one of 19 MAGs was predicted as a spore-former (Christensenellaceae, unknown genus/species). Among the remaining, mostly from poorly characterized families and the Borkfalkiaceae, were classified as non-sporulating. All Clostridiales MAGs (n = 3) were predicted to be spore-forming, with three from the Clostridiaceae family (one identified as Clostridium and two unclassified at the genus level). In contrast, all 24 MAGs from Lachnospirales were classified as non-sporulating, including genera within Lachnospiraceae such as Herbinix, Lachnoclostridium, Butyrivibrio, UBA4285, RUG115, Blautia, Coprococcus, Enterocloster, KLE1615, CAG-194, Faecalimonas, Anaerostipes, Mediterraneibacter, CAG-127, Dorea, Acetatifactor, CAG-603, Bariatricus, and Anaerocolumna, all lacking species-level classification.

The order Monoglobales included one MAG predicted as a spore-former (family UBA1381, genus CAG-41) and one as a non-spore-former (family Monoglobaceae, genus Monoglobus). Among 34 MAGs from Oscillospirales, one was predicted as spore-forming (family Ruminococcaceae, genus UBA3818), potentially representing a novel spore-forming species. All MAGs from Peptostreptococcales, Tissierellales, and UMGS1883 were predicted to be non-spore-forming.

In swine, SpoMAG predicted sporulation in eight Acetivibrionales MAGs, including Acetivibrionaceae (n = 7) and DSM-8532 (n = 1), all unclassified at the species level. Additionally, three Christensenellales MAGs were predicted as spore-formers, including one from Christensenellaceae and two from poorly defined families (CAG-74 and SZUA-584). Two Clostridiales MAGs (Clostridiaceae, genus Clostridium) and one Oscillospirales MAG (Ruminococcaceae, unknown genus and species) were also identified as spore-formers. All MAGs from DUPQ01, Lachnospirales, Peptostreptococcales, Saccharofermentanales, and Tissierellales were predicted as non-sporulating.

In cattle, SpoMAG identified six spore-forming MAGs from Acetivibrionales (three from Acetivibrionaceae, three from DSM-8532), two from Christensenellales (families Christensenellaceae and SZUA-584), and three from Clostridiales (two Clostridiaceae, one Caloramatoraceae, unknown genus). Notably, two MAGs from Lachnospirales (Herbinix genus, family Lachnospiraceae) were also predicted to be spore-forming. One MAG from Peptostreptococcales (Anaerovoracaceae, genus UBA7709) was classified as a spore-former.

In poultry, SpoMAG predicted all seven Acetivibrionales MAGs as spore-formers (six from Acetivibrionaceae and one from DSM-8532). All MAGs from Clostridiales (n = 2) were also classified as spore-forming, along with one Christensenellales (Christensenellaceae, unknown species) and one Peptostreptococcales (Natronincolaceae, genus Alkaliphilus). Two MAGs from Lachnospirales (genus Herbinix) were also predicted as spore-formers.

Gene presence/absence patterns are insufficient to distinguish sporulation in Clostridia MAGs

We next investigated whether presence/absence patterns of sporulation-associated genes could distinguish between predicted spore-forming and non-spore-forming MAGs. Unlike Bacilli, binary gene profiles were insufficient to reliably separate sporulation capability in Clostridia (Fig. 7, Table S7). Although some differences in gene content were observed, there was substantial overlap between phenotypic groups, suggesting that gene presence alone may indicate a latent potential to sporulate rather than active functionality. Importantly, the presence of sporulation-associated genes reflects genetic potential, not active gene expression at the time of sampling.

To evaluate gene content variation in more detail in the Clostridia class, we performed PCA on sporulation gene profiles (Fig. 8). In cattle (Fig. 8A), poultry (Fig. 8B), human (Fig. 8C), and swine (Fig. 8D), there was considerable overlap between predicted spore-forming and non-spore-forming MAGs. Although distinct ellipses were generated for each group, indicating some variation, no clear separation was observed in the multivariate space, further demonstrating that unsupervised ordination methods still insufficient to resolve sporulation status even for the predicted Clostridia MAGs.

We then assessed within-group dispersion using betadisper analysis. Statistically significant differences in dispersion were observed in the spore-forming and non-spore-forming MAGs for all hosts: cattle (F = 7.55, p = 0.0083), swine (F = 43.77, p < 0.001), human (F = 14.19, p = 3 ×10⁻⁴), and poultry (F = 13.42, p = 6 ×10⁻⁴), indicating that the two groups differ in internal variability of sporulation gene content, suggesting distinct levels of genomic variability related to sporulation potential within each group.

SpoMAG reveals potentially sporulating Firmicutes species shared across hosts

Our comprehensive pairwise average nucleotide identity (ANI) analysis of high-quality Firmicutes MAGs revealed distinct distribution patterns among genomes predicted to be sporulating (Table S8). Most MAGs, both predicted sporulating and non-sporulating, exhibited ANI values below the 95% species threshold. A notable proportion clustered at 74% ANI, indicating considerable species-level diversity among the recovered genomes (Fig. 9A).

Despite this diversity, SpoMAG-based predictions revealed seven cases of species-level sharing (ANI value > 95%) among spore-forming MAGs recovered from different hosts (Fig. 9A). These included two shared species between poultry and swine, two between poultry and humans, one between swine and humans, and two between swine and cattle (Fig. 9B). Notably, a species from the Acetivibrionaceae family exhibited ∼99% ANI across poultry, swine, and human hosts.

MAGs exhibit distinct sporulation gene profiles at the Class and Order levels across hosts

To explore patterns in sporulation gene content across hosts, we performed PCA on the binary presence/absence matrices of sporulation-associated genes in MAGs predicted to be spore-forming by SpoMAG. The ordination revealed class- and order-specific clustering, as well as host-associated trends.

Human-derived MAGs (Fig. S1, Table S9) had the first two principal components explaining 42.33% of the total variance. The two Bacilli MAGs were clearly separated from Clostridia. Within Clostridia, seven Acetivibrionales MAGs (all from Acetivibrionaceae family) clustered tightly, suggesting conserved sporulation gene profiles. Notably, two of these MAGs are assigned to Ruminiclostridium genus (bin.488_human and bin.208_human).

In swine MAGs (Fig. S2, Table S9), only Clostridia were predicted as spore-formers. The first two components explained 40.57% of the variance and revealed two distinct clusters: one composed of eight Acetivibrionales MAGs, and another of three Christensenellales MAGs. The remaining MAGs were mixed into a cluster of two MAGs from Clostridiales and one from Oscillospirales.

Cattle-derived MAGs (Fig. S3, Table S9) had PC1 and PC2 explaining 34.88% of the variance. Six Acetivibrionales MAGs formed a cluster, adjacent to two Lachnospirales MAGs, indicating partial similarity in gene profiles. A single Bacilli MAG (order Paenibacillales) appeared well separated from Clostridia. A separate group included MAGs from Christensenellales (n = 2), Clostridiales (n = 3), and Peptostreptococcales (n = 1).

In poultry-derived MAGs (Fig. S4, Table S9), the first two components explained 49.20% of the variance. Two Bacilli MAGs and four from Paenibacillales formed distinct clusters outside the main Clostridia group. Within Clostridia, the mixed cluster included MAGs from Acetivibrionales (n = 7), Christensenellales (n = 1), Clostridiales (n = 2), Lachnospirales (n = 2), and Peptostreptococcales (n = 1), again reflecting compositional variability.

When combining all 63 MAGs predicted to be spore-formers, the PCA revealed a separation between Bacilli and Clostridia, with the first two components explaining 32.21% of the variance (Fig. 10, Table S9). Five Paenibacillales MAGs grouped together, including one MAG from Bacillales. A well-defined cluster of 28 Acetivibrionales MAGs confirmed the conservation of sporulation gene profiles within this order, even among hosts. The remaining Clostridia MAGs formed a broader group, representing the orders Christensenellales (n = 7), Clostridiales (n = 10), Lachnospirales (n = 4), Monoglobales (n = 1), Oscillospirales (n = 2), and Peptostreptococcales (n = 2).

Consistent genes present across the 63 predicted sporulating MAGs

We identified nine genes consistently present in all 63 MAGs predicted as spore-formers by SpoMAG (Fig. 11). These include: pth, yaaT, spoIIAB, spoIIIAE, spoIIIAD, ctpB, ftsW, spoVD, and lgt. Although an additional 13 genes were found in 62 of the 63 genomes, we focused on those shared across all genomes to ensure consistency in gene presence.

Functionally, pth and yaaT are associated with early-stage regulation, while spoIIAB, spoIIIAE, and spoIIIAD contribute to initiation and engulfment; ctpB and ftsW are linked to the SigE regulon; spoVD is involved in spore cortex formation; and lgt plays a role in spore germination.

SHAP analysis identified 16 genes consistently contributing to sporulation prediction in both Random Forest and Support Vector Machine models (Fig. 12). These include: lytH, cotP, spoIIIAG, spoIIR, spoVAD, gerC, yabP, yqfD, gerD, spoVAA, gpr, ytaF, gdh, ypeB, spoVID, and ymfJ. All genes exhibited positive SHAP values, indicating their influence on model classification decisions.

The genes lytH, yabP, and yqfD are associated with the spore cortex; cotP and spoVID are involved in the spore coat; spoIIIAG in engulfment; spoIIR, ymfJ are part of the SigF regulon; spoVAD, spoVAA are related to the SigG regulon; ytaF in SigE regulon; and gerC, gerD, gpr, gdh, ypeB are involved in germination.