Potential Using of Molecular Characterization and Structure of The Gut Bacterial Community for Postmortem Interval Estimation in Sprague Dawley Rats

Once the body dies, the inherent microbes of the host begin to break down from the inside and play a key role thereafter. It is hypothesized that after the death certain rectal microbes would change during the decomposition course in the body. This study aimed to investigate the probable shift in the composition of the rectal ora at different time intervals up to 15 days after death and to explore bacterial taxa important for estimating the time of death. At the phylum level, Proteobacteria and Firmicutes showed major shifts, when checked at 11 different intervals, and emerged at most of the postmortem intervals. At the species level, Enterococcus faecalis and Proteus mirabilis existed at most postmortem intervals; the former showed a downward trend after day 5 postmortem, while the latter showed an upward trend. There were obvious differences in bacterial community structure and richness at the phylum, genus, and species levels during the decomposition of the corpse of rats. The phylum, genus, and species taxa richness decreased initially and then increased signicantly. The turning point came on day 9 when genus, rather than phylum or species, contained the most information for estimating the time of death. We constructed a prediction model using genus taxon data from high-throughput sequencing, which explained 87.2% of the time since the rst sampling within 1 h. Seven bacteria, namely Enterococcus, Proteus, Lactobacillus, unidentied Clostridiales, Vagococcus, unidentied Corynebacteriaceae, and unidentied Enterobacteriaceae, were included in this model. The above-mentioned bacteria showed a promising future for estimating the shortest time of death and results of current study were agreeing with the proposed hypothesis.


Visible decay progression
The decomposition process of the corpses of rats during 15 days was recorded and classi ed into ve stages: The fresh stage began at 0.0 ± 0.0 h with no odor emitted; the bloat stage started along with body expansion giving off odorous gases at 2.6 ± 1.1 days; the active stage started at 5.0 ± 1.0 days with the body being ruptured by accumulated gases and several parts of the tissues were broken down along with plenty of liquid owing out; the advanced stage began at 8.0 ± 1.0 days along with most parts of the tissues being removed; the dry stage happened at 12.5 ± 1.9 days with no soft tissue left.
Relative abundance of gut ora in different groups A total of 7,029,815 raw and 6,674,323 clean reads were obtained by performing the high-throughput sequencing, with an effective rate of 94.97% (the ratio of clean to raw reads). A total of 22,625 OTUs were identi ed based on 97% similarity, with an average of 257 OTUs per sample. The total usable sequences were classi ed into 33 phyla, 49 classes, 108 orders, 203 families, 465 genera, and 306 species. Species accumulation boxplot and rarefaction curves of all the samples were smooth as the number of sequences increased, demonstrating that this sequencing profundity could mirror the complete bacterial species richness among the samples (See Additional File 1). The basic information regarding the number of OTU and alpha diversity indices in individual rats revealed that these parameters in the samples collected at different time points after death reduced sharply in comparison to the living individuals (See Supplementary  Table 1). A Venn diagram was plotted to compare the similarities and variances among the communities obtained in the different groups. The eleven groups showed communities of 36 OTUs in common, with the unique OTUs composed of 84.94%, 90.00%, 74.83%, 87.84%, 21.74%, 14.29%, 5.26%, 10.00%, 16.28%, 23.4%, and 26.53% at time points corresponding to alive, 0 h, 8 h, 16 h, 1 day, 3 days, 5 days, 7 days, 9 days, 13 days, and 15 days, respectively, (Fig. 1).
Ternary plots were constructed using the Ternary plot order of the VCD package in R software. In Fig. 2a, we observed that Corynebacterium amycolatum exhibited the highest abundance and the highest load at 0 h time point among the three-time points and that Bacterium mpn isolate group 2 and Falsiporphyromonas endometrii were the highest in living individuals and at 8 h after death time point, respectively. In Fig. 2b, we observed that the richness of Enterococcus faecalis increased at 16 h time point, and at 1 and 3 days postmortem, followed by Proteus mirabilis, which was highest at 3 days postmortem as compared to the other two time points. In Fig. 2c, we observed that the relative abundance of E. faecalis increased at 5 days postmortem, and exhibited almost the same abundance at 7 and 9 days postmortem. In Fig. 2d, we observed that P. mirabilis exhibited the highest abundance among all the three-time points on day 15 after death, followed by Vagococcus lutrae and E. faecalis that were much higher at 9 and 13 days postmortem than that at 15 days after death, respectively.

Microbial analysis at different levels
The microbial community structure was determined during the succession of decomposition, and all the 16S rRNA sequences were classi ed at the phylum, genus, and species levels. The notable tendencies and uctuations exhibited the relative richness of the diverse bacterial taxa in the rectum of the rat cadavers through the decaying process ( Fig. 3a, b, c; Fig. 4a, b, c; Table 1). Figures 3a, b, c showed the variations in the proportions of bacteria at different levels, and Figs. 4a, b, c showed the relative abundance of the ten topmost bacteria identi ed in the study samples. the living samples was much lower than that at 8 h, 1 day, 3 days, and 9 days postmortems, while the relative abundance of Bacteroidetes was much higher in the living samples than that at 8 h, 5 days, 7 days, 9 days, 13 days, and 15 days postmortem, and the relative abundance of Actinobacteria was signi cantly higher in 0 h than that at days 5, 7, 9, 13, and 15 postmortems (P < 0.05).
At the genus level, Lactobacillus and Enterococcus appeared as the dominant genera at day 1 and 3-13 days postmortem, respectively, Helicobacter disappeared at days 7, 9, and 15 PMIs and Proteus was the most abundant at day 15 of postmortem. The relative abundance of rectum ora in the living samples was signi cantly lower than that at days 3, 5, 7, 9, and 13 PMIs; nevertheless, Helicobacter was much higher in the living samples than those in the postmortems of the above-mentioned days (P < 0.05). The proportion of Proteus was signi cantly lower in the living samples and the 0 h, 8 h, 16 h, and day 1 postmortem samples than those in the 13-and 15-day postmortem samples, while Lactobacillus exhibited an opposite result (P < 0.05).
At the species level, among the ten topmost species that existed at 8 h postmortem, Clostridium sporogenes and F. endometrii disappeared before 1-day postmortem and after 3 days postmortem, respectively. E. faecalis and P. mirabilis appeared during the whole decomposition process of 15 days after death; however, the former showed a downward trend from day 5 of postmortem, while the latter showed an upward trend. Bacterium mpn isolate group 2 disappeared during 5-13 days postmortem or decreased after death. The relative abundance of E. faecalis at days 3, 5, 7, 9, 13, 13, and 15 after death was much higher than that in the living samples, while Bacterium mpn isolate group 2 was distinctly lower (P < 0.05). The relative abundance of C. amycolatum was markedly lower in 0 h, 8 h, 16 h than in 7, 9, 13, and 15 days postmortem samples, while Lactobacillus reuteri and L. intestinalis were higher in 0 h, 16 h than at 13-and 15-days postmortem samples.

Characterization of bacterial diversity and community structure
The complete rectal ora community was evaluated based on diversity and richness, which was calculated at 97% similarity. Alpha diversity indices of the observed species, abundance-based coverage estimators (ACE), and chao1 values for the rectal bacteria in the living samples were signi cantly higher than those in the 5, 7, 9, 13, and 15 days postmortem samples, suggesting that the richness and diversity of the rectal ora declined signi cantly after day 5 of postmortem compared with the alive group ( Table 2). The richness indices (ACE and chao1) went up slightly, however, there was no signi cant difference compared with other time points. All the alpha diversity indices are presented in Table 2, and there were signi cant differences in the overall rectal bacterial community structure among the eleven postmortem intervals.
The similarities in the gut ora communities of rats among the eleven groups were estimated using the beta diversity metrics, such as NMDS and beta diversity heatmap. As presented in Fig. 6a, the differences in coe cients among all the groups were almost higher than 0.5, indicating that the bacterial community in different groups exhibited great diversity. All the samples were clustered into 11 prime clusters. According to the NMDS (stress = 0.152), the bacterial communities of the gut samples were separated into three clusters between the late and early PMI (Fig. 6b). Conspicuously, the 8-h postmortem interval could be signi cantly separated from the other groups, indicating that the gut ora at 8 h after death differed from the other two clusters.
LEfSe is a biomarker detection and descriptive tool for performing high-dimensional statistics. The LEfSe analysis was performed to compare the projected bacterial community among the 11-time intervals at different levels ( Fig. 6c). The results of this analysis suggested that the provision of related taxa was signi cantly diverse among all the groups. The LDA scores revealed that the relative abundances of C. amycolatum, Entero isolate group 2, Bacteroides uniformis, E. faecalis, Streptococcus gallolyticus subsp macedonics, and C. sporogenes were most abundant at postmortem intervals of 0 h, 1 day, 3 days, 5 days, 7 days, and 13 days, respectively, while P. mirabilis and V. lutrae were most abundant on day 15 of postmortem.

Constructing a model for PMI
The best subset selection was combined with phylum, genus, and species indicators to construct the model that best explained aberrance in time of death from 0 h to day 15 (See Additional le 2 and 3; Table 3). Table 3 showed that the poorest model for PMI appeared at the phylum level (16.1% of the variation). The taxon that most contributed to postmortem existed at the genus level. The best subset selection results showed that seven bacteria at the genus level were selected as the best features to develop the model. These seven genera were: Enterococcus, Proteus, Lactobacillus, unidenti ed Clostridiales, Vagococcus, unidenti ed Corynebacteriaceae, and unidenti ed Enterobacteriaceae, and this model contained the most information and explained 87.2% (generalized cross-validation score (GCV) = 0.307) variation of time of death in this study. The species, including four bacteria, were identi ed from best subset selection as the most informative and explained 56.6% (GCV = 0.515). This model was poorer in explaining the variation in time of death than the model constructed by the genus features. The whole models were evaluated by the adjusted R 2 value and the generalized cross-validation score (GCV). A higher R 2 and lower GCV suggest better model. The percent variation of PMI (%) explained by per model.

PICRUSt
The shifts in the probable functions of the gut ora of rats before and after death were inspected by predicting the 16S rRNA genes using PICRUSt (See Additional le 4 and 5). The top four different pathways at the Kyoto Encyclopedia of Genes and Genomes (KEGG) at levels 1 and 2 have been shown in Additional le 4, 5 and Table 4, 5. Among them, the functional pathways associated with metabolism, including carbohydrate and amino acid metabolisms, corresponded to a large number of related genes in all the samples (Table 4 and 5). The pathways related to environmental information processing and organismal systems were signi cantly higher in the rectal bacterial community on days 5, 7, 9, 13, and 15 after death as compared to the living samples (Table 4 and 5). Amino acid metabolism, energy metabolism, and metabolism of cofactors and vitamin pathways were signi cantly different between the living and 5, 7, and 9 days postmortem samples, while energy metabolism showed notable differences between the living and 5 days postmortem samples (Table 4 and 5) (P < 0.05). Although the pathway of carbohydrate metabolism did not differ signi cantly between the living and postmortem samples, the relative abundance increased in the 0 h and 3 days postmortem samples compared to that in the living sample.

Discussion
The microbiome composition in the living body is complex, and multiple signi cant differences have been observed between and within the individuals [21].
Microorganisms exist both inside and outside the dead body and exhibit a distinct and temporal shift during the process of decay [22]. The gut, as an organ contains a large number of bacteria, and some bacteria have been shown to exhibit variations during the decomposition of corpses [11,23,24], which prompted this study. If researchers want to use bacteria to infer the approximate postmortem interval, a database of the regulation (structure and composition) under different conditions during the cadaver degradation process should be developed rst. This study aimed to construct a usable model to estimate time of death and examine the bacterial composition and structure under speci c conditions of temperature (20.63 °C ± 0.93 °C) and humidity (15.37% ± 2.79%) in the process of rat corruption during the 15 days of decomposition, thereby contributing to a basic data accumulation for future use. A rat model was used because a large number of samples provided convenience in evaluating the intra-individual microbial distinction during decomposition, and the changes in microbiota composition existing in living individuals and other mammals are unknown [25].
The bacterial communities in the decomposing rats showed a distinct shift in the categories and relative abundance compared with those in the living rats.
Indices for evaluating bacterial community richness and diversity, such as observed species, Shannon, chao1, and ACE were signi cantly reduced after day 5 postmortem samples as compared to those in the premortem samples, which was in agreement with a previous study [24]. We also found that the richness indices (ACE and chao 1) showed consistent variation before 9 days postmortem, but displayed the opposite trend after 9 days, suggesting that the decomposition of the rats led to an evident change in the internal environment on day 9 and also external environmental factors in uenced the gut bacterial community structure and composition after day 9. In the decomposition process of rats, gas accumulation gave rise to bloating and rupture of the corpses, exchanging the internal condition with the external environment, which resulted in a decrease in the anaerobes such as Lactobacillus that even disappeared, while the facultative anaerobe Enterococcus seized the opportunity to thrive.
This study indicated a notable variation in the gut bacterial diversity and relative abundance during the decomposition course of 15 days at the phyla, genera, and species levels. The results of this study showed that the main phyla in the intestinal samples that existed in high abundance before death were Firmicutes and Bacteroidetes, which were similar to the NIH Human Microbiome Project [26] (Fig. 4a). Although the dominant bacterial community diversity at the phyla level showed no evident uctuation, we found that changes in the phyla Firmicutes and Proteobacteria were in the opposite direction of the relative abundance over time, which was consistent with a previous study [27]. Previous research aiming at the epinecrotic communities on human corpses showed that Firmicutes formed the inherent phylum for corpses, while Proteobacteria was initially identi ed from the environment, and the reason for the increase in Proteobacteria was outer bacterial migration [28][29][30]. This suggests that the increased relative abundance of Proteobacteria in the rectum of dead rats after 8 h may be due to the migration from the external environment. At the genus level, we found that Lactobacillus was the dominant bacteria at the time points before day 1, while after day 1, Enterococcus and Proteus increased as the most abundant bacteria. We suggested that this was because oxygen pressure decreased, and pH value increased because of protein and carbohydrate degradation. The relative abundance of Bacteroidetes declined notably after death, which was similar to the preceding analysis results [10]. As previously reported [31], high levels of Vagococcus might account for the larvae of blow ies, indicating that the carrion insects might participate in the decomposition process (the eggs of larvae might be from the rats themselves; for Vagococcus in the alive samples abundance was 0.001%). LEfSe results suggested that the seven bacteria at the species level were identi ed as seven potential PMI indicators. B. uniformis belongs to Bacteroides spp. and can be regarded as a PMI indicator of day 3 postmortem. It has been reported that [32] Bacteroides spp. could be used as a quantitative indicator of PMI. E. faecalis, Streptococcus gallolyticus subsp. macedonics, and C. sporogenes were also included from the Enterococcus, Streptococcus, and Clostridium spp. that could also be used as PMI indicators on days 5, 7, and 13 postmortems, respectively, which have previously been reported as the most abundant species during decomposition [33]. Our study also showed an evident shift in C. amycolatum at early postmortem, whereas E. faecalis, P. mirabilis, C. sporogenes, and V. lutrae were identi ed to be most abundant at late postmortem intervals. Interestingly, the present ndings on P. mirabilis and V. lutrae indicate almost parallel changes, although the percentages were different. The above-mentioned two species of bacteria after day 15 postmortem showed a downward trend, while E. faecalis and C. sporogenes showed an upward trend. According to a previous study, P. mirabilis was able to attract blow ies, which was the reason for the high percentage of V. lutrae [34].
We constructed a model that could explain 87.2% variation of the time of death within 1 h. Regarding the decomposition process of a corpse as continuous variables for analysis and to develop a particularly ne model to estimate postmortem interval was similar to a previous study [35]. The features contained most information for evaluation death time were selected by Best Subset Selection coming from the genus level and features of the poorest model belong to the phylum. Previous research identi ed that the poorest model for estimating physiological time by epinecrotic bacterial community was composed of phylum features, which was consistent with our study [7], however, features of the best model were from the family taxon. They did not incorporate the genus of the bacteria into the research and the sampling position was different from ours. A study by Hunter et al. reported their focus on human skin useful bacterial community to develop a promising tool for time of death estimation. Their results were similar to those of this study and agree that genus was the most informative taxon [36]. Therefore, the genus taxonomic level was the most promising community for variation in PMI and was worth further research.
The results of PICRUSt suggested that the microbial function shifted signi cantly between pre-and postmortem intervals (Tables 4 and 5). In this study, we found that eight categories of pathways, including preponderant metabolism patterns, amino acid metabolism, and carbohydrate metabolism were associated with Bacteroides, Lactobacillus, Enterococcus, and Proteus, which were recognized for their involvement in proteolysis. These four bacteria participated in amino acid and carbohydrate metabolisms as the main force at different time points, transforming proteins into smelly gases, such as H 2 S, methane, ammonia, sulfur dioxide, and organic acids (e.g., propionic and lactic acids). As reported previously, primary degradation of carbohydrates was performed by Bacteroides, and oligosaccharides were fermented by Lactobacillus to release gaseous byproducts such as hydrogen, carbon dioxide, H 2 S, and methane [32,37]. In our study, Bacteroides and Lactobacillus showed the most evident variations at the early death time points, while Enterococcus and Proteus were at the late time points, suggesting that in the early death, carbohydrates might be the main energy source for the anaerobic bacteria. The above-mentioned processes were considered to result from a reduction in the ability to obtain oxygen [38]. The byproduct gases lead to the rupture of the abdominal cavity, resulting in the shift in bacteria in the cavity to that of aerobic bacteria [2]. Lactobacillus was also identi ed for breaking down the aromatic amino acid (tryptophan) to a smelly indole substance with proper enzyme and protein fermentation occurring at an elevated pH [39,40]. This implies that the pH of the rectum of rats may rise after death, contributing to protein lysis and odor emission. Proteus, which produces urease that is capable of dissolving urea to ammonia and CO 2, usually increases in dysbacteriosis and pathologic conditions [41], indicating that the mass of the protein of corpses was dissolved after day 3 in this study. Therefore, the accumulation of ammonia and other gases in the rectum further leads to an increase in the pH and an increase in the oxygen pressure. Furthermore, following a decrease in the numbers of strict anaerobic bacteria, the facultative anaerobic Enterococcus showed an evident tendency to increase. It has been found that some kinds of proteins can be decomposed into amides and amino acids by Enterococcus. Overall, our results maintained the concept that the four above-mentioned genera might act as imperative donors in the process of decomposition. C. sporogenes mainly occurs before day 1 postmortem as a Gram-positive, obligate anaerobic species that possesses the capacity to decompose carbohydrates and peptones into organic acids and alcohols, and also participates in tryptophan metabolism producing smelly 3-indolepropionic acid [42], suggesting that the emitted odors in the fourtime points (0 h, 8 h, 16 h, and 1 day postmortem) were probably produced by C. sporogenes. In our study, P. mirabilis, which is a Gram-negative facultative anaerobic bacterium, appeared on 13 and 15 days PMIs, which could also produce H 2 S gas, indicating that the smell in the last two time points of this study was associated with P. mirabilis. Nevertheless, further studies will be required to better understand the role of the above-mentioned species in PMI after day 15 of postmortem.
Although the results of this study provided a detailed description of the bacterial composition within a decomposing cadaver system and recommend that the microbial community data can be evolved into a legal medicine means for estimating the PMI, additional research will be required to better comprehend this perception.

Conclusion
Taken together, the bacterial community exhibited a distinct shift during the 15 days of decay in both composition and structure. Proteobacteria and Firmicutes exhibited opposite patterns in the whole decomposition process in this study. The most abundant bacteria at the genus level showed signi cant differences between pre-day 1 and post-day 1. Additionally, the genus taxon bacteria were the most potential features for estimating PMI. Therefore, these ndings offer the foundation for the analysis of the bacterial community at speci c time points after death.

Method of sample collection and selection of experimental time points
Rectal samples were collected after the death of living rats. The postmortem putrefaction of all the dead rats took place in Xi'an City, Shaanxi Province, China (34°15′39.9″°N, 108°56′33.32″°E) in December 2017. Samples of the rectal ora from the dead rats were swabbed with sterile cotton swabs dipped in sterile saline for one minute before and after the execution to obtain more rectal bacterial samples according to the previously described protocols [10][11][12].

DNA extraction
The total genomic DNA of the bacterial samples was extracted from the swabs using the QIAamp DNA Mini Kit (Qiagen, Germany), and the speci c procedures were performed according to the manufacturer's instructions. The DNA concentration was determined using NanoDrop2000 (Thermo Scienti c, Waltham, MA, USA) and the extracted DNA was stored in a refrigerator at -80°C until further use.

High-throughput sequencing
Operation ow Single-end sequencing (SE600) was performed on Thermo sher's IonS5TMXL sequencing platform to obtain high throughput sequences. The total DNA extracted from the bacteria was ampli ed with the 16S rRNA V3+V4 universal primers, 341F (5′-CCTAYGGRBGCASCAG-3′) and 806R (5′-GGACTACNNGGTATCTAAT-3′). The ampli ed products were recovered, puri ed, and quanti ed, and then the corresponding mixing ratio of each sample was adjusted according to the quantitative results. Thereafter, the library was prepared, and sequencing was performed using the sequencing platform.
Processing of the sequencing data The Cutadapt (v 1.9.1) software [13] was used to remove ambiguous bases (N), organize the sequencing data according to the barcodes, remove low-quality bases, barcodes, and primers, and then to obtain the original data. Clean reads were obtained using the UCHIME algorithm and comparing with the species annotation database to detect the chimera sequences and remove them.

OTU clustering and species annotation
The Uparse (v 7.0.1001) software [14] was used to cluster the clean reads with 97% identity to form the operational taxonomic units (OTUs). The OTUs with the highest frequencies were selected as representative sequences. The Mothur method was used to analyze the species annotations using the SSUrRNA [15] database of SILVA132 [16]. The MUSCLE (v 3.8.31) software [17] was used for rapid multi-sequence alignment to obtain the phylogenetic relationships of all OTU sequences.

Sample complexity analysis
Alpha diversity (chao1, ACE, Shannon index, and Simpson) was calculated using the QIIME (v 1.9.1) software [18], and the species accumulation curves were drawn using the Rarefaction Curve tool of the R software (v 2.15.3). Alpha diversity among the different groups was compared, and the differences were determined by performing the non-parametric Kruskal-Wallis test of the Dunn's multiple comparison test, and the results were considered statistically signi cant when P < 0.05.

Diversity comparison analysis
UniFrac distance was calculated using the QIIME (v 1.9.1) software, and the hierarchical cluster of the samples was constructed using the unweighted paired group method with arithmetic mean. The non-metric dimensional scaling (NMDS) diagram was drawn using the vegan package of the R software. The beta diversity index was analyzed using the R software, and the parametric and non-parametric tests were performed subsequently. Linear discriminant analysis effect size (LEfSe) [19] was used to identify the microbial taxa and predict functional genes (PICRUSt) that were abundant in the gut at different successions of time points, which was based on the LDA score > 2.0 and P < 0.05. Analysis of similarities was performed by using the Adonis function of the vegan package of the R software, and the species analysis with signi cant differences between different groups was performed using the R software.

PICRUSt
The phylogenetic investigation of communities by reconstruction of unobserved states (PICRUSt) [13] was used to determine the bacterial function based on the Kyoto Encyclopedia of Genes and Genomes (KEGG) database [20], and the abundance of metabolic pathways and the Bootstrap Mann-Whitney U-test was applied for the detection of gene pathways or OTUs with signi cantly diverse abundance among different groups.

Statistical analysis
All statistical analyses were performed using the GraphPad Prism (v 5.01) software or the R (v 2.11.1) package, and the Kruskal-Wallis test of the Dunn's multiple comparison tests was performed. P-value < 0.05 was considered statistically signi cant. The R statistical package (v 3.6.3) was used to construct the best subset selection and the mgcv (v1.

Availability of data and materials
The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request.                Beta diversity of bacterial populations in eleven time points re ects inter-group differences. Heatmap was drawn by the Unweighted Unifrac distance (a). The microbial diversity in a certain time point increase follow with the size of value. NMDS (non-metric multi-dimensional scaling) coordination plot among the whole samples (b). Stress value less than 2 indicates that NMDS can accurately re ect the difference between samples. Signi cant taxa obtained in the sampling time points using LEfSe analysis. Linear discriminant analysis (LDA) plots of bacteria at different levels as results of all time points (c).