Stochasticity Prevails but Differs: Tissue‐Specific Assembly of Gut Microbiomes Across Seasons in an Amphibian Model

TL;DR

This study explores how gut microbiomes in Black-spotted frogs change across seasons and gut regions, finding that randomness plays a key role in their assembly.

Contribution

The study reveals tissue-specific and seasonal dynamics of amphibian gut microbiomes and identifies differing roles of stochastic and deterministic processes in microbiome assembly.

Findings

Gut microbiomes show tissue-specific succession with temporal fluctuations in structural discrepancies.

Stochastic processes dominate microbiome assembly, with deterministic factors more influential in large-intestine microbiomes.

Small-intestine microbiomes are more randomized yet centralized, with stable taxonomic composition over time.

Abstract

Gut microbiota generally undergoes dynamic remodeling in concert with multifaceted self‐regulation of amphibian hosts during key life stages, such as metamorphosis and hibernation. However, the spatiotemporal dynamics of amphibian gut microbiomes across the lifecycle remain poorly understood. In this study, we applied 16S rRNA gene amplicon sequencing to characterize the gut microbiomes of cultivated Black‐spotted frog (
Pelophylax nigromaculatus
) across seasons. The gut microbiomes exhibited tissue‐specific succession, and structural discrepancies between gut regions fluctuated temporally. Both small‐ and large‐intestine microbiomes showed temporal decay patterns in abundance‐unweighted intercommunity indices, but not in abundance‐weighted indices. Compared with large‐intestine microbiomes, small‐intestine microbiomes were more randomized yet more centralized in terms of amplicon sequence variants, particularly within Proteobacteria (especially Pseudomonas). The alpha diversity of small‐intestine microbiomes was comparatively lower, and their taxonomic composition was more stable over time. We further elucidated the assembly mechanisms of gut microbiomes by systematically analyzing dominant driving factors, ecological processes, phylogenetic traits, source‐sink relationships, and co‐occurrence networks. Stochastic processes played a dominant role in gut microbiome assembly, while deterministic processes (e.g., habitat filtering and microbial interaction) contributed more strongly to large gut microbiomes than to small gut microbiomes. Overall, this study provides insights into the ecological dynamics and assembly mechanisms of amphibian gut microbiomes across the lifecycle and may inform targeted microbiome modification for amphibian breeding and conservation.

Amphibian gut microbiomes exhibit a temporal‐decay pattern in abundance‐unweighted intercommunity indices rather than abundance‐weighted indices. Stochastic processes play a dominant role in the gut microbiome assembly. Deterministic processes contribute more to large‐intestine microbiomes than to small‐intestine microbiomes.

Figures

(a) Photographs of Black‐spotted frog ( Pelophylax nigromaculatus ) at three representative life stages. (b) Experimental field for rice‐frog co‐cultivation. (c) Research workflow in this study.

Structural comparisons across months or life stages using pairwise shared‐ASV ratios, alpha diversity, and beta diversity indices. The shared_abund_rat and shared_rich_rat represent abundance‐weighted and abundance‐unweighted pairwise shared‐ASV ratios, respectively. The significance level for Wilcoxon rank sum tests was set to 0.05.

Spearman correlations between inter‐group indices and their interval days. SG and LG microbiomes were analyzed separately, excluding female samples. Linear regression lines with 95% confidence intervals are shown in gray.

Structural comparisons between tissues and between sexes using pairwise shared‐ASV ratios, alpha diversity, and beta diversity indices. The significance level for Wilcoxon rank sum tests was set to 0.05.

Variance partitioning analysis (VPA) on beta diversity indices of gut microbiomes for three groups of factors. Each group included three dominant factors (Tissue, Month, Life stage) and one subdominant factor (Sex, SVL, or Temperature). Negative values are not shown. *The conditional significance of the deterministic factor < 0.05.

Relative importance of ecological processes deduced with bin‐based null model of R package iCAMP (ver 1.5.12) for: (a) gut microbiome groups; (b) ASV categories; (c) bins; (d) between male SG and LG microbiomes. A, abundant; DL, dispersal limitation; DR, drift; HD, homogenizing dispersal; HeS, heterogeneous selection; HoS, homogeneous selection; M, moderate; R, rare.

(a) Wilcoxon rank sum tests on the CMPD and CMNTD indices of contemporary and temporally adjacent gut microbiomes. The significance level was set to 0.05. (b) The adjusted R‐Squares of three groups of factors in variance partitioning analysis (VPA) across different taxonomic levels. g1, Tissue, Month, Life stage and Sex factors; g2, Tissue, Month, Life stage and SVL factors; g3, Tissue, Month, Life stage and Temperature factors.

Comparisons of co‐occurrence networks across eight gut microbiome groups: (a) overlapping nodes; (b) overlapping edges; (c) taxonomic composition of positive edges. M_AUG, male gut microbiomes in August; M_JAN, male gut microbiomes in January; M_JUL, male gut microbiomes in July; M_LG, male large gut microbiomes; M_NOV, male gut microbiomes in November; M_SEP, male gut microbiomes in September; M_SG, male small gut microbiomes; WG, whole gut microbiomes in June.

Tables

TABLE 1. Sample sizes of 17 gut microbiome groups (n = 110), named using an integrated format of life stage, sex, tissue type, and sampling month.

Group	Sample size
TH‐U‐WG‐JUN	5
TF‐U‐WG‐JUN	6
FL‐U‐WG‐JUN	5
FA‐M‐SG‐JUL	9
FA‐M‐LG‐JUL	8
FA‐F‐SG‐JUL	5
FA‐F‐LG‐JUL	5
FA‐M‐SG‐AUG	8
FA‐M‐LG‐AUG	9
FA‐M‐SG‐SEP	10
FA‐M‐LG‐SEP	10
FA‐F‐SG‐SEP	5
FA‐F‐LG‐SEP	5
PH‐M‐SG‐NOV	6
PH‐M‐LG‐NOV	6
MH‐M‐SG‐JAN	4
MH‐M‐LG‐JAN	4

TABLE 2. Proportions of ASV categories in SG and LG microbiomes across months or life stages.

ASV number	Group	N (%)	R (%)	M (%)	A (%)	p (Χ 2)
5274 (SG microbiomes)	TH‐U‐WG‐JUN	83.52	15.28	0.93	0.27	< 0.001
	TF‐U‐WG‐JUN	89.29	10.60	0.06	0.06
	FL‐U‐WG‐JUN	89.29	10.43	0.23	0.06
	FA‐M‐SG‐JUL	85.93	13.84	0.17	0.06
	FA‐M‐SG‐AUG	90.77	8.84	0.32	0.08
	FA‐M‐SG‐SEP	80.32	19.45	0.19	0.04
	PH‐M‐SG‐NOV	48.84	51.02	0.09	0.04
	MH‐M‐SG‐JAN	60.50	39.40	0.06	0.04
6294 (LG microbiomes)	TH‐U‐WG‐JUN	86.19	12.81	0.78	0.22	< 0.001
	TF‐U‐WG‐JUN	91.02	8.88	0.05	0.05
	FL‐U‐WG‐JUN	91.02	8.74	0.19	0.05
	FA‐M‐LG‐JUL	84.40	14.19	1.10	0.32
	FA‐M‐LG‐AUG	87.75	11.17	0.84	0.24
	FA‐M‐LG‐SEP	77.30	21.58	0.84	0.29
	PH‐M‐LG‐NOV	58.61	40.48	0.75	0.16
	MH‐M‐LG‐JAN	43.96	55.54	0.40	0.10

TABLE 3. Proportions of ASV fluctuation types in SG and LG microbiomes across months or life stages.

ASV number	Group	D (%)	I (%)	U (%)	p (Χ 2)
5274 (SG microbiomes)	TH‐U‐WG‐JUN|TF‐U‐WG‐JUN	4.04	0.66	95.30	< 0.001
	TF‐U‐WG‐JUN|FL‐U‐WG‐JUN	0.28	1.12	98.60
	FL‐U‐WG‐JUN|FA‐M‐SG‐JUL	1.02	0.34	98.63
	FA‐M‐SG‐JUL|FA‐M‐SG‐AUG	0.97	0.25	98.79
	FA‐M‐SG‐AUG|FA‐M‐SG‐SEP	0.21	0.30	99.49
	FA‐M‐SG‐SEP|PH‐M‐SG‐NOV	0.28	2.07	97.65
	PH‐M‐SG‐NOV|MH‐M‐SG‐JAN	0.00	0.78	99.22
6294 (LG microbiomes)	TH‐U‐WG‐JUN|TF‐U‐WG‐JUN	3.38	0.56	96.06	< 0.001
	TF‐U‐WG‐JUN|FL‐U‐WG‐JUN	0.24	0.94	98.82
	FL‐U‐WG‐JUN|FA‐M‐LG‐JUL	1.46	1.67	96.87
	FA‐M‐LG‐JUL|FA‐M‐LG‐AUG	1.26	0.76	97.98
	FA‐M‐LG‐AUG|FA‐M‐LG‐SEP	0.76	2.57	96.66
	FA‐M‐LG‐SEP|PH‐M‐LG‐NOV	1.59	4.53	93.88
	PH‐M‐LG‐NOV|MH‐M‐LG‐JAN	0.29	1.13	98.59

TABLE 4. Proportions of gut microbiomes contributed from three sources: previous month or life stage, contemporary different tissue, and unknown.

Sinks	Sources	Proportions (Mean ± SD)
TF‐U‐WG‐JUN	TH‐U‐WG‐JUN	0.010 ± 0.002
	Unknown	0.990 ± 0.002
FL‐U‐WG‐JUN	TH‐U‐WG‐JUN	0.032 ± 0.006
	TF‐U‐WG‐JUN	0.890 ± 0.058
	Unknown	0.078 ± 0.058
FA‐M‐LG‐JUL	TH‐U‐WG‐JUN	0.015 ± 0.024
	TF‐U‐WG‐JUN	0.021 ± 0.054
	FL‐U‐WG‐JUN	0.042 ± 0.047
	FA‐M‐SG‐JUL	0.491 ± 0.162
	Unknown	0.431 ± 0.179
FA‐M‐SG‐JUL	TH‐U‐WG‐JUN	0.002 ± 0.001
	TF‐U‐WG‐JUN	0.610 ± 0.120
	FL‐U‐WG‐JUN	0.290 ± 0.070
	FA‐M‐LG‐JUL	0.062 ± 0.129
	Unknown	0.036 ± 0.074
FA‐M‐LG‐AUG	FA‐M‐SG‐JUL	0.002 ± 0.002
	FA‐M‐LG‐JUL	0.554 ± 0.160
	FA‐M‐SG‐AUG	0.006 ± 0.003
	Unknown	0.438 ± 0.161
FA‐M‐SG‐AUG	FA‐M‐SG‐JUL	0.673 ± 0.323
	FA‐M‐LG‐JUL	0.001 ± 0.002
	FA‐M‐LG‐AUG	0.012 ± 0.009
	Unknown	0.314 ± 0.324
FA‐M‐LG‐SEP	FA‐M‐SG‐AUG	0.001 ± 0.001
	FA‐M‐LG‐AUG	0.801 ± 0.208
	FA‐M‐SG‐SEP	0.054 ± 0.087
	Unknown	0.144 ± 0.149
FA‐M‐SG‐SEP	FA‐M‐SG‐AUG	0.960 ± 0.029
	FA‐M‐LG‐AUG	0.002 ± 0.002
	FA‐M‐LG‐SEP	0.017 ± 0.015
	Unknown	0.022 ± 0.014
FA‐M‐LG‐NOV	FA‐M‐SG‐SEP	0.094 ± 0.055
	FA‐M‐LG‐SEP	0.430 ± 0.110
	FA‐M‐SG‐NOV	0.262 ± 0.126
	Unknown	0.214 ± 0.176
FA‐M‐SG‐NOV	FA‐M‐SG‐SEP	0.853 ± 0.173
	FA‐M‐LG‐SEP	0.020 ± 0.048
	FA‐M‐LG‐NOV	0.012 ± 0.008
	Unknown	0.115 ± 0.177
FA‐M‐LG‐JAN	FA‐M‐SG‐NOV	0.273 ± 0.074
	FA‐M‐LG‐NOV	0.284 ± 0.170
	FA‐M‐SG‐JAN	0.221 ± 0.116
	Unknown	0.221 ± 0.209
FA‐M‐SG‐JAN	FA‐M‐SG‐NOV	0.774 ± 0.143
	FA‐M‐LG‐NOV	0.003 ± 0.001
	FA‐M‐LG‐JAN	0.089 ± 0.016
	Unknown	0.134 ± 0.131

Full text

Introduction

1

Gut microbiota (i.e., microbial community) and amphibian hosts function as a holobiont through a stable yet dynamic symbiosis (Zilber‐Rosenberg and Rosenberg 2008). These microorganisms play essential roles in amphibian development, ecological adaptation, and health maintenance. For example, gut microbiota can promote metamorphosis and neurodevelopment by mediating the gut–liver and gut–brain axes (Emerson and Woodley 2024; Zhu et al. 2024). The maturation and homeostasis of the immune system are also likely supported by microbial metabolites of gut microbes, such as short‐chain fatty acids and secondary bile acids (Brown and Rudensky 2023; Tran et al. 2025). Conversely, dysbiosis of the gut microbiome (i.e., microorganisms and their theater of activity) may increase disease susceptibility and weaken resistance to parasites, particularly during early developmental stages (Berg et al. 2020; Knutie et al. 2017; Warne et al. 2019).

Amphibians undergo multifaceted self‐regulation across seasons, including changes in body morphology, dietary behavior, and physiological function during key life stages (Wang et al. 2025). For instance, larval organs such as intestines restructure extensively during metamorphosis to facilitate the transition from aquatic to territorial habitats. Adult amphibians enter hibernation to survive harsh environmental conditions, including prolonged food scarcity and cold exposure. Concurrently, gut microbiomes also remodel dramatically (Cao et al. 2023; Chai et al. 2018; Huang and Liao 2021; Park and Do 2024; Scalvenzi et al. 2021; Shi et al. 2023; Tong, Cui, et al. 2020; Tong et al. 2023; Wiebler et al. 2018; Zhang et al. 2019, 2020), driven by synergistic effects of environmental factors (e.g., habitat conditions) and host‐related factors (e.g., physiological status) (Levin et al. 2021; Martinez‐Guryn et al. 2019). However, gut microbiome succession across the amphibian lifecycle remains insufficiently studied (Scalvenzi et al. 2021; Shi et al. 2023; Tong, Cui, et al. 2020). In particular, temporal fluctuations of microbial communities in different intestinal regions have been largely overlooked (Scalvenzi et al. 2021; Zhang et al. 2019), limiting a comprehensive understanding of gut microbiomes across both temporal and spatial scales.

Beyond deterministic forces such as habitat filtering and microbial interactions, gut microbiome assembly can also be shaped by stochastic ecological processes, including dispersal and ecological drift (Kern et al. 2021; Zhou and Ning 2017). Our pilot studies have highlighted the importance of stochastic processes in shaping gut microbiomes of hibernating amphibians (Song et al. 2021, 2023). Nevertheless, most present studies of amphibian gut microbiomes have paid limited attention to stochasticity. This knowledge gap hinders a mechanistic understanding of gut microbiome assembly and its role in amphibian physiological adaptation.

In this study, we aimed to characterize the spatiotemporal dynamics of amphibian gut microbiomes across seasons and to elucidate their underlying assembly mechanisms, with particular emphasis on stochastic processes. Cultivated Black‐spotted frog ( Pelophylax nigromaculatus ) was selected as the model organism because its living conditions and genetic backgrounds can be effectively controlled (Figure 1). Frog individuals were sampled successively over 6 months, and gut microbiomes were profiled using 16S rRNA gene amplicon sequencing. Microbial community structure was assessed using multiple diversity and compositional indices. Furthermore, we investigated gut microbiome assembly by systematically analyzing dominant driving factors, ecological processes, phylogenetic traits, source‐sink relationships, and co‐occurrence networks.

(a) Photographs of Black‐spotted frog ( Pelophylax nigromaculatus ) at three representative life stages. (b) Experimental field for rice‐frog co‐cultivation. (c) Research workflow in this study.

Materials and Methods

2

Animal Cultivation and Field Experiment

2.1

A rectangular field for rice‐frog (RF) co‐cultivation was established in January 2019 (Figure 1). The experimental field consisted of three main components: a ditch for water inflow and outflow, a rice paddy, and a feeding‐table area. The aquaculture water was pumped into the field from a nearby river. Cultivation procedures were manually adjusted to meet the requirements of frogs at different life stages, including hind‐limbed tadpoles (TH), fore‐limbed tadpoles (TF), froglets (FL), active frogs (FA), prehibernation frogs (PH), and mid‐hibernation frogs (MH). All animal cultivation followed local aquaculture regulations, and the experimental procedures were approved by the Animal Care and Use Committee in Xinyang Normal University (XFEC‐2018‐09).

Frog and Tissue Sampling

2.2

A total of 246 frog individuals were collected from the RF field over 6 months: June 1, 2019; July 7, 2019; August 24, 2019; September 27, 2019; November 4, 2019; and January 12, 2020. Randomly sampled frogs were euthanized by immersion in 95% ethanol (tadpoles and froglets) or by double pithing of the brain and spinal cord. Snout–vent length (SVL) and body mass were measured, and sex was determined by visual inspection of the gonads. The intestinal tissues were dissected under aseptic conditions and immediately stored at −20°C. The small gut (SG) and large gut (LG) were separately collected for each individual, except for tadpoles or froglets sampled in June, for which the whole gut (WG) was collected.

DNA Extraction and 16S rRNA Gene Amplicon Sequencing

2.3

To profile gut microbiomes, a subset of 75 frog individuals was selected for DNA extraction. Intestine tissues were rinsed with sterile water and opened using sterilized scissors under aseptic conditions. Total microbial genomic DNA was extracted using a modified phenol‐chloroform method (Song et al. 2018). After assessing DNA purity and concentration, 16S rRNA gene amplicon libraries targeting the V3–V4 hypervariable regions (primers 341F/805R) were constructed and sequenced using paired‐end sequencing on an Illumina MiSeq platform (Biobit Biotech Inc., China) (Song et al. 2021). In total, 110 gut microbiomes were obtained and categorized into 17 sample groups based on life stage, sex, tissue type, and sampling month (Tables 1 and S1). The sequencing data have been deposited in the National Microbiology Data Center under project number NMDC10020041.

ASV Identification and Taxonomic Classification

2.4

Amplicon sequence variants (ASVs) were identified from the 110 gut microbiomes in this study, together with an additional 131 microbiomes that will be reported elsewhere. ASV inference was conducted using a QIIME2‐based pipeline (ver 2019.10) (Bolyen et al. 2019). Quality control included removal of sequences containing ambiguous bases, filtering by PHRED score (≥ 4), trimming to 200 bp, and denoising using the deblur algorithm (Amir et al. 2017). These ASVs were phylogenetically analyzed using MAFFT and FastTree (Katoh and Standley 2013; Price et al. 2010). Taxonomic annotation was performed using the classifier of the SILVA 16S rRNA database (ver 132) with a confidence threshold of 0.8 (Quast et al. 2012). The ASV table was filtered to retain ASVs with a total frequency ≥ 10, occurrence in at least two samples, exclusion of mitochondrial and chloroplast sequences, and inclusion of only bacterial and archaeal taxa. Two samples (IIIAS30 and IIIAS51) with insufficient sequencing depth were excluded from further analyses.

Structural Analyses of Gut Microbiomes

2.5

Structural analyses were performed in R (ver 4.4.0) using both unrarefied and rarefied datasets (R Core Team 2024). Major R packages included microeco, vegan, ape, magrittr, ggpubr, tidyr, and data.table (Bache and Wickham 2022; Barrett et al. 2024; Kassambara 2023; Liu et al. 2021; Oksanen et al. 2022; Paradis and Schliep 2019). For unrarefied data, ASV relative abundances and pairwise shared‐ASV ratios were calculated. Rarefaction was performed to a sequencing depth of 8868 (minimum depth among retained samples), after which alpha diversity (Observed ASVs, Shannon, Simpson, and Faith's phylogenetic diversity), beta diversity (Bray–Curtis, Jaccard, weighted and unweighted UniFrac), and taxonomic relative abundances (phylum, family, genus levels) were computed.

To examine temporal succession patterns, 14 group pairs were compared while controlling for tissue and sex effects. Tissue‐ and sex‐related differences were assessed using 11 pairs of contemporaneous samples. ASVs were classified into four abundance categories: None (0), Rare (0–0.001), Moderate (0.001–0.01), and Abundant (> 0.01). Temporal changes of ASVs were categorized as Decreased, Increased, or Unchanged based on Wilcoxon rank‐sum tests between temporally adjacent samples. Chi‐squared tests were used to evaluate differences in ASV categories and fluctuation types. The pairwise shared‐ASV ratios and alpha diversity were visualized using boxplots and tested with Wilcoxon rank‐sum tests. Beta diversity was visualized using boxplots and principal coordinates analysis (PCoA), and assessed using Wilcoxon tests and PERMANOVA. Taxonomic composition was visualized with alluvial plots, and differential taxa were identified using LEfSe and random forest analyses (Liu et al. 2021; Segata et al. 2011). The filtering thresholds of taxonomic relative abundance, logarithmic LDA score, and mean Gini decrease were set as 0.001, 4, and 0.5. Statistical significance was set at p < 0.05 unless otherwise stated, and p values were adjusted using the Benjamini–Hochberg FDR method.

Temporal decay was evaluated by calculating Spearman correlations between structural dissimilarities and temporal intervals (Shade et al. 2013). The analyses were conducted separately for SG and LG microbiomes, excluding female samples. Structural divergence between SG and LG microbiomes was also assessed across months using male samples only.

Assembly Mechanism Analyses of Gut Microbiomes

2.6

Variance partitioning analysis (VPA) was applied to evaluate the explanatory power of factors on structural variation of gut microbiomes (Legendre 2008; Oksanen et al. 2022). Separate VPAs were performed for three sets of factors, each including three dominant factors (tissue, month, life stage) and one subdominant factor (sex, SVL, or temperature). Factor significance was assessed using the “anova.cca” function of the vegan package.

Ecological processes governing gut microbiome assembly were inferred using phylogenetic bin‐based null model analysis methods, that is, “icamp.cm”, “icamp.bins” and “icamp.cate” functions implemented in the iCAMP package (ver 1.5.12) (Ning et al. 2020). ASVs were grouped into phylogenetic bins with “bin.size.limit = 48”. Ecological processes were classified as homogeneous selection, heterogeneous selection, homogenizing dispersal, dispersal limitation, and drift based on βMNTD_RC and RCbray values.

Phylogenetic structure of gut microbiomes was evaluated using confidence indices (cMPD and cMNTD) (Ning et al. 2020). Phylogenetic conservation in the niche adaptation of gut microbiomes was assessed by examining whether explanatory powers of factors had upward or stable trends across taxonomic levels (Lu et al. 2016). The structural variation data were allocated to three sets of factors in a similar manner as above VPAs.

Source–sink relationships were inferred using FEAST (version 0.1.0) to quantify contributions from temporal transmission, spatial dispersal, and unknown sources (Shenhav et al. 2019).

Co‐occurrence networks were constructed using the microeco package (ver 1.6.0) for eight gut microbiome groups (Liu et al. 2021). Spearman correlations were calculated among ASVs, and network modules were partitioned using the Louvain method. The network topology and robustness were compared using the meconetcomp package (ver 0.6.1) (Liu et al. 2023).

Results

3

The Dynamics of Gut Microbiomes Across Months or Life Stages

3.1

Most ASVs were rare or absent in gut microbiome groups and remained stable in relative abundance during succession (Tables 2 and 3). However, almost all temporally adjacent groups showed different proportions of ASV categories and ASV fluctuation types (Tables S2 and S3). Gut microbiomes exhibited complex succession trends in terms of pairwise shared‐ASV ratios, alpha diversity, and beta diversity indices (Figure 2). The abundance‐unweighted pairwise shared‐ASV ratios and beta diversity indices gradually decreased and increased across months, respectively, whereas abundance‐weighted indices remained comparatively stable. Furthermore, the dynamics of gut microbiomes was characterized by significant tissue specificity (Figures 2 and S1, Table 4). For instance, the abundance‐weighted indices of LG microbiomes presented greater variation than those of SG microbiomes. Significant differences were detected in abundance‐weighted indices of temporally adjacent LG microbiomes, but not in SG microbiomes. LG microbiomes displayed greater temporal fluctuation of taxonomic composition than SG microbiomes (Figure S2). Proteobacteria and Firmicutes were the most dominant taxa in SG and LG microbiomes from July through November, respectively. Several genera markers were consistently identified by LEfSe and random forest methods across four pairs of temporally adjacent groups: (1) TH‐U‐WG‐JUN (Mycobacterium) and TF‐U‐WG‐JUN (Pseudomonas); (2) FA‐M‐LG‐JUL (Pseudomonas) and FA‐M‐LG‐AUG (Cetobacterium); (3) FA‐M‐LG‐AUG (Cetobacterium, Escherichia−Shigella) and FA‐M‐LG‐SEP (Akkermansia); (4) FA‐M‐LG‐SEP (Bacteroides, Erysipelatoclostridium) and PH‐M‐LG‐NOV (Pseudomonas). Temporal‐decay phenomena occurred in the abundance‐unweighted between‐group indices of SG and LG microbiomes (Figure 3), with significant negative Spearman correlations observed between unweighted indices and the interval days.

Structural comparisons across months or life stages using pairwise shared‐ASV ratios, alpha diversity, and beta diversity indices. The shared_abund_rat and shared_rich_rat represent abundance‐weighted and abundance‐unweighted pairwise shared‐ASV ratios, respectively. The significance level for Wilcoxon rank sum tests was set to 0.05.

Spearman correlations between inter‐group indices and their interval days. SG and LG microbiomes were analyzed separately, excluding female samples. Linear regression lines with 95% confidence intervals are shown in gray.

The Structural Divergence Between Small and Large Gut Microbiomes

3.2

Contemporary SG and LG microbiomes, except those in January, showed significant structural differences in ASV categories, ASV fluctuation types, pairwise shared‐ASV ratios, alpha diversity, and beta diversity indices (Figure 4, Tables S5 and S6). By contrast, differences between sexes were weak or absent. Specifically, alpha diversity of SG microbiomes was lower than that of LG microbiomes. The structural differences between tissues were opposite for the abundance‐unweighted and abundance‐weighted intercommunity indices. SG and LG microbiomes initially diverged and later converged during succession (Figures S1 and S3). Taxonomic composition differed significantly between tissues, but not between sexes. Tissue‐associated genera markers consistently identified by LEfSe and random forest analyses included: (1) FA‐F‐SG‐JUL (Pseudomonas) and FA‐F‐LG‐JUL (endosymbiont ‘TC1’ of Trimyema compressum); (2) FA‐M‐SG‐JUL (Pseudomonas) and FA‐M‐LG‐JUL (endosymbiont ‘TC1’ of Trimyema compressum, Erysipelatoclostridium and Ruminococcaceae_UCG−014); (3) FA‐M‐SG‐AUG (Pseudomonas) and FA‐M‐LG‐AUG (Cetobacterium); (4) FA‐M‐SG‐SEP (Pseudomonas) and FA‐M‐LG‐SEP (Breznakia, Akkermansia and Erysipelatoclostridium); (5) PH‐M‐SG‐NOV (Pseudomonas) and PH‐M‐LG‐NOV (endosymbiont ‘TC1’ of Trimyema compressum).

Structural comparisons between tissues and between sexes using pairwise shared‐ASV ratios, alpha diversity, and beta diversity indices. The significance level for Wilcoxon rank sum tests was set to 0.05.

The Dominant Factors for the Structural Variation of Gut Microbiomes

3.3

The main factors explaining structural variation of gut microbiomes were tissue, month, and life stage, rather than sex, SVL, or temperature (Figure 5). The host and environmental factors had greater explanatory power for abundance‐weighted structural variation than for abundance‐unweighted variation. For example, tissue, month, life stage, and sex explained 40.2% of Bray‐Curtis and 43.8% of weighted UniFrac distances, compared to only 9.7% of Jaccard and 17.2% of unweighted UniFrac distances. Nevertheless, most structural variation remained unexplained by these factors.

*Variance partitioning analysis (VPA) on beta diversity indices of gut microbiomes for three groups of factors. Each group included three dominant factors (Tissue, Month, Life stage) and one subdominant factor (Sex, SVL, or Temperature). Negative values are not shown. The conditional significance of the deterministic factor < 0.05.

The Ecological Processes to the Assembly of Gut Microbiomes

3.4

Gut microbiome assembly was dominated by stochastic ecological processes, with the relative importance ranging from 80.43% to 99.71% (Figure 6a). Ecological processes varied spatiotemporally and showed clear heterogeneity between ASV categories and bins (Figure 6b,c). For active frogs, drift accounted for almost all ecological processes in SG microbiomes, whereas dispersal limitation and homogeneous selection contributed substantially to LG microbiomes. During hibernation, dispersal limitation and homogeneous selection decreased dramatically in LG microbiomes, while homogenizing dispersal increased, especially in MH. Ecological processes between SG and LG microbiomes fluctuated significantly from July through January (Figure 6d). Homogeneous selection accounted for the smallest proportion in August (20.99%) and the largest in January (89.07%). Dispersal limitation was detected from July through November but was absent in January.

Relative importance of ecological processes deduced with bin‐based null model of R package iCAMP (ver 1.5.12) for: (a) gut microbiome groups; (b) ASV categories; (c) bins; (d) between male SG and LG microbiomes. A, abundant; DL, dispersal limitation; DR, drift; HD, homogenizing dispersal; HeS, heterogeneous selection; HoS, homogeneous selection; M, moderate; R, rare.

The Phylogenetic Traits of Gut Microbiomes

3.5

Gut microbiomes exhibited variable phylogenetic convergence across space and time (Figure 7). All taxa showed high phylogenetic divergence during metamorphosis, whereas nearest taxa in TF and FL displayed greater convergence. SG and LG microbiomes generally exhibited divergent and convergent phylogenetic traits, respectively. Nearest taxa of SG microbiomes in July retained convergent traits of TF and FL, which differed from subsequent months. During hibernation, LG microbiomes gradually shifted from high convergence to near randomization. In addition, phylogenetic conservation signals were detected for gut microbiomes in the niche adaptation, particularly in tissue‐associated niches (Figure 7).

(a) Wilcoxon rank sum tests on the CMPD and CMNTD indices of contemporary and temporally adjacent gut microbiomes. The significance level was set to 0.05. (b) The adjusted R‐Squares of three groups of factors in variance partitioning analysis (VPA) across different taxonomic levels. g1, Tissue, Month, Life stage and Sex factors; g2, Tissue, Month, Life stage and SVL factors; g3, Tissue, Month, Life stage and Temperature factors.

The Source‐Sink Relationships of Gut Microbiomes

3.6

Source‐sink relationships among gut microbiomes varied spatiotemporally (Table 4). During metamorphosis, TH sources contributed minimally (~1%) to TF sinks, whereas TF sources contributed the largest proportion (~89%) to FL sinks. Temporal transmission contributed far more than spatial dispersal to SG microbiomes. For instance, SG microbiomes in July mainly derived from TF (~61%) and FL (~29%), rather than contemporary LG microbiomes (~6%). In LG microbiomes, contributions from spatial dispersal and temporal transmission fluctuated across months. LG microbiomes in July were primarily contributed by contemporary SG (~49%) and unknown sources (~43%). LG microbiomes in August and September derived mainly from counterparts in the previous month (~55% and ~80%), with minimal spatial dispersal (< 1% and ~5%). During hibernation, the contribution of spatial dispersal to LG microbiomes increased (~26% in November and ~22% in January).

The Co‐Occurrence Networks of Gut Microbiomes

3.7

Co‐occurrence networks exhibited significant temporal dynamics and tissue dependence (Figure 8 and Table S7). Overlapping nodes and edges between SG and LG networks were low (10.8% and 0.8%). Across 6 months, the networks shared 2 nodes (0.4%) and no edges. Many key nodes displayed temporal fluctuations or tissue‐specific differences (Table S8). For example, four key nodes in the WG network corresponded to genus markers of TH‐U‐WG‐JUN and TF‐U‐WG‐JUN, that is, Mycobacterium and Pseudomonas. The July network included six key nodes associated with genus markers of FA‐M‐SG‐JUL and FA‐M‐LG‐JUL, that is, Pseudomonas, endosymbiont ‘TC1’ of Trimyema compressum, Erysipelatoclostridium, and Ruminococcaceae_UCG−014. Positive network edges showed divergent taxonomic composition (Figure 8). For example, the June network had more edges associated with Actinobacteria (~25%) and Planctomycetes (~45%) than Firmicutes (~15%), whereas subsequent months had more edges associated with Firmicutes. LG network possessed more positive edges associated with Bacteroidetes than SG network. The network during metamorphosis showed greater robustness than those in later months (Figure S4), and the LG network was more robust than the SG network.

Comparisons of co‐occurrence networks across eight gut microbiome groups: (a) overlapping nodes; (b) overlapping edges; (c) taxonomic composition of positive edges. M_AUG, male gut microbiomes in August; M_JAN, male gut microbiomes in January; M_JUL, male gut microbiomes in July; M_LG, male large gut microbiomes; M_NOV, male gut microbiomes in November; M_SEP, male gut microbiomes in September; M_SG, male small gut microbiomes; WG, whole gut microbiomes in June.

Discussion

4

Amphibian gut microbiomes undergo adaptive restructuring in response to dramatic changes of environmental and host factors across the lifecycle (Cao et al. 2023; Chai et al. 2018; Huang and Liao 2021; Park and Do 2024; Scalvenzi et al. 2021; Shi et al. 2023; Tong, Cui, et al. 2020; Tong et al. 2023; Wang et al. 2025; Wiebler et al. 2018; Zhang et al. 2019, 2020). Seasonal dynamics of gut microbiomes represent a widespread ecological strategy among animals, including insects (Fluch et al. 2024; Hendrycks et al. 2025), fish (Xu et al. 2024), and mammals (Fan et al. 2022). Here we also depict the complex gut microbiome succession of cultivated Black‐spotted frog across seasons. However, the underlying assembly mechanisms vary depending on host phylogeny, physiological traits, and ecological niche. For example, deterministic and stochastic processes play a considerable role in shaping the gut microbiomes of cucurbit‐feeding fruit flies (Hendrycks et al. 2025), while deterministic processes dominate in gut microbiota assembly of cold‐water fish (Xu et al. 2024). By contrast, stochastic processes play dominant roles in shaping frog gut microbiomes, although deterministic processes, particularly for LG microbiomes, remain significant. Structural divergence between SG and LG microbiomes appears to be driven largely by fluctuations of homogeneous selection and dispersal limitation.

During metamorphosis, amphibian larva experience body reconstruction for the adaptation from aquatic to territorial habitats (Hourdry et al. 1996). Once tadpoles enter the forelimb phase, they cease feeding, and their tails are absorbed. The intestines transform from a fish‐like form (long and undivided) into terrestrial form (shortened with enhanced differentiation) (Ishizuya‐Oka 2007). Concurrently, gut microbiomes shift from fish‐like to amniote‐like composition (Kohl et al. 2013). In line with previous studies (Chai et al. 2018; Shi et al. 2023; Tong, Cui, et al. 2020), species richness, evenness, and phylogenetic diversity decrease from TH to TF, and stabilize in FL. Compared to TH microbiomes, TF and FL microbiomes are more randomized, with core ASVs more centralized. Proteobacteria (especially Pseudomonas) accounts for a greater proportion in TF and FL, while other phyla such as Actinobacteria and Planctomycetes decrease. Similarly, Actinobacteria is enriched in the gut microbiomes of other anurans during early metamorphosis (Chai et al. 2018; Scalvenzi et al. 2021; Shi et al. 2023; Tong, Cui, et al. 2020). However, the dynamics of Proteobacteria vary across amphibian species. Planctomycetes, though widespread in aquatic environments (Kaboré et al. 2020), are rarely detected in amphibian gut microbiomes (Shi et al. 2023; Zhang et al. 2020). The high phylogenetic divergence in TH, TF, and FL microbiomes coincides with the predominance of stochastic processes. The remodeling of gut microbiomes from TH to TF is unlikely to be driven by different species pools, as the habitat remains aquatic. The intestinal microenvironment in the forelimb phase may selectively filter microorganisms that are unsuitable for future terrestrial feeding. Temporal transmission contributes most during the final phase of metamorphosis, suggesting ecological niches remain relatively stable as fore‐limbed tadpoles develop into froglets. Comparisons of co‐occurrence networks indicate that microbial interactions during metamorphosis differ markedly from other periods. Actinobacteria and Planctomycetes probably have critical functions in network maintenance during the hindlimb stage, whereas Firmicutes become functionally dominant in later stages.

After metamorphosis froglets actively feed on land and grow rapidly. The gastrointestinal tract differentiates structurally and functionally: the slender, coiled SG is connected via mesentery, whereas the stubby, straight LG connects to the cloaca. Differences in internal environments (e.g., microvilli length and oxygen content) create distinct ecological niches for microbial colonization (Martinez‐Guryn et al. 2019). Higher oxygen in SG favors Proteobacteria colonization, while oxygen consumption by Proteobacteria may create hypoxic conditions conducive to Firmicutes and Bacteroidetes in LG. Although stochastic processes (drift) predominantly shape both SG and LG microbiomes, deterministic processes (homogeneous selection, habitat filtering, and microbial interactions) are more influential in LG. Temporal transmission contributes far more than spatial dispersal to SG microbiomes, likely due to unidirectional intestinal flow, while spatial dispersal contributes more to LG microbiomes in July, reflecting incompletely differentiated intestines at the subadult stage (Gao et al. 2018). Notably, gut microbiomes in August exhibit maximal structural differences between tissues, weakest homogeneous selection, and least robust co‐occurrence network, possibly due to decreased food intake under high temperatures and attenuated intestine function (Fontaine et al. 2018; Lian et al. 2020).

Adult amphibians typically enter hibernation to survive prolonged fasting and cold exposure, accompanied by dramatic gastrointestinal restructuring to optimize energy allocation (Naya et al. 2009). Gut microbiomes also undergo significant remodeling (Cao et al. 2023; Shi et al. 2023; Tong, Cui, et al. 2020; Tong et al. 2023; Tong, Hu, et al. 2020; Weng et al. 2016). The remodeling patterns vary highly among species. For instance, alpha diversity decreases in hibernating Rana amurensis , Rana dybowskii , Polypedates megacephalus , and Strauchbufo raddei , whereas Odorrana Tormota exhibits decreased species richness but initially increased evenness. In this study, species richness and phylogenetic diversity increase during pro‐ and mid‐hibernation, regardless of tissue. SG microbiomes maintain low evenness, whereas the evenness of LG microbiomes declines gradually. Proteobacteria consistently dominates SG microbiomes, while Proteobacteria and Planctomycetes increase and Firmicutes and Fusobacteria decline in LG microbiomes. Gut symbionts act as an essential regulator in material and energy metabolism during hibernation (Kurtz et al. 2021), with long‐term fasting favoring microbes capable of metabolizing endogenous substances such as mucosal proteins, apoptotic cells, and urea. For example, Pseudomonas in LG microbiomes may enhance urea nitrogen recycling (Regan et al. 2022; Wiebler et al. 2018). Enhanced stochasticity during hibernation also contributes to greater species richness, lower evenness, and increased rare taxa. Increased similarity between SG and LG microbiomes during hibernation likely reflects stronger homogeneous selection and reduced dispersal limitation, consistent with increased SG‐to‐LG dispersal in January.

Conclusions

5

This study demonstrates that amphibian gut microbiomes exhibit tissue‐specific structural succession, with temporal fluctuations in structural discrepancy between SG and LG microbiomes. Stochastic processes dominate the gut microbiome assembly, although deterministic processes—including habitat filtering and microbial interactions—remain significant, particularly for LG microbiomes. Comparisons between sexes, performed only in July and September, revealed minimal sex effects; however, gonadal differentiation likely occurs in the hindlimb stage of metamorphosis (Li et al. 2018), so the finding should be interpreted cautiously. To address small sample size, co‐occurrence networks were constructed by merging groups across tissues and months, potentially limiting detection of tissue‐ or month‐specific microbial interactions. Additionally, gut microbiomes in later and after hibernation were not profiled, though they likely restore prehibernation structure (Park and Do 2024). Further research is warranted to characterize the spatiotemporal dynamics of amphibian gut microbiomes during and after hibernation.

Author Contributions

Xiaowei Song: conceptualization (lead), data curation (lead), formal analysis (lead), funding acquisition (lead), investigation (lead), methodology (lead), project administration (lead), writing – original draft (lead), writing – review and editing (lead). Yuanyuan Zhai: investigation (equal). Mengyang Zhang: investigation (equal). Jingyuan Guo: investigation (equal). Benjun Guo: data curation (equal), software (equal). Chaolong Zhang: data curation (equal), software (equal). Jin Jin: data curation (equal), software (equal). Weiye Wang: data curation (equal), software (equal). Yuanping Xu: data curation (equal), software (equal). Bicheng Zhu: writing – original draft (equal), writing – review and editing (equal). Xiangzhen Li: conceptualization (equal), writing – original draft (equal), writing – review and editing (equal).

Funding

The study was supported by the National Natural Science Foundation of China (NSFC 31600104), Key R&D and Promotion Projects of Henan Province (192102110005), and Research Startup Foundation in Chengdu University of Information Technology (KYTZ202279).

Ethics Statement

The animal cultivation complied with the local aquaculture regulations, and the procedures used in the research were approved by the Animal Care and Use Committee in Xinyang Normal University (XFEC‐2018‐09).

Conflicts of Interest

The authors declare no conflicts of interest.

Supporting information

Data S1: ece373041‐sup‐0001‐Supinfo1.docx.

Table S4: ece373041‐sup‐0002‐TableS4.tsv.

Table S8: ece373041‐sup‐0003‐TableS8.xlsx.