Gut microbiome responds compositionally and functionally to the seasonal diet variations in wild gibbons

Abstract

Wild animals may encounter multiple challenges especially food shortage and altered diet composition in their suboptimal ranges. Yet, how the gut microbiome responds to dietary changes remains poorly understood. Prior studies on wild animal microbiomes have typically leaned upon relatively coarse dietary records and individually unresolved fecal samples. Here, we conducted a longitudinal study integrating 514 time-series individually recognized fecal samples with parallel fine-grained dietary data from two Skywalker hoolock gibbon (Hoolock tianxing) groups populating high-altitude mountainous forests in western Yunnan Province, China. 16S rRNA gene amplicon sequencing showed a remarkable seasonal fluctuation in the gibbons' gut microbial community structure both across individuals and between the social groups, especially driven by the relative abundances of Lanchnospiraceae and Oscillospiraceae associated with fluctuating consumption of leaf. Metagenomic functional profiling revealed that diverse metabolisms associated with cellulose degradation and short-chain fatty acids (SCFAs) production were enriched in the high-leaf periods possibly to compensate for energy intake. Genome-resolved metagenomics further enabled the resolving metabolic capacities associated with carbohydrate breakdown among community members which exhibited a high degree of functional redundancy. Our results highlight a taxonomically and functionally sensitive gut microbiome actively responding to the seasonally shifting diet, facilitating the survival and reproduction of the endangered gibbon species in their suboptimal habitats.

Fig. 1. Microbiome composition of 514 gibbon fecal samples.

a Home ranges of two social groups Nankang (orange, NK) and Banchang (purple, BC) during sample collection. b–e Relative abundances of phyla (b NK; c BC) families (d NK; e BC). Samples were sorted by the sampling time. Each column represented one sample and samples from different months were separated by the blank space. Significance was calculated by the Wilcoxon rank-sum test with samples classified by social groups.

Fig. 2. Covariation of diet and gut microbiota across the social groups.

a and b Hierarchical clustering with the dietary composition of the social groups (a NK; b BC). Detailed information was provided in Supplementary Table 7. c The correlation between the first principal component (PC1) of between-sample dissimilarity and the proportion of leaf in the diet with colors representing different subjects (orange, NK; purple, BC). The p-values were obtained from linear regression analysis. d Gut microbial composition of NK and BC social groups in the high-fruit and high-leaf periods (the dominant families were shown). e and f Loading PC1 scores of each ZOTU within social groups (e, NK; f, BC) with a loading score of >0.05 and <0.05 highlighted.

Fig. 3. The relationship between the five properties (Size, Total link, Average degree, Edge connectivity, and Average path length) of gut microbiome co-occurrence network and the diet.

The points and error bars on the panels represented mean ± standard deviation after subsampling 999 times. Colors denoted different social groups (orange, NK; purple, BC). The adjusted R² values for linear regression were presented. The p-values were obtained from linear regression analysis.

Fig. 4. Covariation of diet and gut microbiota across the individuals.

a NMDS (center panel) of the Bray–Curtis dissimilarity based on the microbial community composition. Different colors corresponded to different gibbon individuals in the top panel. Boxplots (top panel) indicated the distribution of each individual along the first axis of NMDS. The boxplots on the left panel depicted different diet types according to their gut microbial community placement on the second axis of NMDS. Different letters represented significant differences (ANOVA followed by Duncan’s test). b The number of the diet-specific responsive ZOTUs. Colors indicated the types of responsive ZOTUs (yellow, responding to fruit; blue, responding to leaf). c Within-individual correlations between diet and the relative abundances of ZOTUs (Spearman’s correlation r > 0.6; FDR-corrected p-value < 0.01). The colored circles in the middle represented different gibbon individuals. The tree outside was obtained by 16S rRNA gene sequences, and the colors represented diverse microbial families. The lines with different colors indicated the types of responses.

Fig. 5. Overview of the functional profiles depicting the pathways enriched in different diet periods with colors differentiating gibbon subjects (orange, A2; purple, B2).

The diagrams showed metabolic pathways enriched during high-leaf (left) and during high-fruit (right). The colors of different solid lines corresponded to significantly differential metabolic pathways (Wilcoxon rank-sum test, and FDR-corrected p-value < 0.05). The genes with no significant difference between the two diets were marked with gray lines. The detailed information on the genes was summarized in Supplementary Table 9.

Fig. 6. CAZymes profile adapted to polysaccharide metabolism.

a and b Heatmap of the CAZymes profile in the gibbon gut metagenomes (a A2; b B2). CAZymes associated with cellulase, hemicellulase, pectinase, debranching enzymes, amylases, and oligosaccharide degradation were presented. CAZymes families with significant enrichment (Wilcoxon rank-sum test, and FDR-corrected p-value < 0.05) in the leaf- and fruit-dominated periods were marked. The standardized relative abundances of each CAZymes family were shown. c and d Sankey diagram describing the distribution of the top 10 microbial families associated with all selected CAZymes families ranked by the number of genes in gibbons A2 (c) and B2 (d). CAZymes were represented by different colors according to functional types and microbes were colored according to the associated families. The heights of the rectangles indicated the numbers of the CAZymes (left) and microbial species (right). The detailed information on the genes was summarized in Supplementary Table 10.

Fig. 7. Metabolic reconstruction of MAGs enriched in high-leaf diets.

Heatmap indicated the proportions of MAGs with the abilities for major polysaccharides degradation, sugar utilization, and fermentation. The MAGs were summarized into family-level taxonomic clades based on the phylogeny inferred by GTDB-Tk. Functional profiles were reconstructed inferred from the presence of the key marker genes within the metabolic pathway. Detailed information was provided in Supplementary Table 5.

Fig. 8. Prediction performance between different biotic levels in different sampling scales.

The comparison of cross-validation results (a–c) and prediction accuracies (d–f) between different biotic levels in different sampling scales. a–c Scatter plots showing the predicted and observed values of relative abundances of different biotic levels inferred from dietary compositions. Colors represented different biotic levels (dark gray: phylum; light gray: family; black: ZOTU). The diagonal line represented perfect prediction (predicted value = observed value). d–f The accuracy of random forest models in predicting relative abundances in different taxa. Significant differences among different datasets were denoted by lettering (P < 0.05, ANOVA followed Duncan’s test). The detailed information was summarized in Supplementary Table 6.