Multi-omics analysis of miRNA-mediated intestinal microflora changes in crucian carp Carassius auratus infected with Rahnella aquatilis

Infection by an emerging bacterial pathogen Rahnella aquatilis caused enteritis and septicemia in fish. However, the molecular pathogenesis of enteritis induced by R. aquatilis infection and its interacting mechanism of the intestinal microflora associated with microRNA (miRNA) immune regulation in crucian carp Carassius auratus are still unclear. In this study, C. auratus intraperitoneally injected with R. aquatilis KCL-5 was used as an experimental animal model, and the intestinal pathological changes, microflora, and differentially expressed miRNAs (DEMs) were investigated by multi-omics analysis. The significant changes in histopathological features, apoptotic cells, and enzyme activities (e.g., lysozyme (LYS), alkaline phosphatase (AKP), alanine aminotransferase (ALT), aspartate transaminase (AST), and glutathione peroxidase (GSH-Px)) in the intestine were examined after infection. Diversity and composition analysis of the intestinal microflora clearly demonstrated four dominant bacteria: Proteobacteria, Fusobacteria, Bacteroidetes, and Firmicutes. A total of 87 DEMs were significantly screened, and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses revealed that the potential target genes were mainly involved in the regulation of lipid, glutathione, cytosine, and purine metabolism, which participated in the local immune response through the intestinal immune network for IgA production, lysosome, and Toll-like receptor (TLR) pathways. Moreover, the expression levels of 11 target genes (e.g., TLR3, MyD88, NF-κB, TGF-β, TNF-α, MHC II, IL-22, LysC, F2, F5, and C3) related to inflammation and immunity were verified by qRT-PCR detection. The correlation analysis indicated that the abundance of intestinal Firmicutes and Proteobacteria was significantly associated with the high local expression of miR-203/NF-κB, miR-129/TNF-α, and miR-205/TGF-β. These findings will help to elucidate the molecular regulation mechanism of the intestinal microflora, inflammation, and immune response-mediated miRNA–target gene axis in cyprinid fish.


Introduction
The intestine and its normal microbiota play an important role in the physiological and mucosal immune processes in fish (1,2).The changes in gut microbiota influence the occurrence and development of diseases such as bacterial enteritis (3)(4)(5).The gut microbial community is involved in numerous other functions, such as protecting the host from pathogens, modulating immunity, and regulating metabolic processes (6).The intestinal microflora is of functional importance to all aspects of host physiology, and regulation of the microbiota is a feasible strategy to alleviate emerging diseases in teleosts (7).In mammals, the normal development and behavior of the immune system are strongly affected by gut microbial metabolites (8).There is mounting evidence that the intestinal microflora of fish plays an important role in growth and development, nutrition absorption, immunity, and disease resistance (9)(10)(11).Bacterial pathogens from the intestine can migrate through the epithelial mucosa to infect internal organs and injure the host.In previous studies, crucian carp infected with pathogenic Salmonella and Vibrio anguillarum were found to have damaged intestinal cells (12).In fish, the disturbance of the intestinal microflora can lead to impaired physiological function, such as dysfunction of intestinal epithelial cells and weakening of nutrient absorption, metabolism, and immune response (13).Therefore, the intestinal microflora of fish is essential for bacterial enteritis and involved in the regulation of immune function, which is of great significance for preventing diseases and realizing the healthy growth of teleosts.However, there are few reports on the interaction between the intestinal microflora, immunity, and disease in fish, and the characteristics of the ecological structure of the intestinal microflora and bacterial diversity caused by Rahnella aquatilis have not been reported in crucian carp.
MicroRNAs (miRNAs) are small endogenous non-coding RNAs with a length of approximately 22 nucleotides.It guides the degradation of target genes by binding to the 3′ untranslated region (3′-UTR) or inhibiting translation to negatively regulate gene expression (14,15).Accumulating evidence demonstrated that miRNAs participate in the regulation of cell proliferation, differentiation, apoptosis, and immune response (16,17).In higher vertebrates, a subset of miRNAs has been identified as an important regulator of numerous key genes in the immune system gene network.Mature miRNAs and their corresponding genes were characterized in a number of teleost species (18).As an important molecular regulator of fish, the immune system protects the host from infections by pathogens.MiRNAs regulate systemic and mucosal immunity by targeting pattern recognition receptors and downstream signaling factors (19,20).Several kinds of miRNAs (e.g., miR-21b-3p, miR-146a-3p, and miR-155-5p) were found to effectively inhibit harmful inflammation or promote early immune response in Atlantic salmon Salmo salar (21).Research reported that the changes in Toll-like receptor (TLR) 2 in fish may be related to miR-155 and its target SOCS1.The TLR2/MyD88 signaling pathway was partially involved in the influence of miR-155 on the expression of immunocytokines (22).In miiuy croaker Miichthys miiuy, miR-203 and miR-148-1-5p were involved as negative regulatory factors in the innate immune response (23,24).In addition, miR-205-5p was an important regulatory factor in genetically improved farmed tilapia (Oreochromis niloticus), controlling lipid metabolism (25).Very recently, the integration of miRNA and mRNA analyses revealed the potential roles of the identified known miRNAs that are involved in the regulation of glycometabolism-related pathways in common carp Cyprinus carpio (26).Previous studies revealed that mucosal immune responses in crucian carp include the TLR, NOD-like receptor (NLR), and complement and coagulation pathways (20,27,28).However, there are few transcriptome analyses of miRNAs (e.g., miR-17, miR-26a, miR-144, and miR-146a) related to the involvement in the mucosal immune response of miRNA and the intestinal microflora regulation function in crucian carp (20).
Recently, intensive culture of commercial fish such as crucian carp has increased the severity and frequency of diseases in China (27,(29)(30)(31)(32). Crucian carp have fast growth rates, high nutritional quality, and good palatability and have become an important freshwater economically farmed fish in China (33, 34).In recent years, crucian carp were infected with bacteria and caused diseases more frequently due to the high-density culture (27,(30)(31)(32).More recently, R. aquatilis has been recognized as an emerging bacterial pathogen, and more and more infections were detected in fish, which have demonstrated that R. aquatilis caused clinical enteritis and septicemia (30-32, 35, 36), and frequent outbreaks of its infections have caused severe economic losses to fish farmers in China (27, 30-32, 36, 37).The innate immune system involved different signaling pathways that were regulated by complex mechanisms, of which the mucosal immunity (i.e., intestine, gills, and skin) was the host's first line of defense against pathogen infection (38-41).However, the molecular mechanism of the miRNAs in the local intestinal microflora and mucosal immune response is still unclear in fish.
Bacterial diseases, especially bacterial enteritis, are a major problem in crucian carp Carassius auratus, which has caused large economic losses (3).In previous studies by our group, a new pathogenic R. aquatilis KCL-5 strain was isolated from diseased crucian carp (27,(30)(31)(32).Subsequently, we detailedly studied the mucosal immunity (skin and gill) response of bacterial septicemia of crucian carp through multi-omics analysis (20,28,41).However, the molecular mechanism of the intestinal microflora, inflammation, and immune response-mediated miRNA-target gene axis in the mucosal immune regulation infected with R. aquatilis is still unclear.In the present study, bacterial 16S rRNA high-throughput and small RNA (sRNA) high-throughput transcriptome sequencing techniques were used for the first time, and the molecular regulation mechanism of the intestinal microflora and miRNAs in intestinal mucosal immunity of crucian carp infected with R. aquatilis was first revealed by us, which lays a foundation for further study of systemic immunity in cyprinid fish.

Sample and challenge
The pathogenic strain R. aquatilis strain KCL-5 was previously isolated from diseased crucian carp C. auratus according to our previous report (30).All the crucian carp (average weight 250 ± 10 g, average length 21 ± 1 cm) used in this experiment were purchased from Tianjin Aquatic Science and Technology Development Co., Ltd.(Tianjin, China) and acclimatized at a temperature of approximately 24°C ± 0.5°C in freshwater tanks.In the infection experiment, the fish were intraperitoneally injected with R. aquatilis with a final concentration of 6 × 10 8 CFU/mL based on the concentration of 50% lethal dose (LD 50 ) by our group report (30).The fish were randomly selected, and the intestines were sampled at 0, 36, and 72 h post-infection (hpi), according to the clinical signs and our laboratory findings (30,32), such as skin bleeding, ascites, and enteritis; the pathological changes were observed during early infection.The samples were mixed, with three replicates per group, and then stored in a −80°C refrigerator for later experiments.

Pathological observation and enzyme activity changes in intestine
The crucian carp were randomly selected and dissected for histopathological examination.After deep anesthesia using MS-222 (200 mg/L), the intestine samples were taken from fish, treated with 0.8% sterile saline, and fixed with 4% paraformaldehyde.The embedded paraffin tissue was placed on a tissue slicer for continuous slicing (thickness 5 mm) and thereafter stained with hematoxylin and eosin (H&E) and Alcian blue-Periodic acid Schiff (AB-PAS).Detection and analysis of TUNEL apoptosis were performed according to the recent report by Huo et al. (27).Alkaline phosphatase staining was tested according to the instructions.The positive result of alkaline phosphatase staining manifested as the appearance of gray-brown to dark-black granular or sheet-like precipitates in the cytoplasm.Sterile physiological saline was added to the intestinal samples, which were then homogenized in an ice water bath and centrifuged at 2,500 rpm for 10 minutes.The supernatant was taken and diluted into 2% tissue homogenate.The enzyme activities of lysozyme (LYS), aspartate transaminase (AST), alanine aminotransferase (ALT), alkaline phosphatase (AKP), and glutathione peroxidase (GSH-Px) were then detected according to the instructions of the kit (Nanjing Jiancheng, Nanjing, China).Each experiment was repeated three times.

Genomic DNA extraction and highthroughput sequencing
The total DNA of the intestinal content samples was extracted from the fish according to a previous method, with slight modifications (32).The universal primers were designed in the conserved region, and the necessary connectors were added for sequencing at the end of the primers for PCR amplification.Purification, quantification, and homogenization of PCR products were performed to construct three control and infected libraries by 16S rRNA high-throughput sequencing.Further quality inspection of the library was carried out, and sequencing analysis was conducted based on the Illumina HiSeq 2500 platform to obtain the original image data files, which were then converted into the original sequenced sequences through base recognition analysis.

Data statistics and analysis of intestinal microflora
Statistical analysis of the intestinal microflora was used to filter the original data using Trimmomatic (V0.33) software, and partial low-quality cropping was performed on all reads obtained from sequencing using Cutadapt (V1.9.1).The data of each sample were split from the obtained reads, the barcode and primer sequences were truncated, and the data were obtained through preliminary quality control.Subsequently, USEARCH (V10) was used to splice the dual-ended reads, and chimeras were removed using UCHIME (V8.1), ultimately obtaining high-quality sequences for subsequent analysis.
At a similarity level of 97%, Usearch was used to cluster the effective reads of multiple intestinal content samples to obtain operational taxonomic units (OTUs).Using a naive Bayesian classifier with the SILVA as a reference database, taxonomic annotation was performed on sequences with features to obtain species classification information for each feature.Further, QIIME software was used to generate species abundance at various taxonomic levels (i.e., phylum, class, order, family, genus, and species), and R language was used to draw species composition structure maps at different taxonomic levels.The sample with the least amount of data was the standard for homogenization processing.Alpha-diversity index (i.e., Chao1, Ace, Shannon, and Simpson) and beta diversity of homogenized sample data were evaluated using PAST, Mothur, and QIIME software.

KEGG function prediction of intestinal microflora
The 16S rRNA high-throughput sequencing data with Tax4Fun2 software was used to classify species based on the SILVA database using the QIIME or SILVAngs platform.The 16S rRNA copy number was standardized based on the classification results and the genome annotation of the National Center for Biotechnology Information (NCBI).Finally, a linear relationship between SILVA classification and prokaryotic classification was constructed in the Kyoto Encyclopedia of Genes and Genomes (KEGG) database, thereby achieving KEGG functional prediction of prokaryotic microbial communities.

Extraction of total RNA and construction of sRNA library
For qualified intestinal samples from the crucian carp, total RNA was used as the starting sample, and a library was constructed using a sRNA Sample Pre Kit.T4 RNA Ligase 1 and T4 RNA Ligase 2 were used to connect connectors at the 5′ and 3′ ends of sRNA and then reversed transcribed to synthesize cDNA for PCR amplification.The target fragment was screened via polyacrylamide gel electrophoresis (PAGE) gel electrophoresis, the gel was cut for gel recovery, and the resulting fragment was the sRNA library.After the library construction was completed, Qubit 2.0 and Agilent 2100 bioanalyzer were used to detect the concentration of the library and sequence insert size.In addition, the qRT-PCR method was used to accurately quantify the effective concentration of the library, and after passing the quality inspection, the NovaSeq 6000 platform was used for high-throughput sequencing analysis.

Identification and predictive analysis of miRNAs
According to the recent report by Bai et al. (42), the reads from the reference genome of crucian carp were compared with known miRNAs in the miRBase (V22) database, where the identified reads were considered known miRNAs.Based on the biological characteristics of miRNAs, miRDeep2 software was used for the unidentified miRNA sequences, and the possible precursor sequences were obtained by comparing the position information on the genome through reads.Based on the distribution information of reads on the precursor sequence and the energy information of the precursor structure (RNAfold randfold), the Bayesian model was used to score and finally predict novel miRNAs.In addition, typical miRNA base ratios were obtained by analyzing the base preference of miRNAs.

Prediction and functional enrichment analysis of DEM target genes
As previously reported (20), the differential expression level of the miRNAs was analyzed.The candidate gene functional enrichment of the differentially expressed miRNAs (DEMs) was analyzed using Gene Ontology (GO) and KEGG.The RNAhybrid method was perfected to predict the binding site of key miRNAs to target genes, and SWISS-MODEL and PyMOL online software were used to generate the tertiary structure of the protein analysis.The miRBase and MEGA 7.0 software were used to select sequences and multiple sequence alignment of pre-miR-205, and maximum likelihood estimation (MLE) was utilized for the construction of pre-miR-205 phylogenetic tree (bootstrap = 1,000) (Supplementary Table 1).

RNA isolation and quantitative real-time PCR analysis
RNAiso Plus (Takara, Dalian, China) and RNAeasy ™ Animal RNA Isolation Kit with Spin Column (Beyotime, Shanghai, China) were used to extract total RNA from intestinal samples according to the instructions of the kits.The quality was detected using a Nanodrop ND-2000 ultra-micro spectrophotometer, and the purified total RNA was stored in a cryogenic refrigerator at −80°C.Then, the complementary DNAs (cDNAs) were obtained by reverse transcription using PrimeScript ™ RT reagent Kit with gDNA Eraser (Takara, Dalian, China), with each sample containing three repeats.In addition, miRanda and TargetScan software were used to predict the candidate target genes of key miRNAs (e.g., miR-129, miR-203, and miR-205).Eleven mucosal immune-related target genes (e.g., TGF-b, TNF-a, MyD88, TLR3, and NF-kB) were selected, and b-actin was used as the internal reference to analyze the changes in gene expression level after infection by qRT-PCR.Supplementary Table 2 lists the primer sequences used for the qRT-PCR analysis.

Data processing and analysis
The data were processed using SPSS 25.0 software and expressed as mean ± standard deviation.Single-factor analysis of variance was used to compare the significance of the data, with p < 0.05 as the significant difference.

Histopathological changes in the intestine of C. auratus infected with R. aquatilis
The crucian carp infected with R. aquatilis strain KCL-5 caused symptoms of enteritis and septicemia.Histopathological observation showed necrosis, degeneration and exfoliation of intestinal mucosal epithelial cells, severe injury and shortening of intestinal villi, bleeding of the submucosa, mucosal layer and muscular layer thinning, intestinal hemorrhage, necrosis, and vacuolation (Figure 1A); intestinal mucous cells increased and stained blue after AB-PAS staining, distributed with goblet cells (Figure 1B).Gomori calcium-cobalt staining showed that the positive reaction of alkaline phosphatase was brown-black granules, and the reaction was enhanced in the intestinal mucosa of the infected group (Figure 1C).In addition, the apoptotic cells were fluorescently demonstrated in the intestinal mucosal layer through the TUNEL fluorescence assay (Figure 1D).
The changes in enzyme activity, including ALT, AST, GSH-Px, AKP, and LYS, were determined in the intestines of infected crucian carp (Figure 1E).The activity of LYS decreased by 50% at 36 hpi after infection and then increased by 66.67% at 72 hpi.In contrast, the activities of AKP, ALT, AST, and GSH-Px increased at first and then decreased.It is suggested that five enzymes play a key role in regulating intestinal inflammation, which can effectively regulate the mechanism of intestinal immune response.

Evaluation of intestinal microflora in C. auratus
The intestinal contents of crucian carp were sequenced by 16S rRNA high-throughput sequencing, a total of 720,812 pairs of reads were obtained, and a total of 718,785 clean reads were produced after double-terminal read quality control and splicing (Supplementary Figure 1A).Based on the preliminary analysis of the sequencing results of intestinal bacteria, the dilution curve and species accumulation curve were obtained as shown in Supplementary Figures 1B, C. With the increase in the number of sequencing, the changes in the OTUs and species accumulation curve showed no significance and gradually tended to enter the sequencing stable phase at 30,000 reads.It turned out that the intestinal contents were sufficient and the quality of 16S rRNA sequencing was reliable.The results truly reflect the changes in the intestinal microflora diversity induced by intraperitoneal injection of R. aquatilis in crucian carp.

Diversity analysis of intestinal microflora in C. auratus
The changes in species richness and diversity of the intestinal microflora in crucian carp infected with R. aquatilis were investigated by alpha analysis.The results showed that the Chao1 and ACE index increased at first and then decreased, while the Shannon index and Simpson index decreased at first and then increased.Compared with that at 0 hpi, the Chao1 index of the microflora increased significantly at 72 hpi (p = 0.019), indicating that the species richness of the intestinal microflora increased after 72 hpi, but the ACE, Shannon, and Simpson indexes of the microflora did not significantly change (ACE, p = 0.064; Shannon, p = 0.25; Simpson, p = 0.19) (Figure 2A).Principal coordinates analysis (PCoA) and unweighted pair group method  with arithmetic mean (UPGMA) were used to further analyze the species beta diversity.The results showed that there were obvious groups at 0 and 72 hpi (Figures 2B, C and Supplementary Figures 2A, B), illustrating that R. aquatilis had an effect on the composition and structure of the intestinal microflora in crucian carp.
Component analysis in the intestinal microflora of C. auratus Due to the changes in the composition and structure of the intestinal microflora in crucian carp, clustering analysis was carried out to reflect the similarity among the three groups of intestinal samples.According to the Venn diagram (Figure 3A), the number of public OTUs was 634, of which the unique number of OTUs at 0 and 36 hpi was 1, while the unique number of OTUs at 72 hpi was 12, indicating that there were some important differences in OTUs among groups after R. aquatilis was injected into crucian carp at different infection time.According to the results of species annotation, the composition at the phylum classification level is shown in Figure 3B.The top 10 bacteria with the largest species abundance were selected, among which the relative abundance of phylum Firmicutes was the highest, followed by Proteobacteria.The relative abundance of Proteobacteria decreased significantly at 36 hpi after infection, while that of Fusobacteria increased significantly (p < 0.05).In addition, the largest species abundance of bacteria at the genus level is shown in Figure 3C.The relative abundance of the Cetobacterium increased significantly at 36 hpi.In contrast, the relative abundance of Escherichia-Shigella decreased continuously and especially decreased significantly at 72 hpi (p < 0.05).Linear discriminant analysis effect size (LEfSe) analysis was carried out to explore the key region taxa.The results showed that there were 13 groups with significant differences between 0 and 36 hpi, among which mainly included Fusobacteriales, Cetobacterium, and Bacteroidaceae at 36 hpi.However, at 72 hpi, only phylum Acetoanaerobium showed a significant difference, which indicated that the abundance of Proteobacteria, Actinobacteria, and Firmicutes bacteria in the crucian carp intestine was higher (Figure 3D).

KEGG functional and relevance analysis in intestinal microflora of C. auratus
KEGG function of the intestinal microflora was analyzed to explore the changes in metabolic pathways related to crucian carp infected with R. aquatilis.The results showed that 14 different metabolic-related pathways were annotated (Figure 4A).The main KEGG pathways enriched by functional genes at 36 hpi included glycolysis/gluconeogenesis, biosynthesis of amino acids, pyruvate metabolism, and purine metabolism, while the main enriched KEGG pathways at 72 hpi included biosynthesis of secondary metabolites and antibiotics.Furthermore, according to the abundance and changes in the intestinal microflora at the phylum level, the correlation analysis was carried out by drawing a coexpression analysis network map based on python by Spearman.The results showed that phylum Firmicutes was negatively correlated with Verrucomicrobia and that Proteobacteria was negatively correlated with Fusobacteria (Figure 4B), which were the same as the changing trend of species composition at the phylum level.

Identification and predictive analysis of miRNAs in the intestine
In order to identify the immune-related miRNA expression profile, six sRNA libraries of fishes in the control group (IC1, IC2, and IC3) and infected group (IT1, IT2, and IT3) were constructed based on the NovaSeq 6000 platform at 0 and 36 hpi.The results showed that a total of 86.84M clean reads was obtained, and each sample was more than 12.15M clean reads, with a base mass value of Q30 in each group greater than 95% (Supplementary Table 3).The DEMs were further screened and identified, and a total of 4,422 miRNAs were obtained, of which 2,873 were known miRNAs and 1,549 were novel miRNAs (Table 1).The length distribution of the identified known miRNAs and novel miRNAs was mainly in the range of 20 nt to 23 nt (Supplementary Figure 3A).In addition, the nucleotide bias of the DEMs was annotated and analyzed.The results showed that the composition of miRNA in the intestine of crucian carp had the uracil (U) bias.At the same time, it was also found that the ninth base has a higher U bias compared to others (Supplementary Figures 3B, C), indicating that the first and ninth bases may be closely related to the regulatory function of miRNAs.

Differentially expressed miRNAs in the intestine of C. auratus
DEMs were further screened for intestinal inflammation and immune-related genes in the intestine of crucian carp after the infection.The results showed that 4,422 miRNAs showed differential expression at three levels: upregulated, equal, and downregulated (Figure 5A).Particularly, 87 miRNAs showed significant differences in expression (p < 0.05), among which 49 miRNAs were significantly upregulated, including miR-133-3p, miR-206, miR-499, novel-miR-234, and novel-miR-485, and 38 miRNAs were significantly downregulated, such as miR-10b-3p, miR-203, miR-205, novel-miR-109, and novel-miR-129 (Figure 5B).In addition, 22 representative known miRNAs and novel miRNAs were screened (Supplementary Table 4), and the significant differential expression of these miRNAs indicates the involvement in the immune response process of the intestine.Furthermore, in order to analyze the different change patterns of miRNA expression abundance among different samples, a total of 87 miRNAs with significant differences in expression were analyzed for co-expression trends.The results showed that K-means clustering ultimately divided the data into two categories, with six DEMs such as miR-206, miR-206-3p, and novel-miR-1631 clustered into one category, showing an increasing co-expression trend.The other 81 DEMs clustered into one category, showing an overall decreasing co-expression trend (Supplementary Figure 4).

GO and KEGG enrichment analyses of the target genes
According to the DEM target genes annotated between the control group and infected group fishes, a total of 1,024 DEM target genes were annotated, among which GO analysis annotated the distribution of DEM target genes in three terms, including biological process, cellular component, and molecular function.in biological processes; the cells and cell parts had the highest enrichment of differentially expressed genes in cell components; the target genes that bind molecular functions were enriched the most in the molecular functional categories (Figure 6A).The result of KEGG pathway enrichment analysis showed that 921 DEM target genes were enriched into the metabolic pathway, mainly including glycerophospholipid metabolism, carbon metabolism, nitrogen metabolism, ether lipid metabolism, glutathione metabolism, arginine and profile metabolism, cytosine and metabolism, purine metabolism, and other metabolic pathways.However, the DEM target genes enriched to endocytosis were the most enriched, followed by melanogenesis (Figure 6B).In addition, DEM target

The cellular process was the most enriched target gene for DEMs
Analysis of Kyoto Encyclopedia of Genes and Genomes (KEGG) function (A) and correlation (B) of intestinal microflora in Carassius auratus infected with Rahnella aquatilis.genes were also enriched in the TGF-b signaling pathway, intestinal immune network for IgA production, lysosome, and other KEGG pathways.Analysis of the GO and KEGG functional enrichment preliminarily revealed the expression patterns of immune-related genes of crucian carp intestine infected with R. aquatilis.

Expression analysis of miRNAs and immune-related target genes
To understand the function of DEMs in intestinal mucosal immunomodulatory response, the key DEMs were then screened, including miR-203, miR-205, miR-206, miR-129, and miR-499, from the intestines of crucian carp.Further analyses of the interaction network between miRNAs and their target genes were performed, and the results showed that miR-203 and miR-205 jointly target transforming growth factor-b (TGF-b), miR-129 targets tumor necrosis factor a (TNF-a), miR-205 targets myeloid differentiation factor 88 (MyD88), and type C lysozyme (LysC), interleukin 22 (IL-22), miR-203, and miR-206 jointly target allograft inflammatory factor 1 like (AIF1L) (Figure 7A).Additionally, immune-related gene expression in the intestine was induced by R. aquatilis intraperitoneally injected in crucian carp.The results showed that the expression levels of TLR3 and TNF-a gene were upregulated until reaching a maximum at 72 hpi.However, the expression of C3 and LysC decreased significantly at 36 hpi (p < 0.001).F2 and F5 were significantly downregulated at 36 hpi while significantly upregulated at 72 hpi (p < 0.05).On the contrary, the expression levels of MHC II, NF-kB, and TGF-b were significantly upregulated and then decreased (Figure 7B).In addition, based on the immune-related pathway of the KEGG database, the protein interaction of intestinal immune-related genes was analyzed.The results showed that there were high-order interactions between 32 network edges and 18 nodes (p < 1.0e−16).It mainly includes TGFb and TNF-a in the TGF-b signaling pathway; F2, F5, and C3 in the complement and coagulation cascade signaling pathway; and Tolllike receptor (TLR3), MyD88, and NF-kB in the TLR signaling pathway (Figure 7C).Pathological observation results showed that the cell apoptosis was induced by R. aquatilis infection through apoptosis-related genes, such as expression of TAP1 and TNF-a in the kidney (data not shown).Furthermore, it was found that miR-129, miR-203, and miR-205 could bind to the 3′-UTR of the target genes TNF-a, NF-kB (nfkb2) and TGF-b, respectively, suggesting potential effects on the structure and function of three inflammation and immune-related proteins (Figures 8A, B).The phylogenetic tree results also showed that miR-205 of crucian carp was closely clustered with zebrafish Danio rerio and common carp C. carpio, which revealed the high conservation of its precursor sequences from the same family in different species of cyprinid fish (Figures 8C, D).Thus, they might participate in the local intestinal immune regulation of the TLR signaling pathway of crucian carp infected with R. aquatilis.

Analysis of the interaction among intestinal microflora and miRNAs
Spearman's correlation analysis was used to explore the potential relationship between the key DEMs of intestinal mucosal immunity and the intestinal microflora in the crucian carp intestine infected with R. aquatilis.The results showed that the expression of miR-203 and miR-205 in intestinal was positively correlated with the abundance of Proteobacteria but negatively correlated with that of Firmicutes (Figure 9A).Further analysis of the interaction between the intestinal microflora, miRNAs, and immune genes showed that miR-203 could inhibit NF-kB and that butyric acid may be produced by Firmicutes as short-chain fatty acids (SCFAs), which also inhibit NF-kB to inhibit intestinal mucosal inflammation.Proteobacteria may produce propionic acid, which acts as a histone deacetylase (HDAC) inhibitor, and then promote the release of the anti-inflammatory factor TGF-b.At the same time, miR-205 also directly inhibited TGF-b, thereby regulating the local immune response (Figure 9B), indicating that they were involved in local inflammation response related to target genes mediated via three key DEMs (i.e., miR-129, miR-203, and miR-205).It is suggested that miRNA may cooperate with the intestinal microflora to regulate intestinal inflammation and immune response and elucidate the potential molecular regulatory mechanism of apoptosis, intestinal microflora, metabolism, inflammation, and immune response-mediated miRNA-target gene axis in cyprinid fish.

Discussion
Crucian carp is one of the most common freshwater-farmed fish in China (29,43,44).In previous studies, a pathogenic R. aquatilis strain KCL-5 was isolated from the diseased crucian carp, which was recognized as a new aquatic pathogen that causes high lethality in fish (30)(31)(32).R. aquatilis is an emerging pathogen affecting both fish and humans, causing enteritis, septicemia, and death in severe cases (30,31,36,45,46), and has resulted in significant economic losses to fish farmers in China.Recent studies indicated that the thickness of the mucosal layer, muscular layer, and intestinal villus length was decreased after bacterial infection (47,48), and the pathological changes were basically consistent with those of the presently studied crucian carp infected with R. aquatilis.The enzyme activities of ALT, AST, and AKP were usually used in the disease-resistance immune and pathological damage (49).In this study, AKP was significantly increased after infection with R. aquatilis, indicating that it plays a crucial role in regulating the host immune response to bacterial infection (50).However, there is still no report on the molecular pathogenesis of enteritis and the interacting mechanism of the intestinal microflora associated with miRNA immune regulation induced by R. aquatilis infection in fish such as crucian carp.
In fish, the intestinal microflora is regarded as an organ that helps the host maintain a key physiological function of health (51).The relationship between fish and surrounding microorganisms may be mutually beneficial or pathogenic.Fish-related symbiotic intestinal microflora plays a role in nutritional supply, metabolic homeostasis, and immune defense (52).The intestinal microflora influences the physiology, development, and immunity of fish and serves as a barrier against pathogens.It has been found that the intestinal microflora can help the host to perform a variety of physiological and biochemical functions, thereby playing an important role in digestion, nutrient absorption, osmotic regulation, and the dynamic balance of the immune system (10,53).In fact, the diversity of intestinal microbial reflects the health and metabolic capacity of fish (54).Low diversity and stability of microbiota are closely associated with fish diseases (12,13,55).A great deal of evidence showed that improving the diversity of fish intestinal microbiota can more effectively protect the host against pathogen invasion (56)(57)(58).In particular, the genus Cetobacterium played an important role in the intestinal biochemical process of fish (59).In this study, the composition of the intestinal microflora of crucian carp changed significantly after infection with R. aquatilis, with the relative abundance of Cetobacterium first increasing significantly and then decreasing.The study of the intestinal microflora of the grass carp Ctenopharyngodon idella showed that the relative abundance of Proteobacteria and Bacteroidetes increased significantly, indicating that Bacillus subtilis plays a role in regulating the intestinal immune response, fat metabolism, and bacterial composition (60).By stimulating innate immunity and regulating the intestinal microflora, Bdellovibrio bacteriovorus reduced the abundance of pathogenic Aeromonas, decreased intestinal mucosal and epithelial cell damage, and improved the growth performance and survival rate of the mandarin fish Siniperca chuatsi (61).In this study, the Proteobacteria, Fusobacteria, Bacteroidetes, and Firmicutes were found to be predominant in the intestine of crucian carp, suggesting that R. aquatilis can cause intestinal flora disorder by affecting the changes in intestinal microecology and physiological metabolism, which then lead to intestinal pathology and immune response, but its immunomodulatory mechanism still requires further research.
In addition, the intestinal microflora is increasingly recognized as an endocrine system that can regulate host lipid metabolism (62, 63).In this study, several pyruvate and purine metabolism and other associated metabolic pathways were detected in crucian carp after infection with R. aquatilis.The host metabolism of shrimp is closely related to the intestinal microflora (64), and Proteus was found to play a crucial role in carbohydrate metabolism and amino acid biosynthesis in Litopenaeus vannamei (65).A previous study confirmed that the microflora in zebrafish D. rerio was involved in intestinal nutrition absorption and fatty acid metabolism (66).It is worth mentioning that the intestinal microflora is also involved in the regulation of immune responses (67).The intestinal mucosa recognizes bacterial antigens through pattern recognition receptors (PRRs; e.g., Toll-like receptors, RIG-I-like receptors, and NOD-like receptors) and activated signaling cascades to regulate the immune response (68).A study in humans had shown that the intestinal microflora regulated TLR4/NF-kB signaling pathways that trigger inflammation (69).In D. rerio, the intestinal expression of proinflammatory cytokines and antiviral mediators of the TLR/MyD88/NF-kB signaling pathway was induced by gut bacteria Vibrio and Aeromonas (70).Recently, a study reported that the intestinal bacteria Vibrio and Aeromonas activated the MyD88dependent signaling and induced innate immune response in D. rerio (71).The intestinal microflora of fish has received increasing attention due to its important role in maintaining host health.Most previous studies have focused on factors that control healthy intestinal microflora, such as diet, feeding conditions, and fish genes (72).However, there are few studies on the interaction between intestinal microflora, fish immunity, and disease.
Recently, miRNA has been shown to be an important regulator of gene expression by pairing with the complementary site of the target mRNA (14, 19, 73).The miRNAs that hold clues to the regulation of gene networks in the basic metabolic process and diseases caused by impaired miRNA function have been extensively studied in humans and model organisms (18).Although there have been many studies on immunomodulatory miRNA networks in higher vertebrates, the regulation of these responses in fish is not perfect (23).In common carp C. carpio, miR-155 targeted TNF-a, IFN-g, IL-1b, IL-6, and IL-10 and is further involved in regulating the expression of cytokines of the SOSC1 signaling pathway (22).MiR-122 and miR-200a-3p targeted TLR14 and TLR1, respectively, in the TLR signaling pathway and were thus involved in the regulation of gene expression in M. miiuy (74, 75), while miR-144 negatively regulated NF-kB signal transduction by targeting IL-1b (76).MiR-203 was a direct negative regulator of MDA5 in M. miiuy, which inhibited cytokines by regulating the RLR signaling pathway in teleosts (77).In addition, miR-203 and miR-148-1-5p inhibited the expression of the target genes IRAK4 and IRAK1, thereby regulating the NF-kB signaling pathway (23,24).Recently, a study reported that the miR-29a and miR-143 ameliorated the insulin resistance pathways in the common carp by targeting genes pik3r1 and pik3r3 (26).MiRNA is also involved in glucose and lipid metabolism in fish (19).Mennigen et al. detected that miR-122 was involved in the regulation of insulin metabolism and adipogenesis in rainbow trout Oncorhynchus mykiss (78).MiR-205 targeted and regulated acetyl-CoA carboxylase b (ACAC b), thus inhibiting lipid synthesis in the liver of tilapia Oreochromis mossambicus (25).In this study, the expression of miR-206 in the intestine of crucian carp was upregulated after infection with R. aquatilis.The expression of miR-206 was also upregulated in D. rerio by Mycobacterium balnei infection, resulting in a reduction of neutrophil response in the host (79).Therefore, inhibiting the expression of miR-206 may have a therapeutic influence on the development of the disease.Furthermore, this study found that the major DEMs (i.e., miR-129, miR-203, and miR-205) in the intestine of crucian carp were involved in the local immunomodulatory mechanism of the TLR signaling pathway by targeting TNF-a, NF-kB, and TGF-b.In addition, we speculated that butyric acid may be produced by Firmicutes inhibiting NF-kB.Propionic acid may be produced by Proteobacteria, which could act as an HDAC inhibitor and thereby promote the release of the anti-inflammatory factor TGF-b, which plays an important role in the regulation of intestinal inflammation (80,81).Previously, the isobaric tags for relative and absolute quantitation (iTRAQ) analysis of differentially expressed proteins of crucian carp infected with Aeromonas hydrophila was performed by our group, and enrichment analysis obtained TLR, TGF-beta, and NF-kB signaling pathways (20,28).However, the functional validation and analysis of target gene proteins were of great significance in fish immunology.These studies focused on multi-omics analysis, laying the foundation for investigations of the roles of miRNA and its target gene functions in crucian carp, preliminarily elucidating the potential molecular regulatory mechanism that miRNA may mediate the intestinal microflora, inflammation, and immune response in fish, but its specific regulatory mechanism needs further research.

A B
Analysis of the interaction network of regulatory mechanisms of key miRNAs/intestinal microflora/target genes in Carassius auratus intestinal infected with Rahnella aquatilis.(A) Heatmap of Spearman's correlation between the intestinal microflora and key miRNAs.(B) Potential molecular regulatory mechanism of intestinal microflora-mediated miRNA-target gene axis.

Conclusion
In general, we combined the intestinal microflora with miRNA transcriptome for the first time to explore the molecular regulation mechanism of miRNAs in mediating intestinal mucosal immunity.Through the multi-omics analysis of the intestinal microflora diversity in crucian carp induced by R. aquatilis infection, it was revealed that phyla Proteobacteria, Fusobacteria, Bacteroidetes, and Firmicutes were dominant bacterial communities.Further, miRNA transcriptome analysis showed that the key DEMs (i.e., miR-129, miR-203, and miR-205) participated in the local mucosal immune regulation of the TLR signaling pathway by targeting genes TNF-a, NF-kB, and TGF-b.In addition, the short-chain fatty acids may be produced by Firmicutes and Proteobacteria in the intestinal microflora, thereby inhibiting or promoting inflammatory factors TNF-a and TGF-b to regulate intestinal inflammation.MiRNA and the intestinal microflora may interact with the immune-related genes to play an important role in the mucosal immune response of crucian carp infected with R. aquatilis.This provides a new idea for elucidating the potential mechanism of intestinal inflammation and immune response-mediated miRNA-target gene axis in cyprinid fish.

2
FIGURE 2Analysis of intestinal microflora diversity of Carassius auratus infected with Rahnella aquatilis.(A) Box diagram of alpha-diversity index among groups.(B) Principal coordinates analysis (PCoA) of the intestinal microflora.(C) Unweighted pair group method with arithmetic mean (UPGMA) cluster tree analysis of intestinal microflora.Note: Groups a, b, and c were infected for 0, 36 and 72 hpi, respectively, the same as below.

3
FIGURE 3 Analysis of the composition of intestinal microflora in Carassius auratus infected with Rahnella aquatilis.(A) Operational taxonomic unit (OTU) Venn diagram of intestinal microflora.(B) Species distribution map of intestinal microflora at phylum level.(C) Species distribution map of intestinal microflora at genus level.(D) Histogram of linear discriminant analysis (LDA) value distribution of intestinal microflora.*P<0.05.

5
FIGURE 5 Analysis of differentially expressed miRNAs (DEMs) in the intestine of Carassius auratus infected with Rahnella aquatilis.(A) Volcanic map of DEMs.Each point in the figure represents a miRNA.(B) Cluster diagram of DEMs.The column represents different samples, and the row represents different miRNAs, which were clustered by log10 (TPM + 1e−6) value.Red stands for high expression of miRNAs and blue for low expression of miRNAs.

6
FIGURE 6Functional analysis of differentially expressed miRNA (DEM) target genes in Carassius auratus intestine infected with Rahnella aquatilis.(A) Gene Ontology (GO)-annotated classification chart.(B) Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment and scattered point graph.

7
FIGURE 7 Expression analysis of miRNA immune-related genes in Carassius auratus intestine infected with Rahnella aquatilis.(A) Network diagram of the interaction between key differentially expressed miRNAs and target genes.Blue represents downregulation, and red represents upregulation.(B) Temporal changes in the expression of target genes.*P<0.05,**P<0.01,***P<0.001.(C) Analysis of the interaction between genes and proteins related to intestinal immunity.Dark blue lines represent known interactions from curated databases.Purple lines represent experimentally determined known interactions.Green lines represent predicted interactions with neighboring genes.Red lines represent gene fusions.Blue lines represent gene clusters.Yellow lines represent text-mining evidence.Black lines represent co-expression.Light blue lines represent protein homology-based interactions.

8
FIGURE 8 Molecular characteristics of potential interaction between miRNAs and target genes proteins.(A) Three binding sites in miR-129, miR-203, and miR-205 to their target genes.(B) Tertiary structure of TNF-a, NF-kB, and TGF-b proteins.(C) Species genetic relationship of miR-205 sequences.(D) Phylogenetic tree analysis of pre-miR-205 from the 26 species.

TABLE 1
MiRNA statistics of intestinal in Carassius auratus infected with Rahnella aquatilis.