Analysis of postmortem intestinal microbiota successional patterns with application in postmortem interval estimation

Microorganism plays a vital role in the decomposition of vertebrate remains in nature nutrient cycling, where the postmortem microbial succession patterns during decomposition remain unclear. Here, hierarchical clustering based on Manhattan distances analyzed the similarity and difference among postmortem intestine microbial succession patterns from microbial 16S rDNA sequences in a mouse decomposition model. Based on the similarity, seven different classes of succession patterns are obtained. Generally, the normal intestinal ora in the cecum gradually decreases with changes in living conditions after death, while some facultative anaerobes and obligate anaerobes grew and multiplied with oxygen consumption. Furthermore, a random forest regression model is developed to predict postmortem interval with the microbial succession trend data set. The model demonstrates a mean absolute error of 20.01 hours and a squared correlation coecient of 0.95 during 15-day decomposition. Lactobacillus, Dubosiella, Enterococcus, and Lachnospiraceae NK4A136 group are considered as signicant biomarkers for this model according to the ranked list. This study explores microbial succession patterns in terms of relative abundances and variety, aiding in predicting postmortem intervals, offering some knowledge of microbial behaviors in decomposition ecology.


Introduction
As a key aspect of the material cycle and energy ow of the ecosystem, the decomposition of vertebrate remains is one of the main processes involved in maintaining ecosystem functions. This process begins immediately after death. A lack of oxygen results in cell autolysis, and tissue breakdown, releasing cellular components and macromolecules. These nutrients are subsequently metabolized by predominately microbial communities, particularly those concentrated in the gastrointestinal tract. These microbes begin to exponentially proliferate, transmigrate, and create specialized proteins that digest host tissues during putrefaction [1]. That the intestinal microbiome plays an essential role in the remaining decomposition [2,3]. However, the number of studies of the dynamic succession patterns in these microbial communities is relatively small, but these studies may facilitate estimation of the postmortem interval (PMI) in forensic investigations and improve understanding of the mechanisms governing succession processes during decomposition to accelerate decomposition, facilitate nutrition recycling, and reduce environmental pollution.
Recently, some studies reported microbial succession trends during decomposition. For example, Liu R et al. [4] established a microbial clock providing an accurate estimate of PMI in a mouse model system and detected Enterococcus faecalis, Anaerosalibacter bizertensis, Lactobacillus reuteri, and others as the most informative species in the decomposition process for PMI estimation. Zhang J. et al [5]. also established a random forest model to estimate PMI with burial microbial sequence data set. Additionally, Debruyn et al. [6] documented that Bacteroidales (Bacteroides, Parabacteroides) signi cantly declined while Clostridiales (Clostridium, Anaerosphaera) and the y-associated Gammaproteobacteria Ignatzschineria and Wohlfahrtiimonas increased in human gut bacterial communities. Adserias et al. [7] monitored the oral microbiota of decaying human bodies and found that Firmicutes and Actinobacteria were the predominant phyla in the fresh stage, while the presence of Tenericutes corresponded to the bloat stage, and Firmicutes was the predominant phylum in advanced decay. As decomposition progressed, Firmicutes replaced Proteobacteria as the dominant phylum in the late decay stages of swine models [8]. Additionally, researchers found that the relative abundances of Clostridiaceae and Enterococcaceae increased during decomposition in murine models [2,9]. These similar results regarding microbial succession between animal decomposition models and human remains showed the possibility to use animal decomposition models as alternatives for microbial succession in human corpses to mathematically establish microbial community change models during decomposition.
However, these studies only partially focused on microbial succession during decomposition. Few studies have tried to mathematically describe microbial succession patterns. Fortunately, recent studies have shown that utilizing mathematical models to describe temporal patterns of microbial diversity is a developing strategy that facilitates the inference of relative growth rates, mechanisms of interaction, and responses to perturbations [10][11][12]. A recent increase in the number of microbial time series studies with mathematical methods offers some insights into the stability and dynamics of microbial communities in a wide range of environments [13][14][15][16]. In a meta-analysis of temporal dynamics in microbial communities (including 76 sites representing air, aquatic, soil, brewery wastewater treatment, human-and plantassociated microbial biomes), approximately half had time-decay curves with negative slopes [15].
Besides, the temporal variability in microbial community diversity was comparable across studies within the same environment but varied across environments, being lowest in soil and brewery wastewater and highest in the human palm and the infant gut [15,16]. Some microbial communities go through a series of predictable states after colonization (primary succession), which occurs for instance during the formation of dental plaque, where oxygen-tolerant early colonizers prepare the ground for later oxygensensitive colonizers [17]. In some cases, such as infant gut microbiota colonization [18], communities vary considerably in the initial stages of succession, but stabilize at similar states. These studies indicated the potential of using mathematical methods in microbial dynamic analysis.
Thus, this study aimed to explore intestinal microbial succession patterns during decomposition by a mathematical algorithm and establish a PMI estimation model based on successional trend data set. The statistically signi cant changes in microbial relative abundance during decomposition were analyzed by the Kruskal-Wallis test (P< 0.001). Utilizing hierarchical classi cation analysis to classify the successional patterns based on the Manhattan distance, seven classes of successional patterns were obtained. Moreover, A random forest model combined with the postmortem intestinal microbial successional trend data set yielded a mean absolute error of 20.01 hours within 15-day decomposition.
Some initial colonizers of the mammalian gastrointestinal microbiota and some obligately anaerobic microbes were regarded as important taxa in this model according to the ranked list of all taxa.

Materials And Methods
Sample collection and experimental setup All experiments were approved by the Institutional Animal Use and Care Committee of Xi'an Jiaotong University (approval No: 2017-288). Adult male C57BL/6J mice (18-25 g, 6-10 weeks, n = 240) were obtained from the Experimental Animal Center of Xi'an Jiaotong University. All mice were housed in micro isolator cages under standard lighting (light/dark periods of 12h) and were given autoclaved food and water ad libitum. After a sacri ce by cervical dislocation, the mice were placed on sterile plates with ambient temperature (Ta, 25.0±1.5℃) and moderate relative humidity (RH, 50.0±7%) in an arti cial climate box (PRX-35B, Zuole, Shanghai, China). The arti cial climate box and clean bench were treated with 75% alcohol solution and an hour of UV radiation to minimize and fragment environmental contaminant DNA. The vent of the arti cial climate box was wrapped by sterile gauze at 200 mesh to resist insect access. Every surgical instrument was autoclaved and burned with the outer ame of an alcohol lamp. Ten-time points (0 hours, 8 hours, 12 hours, 1 day, 2 days, 4 days, 7 days, 10 days, 13 days, 15 days after death) were set spanning 15-day decomposition. At each time point, ceca were obtained from 20 mice under strict aseptic operation by multiple people in a short period, to avoid the in uence of sampling on microbiota variability. This 20*10 sample set was used for the following model training and validation. Besides, we set another four mice to every time point as an external sample set (4*10) according to the single-blind principle. Concretely, another experimenter collected postmortem ceca without telling the PMI of the samples. Samples were stored at −80°C immediately after harvest until further processing.
Following bead beating, DNA was isolated from each cecum sample using the QIAamp® DNA Mini Kit (Qiagen, Germany) according to the manufacturer's speci cations on a super clean bench. The variable region 3-4 (V3, V4) of the 16S rRNA gene was then ampli ed from each DNA sample by the bacterial primers 341F (5'-CCTAYGGGRBGCASCAG-3') and 806R (5'-GGACTACNNGGGTATCTAAT-3'). Each amplicon was sequenced on an IonS5 TM XL sequencer. DNA sequences were analyzed using the mothur (version 1.33.3) pipeline. Potentially chimeric sequences were removed using the UCHIME algorithm [19]. Sequences were assigned to OTUs based on a 3% dissimilarity cutoff using the average neighbor algorithm in the UPARSE pipeline (version 7.0.1001) [20]. Each representative sequence was matched against the SSUrRNA database of SILVA reference alignments using mothur for taxonomic annotation [21]. OTU abundance information from each sample was homogenized, with the standard of the fewest sequences in the samples. Rarefaction and rank-abundance curves were drawn using R software (version 2.15.3). Alpha diversity was measured by the Shannon, Chao1, Simpson, and ACE indexes in QIIME software [22][23][24]. Additionally, phylogenetic beta diversity was measured using unweighted UniFrac and weighted UniFrac distances in QIIME software [22,25,26].

Data preprocessing
Every individual sample contained the relative abundance values of 5099 different taxa. First, the data were ltered to remove the noninformative taxa because the relative abundance values of many taxa were mostly zero. Speci cally, samples were grouped according to the time point of death. For the data set, if the relative abundance of a taxon was zero in more than 60% of samples in all groups, then this taxon was discarded. Thus, the number of taxa was reduced from 5099 to 1810. Second, the dataset was normalized by Z-score standardization. After normalization, the relative abundance values for each taxon were scaled to be between 0 and 1 so that the change tendencies could be explored.

Microbial community dynamic analysis
The Kruskal-Wallis test was utilized to analyze the signi cance of the differences in the relative abundances between adjacent time points. We counted a signi cant increase as +1 and a signi cant decrease as -1. The trend values were accumulated in sequence according to the PMI time series. Finally, we obtained a microbial community succession trend dataset. In gure 2B, the red lines indicate that the changing trend of that taxon is signi cant, while the blue lines indicate that the changing trend of that taxon is not signi cant (P < 0.001, Kruskal-Wallis test). The decomposition patterns were explored by an HCA method based on the Manhattan distances between the trend values of different taxa. HCA is a commonly used unsupervised pattern recognition method that seeks to build a hierarchy of classes based on similarity measurement. The results of HCA are presented in a dendrogram. Different arrangements of classes were obtained when choosing different threshold values. Besides, the class results were signi cantly affected by the choice of dissimilarity or distance measurement method. In this study, the Manhattan distance, or taxicab distance, was used, in which the distance between two points is the sum of the absolute differences of their Cartesian coordinates [27], which was de ned as Evaluation of microbial succession patterns Random forest, proposed by Breiman [28], is a popular approach in applied statistics due to its easy applicability to classi cation and regression problems. It creates random decision trees based on subsets of features (e.g., taxa) in the data and then chooses the subsets of features that are best able to classify samples into prede ned groups using only part of the data. This technique shows high predictive accuracy and applicability even in high-dimensional problems with highly correlated variables, which often occurs in bioinformatics [29]. Recently, some researchers used RF in PMI prediction studies [8,30] and other microbial ecology studies [31], providing high classi er accuracy. In this study, a PMI prediction model was established based on the microbial community succession trend dataset mentioned above with the RF algorithm to evaluate the microbial succession patterns. The X-variable corresponded to the matrix of microbial community succession trends and the responding Y-variable was associated with PMI values. The data set was randomly divided into two parts: 80% as the training set, 20% as the validation set. Training refers to tting or building the model while validation is equivalent to predicting. For murine body decomposition, each subset was partitioned based on an individual for cross-validation so that the samples from the same individuals are either in the training set or validation set, but not both. In the current study, the best RF regression model after hyperparameter tuning through cross-validation was selected to represent the nal model. To evaluate the stability and predictive ability of the RF model, the external set was preprocessed and analyzed following the above method and then formed the testing data set. This external testing set was performed using constructed models to predict the PMI. The MAE, RMSE, and , as the three main parameters of the model calibrated and predicted results, were used to evaluate the regression model reliabilities. The MAE was calculated as the deviation of the predicted from observed values and representing the average prediction error in the same unit of the original data. The RMSE represents the quadratic mean of these differences between predicted values and observed values. The can evaluate the goodness of tting between actual and predicted values. High values of and low values of RMSE and MAE demonstrate a well-established RF regression model. In particular, the model was run and evaluated fteen times in the validation and test sets, and the mean value of the MAE and as shown in Table 1. The modeling was done with Matlab software version R2019b (MathWorks, Natick, MA, USA) equipped with Windows-Precompiled-RF_Mexstandalone-v0.02toolbox (https://randomforestmatlab.googlecode.com/ les/Windows-Precompiled RF_Mexstandalone-v0.02-.zip), and the number of trees and predictors sampled for splitting at each node were both 500.
Additionally, mean decrease accuracy (MDA) [32] values were used to evaluate variable importance in this PMI prediction model. The MDA takes the importance of out-of-bag error and the impurity of the variable into consideration. The high value of the MDA showed the great importance of the variables to this model.

Statistics and reproducibility
Graphical illustrations are plotted using Prism software (GraphPad, San Diego, CA). The P-value cut-off used in this study is 0.05. *, **, and *** in the gures refer to P-values ≤0.05, ≤0.01, and ≤0.001, respectively. For all comparisons between two groups, the unpaired Student's t-test or the Kruskal-Wallis test was performed depending on the data normality and homogeneity of variance. We considered biological replicates the sequencing results of biologically distinct samples.

Overview of microbial composition during decomposition
The average raw library size after paired-end assembly was 78,545 ± 10,670 (SD). After sequence trimming, quality ltering, and removal of chimeras, 17,696,903 high-quality sequences remained, with an average of 73,771 ± 9,688 (SD) reads per sample. Samples with library sizes close to, or larger than, 3000 reads showed convergence in their sampling effort curves. These high-quality sequences were assigned to a total of 5,099 operational taxonomic units (OTUs) by the UPARSE pipeline using a threshold of 97% identity. The OTUs belonged to 29 phyla, 67 classes,160 orders, 250 families, 473 genera, and 835 species. Among the 29 phyla, Firmicutes, Bacteroidetes, and Actinobacteria were dominant in all samples. (Fig. 1a). Lactobacillus was dominant in all groups at the genus level (Fig. 1b). The Lachnospiraceae NK4A136 group accounted for more than 15% of the total sequences before day 1. At the species level, Lactobacillus reuteri, Lactobacillus johnsonii, Enterococcus faecalis, Firmicutes bacterium M10-2, and Clostridium tetani E88 were the most abundant species (Fig. 1c). In general, the microbial community was dominated by OTUs with the Firmicutes (Fig. 1d) a liated within the Lactobacillaceae.

Dynamics of microbial communities during decomposition
Although there appeared to be a signi cant amount of stochastic variation across the duration of this experiment, the overall reproducibility of these patterns between animals suggested that deterministic microbial community succession decreased with time (Fig. 1). As a result, only the following genera remained at day 15 postmortem: Gordonibacter, Bi dobacterium, Enterorhabdus, Lactococcus, Clostridium sensu stricto 18, Clostridium sensu stricto 15, Anaerosalibacter, Enterococcus, Dubosiella, and Lactobacillus (Fig. 1b). At the species level (Fig. 1c), the relative abundance of Lactobacillus reuteri, Clostridium tetani E88, Lactobacillus johnsonii, and Enterococcus faecalis showed notable variation during decomposition, which may aid in PMI estimation. Lactobacillus reuteri increased and exhibited peaks in abundances at both day 7 and day 15. Clostridium tetani E88 appeared on day 7 and decreased until day 15. Lactobacillus johnsonii showed a higher relative abundance one week after death than in the advanced stage of decomposition. Besides, Firmicutes bacterium M10-2 appeared on day 2 and showed an immediate increase from day 2 to day 4. Enterococcus faecalis appeared on day 2 and increased until day 10. Clostridium tetani E88 appeared on day 4 and increased and then decreased from then on. At the OTU level, OTU 1613 gradually increased in relative abundance during the 15-day decomposition and nally became the dominant OUT at day 15 (Fig. 1d). Among the Bacteroidetes and Actinobacteria, the OTUs had a similar change trend, increasing up to their peak abundances at day 2 and decreasing until day 4. However, OTU 2202, annotated to Helicobacter, gradually decreased during the 15day decomposition (Fig. 1d). A heat map also showed discrepancies in microbial community abundance change trends during decomposition, with a color gradient from deep blue (low abundance) to deep red (high abundance) (Fig. 1e).

Classi cation of microbial succession patterns during decomposition
To analyze the microbial succession patterns, we applied hierarchical cluster analysis (HCA) based on the Manhattan distances between different taxa (Fig. 2a). Seven classes of taxon change patterns were obtained (Fig. 2b). Classes 1, 2, and 5showed a similar decreasing trend, which appeared with different PMIs. The taxa in lass 1 and class 5exhibited a fast decline before four days after death and a slow decline from then on. However, the declination range was greater in class 5 than in class 1. The taxa in Class 2 exhibited a continuous decrease during 15-day decomposition. In contrast, classes 3 and 4 showed an increase during the early PMI points. Class 3 showed a small increase two days after death, while class 4 showed a fast increase, a slow increase during 7-day decomposition, and a decrease from day 7 to day 13. Class 7 and class 6 assembled in clustering tree as both contained a lot of taxa with non-signi cant change trend during decomposition. Class 6 showed a relatively continuous pro le, while class 7 showed a slow decline during the whole 15-day decomposition. After species annotation, class 1 mainly contained Lachnospiraceae, Ruminococcaceae, and Muribaculaceae at the family level and Lachnospiraceae NK4A136 group, Lactobacillus, Bacteroides, Lachnoclostridium, and Prevotellaceae UCG-001 at the genus level. Class 2 mainly contained Muribaculaceae, Lachnospiraceae, Ruminococcaceae, and Lachnospiraceae at the family level, Lactobacillus, Lachnospiraceae NK4A136 group, Ruminococcaceae UCG-014, Ruminiclostridium 5, Helicobacter, Ruminiclostridium 9, Alistipes, and Lachnoclostridium at the genus level. The microbial composition in class 1 and class 2 were similar. Class 5 also mainly contained Muribaculaceae, Lachnospiraceae, and Ruminococcaceae at the family level, and Lachnospiraceae NK4A136 group, Ruminococcaceae UCG-014, Lachnoclostridium, Bacteroides, Prevotellaceae UCG-001, Lactobacillus, Helicobacter, and Ruminiclostridium 9 at the genus level. Class 7 contained Enterorhabdus and Dubosiella. Dubosiella, Lactobacillus, Lactococcus, and Enterococcus were dominant genera in class 3 and class 4, increasing during early decomposition. We also counted the number of taxa in each class (Fig. 2b). Classes 6 and 7 contained 58.95% of the microbial taxa in all patterns. However, class 4 contained the fewest microbial taxa among all classes and was unimodal with an increasing pro le. The blue lines in the plot represent relatively invariable taxa based on the Kruskal-Wallis test, which were mostly included in class 6. The other classes are mainly represented by the red lines, with signi cant variation during decomposition (P < 0.001, Kruskal-Wallis test). To display the changes in the patterns more obviously, the relative abundance changes of microbes in every class are shown in Fig. 3 (at the phylum level) and Supplementary Figure 1 (at the genus level). Speci cally, in class 1, the relative abundance of Firmicutes gradually decreased, while Proteobacteria increased with PMI. The relative abundance of Bacteroidetes decreased during decomposition in class 2. Firmicutes in class 3 gradually increased. Latescibacteria in class 4 increased during 15-day decomposition. Bacteroidetes in class 5 gradually decreased during decomposition. The phyla in class 6 and class 7 did not show obvious variation during decomposition. In conclusion, the Lachnospiraceae NK4A136 group, Lactobacillus, Muribaculaceae, and Lachnospiraceae showed a declining pro le, while Dubosiella, Enterococcus, Lactococcus, Clostridium_sensu_stricto_15, and Enterorhabdus were mainly increased during the 15-day decomposition. Then, we calculated the percentage of OTU annotation for the top 20 genera and every phylum. Firmicutes were the most dominant phylum in all classes (Fig. 4a). Class 6 contained the greatest number of phyla among the classes. Class 4 contained the fewest phyla, including Firmicutes, Actinobacteria, Bacteroidetes, Deferribacteres, Epsilonbacteraeota, and Proteobacteria. At the genus level (Fig. 4b), Lachnospiraceae NK4A136 group was the dominant genus in class 1. Lactobacillus was dominant in class 2. Dubosiella, Lactobacillus, and Enterococcus were dominant genera in class 3. Lactobacillus was the most dominant genus in class 4. Class 6 and class 7 contained many genera, while the distribution of the taxa in class 6 was discrete. Class 4 contained the fewest genera of the classes.
A PMI estimation model based on the microbial community succession trend data set A random forests (RF) regression model was established to predict the PMI based on the microbial succession trend data set. The root mean square error (RMSE), the mean absolute error (MAE), and the squared correlation coe cient (R-squared, ) values in the validation set were 30.03 hours, 20.42 hours, 0.945, respectively. Those values in the testing set were 28.38 hours, 20.01 hours, and 0.947, respectively. Besides, the importance of variables was evaluated using the mean decrease accuracy (MDA) in the prediction illustrated in Fig. 5. In this model, the taxa with MDA values above 17.93, as 1.5 times the mean value of all MDA values, were regarded as signi cant strains, which was represented by the orange line in Fig. 5. The microbial taxa with MDA values over 17.93 were from all classes. According to the species annotation, signi cant taxa were mainly related to Lactobacillus, Dubosiella, Enterococcus, Lachnospiraceae NK4A136 group, Prevotellaceae UCG-001 and Clostridium sensu stricto 15 at the genus level.

Discussion
Decay and decomposition in mammalian and other vertebrate taxa are a key process in biological nutrient cycling. Cooperation between vertebrate and invertebrate scavengers (including bacteria, archaea, fungi, and protists) plays an important role in the chemical decomposition of animal waste. In particular, microbial decomposers are expected to lead to the conservation of key biochemical metabolic pathways and cross-kingdom ecological interactions for the e cient recycling of nutrient reserves. When a mammalian body is decomposing, microbial and biochemical activities result in a series of decomposition stages [9] that are associated with reproducible microbial succession across mice [4], swine[8], and human cadavers. Thus, a deep investigation of the behavior of this group during decomposition could be very helpful.
In the current study, we utilized hierarchical classi cation analysis with Manhattan distances to classify the relative abundance changes in the microbiota during decomposition. Seven classes of succession patterns were detected in the current study. The seven changing patterns could be further divided into unimodal, continuous decline, and a single peak followed by an ascending phase. The declining pattern was the most dominant pattern of microbial succession. We focused on intestinal microbiome changes during decomposition. At the very beginning of death, the microbes in the cecum were almost all normal intestinal tract bacterial ora. These microbes play a vital role before death in functions such as the decomposition of food, degradation of drugs, promotion of fat accumulation, and synthesis of vitamins and amino acids. Once the body dies, the living environment of the microbial organism changes, inducing a decline in microbial diversity and richness. The Chao1 estimator showed a decline in microbial community richness within 4 days of death and a subsequent gradual increase until 13 days (Supplementary Figure 2a). The Shannon index demonstrated that the diversity of the microbial community also declined during 15-day decomposition (Supplementary Figure 2b). These results coincided with the dominant microbial succession pattern pro les (Fig. 2b).
However, there remained discrepancies among 7 classes of postmortem microbial successional patterns. Class 1, 2, and 5 showed similar decrease patterns. Speci cally, class 1 showed a fast decline four days after death and a slow decline from then on. The microbiome in classes 2 and 5 exhibited a continuous decrease during 15-day decomposition, while the range of decrease was higher in class 5 than in class 2. The reasons for these changes may be that the microbial compositions in these classes differ. After species annotation, the microbial compositions in classes 1 and 2 were similar. A hierarchical cluster plot ( Fig. 2a) also showed the similarity between class 1 and class 2. However, the relative abundances of Muribaculaceae and Lactobacillus were higher in class 2, while the relative abundance of Lachnospiraceae NK4A136 group was higher in class 1. Lachnospiraceae NK4A136 group, a genus of the gut microbiota, is involved in the pathogenesis of age-related cognitive impairment in mice [33]. Classes 1, 2, and 5 all contained Lachnospiraceae, Ruminococcaceae, and Muribaculaceae at the family level. The bacteria in Ruminococcaceae are obligate anaerobes, while Lachnospiraceae is a family of anaerobic, spore-forming bacteria in the order Clostridiales that ferment diverse plant polysaccharides to short-chain fatty acids (butyrate, acetate) and alcohols (ethanol). Based on NCBI taxonomy, Muribaculaceae, Lachnospiraceae, and Deferribacteraceae were the predominant families identi ed in the healthy mouse gut microbiome. Researchers [34] also found that the genes for carbohydrate metabolism were upregulated in Muribaculaceae, while genes for cofactors and vitamin metabolism and amino acid metabolism were upregulated in Deferribacteraceae. These bacteria were reported in the mammalian gut microbiota [35,36]. Class 7 showed a slow decline during 15-day decomposition.The microbial composition of class 7 was slightly different from those of classes 1, 2, and 5, which was also demonstrated in the hierarchical cluster results. Class 7 included Enterorhabdus and Dubosiella, while classes 1, 2, and 5 did not. Some studies reported that the relative abundance of the genus Enterorhabdus increased in subjects with prediabetes compared to healthy subjects [37]. In conclusion, classes 1, 2, and 5 showed similar decreasing pro les, and they were mainly annotated to Ruminococcaceae, Lachnospiraceae,Deferribacteraceae, and Muribaculaceae at the family level. The transformations listed above might be caused by changes in the internal conditions of the mouse after death; in particular, the oxygen content decreased, which was not suitable for obligate anaerobes. Conversely, classes 3 and 4 showed an increased pro le during the early PMIs. Class 6 showed a relatively invariable pro le during the whole 15-day decomposition. Dubosiella, Lactobacillus, Lactococcus, Clostridium sensu stricto 15, and Enterococcus were dominant genera in class 3 and class 4, with increasing patterns during early decomposition. Dubosiella, a gram-positive genus from the family Erysipelotrichidae, has been isolated from the intestinal content of mice. Lactococcus species are capable of proteolysis, a process that converts amino acids into avor complexes, which may aid in decomposing the protein of corpses. Lactobacillus is a genus of gram-positive, facultative anaerobic or microaerophilic, rod-shaped, non-spore-forming bacteria. The relative abundances of the genus Enterococcus, as facultative anaerobic and gram-positive bacteria, suggested that as time continued after death, the proportion of normal intestinal ora became maladjusted, and opportunistic pathogens rose immediately. Classes 3 and 4 also contained Enterorhabdus at the genus level, which are grampositive, facultative anaerobe, non-spore-forming rods that exist in mammalian intestines [38]. The possible reason for this facultative anaerobe increase is that body rupture permits oxygen entry. At the species level, Lactobacillus reuteri and Lactobacillus animals were dominant in class 4. These species were reported as commensal probiotic bacteria that produce antimicrobials to ght off harmful infections and mediate the immune system [39,40]. They both showed the highest relative abundances on day 7 after death. This nding provides the possibility of extracting useful microbes from the decomposing microbiota. Thus, local microbes are most dominant in the intestinal tract at the early stage. With corpse bloating, the oxygen inside the body is consumed by the decomposition of nutrient substances. Meanwhile, the bacteria of the natural gut microbiome spread as tissue lique es. When the body enters advanced decay, there is a proliferation of obligate anaerobes. After body rupture, facultative anaerobes proliferated inside the gut track. Moreover, we observed the relative abundance change patterns of each OTU in every class. The original curves of taxa were aligned with the tting curves of every class. Thus, it is obvious that the classifying and tting methods used in this study are suitable and reliable. Interestingly, the succession patterns showed similar trends in general, and the change rates for 4-day decomposition were higher than those for the subsequent processes. From this point, the microbial decomposers may have similar behavior during the degradation of bodies, which improves understanding of the ecological cycling rules in the microbial kingdom.
A RF model was established with high at 0.95 and low MAE at 20.01 hours during 15-day decomposition based on postmortem microbial successional trend data (Table 1). To evaluate the PMI prediction ability of each class, the MDA value was calculated for every microbial taxon. Classes 2 and 7 contained the largest number of important variables for predicting PMI based on the MDA value of each taxon (Fig. 5). However, class 6 contained the fewest among the 7 classes. This result is explained in Fig. 2b that class 6 contained many invariable taxa, while the abundances of members of classes 2 and 7 were constantly changing during the 15-day decomposition. The important taxa in class 2 were mainly initial colonizers of the mammalian gastrointestinal microbiota, most of which are anaerobic or aerotolerant anaerobic bacteria, such as Helicobacter, Bacteroides, Alistipes, Muribaculum, Lactobacillus, Enterorhabdus, Adlercreutzia, Anaerotruncus, and Mucispirillum. When the bodies died, the living conditions of the intestinal microbiota gradually changed, resulting in a decline in the relative abundances of these initial microbes. The signi cant taxa in class 7 were mainly obligately anaerobic microbes, for example, Desulfovibrio (found in a variety of habitats, including soil and the intestines, and feces of animals), Roseburia (gram-positive anaerobic bacteria that inhabit the human colon), Anaerovorax (strictly anaerobic), and Allobaculum (strictly anaerobic). The relative abundances of these anaerobic microbes decreased immediately after body rupture, which may facilitate the prediction of the early PMI. As a result, Lactobacillus, Dubosiella, Enterococcus, Lachnospiraceae NK4A136 group, and Clostridium sensu stricto 15 were the most important genera for the PMI estimation model, which also coincides with results of previous studies.
This study provides basic knowledge for future studies on microbial decomposers. However, there are still knowledge gaps existing in the natural decomposing kingdom. Future work needs to focus on the succession patterns of fungi, heterotrophic organisms (scavengers), and their relationship with PMI estimation. Our group is also dedicated to expanding the application of this work to human remains in forensic practice.

Conclusion
The present study contributes to our knowledge about postmortem microbiology as the thorough characterization of microbial succession patterns. We demonstrated that postmortem microbial communities differ over time. However, their changing trends have some similarities, as the successional rate changes faster at early decomposition and slower later. A random forest model was established to predict the PMI with an MAE of 20.01 hours during 15-day decomposition based on the trend data set.
Taken together, these succession patterns are con rmed as critical information to consider when developing the use of microorganisms as physical evidence for PMI estimation.

Declarations
Funding This work has been supported by the funders of the National Natural Science Foundation of China (No. 81730056).

Con icts of interest
The authors declare that there is no con ict of interest.

Availability of data and material
The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.

Code availability
The code used in this work is available on request from the corresponding author.   Figure 1 Taxonomic pro les of postmortem microbial communities in the intestine at different levels (a: at the phylum level, b: at the genus level, c: at the species level, d: at the OTU level) and e: a heat map of the microbial 16S rRNA gene sequences abundances (OTUs) during decomposition, indicating that there might be discrepancies in OTU-relative abundance of different PMIs with a color gradient from deep blue (low abundance) to deep red (high abundance). Classi cation of microbial succession patterns during decomposition. a: Hierarchical classifying tree plot based on Manhattan distances between taxa. b: Postmortem microbial succession pattern pro les. The X-axis represents the postmortem interval (PMI) hour. The Y-axis represents the accumulated value of the changing trend in every taxon at every PMI. The red lines indicate that the changing trend of that taxon is signi cant, while the blue lines indicate that the changing trend of that taxon is not signi cant (P < 0.001, Kruskal-Wallis test). The black line in every class shows the median trend value of taxa in the class. shows the results at the genus level.)