The animal study was reviewed and approved by the Animal Care and Use Committee of Huazhong Agricultural University. The overall design of the experiment is illustrated in Figure 1. This experiment was conducted at the Wuli Town pig farm, Xin Yang of Tecon Biology Co. A total of 1,000 healthy sows of comparable body weight (DanBred Landrace × DanBred Yorkshire, parities 1–2) were selected and randomly assigned to two dietary treatment groups after artificial insemination. Sows in the control (CON) group received a diet primarily composed of beet pulp and barley as fiber sources, while those in the experimental group were provided a diet supplemented with 1% functional fiber [DF, consisted of 85.7% resistant starch (Hangzhou, China) and 14.3% guar gum (Yunzhou, China)]. The functional fiber was provided by Wuhan Fanbo Biotechnology Co., Ltd. (Wuhan, China). This was achieved by partially replacing the energy-yielding ingredients in the control diet while ensuring that the net energy levels in the CON and DF diets was comparable. All sows were housed in individual gestation stalls (2.2 m length × 0.65 m width × 0.6 m height). Sow culling data were recorded during gestation (Supplementary Table S1).

All experimental diets were formulated to meet or exceed the nutritional requirements for gestating and lactating sows as outlined in National Research Council (NRC) (2012) (Supplementary Table S2). The gestation diet was provided at 2.6 kg/d from days 0 to 30, 2.2 kg/d from days 31 to 90, and 2.8 kg/d from days 91 to 109. From day 110 until the day before parturition, sows received 2.0 kg/d. Post-farrowing, all sows received the same lactation diet throughout the lactation period, and sows were fed 1.0 kg/d on the first day postpartum, with the amount increased by 1.0 kg/d over the next 7 days. Thereafter, sows had ad libitum access to feed until weaning.

Sow performance and reproductive parameters were recorded during gestation and lactation. These included total piglets born, live-born piglets, birth weight, stillbirths, mummified fetuses, number of piglets weaned, weaning weight, and the interval to return to estrus. Reproductive performance data were collected, encompassing the number of total piglets born, piglets born alive, healthy piglets (birth weight ≥ 0.8 kg), piglets with intrauterine growth retardation (IUGR) (birth weight < 0.8 kg), piglets born dead (stillborn, mummified, crushed, or abnormal), litter birth weight, and individual piglet birth weight.

Fecal specimens were systematically collected from 20 randomly selected sows per group at critical developmental stages. In detail, fresh fecal samples (~5 g/sow) were collected at gestation days 0 (G0 d), 30 (G30 d), 90 (G90 d), and 109 (G109 d), as well as lactation days 3 (L3 d) and 14 (L14 d). All samples were immediately frozen in liquid nitrogen and stored at −80 °C until further analysis.

Microbial co-occurrence networks were constructed at the genus level using SparCC with 100 bootstrap iterations and a significance threshold of p < 0.05, as implemented in the SpiecEasi R package (v1.1.2) (Friedman and Alm, 2012). Only correlations with |r| > 0.5 were retained to build the final networks (Chen et al., 2023). Network topological properties were calculated using Cytoscape (v.3.10.2) (Shannon et al., 2003). Network visualizations were generated in Gephi (v.0.9.2) (Bastian et al., 2009). Natural connectivity analysis was performed following the method described by Pan et al. (2021), and the visualization of natural connectivity dynamics was generated using the ggplot2 package (version 3.5.1) in R. Natural connectivity (Wu et al., 2010) was introduced to describe network stability differences, the estimation of natural connectivity was based on the following algorithm:

where ave.(k) represents the natural connectivity, N is the number of nodes in the network, and kᵢ denotes the eigenvalue of the adjacency matrix.

The performance and reproductive parameters of sows were analyzed separately for the full dataset (Table 1) and the microbiome subset (Table 2). For the Table 1 (CON: n = 405; DF: n = 411), a one-way analysis of variance (ANOVA) was performed using the GLM procedure in SAS (version 9.4; SAS Institute Inc., Cary, NC, United States), with dietary treatment (CON vs. DF) as the fixed effect and individual sows as the experimental units. Parity was not included as a covariate due to the large sample size and the relatively balanced distribution between groups (CON: 1.59 ± 0.68; DF: 1.79 ± 0.66).

For exploratory microbiome and metabolome analyses, a subset of 32 sows (n = 16 per group) was selected based on extreme reproductive performance. Following an initial random sampling of 20 sows per group, the 16 sows with the lowest total number of piglets born (CON group) and the 16 sows with the highest total number of piglets born (DF group) were chosen to amplify potential differences associated with treatment. Comparisons between CON and DF groups in this subset were performed using an unpaired two-tailed Student’s t-test (GraphPad Prism, version 8.3.0). It should be noted that selecting sows based on extreme reproductive performance may introduce selection bias. This exploratory subset analysis aimed to characterize diet-associated microbial differences, and the results should be interpreted as supportive, hypothesis-generating evidence rather than as definitive inferences at the population level.

Alpha and beta diversity of microbial data were calculated using the vegan package (version 2.6–8) in R software (version 4.1.3). PCoA was performed using Bray-Curtis distance metrics. PERMANOVA was used to evaluate factors shaping microbiota by using the Adonis function of the vegan package (999 permutations). In addition, the significantly different genera were selected by the linear discriminant analysis effect size (LEfSe) method (Segata et al., 2011), performed via Wekemo Bioincloud¹ (Gao et al., 2024).

Orthogonal partial least squares discriminant analysis (OPLS-DA) was also introduced as a supervised model to reveal the metabolome variation between groups, using the “plsda” function in R package “ropls” (version 1.26.4). The analysis was conducted on untargeted metabolomics. The differential metabolites with false discovery rate (FDR) < 0.05 and the absolute value of |log2 (fold change)| > 1 were visualized by a volcano plot. The results were visualized using the ggplot2 (version 3.5.1) package for R. Statistical significance was determined when p value was < 0.05.

As shown in Table 1, compared with the CON group, the numbers of total born and healthy piglets were significantly improved (p < 0.01) in the DF group, accompanied by an increase in litter birth weight (p < 0.05). The number of IUGR piglets was significantly reduced (p < 0.05) in DF group. However, the numbers of white and black stillbirths were also higher in the DF group (p < 0.05). No significant differences were observed between the CON and DF groups in the numbers of piglets born alive, mummified fetuses, abnormal piglets, or average birth weight per piglet (p > 0.05).

In order to study the composition of gut microbiota in sows with different reproductive performance, fecal samples from 32 sows with significantly different reproductive outcomes (Table 2) were subjected to 16S rRNA gene sequencing. Global changes in the gut microbial community were observed by assessing diversity. Microbial diversity (Shannon index) and richness (the Chao1 estimator, Chao1) initially increased, subsequently declined, and later exhibited a recovery (Figures 2A,D). To assess the effects of dietary treatment and time on microbial α-diversity, linear mixed-effects models were employed using the nlme package (version 3.1–162) in R. These models were chosen because they appropriately account for the longitudinal and repeated-measures design of the study. The analysis revealed microbial α-diversity exhibited significant changes over time (p < 0.0001). Dietary treatment also influenced the microbiota, showing a trend toward significance for Shannon diversity (p = 0.051; Figures 2A,B) and a significant effect on richness as estimated by Chao1 (p = 0.038; Figures 2D,E). The alpha diversity of DF microbiota was higher than CON on G30 d and L14 d (p < 0.05, Figures 2B,E). Incorporating the changes in diversity over time, the diversity in the DF group significantly increased at G30 d and remained stable, while the CON group gradually increased with the progression of pregnancy (p < 0.05, Figure 2C). Regarding the Chao1 index, the DF group exhibited an increase in the median value from G0 d to G30 d, although the change was not statistically significant, while the CON group remained relatively stable. Consequently, the results showed a significant difference between the two treatments at G30 d, with no significant difference observed at G109 d. These findings suggest that sows in the DF group rapidly achieved a stable gut microbial community, indicating that DF supplementation enhanced the stability of the gut microbiota. On the other hand, microbial diversity significantly decreased in both treatment groups after farrowing in sows (p < 0.05, Figures 2C,F); however, the diversity in the DF group significantly increased from L3 d to L14 d and showed no significant difference compared to G0 d, whereas the CON group, despite a significant increase at L14 d compared to L3 d, remained significantly lower than G0 d (Figure 2C). These results indicate that DF treatment facilitated, to a certain degree, the recovery of the gut microbiota during lactation.

For beta diversity, Principal coordinate analysis (PCoA) based on Bray–Curtis dissimilarity revealed that the structure of the gut microbiota was distinctly differentiated by different stages of gestation (Figure 2G), particularly with the gestation and lactation periods clearly clustering into two separate groups (Supplementary Figure S1A). And PCoA exhibited significant segregation except on L3 d (PERMANOVA, p < 0.05, Supplementary Figure S1B). PC1 and PC2 analyses revealed similar microbial trajectory patterns between sows in the CON and DF group (Figures 2H,I). Consistent with the changes in α-diversity, at G30 d, both PC1 and PC2 values in the DF group significantly decreased and remained stable at subsequent gestation time points, whereas the CON group exhibited a continuous decline throughout the gestation period (Figures 2H,I).

The relative abundances of the top 10 phyla (Figure 3A) and top 20 genera (Figure 3B) were calculated as the overall mean across all samples, providing a comprehensive overview of fecal microbial composition. The gut microbial community was dominated by Firmicutes (68.33%) and Bacteroidota (23.65%), followed by Spirochaetota (3.78%), Proteobacteria (1.94%) and Actinobacteriota (0.6%) (Figure 3A). Likewise, Christensenellaceae_R-7_group (9.58%), UCG-002 (5.72%), Clostridium_sensu_stricto_1 (4.86%), un_f__Lachnospiraceae (4.51%), Prevotellaceae_UCG-001 (4.08%), Ruminococcus (3.74%), Rikenellaceae_RC9_gut_group (3.74%), Treponema (3.61%), norank_f__Eubacterium_coprostanoligenes_group (3.32%) and norank_f__Muribaculaceae (3.31%) were the 10 most abundant genera (Figure 3B). In order to evaluate the impact of treatment on the composition of the gut microbiota at each gestation stage, we compared the relative abundance of the Firmicutes and Bacteroidetes phyla, as well as the Firmicutes-to-Bacteroidetes ratio, in the CON and DF groups. Significant alterations were observed only on gestation day 109 (Supplementary Figure S2). Specifically, compared with the CON group, the abundance of Firmicutes was markedly increased in the DF group (p < 0.001, Supplementary Figure S2A), whereas the abundance of Bacteroidetes was significantly decreased (p < 0.05, Supplementary Figure S2B). Consequently, the Firmicutes-to-Bacteroidetes ratio was significantly elevated in the DF group (p < 0.01, Supplementary Figure S2C), indicating that the gut microbial composition was markedly altered.

To evaluate the temporal dynamics of the microbial community across reproductive stages and dietary treatments, we applied a linear mixed-effects model at the genus level (Supplementary Table S3). A total of 43 genera exhibited significant differences in response to dietary treatment, 130 genera showed significant temporal variation across reproductive stages, and 23 genera were significantly influenced by the treatment × stage interaction. Of these, 34 genera were simultaneously affected by both treatment and stage (Figure 3C; Table 3), and were considered key genera.

Microbial differences between dietary treatments during gestation were identified using linear discriminant analysis (LDA) coupled with effect size estimation (LEfSe). As shown in Supplementary Table S4, using an LDA score > 2 as the threshold, 28 genera were significantly different at G0 d, 23 genera at G30 d, 24 genera at G90 d, and 55 genera at G109 d. These results indicated that G109d had the largest number of differentially abundant genera, reflecting the greatest compositional differences between the CON and DF groups. At G109d specifically, these 55 genera comprised 18 taxa significantly enriched in the CON group and 37 taxa enriched in the DF group. Therefore, the differential genera identified by LEfSe at G109 d were further integrated with the key genera revealed by the linear mixed-effects model (Figure 3D). In total, 17 genera were identified as differential taxa by both approaches, and their relative abundances were shown in Figure 3E and Supplementary Figure S3. The genera with relative abundance > 0.1% were shown in Figure 3E. During early gestation (G0-G30d), distinct microbial differences between the CON and DF groups were already evident. Treponema showed a significant increase in the CON group starting from G30d (p < 0.05). Streptococcus showed a tendency toward higher abundance at G0d (p < 0.1) in DF group. Compared with the CON group, Family_VIII_AD3011_group exhibited significantly higher abundance at G0d (p < 0.05) and G30d (p < 0.05), while norank_f_Eubacterium_coprostanoligenes_group showed a tendency toward higher abundance at G0d (p < 0.1) and a significant increase at G30d (p < 0.05). Moreover, NK4A214_group was significantly increased in the DF group at G30d (p < 0.05), indicating an early microbial response to DF diet.

As gestation progressed into late pregnancy (G90-G109d) and lactation (L3-L14d), these microbial differences became more pronounced. Treponema remained significantly enriched in the CON group during late gestation (p < 0.05). In the DF group, Streptococcus was significantly enriched during late pregnancy (p < 0.05). And Oscillospira showed a tendency toward higher abundance at G90d (p < 0.1) and a significant increase at G109d (p < 0.05) compared to the CON group. Family_VIII_AD3011_group and norank_f_Eubacterium_coprostanoligenes_group continued to show higher abundances in the DF group and were both significantly enriched at G109d (p < 0.01 and p < 0.001, respectively), with norank_f_Eubacterium_coprostanoligenes_group also exhibiting significantly higher abundance at L14d (p < 0.01). Furthermore, NK4A214_group maintained a relatively high abundance throughout the entire reproductive period, showing a trend toward higher abundance at G90d (p < 0.1), a significant increase at G109d (p < 0.05), and a further increase at L14d (p < 0.01) compared with the CON group. Overall, the temporal dynamics of these genera suggest that Treponema and NK4A214_group may serve as key differential genera between the CON and DF groups, reflecting distinct microbial trajectories under different dietary treatments.

The G109d time point was selected for association analysis with reproductive performance because previous study has indicated that the gut microbiome undergoes the most profound remodeling during late gestation, which supports metabolic adaptations beneficial for fetal development and litter size (Koren et al., 2012; Chen et al., 2023). Spearman correlation analysis was performed to investigate the relationship between the bacterial genera abundance at G109 d and reproduction performance of sows (Figure 3F). Among the genera with high abundance, the relative abundance of Treponema was negatively correlated with the numbers of total born (p < 0.01) and healthy piglets (p < 0.05, Figure 3F). The relative abundance of Streptococcus showed a highly significant positive correlated with the numbers of total born, born alive and healthy piglets (p < 0.01). The relative abundance of Family_VIII_AD3011_group was highly positively correlated with the numbers of total born and healthy piglets (p < 0.01), and positively correlated with the numbers of born alive (p < 0.05, Figure 3F). In addition, NK4A214_group showed a positive correlated with the numbers of total born and born alive (p < 0.05), and a highly significant positive correlated with the numbers of healthy piglets (p < 0.01, Figure 3F).

Moreover, we applied the same analytical strategy at the family level (Supplementary Figure S4; Supplementary Tables S5, S6). A total of 13 families were significantly influenced by both treatment and stage (Supplementary Figure S4A; Supplementary Table S5). And the Venn diagram showed 7 families were identified as differential taxa by both the linear mixed-effects model and LEfSe analysis (Supplementary Figure S4B), and their relative abundances were shown in Supplementary Figure S4C. The correlation heatmap indicated that Spirochaetaceae was significantly negatively correlated with the numbers of total born (p < 0.01) and healthy piglets (p < 0.05, Supplementary Figure S4D), which was enriched in the CON group (Supplementary Figure S4C). In addition, among the 6 families that were enriched in the DF group (Supplementary Figure S4C), Oscillospiraceae, the highest abundance with a median relative abundance of 14.74%, was markedly positively correlated with the numbers of total born and healthy piglets (p < 0.05, Supplementary Figure S4D). Taken together, these findings suggest that Spirochaetaceae, especially Treponema, are negatively associated with reproductive performance, whereas Oscillospiraceae and its member NK4A214_group enriched in DF group may play a positive role in improving reproductive performance of sows.

Furthermore, based on the patterns of diversity dynamics (Figure 2) and differential taxa analysis (Figure 3), we consider G30 d and G109 d to be the most important stages during gestation. Therefore, Co-occurrence networks were constructed using SparCC to evaluate the network properties and stability at these two time points (Figure 4). The network stability at G30 d and G109 d was shown in Figures 5A,B, respectively, where the DF group consistently exhibited higher stability than the CON group, with an obvious increase from G30 d to G109 d. Co-occurrence networks constructed for the CON and DF groups at G30 d and G109 d were illustrated in Figures 4C–F, and the topological properties of the networks were presented in Table 4. We found that the number of nodes among the four networks was comparable, but the numbers of edges and average number of neighbors were obviously different. To be more specific, compared with CON group, both metrics were higher in DF group, with more edges (521 vs. 425 at G30; 737 vs. 542 at G109) and a greater average number of neighbors (6.83 vs. 5.65 at G30; 10.5 vs. 7.13 at G109), which signify the connectivity of network in the DF group was stronger than CON group. Moreover, the clustering coefficient, network density, and network centralization exhibited similar patterns, showing comparable values between the groups at G30 d, whereas all three metrics were higher in the DF group than in the CON group at G109 d (clustering coefficient: 0.358 vs. 0.273; network density: 0.076 vs. 0.047; network centralization: 0.186 vs. 0.106). Collectively, these results indicate that DF treatment enhanced microbial network connectivity and clustering, leading to tighter internal organization, which may partly underlie the observed differences in network stability.

To further characterize the core network structure, the importance of each node was evaluated using MCC scores (Table 5; Supplementary Table S7). In the CON group, g__NK4A214_group exhibited very high MCC scores at both G30 d and G109 d, indicating a consistently central role in the microbial network. And g__Treponema also showed relatively high MCC rankings in the CON networks, particularly at G109 d. By contrast, in the DF group, both the g__NK4A214_group and the g__Treponema displayed substantially lower MCC rankings at G30 d (ranked 39/147 and 41/147, respectively), suggesting reduced centrality during the early stages of gestation. By G109 d, g__Treponema showed an increased MCC ranking (15/144), whereas g__NK4A214_group remained low in centrality, indicating treatment-specific differences in core network organization.

To characterize the metabolome signatures of sows on G109 d, the fecal metabolic profiles were analyzed by non-targeted LC–MS/MS metabolomics. Overall, stool metabolite profiles were found to significantly differ between the CON and DF group by OPLS-DA (p = 0.001; Figure 5A). Compared with CON group sows, DF group sows displayed increased concentrations of 68 metabolites and decreased concentrations of 83 metabolites (Figure 5B). KEGG analysis revealed that most of the differential metabolites in DF group compared with CON group were enriched in bile secretion, prolactin signaling pathway, sphingolipid signaling pathway, breast cancer, alcoholism, prostate cancer and cocaine addiction (Figure 5C; Supplementary Table S8). Following adjustment for multiple testing, two pathways remained significant (FDR < 0.05). Heatmaps of the metabolites in the bile secretion and prolactin signaling pathway were presented in Figures 5D,E, respectively. Enriched metabolites in the DF group sows included lithocholic acid, chenodeoxycholate, cholic acid, glycochenodeoxycholate, and progesterone (Figures 5D,E). Notably, bile secretion was the most significantly altered pathway, with subsequent findings revealing that a series of bile acids (BAs), including HDCA, CDCA, GCDCA, CA, LCA, UDCA, and γ-MCA, were significantly increased in the DF group, whereas TUDCA was significantly reduced in the DF group (Figure 5F). The results of correlation between sow reproduction performance and BAs showed that TUDCA, which was abundant in the CON group, was negatively correlated with the numbers of total born (p < 0.01) and healthy piglets (p < 0.05), GCDCA was positively correlated with the numbers of total born (p < 0.01), born alive (p < 0.05) and healthy piglets (p < 0.01), CDCA and HDCA were positively correlated with the numbers of total born (p < 0.05, Figure 5G).

To further elucidate the potential interactions between gut microbes and metabolites associated with reproductive performance, spearman correlation analysis was conducted using the bacterial taxa and metabolites that showed significant associations with reproductive performance (Figure 6). Among these associations, particular attention was paid to Treponema and NK4A214_group, which were previously identified as key differential genera in the CON and DF groups, respectively. Treponema showed negative correlations with fecal progesterone (p < 0.01), carnitine (p < 0.001) and GCDCA (p < 0.001), but positive correlations with hydrocortisone (p < 0.05), fluvastatin (p < 0.05), levofloxacin (p < 0.05) and Beta-estradiol (p < 0.01) (Figure 6A). In addition, Treponema was also negatively correlated with fecal HDCA (p < 0.05), whereas it was positively correlated with TUDCA (p < 0.01) (Figure 6B). By contrast, NK4A214_group was negatively correlated with hydrocortisone and fluvastatin (p < 0.01) but positively correlated with GCDCA (p < 0.05) and progesterone (p < 0.001) (Figures 6A,B).