Microbiome-Metabolome Analyze The Immune Microenvironment of Cecal Contents, Soft and Hard Feces of Hyplus Rabbits

Background: Intestinal microbiota and its metabolites play a vital role in host growth, development and immune regulation. Methods: This study analyzed the microbial community distribution, cytokines and short chain fatty acids (SCFAs) content of cecal contents (Con) group, soft feces (SF) group and hard feces (HF) group of 140-day-old Hyplus rabbits, and veried the effect of soft feces on cecal immune microenvironment by fasting soft feces. Results: The results showed that there were signicant differences in the levels of phylum and genus composition, cytokines and SCFAs among Con group, SF group and HF group. In addition, metabolic pathway enrichment analysis found that there were two signicantly up-regulated differential metabolic pathways in SF group and HF group compared with con group, P125-PWY, namely the super pathway of (R, R) - butanediol biosynthesis and P341-PWY, namely the glycolysis V pathway. At the same time, Christensenellaceae_R-7_Group and Lachnoclostridium are signicantly enriched in the above two pathways. The correlation analysis of cytokines and SCFAs with differential microbial communities showed that Muribaculaceae and Ruminococcaceae_UCG-014, Ruminococcaceae_NK4A214_group and Christensenellaceae_R-7_Group are closely related to cytokines and SCFAs. After fasting soft feces (CP), the contents of SCFAs and cytokines (IL-4, IL-10) in cecum decreased signicantly, while cytokines (TNF-a, IL-1β) increased signicantly. The results of multiple immunouorescence showed that the expression of claudin-1, occludin and ZO-1 related to intestinal immune barrier increased signicantly in CP Group. Conclusions: In conclusion, soft feces are not only rich in probiotics and SCFAs, but also play a very important role in improving the immune microenvironment of cecum. This study may provide a valuable reference for the treatment of inammatory intestinal diseases. analysis effect size; PCA: Principal component analysis; PCoA: principal coordinate analysis; IBD: Inammatory bowel disease; TJ: tight junctional; ZO-1: Zonula occludens 1.


Introduction
Animal feces contain a large number of intestinal microorganisms. Its main function is to digest and absorb nutrients in food and macromolecular substances that the host itself can not completely decompose, so as to provide the host with energy and nutrients necessary for growth and development [1]. Carbohydrates are an important energy source in animals. There are a large number of polysaccharides that cannot be decomposed and utilized in the food eaten by animals, and microorganisms in the intestine provide the host with enzymes and biochemical metabolic pathways that the host does not have, so that the host can produce energy by fermenting the degraded polysaccharides (starch, cellulose, hemicellulose, colloid) and absorbed oligosaccharides, at the same time, a large number of short chain fatty acids (acetic acid, propionic acid and butyric acid) are produced [2][3][4][5].
Microorganisms can ferment and decompose polysaccharides in food into short chain fatty acids, which are then used by other organs and cells [2,6]. Intestinal microbiota helps to maintain the integrity of intestinal epithelial mucosal barrier function and causes immune system mediated mucosal protection.
Intestinal microbiota regulates the tight junction of intestinal epithelial cells and protects intestinal epithelial cells from damage by controlling the proliferation rate of intestinal epithelial cells and inducing cytoprotective proteins [7,8]. Mucosal epithelial barrier is the rst line of defense against the invasion of intestinal pathogenic microorganisms and toxins. As an important part of intestinal mucosal barrier, the change of tight junction protein can cause abnormal intestinal barrier function and affect the intestinal health of piglets [9], so it plays an important role in the occurrence of intestinal in ammation and diseases [10]. Intestinal microbiota can also induce the intestinal cavity to produce secretory immunoglobulin LGA (SGA) and antimicrobial peptides (AMPS) [11,12]. Microbial colonization can stimulate intestinal goblet cells to secrete mucin, which, plays a vital role in maintaining intestinal health [13]. In addition, more and more studies have shown that host behavior can affect intestinal microbiota, and intestinal microbiota can also affect animal behavior [14,15].
Coprophagy refers to the behavior of animals feeding on feces, which is very common among small-and medium-sized mammalian. Animal feeding on feces includes feeding on their own feces and the feces of other animals to meet their nutritional needs [16,17]. Dung eating behavior is of great signi cance to small and medium-sized herbivores. Banning dung eating can cause weight loss or growth retardation [18]. Rabbit soft feces contain a large number of bacterial proteins and nutrients. Eating soft feces can prolong the time of feed passing through the digestive tract, improve the digestion and absorption e ciency of feed, help to maintain the normal microbiota of the digestive tract and alleviate nutritional diseases [19,20]. Prevent rabbits from eating dung, reduce the content of vitamin A in the blood, slowly lose weight, increase the content of acetic acid and propionic acid in the blood, reduce the content of butyric acid, and skin damage at eyes and ears. It can be seen that dung eating behavior plays a very important role in the nutrition and health of herbivores. In addition, fecal eating behavior can also help herbivores obtain the necessary intestinal microbiota and maintain the diversity and function of intestinal microbiota [21,22].
In recent years, more and more research reports show that fecal eating behavior is of great signi cance to rabbit growth, development and immune regulation [23]. However, the microbiota and metabolites that play a key role in rabbit soft feces are still unknown. The rapid development of high-throughput sequencing technology has helped researchers better understand the composition of host intestinal microbiota and excavate the key microbiota and metabolites, which has played an important role in improving human and animal health [24,25]. Combined analysis of microbiome and metabolome can better nd the key regulatory metabolic pathways, nd out the key microbiota, and explain the molecular mechanisms of biological growth and development, environmental response, physiological state and pathological response. Therefore, in this study, high-throughput sequencing technology and GC-MS were used to analyze the cecal contents of Hyplus rabbit, microbial diversity, SCFAs and cytokines of soft and hard feces [26,27], and then the effect of soft feces on cecal immune microenvironment was veri ed by fasting soft feces, in order to screen the key microbiota and its metabolites that play an important role in rabbit growth, development and immune regulation.

Experimental design, animals and management
In this study, 40 Hyplus rabbits with similar weight provided by Jiyuan sunshine rabbit Technology Co., Ltd. were used as experimental materials. All Hyplus rabbits were fed according to the feeding standard of Jiyuan sunshine rabbit Technology Co., Ltd. After 20 days of adaptive feeding, the soft feces and hard feces of Hyplus rabbits were collected respectively. Six rabbits were randomly selected to collect the cecal contents after euthanasia. The collected samples were divided into three groups: cecal contents group (CON), soft feces group (SF) and hard feces group (HF). All samples were stored in the refrigerator at -80℃ for a short time, and DNA was extracted as soon as possible. 6-8 bodies were randomly selected from each group for detection. Then, 10 remaining Hyplus rabbits were randomly divided into two groups, wearing collars for fasting soft feces (CP) test. These Hyplus rabbits are kept in cages alone and have free access to feed and water. The same rabbit species (same age), basic dietary composition, feeding method and environment were used. All rabbits are raised according to the appropriate guidelines for raising rabbits.
Extraction and quality assessment of total DNA from soft, hard feces and cecal contents and library construction The fecal genomic DNA Extraction Kit (Tian gen, Beijing, China) was used to extract the fecal genomic DNA from all samples. Meanwhile, Nanodrop was used to quantify DNA and the quality of DNA was detected by 1.2% agarose gel electrophoresis. Then, the V3 + V4 region of 16S rRNA was ampli ed by PCR, gel recovery and uorescence quanti cation. Subsequently, the sequencing library was prepared using TruSeq Nano DNA LT Library Prep Kit, and the library was inspected using Agilent Bioanalyzer.
Illumina Novaseq6000 high-throughput sequencing The quali ed library was quanti ed by the quant-iT PicoGreen dsDNA Assay Kit on the Promega quanti uor uorescence quantitative system. After gradient dilution of quali ed on-line sequencing libraries (index sequence cannot be repeated), they are mixed according to the required sequencing amount in corresponding proportion, and denatured into single strand by NaOH for on-line sequencing. PE250 sequencing was performed with Illumina novaseq6000 sequencer, and the corresponding reagent was novaseq 6000 S2 reagent kit v1.5 (300 cycles). Sequencing service were performed by Bioyigene Biotechnology Co., Ltd. (Wuhan, China).
Bioinformatics analysis of 16S rRNA gene sequencing data Firstly, the sequenced original data are subjected to quality ltering, denoising, splicing and de chimerism to ensure that the data can be used for subsequent analysis. 1) By making statistics on the ASV/OTU Table after attening, the speci c composition table of microbial community in each sample at each classi cation level can be obtained, and the column diagram of phylum and genus level can be drawn by QIIME2 software) ggtree is used as an evolutionary tree to show the position of each ASV / OTU in the evolutionary tree and the evolutionary distance between them, and re ect their composition, abundance, taxonomy and other information through heat map and histogram. 3) In this study, Chao1 and observed species indexes are used to characterize richness, Shannon and Simpson indexes are used to characterize diversity, faith's PD index is used to characterize evolution based diversity, Pielou's evenness index is used to characterize evenness, and good's coverage index is used to characterize coverage. On the basis of alpha diversity, QIIME2 is used to draw rarefaction curve. Reference for calculation method of alpha diversity index (http://scikit-bio.org/docs/latest/generated/skbio.diversity.alpha.html#module-Skbio.diversity.alpha). 4) Using the uclust function of stat package in R language, UPGMA algorithm (i.e. average clustering method) is used for clustering analysis of Bray Curtis distance matrix by default, and the R script ggtree package is used for visualization. 5) The LEfSe package in Python software is used to analyze the histogram of LDA value distribution of signi cantly different species to show the signi cantly enriched species in each group (note that there is no signi cant down-regulation) and their importance; Using VennDiagram package in R software and make Wayne diagram according to ASV/OTU abundance table, and count the number of members of each set according to their presence or absence in each sample (Group), that is, the number of unique ASV/OTU in each group and the number of common ASV/OTU among groups (note that it is not abundance value); The principal component coordinate scores of each sample and each taxon were calculated by R script, combined with the contents of cytokines and short chain fatty acids, and presented in the form of interaction diagram (RDA Figure). (RDA analysis is conducted by using Shanghai Personal gene cloud platform, https://www.genescloud.cn/chart/) Functional potential prediction PICRUST2 (Phylogenetic Investigation of Communities by Reconstruction of Unobserved States) is a software for predicting the functional abundance of samples based on the abundance of marker gene sequences in samples (Gavin M. Douglas, et al., preprint). Its important function is to count the metabolic pathway of microbiota, analyze the difference of metabolic pathway and the species composition of metabolic pathway. The main analysis process is as follows; 1) Firstly, the 16S rRNA gene sequence of the known microbial genome was aligned to construct a new evolutionary tree. 2) Using castor hidden state prediction algorithm, according to the copy number of gene family corresponding to the reference sequence in the evolutionary tree, infer the nearest sequence species of the characteristic sequence, and then obtain the copy number of gene family. Note that when calculating the nearest sequence species index (NTSI) of each sequence, if the NTSI of the sequence is > 2 by default, it will be excluded from subsequent analysis. 3) Combined with the abundance of characteristic sequences of each sample, the copy number of gene family of each sample was calculated. 4) Finally, the gene family is "mapped" to various databases, and MinPath is used to infer the existence of metabolic pathways, so as to obtain the abundance data of metabolic pathways in each sample. 5) After obtaining the abundance data of metabolic pathways, Metagenomeseq package in R software was used to nd out the metabolic pathways with signi cant differences between groups. 6) The functional composition of samples/groups was obtained, and the species composition of pathways was analyzed by using the hierarchical sample metabolic pathway abundance table (Supplementary Table S2).

Correlation analysis of cytokines, SCFAs and microbiota
In order to analyze the correlation (Pearson correlation) between the microbiota and cytokines and SCFAs in the soft stool group, Use OmicShare Tools (Pearson correlation analysis was performed using the OmicShare tools, a free online platform for data analysis (https: //www.omicshare.com/tools)) online analysis software to analyze the correlation between the rst 15 genera and cytokines and SCFAs, respectively in the Con group, HF group and SF group, and select the cytokines There are TNF-α, IL-1β, IL-6, TGF-β, IL-4, IL-10 and IL-13; short-chain fatty acids include Formic acid, Isocaproic acid, Heptanoic acid, Propanoic acid, Isobutyric acid, Butyric acid, Valeric acid, Acetic acid, Isovaleric acid and Hexanoic acid.
Fasting soft feces to verify cecal immune microenvironment ELISA kits were used to detect the contents of TNF-α, IL-1β, IL-4 and IL-10 in the cecum content of CN and CP groups. Subsequently, the cecal tissues in the CN and CP groups were subjected to multiple immuno uorescence staining to detect the protein expression levels of Claudin, Occludin and ZO-1. The preparation of para n sections refers to Miao et al. [46]. After the sections are prepared, they are placed in a pathology section scanner (Hungary, 3DHISTECH) to collect images: DAPI UV excitation wavelength 330-380nm, emission wavelength 420nm, blue light; FITC excitation wavelength 465-495nm, emission wavelength 515-555 nm, emitting green light; CY3 excitation wavelength 510-560, emission wavelength 590nm, emitting red light. The excitation wavelength of CY5 is 608-648nm, and the emission wavelength is 672-712. CY5 is originally a positive red color. In order to distinguish it from CY3, it is set to pink light. The nucleus stained by DAPI is blue under ultraviolet excitation, and the positive expression is red, pink and green light labeled with uorescein.

Statistical analysis
A completely randomized trial design was used in this study. All data were analyzed by independent sample t-test in SPSS 24.0, * indicates signi cant difference, P-value < 0.05; ** Indicates a very signi cant difference, P-value < 0.01; NS indicates that there is no signi cant difference between the data, i.e. Pvalue > 0.05.

Results
Detection results of cytokines in soft, hard feces and cecal contents The levels of cytokines TNF-α, IL-1β, IL-6, TGF-β, IL-4, IL-10 and IL-13 in cecum contents, soft feces and hard feces were measured using ELISA kits, respectively, and the results are shown in Fig. 1. TNF-α, IL-1β, IL-6, TGF-β, IL-4, IL-10 and IL-13 were signi cantly higher in both soft and hard feces groups than in the cecum contents group; IL-1β, TGF-β, IL-4 and IL-13 contents were higher in the hard feces group than in the soft feces group, and the data were statistically signi cant between the two groups; however, TNF-α, IL-6 and IL-10 levels were higher in the soft feces group than in the hard feces group, but the difference between them was not signi cant (P-value>0.05) (Fig. 1A-G). In addition, three-dimensional PCA plots (Fig. 1H) were normalized according to the cytokines content of each sample and plotted by R software (v3.6.0), and principal component analysis (PCA) was used to perform a natural decomposition of the community data structure and to rank (sort) the samples to observe the differences between samples. The results showed that there was less dispersion between samples within groups in the cecum contents group, and the dispersion between samples within groups in the soft feces and hard feces groups was greater than that in the cecum contents group; while the samples between the cecum contents group, soft feces group and hard feces group were each independently distributed, with the samples in the soft and hard feces groups being closer together and farther away from the samples in the cecum contents group.

Detection results of SCFAs in soft, hard feces and cecal contents
The GC/MS method was used to determine the content of volatile short-chain fatty acids in the three samples, and the results are shown in Fig. 2. The contents of isobutyric acid, butyric acid, valeric acid and isohexanoic acid in the cecal content group were signi cantly higher than those in the soft stool group and hard stool group (P-value <0.05). The contents of acetic acid, isobutyric acid, butyric acid, isovaleric acid, valeric acid and caproic acid in the soft feces group were signi cantly higher than those in the hard feces group (P-value <0.05). The content of formic acid in the soft feces group and the hard feces group was signi cantly higher than that of the cecal content group (P-value <0.05), and the content of the hard feces group was higher than that of the soft feces group, the difference was statistically signi cant (Pvalue <0.05). There was no signi cant difference in the content of enanthic acid in the cecal contents group, soft stool group and hard stool group (P-value >0.05).

Species composition among groups soft, hard feces and cecal contents
Using QIIME software, the composition and abundance distribution of each sample at two levels of phylum and genera were obtained, and the analysis results were expressed in the form of a histogram. At the phylum classi cation level (Fig. 3A), Firmicutes and Bacteroides are the main phyla of the three groups of samples. The relative abundance of Firmicutes (84.48%) in the soft feces group is signi cantly higher than that in the hard feces group. (69.66%) and cecal contents group (56.82%), but the relative abundance of Bacteroides in the soft stool group (8.30%) was signi cantly lower than that in the hard stool group (26.67%) and cecal contents group (38.31%). At the genus level (Fig. 3B), the relative abundance of Muribculaceae in the cecal content group (31.17%) and the hard feces group (23.74%) was signi cantly higher than that of the soft feces group (3.78%); Ruminococcaceae_NK4A214_group in the soft feces group, The relative abundances of Lachnospiraceae_NK4A136_group, Christensenellaceae_R-7_group and Subdoligranulum (21.73%, 8.18%, 8.13%, and 7.93%, respectively) were signi cantly higher than those of the cecal contents group (5.71%, 2.30%, 1.39%, and 1.56%, respectively) and hard Fecal group (5.07%, 4.76%, 2.92% and 1.34%, respectively); the relative abundance of Ruminococcaceae_UCG-014 (12.84%) in the hard feces group was signi cantly higher than that of the cecal contents group (2.80%) and soft feces group (3.39) %). Then use ggtree in R software as evolutionary tree to display the position of each ASV/OTU in the evolutionary tree, and the evolutionary distance between each other, and re ect their composition, abundance, taxonomy and other information through heat maps and histograms. The result is shown in Fig. 3C. Among them, Ruminococcaceae_NK4A214_group, Ruminococcaceae_UCG-014, Ruminococcaceae_1, Ruminiclostridium_6, Muribaculaceae, Clostridiales_vadinBB60_group, Christensenellaceae_R-7_group and other bacterial genera are closely related. Ruminococcaceae_UCG-013 and Subdoligranulum are closely related, but farther than the other genus mentioned above. It can be seen from the heat map that Ruminococcaceae_NK4A214_group and Christensenellaceae_R-7_group have higher relative abundance in the soft feces group. The relative abundance of bacteria such as Ruminiclostridium_6 and Ruminococcaceae_UCG-014 in the HF group was higher; the relative abundance of bacteria such as Muribaculaceae and Ruminococcaceae_1 in the cecal content group was higher. In addition, the abundance histogram results show that Firmicutes and Bacteroides are the main two phyla.
Alpha and beta diversity in microbiota of groups soft, hard feces and cecal contents Alpha diversity and beta diversity index were used to characterize the diversity of species within and between habitats respectively, so as to comprehensively evaluate their overall diversity. In this study, Chao1 and observed species indexes are used to represent richness, Shannon and Simpson indexes are used to represent diversity, faith's PD index is used to represent evolution based diversity, Pielou's evenness index is used to represent evenness, and good's coverage index is used to represent coverage. Alpha diversity results are shown in Fig. 4A. The indexes of Chao1, faith's PD, Shannon, Pielou's evenness and observed species in cecal content group are signi cantly higher than those in soft feces group and hard feces group. There was no signi cant difference in all the above indexes between soft feces group and hard feces group (P-value > 0.05). The sparse curve results of good's coverage are shown in Fig. 4B. The smoothness of the curves of the three groups of samples re ects the impact of sequencing depth on the diversity of observed samples. The curves of the three groups of samples tend to be gentle, indicating that the sequencing results have enough to re ect the diversity contained in the current samples. If the sequencing depth continues to increase, a large number of new ASV/OTUs that have not been found. Then, the dimension of multidimensional microbial data is reduced by Principal Coordinate Analysis (PCoA), and the main trend of data change is displayed by the distribution of samples on the continuous sorting axis. At the same time, permutational multivariate analysis of variance (PERMANOVA) was used to identify whether there were signi cant differences among the three groups of samples. The beta diversity of the three groups of samples was analyzed by the above methods. The results showed that in the PCoA diagram (Fig. 4C), there was an obvious separation between the samples of cecal content group, soft feces group and hard feces group, and the dispersion between the samples of cecal content group was less than that between the samples of soft feces group and hard feces group. The results of difference analysis between groups are shown in Fig. 4C. There is a very signi cant difference between cecal contents group and soft and hard feces group (P-value < 0.05), but there is no signi cant difference between soft and hard feces groups (P-value > 0.05).

Species differences and marker analysis
The results of LEfSe analysis showed that p_Bacteroidetes, c_Bacteroidia, o_Bacteroidales, f_Muribaculaceae and g_Muribaculaceae were the most signi cant bacterial genera in the Con group. g_Ruminococcaceae_UCG-014, g_Lachnospiraceae_UCG_001 and g_Izimaplasmatales are the most signi cant bacterial genera in the HF group. p_Firmicutes, c_Clostridia and o_Clostridiales are the most signi cant differences in the SF group (Fig. 5A). Use the ASV/OTU abundance table (Supplementary  Table S1) to make a Venn diagram, the Con group, the soft feces group and the HF group. The three groups of samples have 13765, 5442 and 7344 ASV/OTUs respectively. Three groups of samples There are a total of 623 ASV/OUT (Fig. 5B). RDA analysis results show that SCFAs are mainly related to the Con group and SF group, and Cytokines were mainly correlated with fecal microbiome (HF and SF groups) (Fig. 5C). The relative abundance of the top 15 bacterial genera in the three groups of samples was drawn using software as shown in Fig. 5D. The relative abundance of Muribaculaceae and Ruminococcus_1 in the Muribaculaceae Con group and the HF group was higher than that of the SF group; Ruminococcaceae_NK4A214_group in the SF group The relative abundances of, Christensenellaceae_R-7_group, Subdoligranulum and Candidatus_Saccharimonas were higher than the Con group and the HF group.

Functional potential predictive analysis
The relative abundance results of primary and secondary metabolic pathways in metacyc database are shown in Fig. 6A. The relative abundance of pathways related to biosynthesis is high, mainly including amino acid biosynthesis, carbohydrate biosynthesis, coenzyme factor and vitamin biosynthesis. The relative abundances of enzymes related to the synthesis of acetic acid, propionic acid and butyrate in the three groups were further analyzed. It was found that the relative abundances of acetic acid synthase such as EC: 6.3.4.3, EC: 1.5.1.20, EC: 2.3.1.54 and EC: 2.3.1.169 in the soft feces group were signi cantly higher than those in the Con group (Fig. 6B). In propionate synthase, the relative abundance of EC: 2.7.1.1 enzyme in soft feces group was signi cantly higher than that in Con group (Fig. 6C). Among butyrate synthase, the relative abundances of EC: 2.3.1.9, EC: 4.2.1.17, EC: 1.3.8.1 and EC: 2.8.3.8 in soft feces group were signi cantly higher than those in Con group (Fig. 6D). After obtaining the abundance data of metabolic pathways, Metagenome-seq was used to nd out the metabolic pathways with signi cant differences among the three groups. As shown in Fig. 6E and F, the up-regulated metabolic pathways in the hard fecal group were PWY-7198 and PWY-7210 compared with the cecal content group, and the upregulated metabolic pathways in the soft fecal group were PWY-7198, PWY -7210, p341-PWY and P123-PWY compared with the cecal content group, that is, P341-PWY and P123-PWY were up-regulated in soft feces group compared with hard feces group. Finally, according to the signi cantly different metabolic pathways, the hierarchical sample metabolic pathway abundance table (Supplementary Table S2) is used to analyze the species composition of different pathways. The results are shown in Fig. 6G and H. The P125-PWY pathway is superpathway of (R, R)-butanediol biosynthesis, in which the relative abundance of Christensenellaceae_R-7_group in the soft stool group is higher than that of the cecal content group and the hard stool group; the P341-PWY pathway is glycolysis V (Pyrococcus), of which the soft stool group and the hard stool group The relative abundance of Lachnoclostridium in the group was higher than that of the cecal content group.

Correlation analysis of cytokines, SCFAs and microbiota
The correlation analysis results of SCFAs and the microbiota are shown in Fig. 7A. Among them, Ruminococcaceae_NK4A214_group has signi cant positive correlation with acetic acid, butyric acid, isobutyric acid, isovaleric acid and hexanoic acid; Christensenellaceae_R-7_group has positive correlation with acetic acid, butyric acid, isobutyric acid, isovaleric acid and hexanoic acid; Muribaculaceae has a signi cant positive correlation with isocaproic acid. Ruminococcaceae_UCG-014 has a signi cant negative correlation with acetic acid, propionic acid, isobutyric acid, butyric acid, valeric acid and caproic acid. The results of correlation analysis between cytokines and The top 15 relative abundance bacteria in the three groups of samples are shown in Fig. 7B. Muribaculaceae has a signi cant positive correlation with IL-1β, IL-4 and TGF-β, respectively; On the contrary, Candidatus_Saccharimonas has a signi cant negative correlation with IL-1β, IL-4 and TGF-β.

Effect of fasting soft feces on cecal immune microenvironment
After fasting soft feces, the cytokines TNF-α and IL-1β Signi cantly increased, while cytokines IL-4 and IL-10 decreased signi cantly (Fig. 8A). In addition, the expression levels of tight junction related proteins in CN group and CP group were measured by Western blotting and multiple immuno uorescence staining respectively. The results of Western blotting showed that the expression levels of tight junction proteins such as ZO-1, Claudin-1 and Occludin increased signi cantly after fasting soft feces (Fig. 8B). The results of immuno uorescence staining were consistent with those of Western blotting. The staining degree of ZO-1 protein (purple), Claudin-1 protein (green) and Occludin (red) protein in CP Group was deeper (compared with CN group), indicating that the contents of Claudin-1 and Occludin protein in CP group were higher (Fig. 8C). Meanwhile, after fasting soft feces, the cecal mucosa was damaged (Fig. 8D, the purple part indicated by the arrow), while the integrity of cecal mucosa in the CN group was better.

Discussion
Because rabbits have a special colon separation mechanism, they can produce two different forms of feces (soft feces and hard feces) [28,29]. There are a lot of nutrients in soft feces, and rabbits have the behavior of eating soft feces from the beginning. This behavior is equivalent to prolonging the retention time of food in the digestive tract and improving the utilization rate of nutrients in food [30][31][32]. It is found that the chemical composition of soft feces is similar to that of cecal contents. There are a large number of microorganisms in both cecal contents and soft feces. These microorganisms can a variety of enzymes and short chain fatty acids, and play an important role in improving rabbit growth performance and immune function [5,33]. However, the microbiota or related metabolites that play a key role are unknown. Therefore, in this study, the combination of microbiome and metabolome to analyze the cecal contents of Hyplus rabbit, microbial diversity of soft and hard feces, SCFAs and cytokines can better help to nd the key regulatory metabolic pathways, key ora and their metabolites, and play an important role in clarifying the molecular mechanism of soft feces in rabbit growth, development and immune regulation.
The higher the value of Firmicutes / Bacteroides, it may be related to obesity. On the contrary,it may be related to emaciation [34][35][36]. Regardless of fat and thin, Firmicutes may play an important role in the intestine, more effectively absorb nutrients in food and promote the growth and development of the host.
Further analysis found that the relative abundance of Muribaculaceae in the Con group and the HF group is higher than that of the SF group at the level of genus classi cation. The relative abundance of this genus has a very important impact on host development and health. The increase of Akkermansiaceae and the decrease of Muribaculaceae can promote the healthy development of bones [37]. Additionally, Akkermansiaceae, especially Akkermansia muciniphila, plays a protective role in diet induced obesity and other diseases [38]. What's more, the relative abundance of speci c genus Ruminococcaceae_NK4A214_group was reduced in HF group, and it was signi cantly correlated with bile acid and Vitamin A levels, which plays a vital role in the treatment of rats with ulcerative colitis [39]. Christensenellaceae_R-7_group as a probiotic is signi cantly negatively correlated with metabolic diseases such as BMI and in ammation, such as IBD and metabolic syndrome [40]. Through association analysis (microbiome & metabolomics analysis), it was con rmed that Muribaculaceae, Christensenellaceae_R-7_group and Ruminococcaceae_NK4A214_group have a strong positive correlation with SCFAs and IL-10, and they have a strong negative correlation with TNF-α, IL-1β and TGFβ. Candidatus_Saccharimonas have a strong negative correlation with IL-1β, TGF-β and IL-4. These results show that there are abundant genes (microorganisms) that can synthesize SCFAs (acetic acid, propionic acid and butyric acid) in soft feces, which may be the key to promote growth and development and resist external environmental interference [41,42]. These may also explain why the host can improve the digestion and absorption e ciency of feed, maintain the normal microbiota of digestive tract and alleviate nutritional diseases after eating soft feces.
After fasting soft feces, the content of IL-4 and IL-10 in cecal contents decreased and the content of TNFα and IL-1β increased, which further veri ed that some microorganisms in soft feces can increase the expression level of cytokines associated with anti-in ammation and reduce the relative expression level of cytokines associated with in ammatory. Intestinal mucosa is a physical barrier against invading intraluminal pathogens and toxins [43,44]. The results of immuno uorescence staining further con rmed that fasting soft feces caused damage to intestinal mucosa, indicating the importance of soft feces in maintaining intestinal immune microenvironment. Nevertheless, western blot and immuno uorescence staining showed that the expression of tight junction associated proteins increased signi cantly after fasting soft feces. The tight junctional (TJ) mainly maintain the integrity of the epithelial barrier by regulating the paracellular permeability [45]. In the process from never fasting soft feces to fasting soft feces, the body suddenly interrupted the intake of probiotics, SCFAs and cytokines from the external environment. In order to deal with this sudden situation (Similar to stress), the body improved the integrity of intestinal barrier by increasing the expression of tight junction protein in the intestine, so as to avoid the impact of the external environment on intestinal homeostasis and even body health.

Conclusion
In this study, the soft feces are rich in bene cial microorganisms (Muribaculaceae, Ruminococcaceae_NK4A214_group and Christensenellaceae_R-7_group) and their metabolites (SCFAs and cytokines). Feeding soft feces plays a very important role in improving the immune microenvironment of cecum. This study may provide a valuable reference for the treatment of    Content of SCFAs in cecal contents and soft and hard feces. A-J. The contents of short chain fatty acids such as formic acid, isocaproic acid, heptanoic acid, propanoic acid, isobutyric acid, butyric acid, valeric acid, acetic acid, isovaleric acid and hexanoic acid in cecal contents and soft and hard feces, respectively; K. The relative abundance of short chain fatty acids such as formic acid, isocaproic acid, heptanoic acid, propanoic acid, isobutyric acid, butyric acid, valeric acid, acetic acid, isovaleric acid and hexanoic acid in cecal contents and soft and hard feces in the form of heat map.    Functional potential prediction of microbiota in cecal contents and soft and hard feces. A. The analysis results of metacyc metabolic pathway including various pathways involved in primary and secondary metabolism and related metabolites, biochemical reactions, enzymes and genes; B. The relative abundance of enzymes related to acetic acid synthase in cecal contents group and soft and hard feces group; C. The relative abundance of propionate synthase related enzymes in cecal contents group and soft and hard feces group; D. The relative abundance of butyrate synthase related enzymes in cecal contents group and soft and hard feces group; E. Different metabolic pathway between cecal content group and hard fecal group; F. Different metabolic pathway between cecal contents group and soft feces group; G and H. The species composition of up-regulated metabolic pathways such as P341-PWY and P123-PWY in soft feces group compared with hard feces group.