Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Expression characteristics and interaction networks of microRNAs in spleen tissues of grass carp (Ctenopharyngodon idella)

  • Yinli Zhao,

    Roles Funding acquisition, Writing – original draft

    Affiliation College of Biological Engineering, Henan University of Technology, Zhengzhou, Henan Province, P.R. China

  • Shengxin Fan,

    Roles Data curation

    Affiliation College of Animal Science and Technology, Henan Agricultural University, Zhengzhou, Henan Province, P.R. China

  • Pengtao Yuan,

    Roles Validation

    Affiliation College of Animal Science and Technology, Henan Agricultural University, Zhengzhou, Henan Province, P.R. China

  • Guoxi Li

    Roles Project administration, Writing – review & editing

    liguoxi0914@126.com

    Affiliation College of Animal Science and Technology, Henan Agricultural University, Zhengzhou, Henan Province, P.R. China

Abstract

The spleen is an important immune organ in fish. MicroRNAs (miRNAs) have been shown to play an important role in the regulation of immune function. However, miRNA expression profiles and their interaction networks associated with the postnatal late development of spleen tissue are still poorly understood in fish. The grass carp (Ctenopharyngodon idella) is an important economic aquaculture species in China. Here, two small RNA libraries were constructed from the spleen tissue of healthy grass carp at one-year-old and three-year-old. A total of 324 known conserved miRNAs and 9 novel miRNAs were identified by using bioinformatic analysis. Family analysis showed that 23 families such as let-7, mir-1, mir-10, mir-124, mir-8, mir-7, mir-9, and mir-153 were highly conserved between vertebrates and invertebrates. In addition, 14 families such as mir-459, mir-430, mir-462, mir-7147, mir-2187, and mir-722 were present only in fish. Expression analysis showed that the expression patterns of miRNAs in the spleen of one-year-old and three-year-old grass carp were highly consistent, and the percentage of miRNAs with TPM > 100 was above 39%. Twenty significant differentially expressed (SDE) miRNAs were identified. Gene ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis showed that these SDE miRNAs were primarily involved in erythrocyte differentiation, lymphoid organ development, immune response, lipid metabolic process, the B cell receptor signaling pathway, the T cell receptor signaling pathway, and the PPAR signaling pathway. In addition, the following miRNA-mRNA interaction networks were constructed: immune and hematopoietic, cell proliferation and differentiation, and lipid metabolism. This study determined the miRNA transcriptome as well as miRNA-mRNA interaction networks in normal spleen tissue during the late development stages of grass carp. The results expand the number of known miRNAs in grass carp and are a valuable resource for better understanding the molecular biology of the spleen development in grass carp.

Introduction

The spleen is an exclusively vertebrate organ and has important functions in immunity and hematopoiesis [1]. Similar to mammals, the spleen of fish also stores and produces erythrocytes, removes the aged and affected red blood cells [2], captures and destroys blood-borne pathogens [3, 4], and induces the adaptive immune response [5, 6]. Therefore, elucidating the regulatory mechanism of spleen function may assist in understanding the immune system thereby controlling infectious diseases in fish. It is known that the spleen’s basic histological architecture and related roles have been conserved during evolution throughout vertebrate species from fish to mammals [7]. Previous studies have shown that the development of the spleen begins with the condensation of the mesodermal mesenchyme underlying the dorsal mesogastrium epithelium. Next, there is the formation of a distinct splenic primordium on the dorsal stomach. The final organization of the spleen is established with the formation of the white pulp and the red pulp [8]. Despite differences in the time nodes of spleen embryogenesis, as well as the late stage development among vertebrates, the development of the spleen involves many biological processes including the growth of vascular channels, the influx of reticular cells to form a filtering network, migration of cells from other organs, and final maturation of transient and resident cell populations [9]. With the outbreak of infectious diseases in aquaculture environments, the research around the spleen’s physiological function in fish has been increasingly valued [1012]. However, understanding of fish spleen physiology is poor at the molecular level and is still largely restricted to a few species such as zebrafish [10], catfish [11], and sea bass [12]. Therefore, it is necessary to elucidate underlying molecular regulatory mechanisms of the spleen’s physiological function in fish.

MicroRNAs (miRNAs) are endogenous, small, non-coding RNAs approximately 22 nucleotides in length. They have an essential function in gene regulation by binding to the 3′-untranslated region (UTR) of target mRNAs, leading to mRNA degradation or the inhibition of mRNA translation [13]. Accumulating evidence showed that miRNAs participate in the control of many biological processes, such as cell proliferation and differentiation, growth, development, organogenesis, immunity, and metabolism [1416]. In particular, some studies have demonstrated that miRNAs are also involved in the regulation of spleen function in fish. For instance, 509 known miRNAs and 1157 novel miRNAs have been identified from four tissues of Micropterus salmoides, including the spleen [17]. In addition, 116 miRNAs from grouper spleen-derived cells [18], 194 conserved miRNAs and 12 novel miRNAs from the common carp spleen [19], and 452 miRNAs from C. semilaevis immune tissues (liver, head kidney, spleen, and intestine) [20] were also identified separately. It has been demonstrated that miR-182-3p was up-regulated after Vibrio parahaemolyticus flagellin stimulation in grouper spleen cells, and could directly target TLR5M [21]. In the spleen of Japanese flounder infected by Edwardsiella tarda, miR-194a was upregulated to a significant extent suppressing the type I interferon response [22]. These studies indicate that miRNA plays an important role in the regulation of spleen function and also provide new insight into the molecular biology of the spleen’s physiological function in fish. However, these studies largely focused on the role of miRNAs in the spleen upon exposure to antigens and chemicals, as well as other stress conditions. At present, research on miRNA expression and function in normal spleen tissues is relatively lacking.

The grass carp (Ctenopharyngodon idella) is an important economic species for freshwater aquaculture in China. The recent growth in the grass carp aquaculture industry is accompanied by increasingly severe infections induced by viruses, bacteria and parasites [23]. Therefore, the prevention and treatment of infectious diseases have become the most prominent problem in grass carp aquaculture. Considering the economic importance of this fish species, the research regarding the immune system and its function has received considerable attention. Recent studies have focused on miRNAs in the regulation of grass carp spleen function. For example, a total of 160 conserved miRNAs and 18 novel miRNAs were identified from grass carp in eight different tissues including the spleen [24]. There were 169 known miRNAs and 380 novel miRNAs obtained in 11 grass carp tissues including the spleen with high-throughput sequencing [25]. A total of 1208 miRNAs were identified from grass carp spleen at 0, 1, 3, 5, 7, and 9 days post-infection with a grass carp reovirus, of which 36 miRNAs exhibited differential expression [26]. In addition, 185 miRNAs associated with aeromonad septicemia in grass carp were identified from five grass carp tissues including the spleen at 4 h, 1 day, 3 days and 7 days post-infection with motile Aeromonas hydrophila [27]. While providing a basis for understanding the molecular biology of grass carp spleen, these studies investigated the spleen under host-pathogen conditions and were not specific to the spleen. Hence, no studies have reported miRNA transcriptomic profiling, and interaction networks between miRNAs in normal spleen tissue of grass carp.

The production of commercial culinary fish is an important aspect of grass carp aquaculture. During the production stage, grass carp grow rapidly and they are often susceptible to disease outbreaks. This is primarily due to factors such as aquaculture water quality and pathogenic microorganisms. Therefore, discovering the molecular biological characteristics of immune defense functions in the grass carp immune system will help to improve the control and prevention of disease. Considering the importance of the spleen in immune defense, and the important role of miRNAs in the regulation of spleen function, two miRNA libraries were constructed from the spleen tissue of one-year-old and three-year-old healthy grass carp. Through high throughput sequencing and subsequent bioinformatic analysis, the characterization of miRNA transcriptome profiling in normal spleen tissues of grass carp was described. Likewise, the miRNA regulatory network related to the spleen’s physiological function was established. The main aim of this work is to provide new insight into the molecular biological characteristics of the spleen’s physiological function. This will prove to be a valuable resource for future research on the role of miRNAs in the defense function of spleen in grass carp.

Materials and methods

Ethics statement

Experimental and animal care was executed in accordance with the program approved by the Experimental Animal Management Ordinance (Ministry of Science and Technology, China, 2004), and was approved by the Institutional Animal Care and Use Committee of Henan Agricultural University in China.

Sample collection and RNA isolation

The experimental animals in this study were the cultivated grass carp Ctenopharyngodon idella from the Animal Center of Henan Agricultural University. A total of 120 experimental fish at one-year (cis1) and three-year of age (cis3) were taken to the laboratory and temporarily raised under standard conditions for one week. Nine healthy individuals from each developmental stage were randomly selected for spleen tissue collection. The experimental fish were euthanized by tricaine methanesulfonate (MS-222) immersion according to the recommendations on the Operational Specification for the Euthanasia of Laboratory Fishes (local standard of Shandong Province, China; DB37/T 4134–2020). Firstly, 2000 mg/L concentrated stock solutions of MS-222 (Sigma Chemical Company, St Louis, MO, USA) were prepared and stored in amber jars at 4°C. Then, a suitable amount of concentrated stock solutions were taken and further diluted with water from the fish’s environment to a final concentration of 350 mg/L. Subsequently, anesthetic solution (350 mg/L MS-222) was placed in a tank, and fish were completely immersed in the anesthetic solution for 20 mins. Once swimming and opercular movement had completely stopped, fish were removed from the anesthetic tank and transferred to the recovery tank filled with water from the fish’s environment. From euthanized fish, samples of spleen tissues were immediately collected and frozen in liquid nitrogen. Total RNA was isolated from spleen tissues using Trizol reagent (TaKaRa, Dalian, China) according to the manufacturer’s instructions. The spleen tissues were homogenized by grinding in a mortar. The RNA concentrations were determined by NanoDrop 2000 spectrophotometry (Thermo Scientific, Wilmington, DE, USA), the RNA integrity was assayed by agarose gel electrophoresis. The 260/280 ratios and RIN (RNA Integrity Number) values of total RNA samples were 1.701~2.016 and 7.8~8.4, respectively, meeting the requirements of subsequent experiments.

MiRNA library construction and post-sequencing analysis

To characterize a general overview of miRNAs expressed in grass carp spleen tissue during normal development, the high-quality RNA samples from spleen tissues of three individuals were randomly mixed equally at each developmental stage, and the mixed RNA sample was used for cDNA library construction. Two miRNA libraries were constructed from the spleen tissue of grass carp at one year and three years of age using the TruseqTM Small RNA sample prep Kit (Illumina, San Diego, USA). The initial total RNA was 1μg for each miRNA library construction. The cDNA library was purified by 6% Novex TBE PAGE gel (Invitrogen, CA, USA) and enriched by the bridge PCR amplification using cBot Truseq SR Cluster Kit v3-cBot-HS (Illumina, San Diego, USA) and quantified using the TBS380 Picogreen (Invitrogen, CA, USA). Finally, the qualified cDNA library was sequenced using Hiseq4000 Truseq SBS Kit v3-HS (Illumina, San Diego, USA) on the Illumina HiSeq 4000 sequencing platform.

Raw reads were obtained by transforming the Illumina sequencing image data and cleaned using SeqPrep (Releases 34, https://github.com/jstjohn/SeqPrep) and Sickle (Releases 4, https://github.com/najoshi/sickle) software by removing low quality reads, and adaptor sequences. After elimination of redundancy, the 18–32 nt clean reads, called small RNA (sRNA), were used for subsequent analyses. For these clean reads, the reads with identical sequences were merged to obtain the unique sequence (that is, unique small RNA), which were used to count the sRNA species in the samples as well as the common and unique sequences between the samples. The length distribution of clean reads was calculated using the Fastx-Toolkit software (Version 0.0.13, http://hannonlab.cshl.edu/fastx_toolkit/). The clean reads were blasted against RNA families in the Rfam database (11.0, http://Rfam.sanger.ac.uk/) to discard ribosomal RNA (rRNA), small cytoplasmic RNA (scRNA), small nucleolar RNA (snoRNA), small nuclear RNA (snRNA), and transfer RNA (tRNA). Subsequently, the remaining sRNA sequences were mapped to the reference genome using Bowtie software (Releases 1.2.3) [28] to analyze the distribution of sRNAs on the reference sequence. Due to the lack of available genomic information in grass carp, the genome of zebrafish with close relationships was used as a reference in this study.

Identification of miRNAs

The conserved miRNAs in the spleen were identified in grass carp. All clean sequences after removing non-miRNA sequences were searched against the sequences of mature metazoan miRNAs in the miRBase (Release 22.1) [29]. The sRNA sequences with identical or related sequences (no more than one mismatch in the seed region or a few end nucleotide variations in the entire length) to the reference mature miRNAs were mapped to the reference genome using the Bowtie software (Releases 1.2.3) [28]. The 300 base sequence flanking each sRNA sequence was extracted as the potential precursor of miRNAs, and their folding secondary structures were evaluated using RNAfold (Version 2.4.13) [30] and miRDeep2 (Version 0.1.3) [31]. Based on the above analysis, the sequences that matched a miRNA in the miRBase, and whose precursor could form a typical stem-loop hairpin structure, were identified as a known grass carp conserved miRNA. In addition, all unannotated sequences were put into the miRDeep2 modules program and used to identify potential novel miRNAs by miRDeep2 under the system default parameter. These novel miRNAs were named as described by Ambros et al.(2003) [32].

Identification of significantly differentially expressed miRNAs

The expression levels of the identified miRNAs were estimated using transcripts per million clean reads (TPM). Based on the normalized TPM values, the differentially expressed miRNAs were analyzed by DEGseq (Releases 3.12) [33]. The fold-change for each miRNA between one-year and three years of age was calculated. The p-value was determined by Fisher’s exact test and adjusted using q-value [34]. The q-value <0.05 and |log2 (fold change)| >1 were set as the threshold for SDE miRNAs.

MiRNA target prediction, functional enrichment analysis and interaction analysis between miRNAs and mRNAs

To predict the potential target genes of miRNAs, transcriptome profiles of spleen tissue from grass carp one year and three years of age were also constructed using the same tissue samples and the miRNA libraries [35]. The raw sequencing data had been deposited into the NCBI database Sequence Read Archive under number SRP078553. Based on the above mentioned transcriptome profile data, the potential target genes of the SDE miRNAs were predicted using the miRanda (Releases 3.3a) [36] and TargetScan (Releases 6.2) [37] software. When the TargetScan software was used to predict the target genes, only the sites belonging to the top quartile of ranked predictions, present in zebrafish, mouse, and human, were selected as the true binding sites. For the miRanda software, all detected targets with the threshold parameters of S >140 and ΔG <-10 kcal/mol and strict 5’ seed pairing were considered as potential targets. In order to make the identification of target genes robust, the target genes that were simultaneously predicted by two programs were selected for subsequent functional enrichment analysis. GO and KEGG pathway enrichment analyses were performed using the Goatools (Version v0.7.6) [38] and KOBAS (Version 3.0) [39] software, respectively. The corrected p-value ≤0.05 was used as the threshold for significant enrichment of a GO term or pathway. Based on the results of the KEGG pathway enrichment analysis, the negatively correlated SDE miRNA-mRNA pairs expressed in the spleen of grass carp were selected from the pathways associated with immunity, lipid metabolism, and cell proliferation and differentiation, and their interaction networks were constructed using Cytoscape software (Version 3.2) (http://www.cytoscape.org/).

Quantitative real-time PCR (qRT-PCR) analysis

Using the qRT-PCR analysis, six SDE miRNAs including miR-18c, miR-451, miR-144-5p, miR-100-3p, miR-363-3p and miR-2188-5p were validated in spleen from one year and three years of grass carp. The qRT-PCR was performed using the customized kits containing the stem-loop reverse transcription primer and miR-specific molecular beacon probe (RiboBio, Guangzhou, China) on the LightCycler® 96 instrument (Roche, Basel, Switzerland). The total RNA from spleen tissues were treated with DNase I (RNase-free) (TaKaRa, Dalian, China), and then reverse-transcribed with M-MLV Reverse Transcription Reagents (Invitrogen, CA, USA) and the stem-loop reverse transcription primer according to the manufacturer’s instructions. The PCR amplification reaction system was as follows: 1 μL cDNA product, 5 μL 2×SYBR Premix Ex Taq II, 1 μL specific primer (10 μmol/L), and 3 μL ddH2O. The PCR amplification program for these miRNAs was as follows: 95°C for 3 mins; 40 cycles of 95°C for 12 s, 60°C for 40 s, and 72°C for 30 s; and 10 mins extension at 72°C. All qRT-PCR analyses included three biological replicates, and each biological replicate contained three technical replicates. U6 small nuclear RNA was used as endogenous control, and the relative expression levels were calculated by the 2−ΔΔct method [40]. The t-test was used to determine the p-value, and p-value<0.05 was considered to be a significant difference.

Statistical analysis

The qRT-PCR data and graphs were generated in GraphPad Prism (version 5.0) software (San Diego, CA, USA). Statistical significance of the qRT-PCR results was tested by performing two-tailed unpaired t-tests. Data are presented as the means containing three replicates.

Results

MiRNA library sequencing and sequence analysis

Through Illumina HiSeq sequencing, 15,839,580 and 16,587,512 high-quality clean reads with a length range of 18–32 nt were obtained from spleen tissue of grass carp at one year (cis1) and three years (cis3) of age, respectively. The length of these cleaned reads consisted of 20–23 nucleotides, and peaked at 22 nt. This represented approximately 35% of the total number of clean reads in the libraries (S1A Fig). The statistical results showed that a total of 31,811,397 sRNAs representing 100,557 unique sRNAs were shared in the cis1 and cis3 libraries (S1B and S1C Fig). In addition, 293,603 sRNAs were annotated with the RFam 11.0 database, accounting for 86.49% of the total number of clean reads. In particular, the number of matched miRNA sequences reached 80% (S1D Fig), indicating a high miRNA abundance in the library. However, the matched unique sRNAs in the cis1 and cis3 libraries only accounted for 23.01% and 19.14%, respectively, indicating that there were many specific sRNAs in the spleen of grass carp. All clean sRNA sequences, after removing rRNA, tRNA, snoRNA, and other snRNAs, were aligned with the reference genome sequence, and a total of 26,731,633 (87.92%) reads were perfectly mapped. Among them, the number of mapped sequences on chromosome 14 (25.02%) was the largest, followed by chromosome 8 (13.37%), chromosome 2 (12.14%), chromosome 23 (11.25%), and chromosome 15 (10.98%) (S1E Fig).

MiRNAs expressed in grass carp spleen tissue

To identify conserved miRNAs in the spleen of grass carp, all clean sequences mapped to reference genomes were compared to the mature sequences of miRNAs from all animals in the miRBase (release 22.1). A total of 26,022,975 (80.25%) reads with homology to known miRNAs were annotated as conserved miRNAs. When combined with the structural predictions of precursor sequences, 324 known conserved miRNAs were identified in grass carp spleen tissue from two developmental stages (S1 Table). Of the 324 conserved miRNAs, 132 miRNAs were perfectly matched with known miRNAs of zebrafish. However, 157 miRNAs were matched with known zebrafish miRNAs but showed variation at the end of mature sequences or differences in genomic localization of precursor sequences. The remaining 35 miRNAs were homologous with known miRNAs from 24 animals, including Cyprinus carpio, Ictalurus punctatus, Oryzias latipes, Salmo salar, Ciona intestinalis, Ornithorhynchus anatinus, Taeniopygia guttata, Ascaris suum, Bombyx mori, Caenorhabditis elegans, Cricetulus griseus, Drosophila melanogaster, Drosophila willistoni, Eptesicus fuscus, Equus caballus, Homo sapiens, Macaca mulatta, Monodelphis domestica, Mus musculus, Pan troglodytes, Pristionchus pacificus, Rattus norvegicus, Schistosoma mansoni, and Tupaia chinensis. In addition, sRNAs sequences that did not match known miRNAs were compared with the reference genome sequence, and novel miRNAs were predicted by miRDeep2 according to the criteria of potential miRNAs identity mentioned above. Ultimately, nine novel miRNAs were identified in the spleen of grass carp (S2 Table).

In addition, a family analysis of 324 known miRNAs was carried out using the sequence similarity searches. The results showed that 300 miRNAs belonged to 105 miRNA gene families. Among these miRNA families, some families such as let-7, mir-1, mir-10, mir-124, mir-8, mir-7, mir-9, mir-153, mir-190, mir-137, and mir-25 were present both in invertebrates and vertebrates, while fourteen miRNA families existed only in fish. These were possibly fish-specific miRNAs that included mir-459, mir-430, mir-462, mir-7147, mir-2187, mir-722, mir-724, mir-737, mir-725, mir-728, mir-729, mir-731, mir-734, and mir-738 (Fig 1). More than two family members in 64 miRNA families were expressed in the spleen of grass carp. In particular, the miRNA families mir-10, let-7, mir-15, mir-17, mir-130, mir-30 and mir-8 had twenty, thirteen, twelve, eleven, nine, nine, and nine members, respectively.

thumbnail
Fig 1. Highly conserved and fish-specific miRNA family in grass carp spleen.

The presence of miRNA is indicated by a plus (+); the absence of miRNA is indicated by a minus (-). Abbreviations: hsa, Homo sapiens; ptr, Pan troglodytes; mmu, Mus musculus; eca, Equus caballus; bta, Bos taurus; ssc, Sus scrofa; cfa, Canis familiaris; oan, Ornithorhynchus anatinus; gga, Gallus gallus; aca, Anolis carolinensis; xtr, Xenopus tropicalis; ccr, Cyprinus carpio; dre, Danio rerio; hhi, Hippoglossus hippoglossus; ipu, Ictalurus punctatus; ola, Oryzias latipes; ssa, Salmo salar; pma, Petromyzon marinus; bbe, Branchiostoma belcheri; cin, Ciona intestinalis; lva, Lytechinus variegatus; bmo, Bombyx mori; dme, Drosophila melanogaster; dpu, Daphnia pulex; lgi, Lottia gigantea; asu, Ascaris suum; bma, Brugia malayi; egr, Echinococcus granulosus; xbo, Xenoturbella bocki. P, protostomia; D, deutostomia; V, vertebrata; Pi, pisces.

https://doi.org/10.1371/journal.pone.0266189.g001

The expression characteristic of miRNAs in grass carp spleen tissue

Among the 324 known miRNAs identified in this study, 308 miRNAs were shared in two development stages of the spleen, seven miRNAs were expressed only in the spleen of one-year-old grass carp, and nine miRNAs were expressed specifically in the spleen of three-year-old grass carp. The distribution of miRNA expression levels showed that the expression pattern of miRNAs was highly consistent between one-year-old and three-year-old spleen (S2A–S2C Fig). In addition, the abundance of these identified miRNAs was different in the grass carp spleen. The percentages of miRNAs with TPM > 100 in the cis1 and cis3 were 40.7% and 39.2%, respectively. These abundant miRNAs were primarily concentrated in families such as let-7, mir-10, mir-30, and mir-15, which are known to be involved in immunity (Table 1).

thumbnail
Table 1. The abundance and function of miRNAs with TPM > 1000 in grass carp spleen tissue.

https://doi.org/10.1371/journal.pone.0266189.t001

The expression levels of spleen miRNAs in grass carp at one year and three years of age were compared (S2D Fig). Under the criteria of |log2(foldchange)|>1 and q-value <0.05, 20 SDE miRNAs were identified compared with one-year old spleen (S2E Fig and S3 Table). Surprisingly, only cid-miR-212-5p and cid-miR-100-3p among these SDE miRNAs were upregulated. However, most miRNAs with high abundance were not differentially expressed between the two development stages in the spleen. In addition, six SDE miRNAs including cid-miR-144-5p, cid-miR-100-3p, cid-miR-451, cid-miR-18c, cid-miR-2188-5p, and cid-miR-363-3p were verified by qRT-PCR analysis. Compared with the RNA-seq data, the expression levels determined by qRT-PCR showed similar expression patterns (Fig 2), indicating that the sequencing analysis results were reliable in this study.

thumbnail
Fig 2. The qRT-PCR verification of the selected SDE miRNAs in the spleen of grass carp at one year and three years of age.

The A, B, C, D, E and F in figure are the qRT-PCR results and TPM value of RNA-seq corresponding to cid-miR-144-5p, cid-miR-100-3p, cid-miR-451, cid-miR-18c, cid-miR-2188-5p and cid-miR-363-3p, respectively. Data are represented by the mean ± SEM (n = 9) in qRT-PCR results. **p-value < 0.01 represent significant difference.

https://doi.org/10.1371/journal.pone.0266189.g002

The potential function of SDE miRNAs

Based on the spleen transcriptome data, a total of 6,887 potential target genes of the SDE miRNAs were predicted in grass carp spleen tissue. Subsequently, the KEGG pathway enrichment analyses for the predictive target genes were performed, and 37 significantly enriched pathways were identified. These pathways were clustered into seven classes including environmental information processing, genetic information processing, cellular processes, organismal systems, drug development, human diseases and metabolism (Fig 3). Additionally, some of these pathways were involved in immunity and lipid metabolism, such as the B cell receptor signaling pathway, T cell receptor signaling pathway, insulin signaling pathway, PPAR signaling pathway, and fatty acid biosynthesis. Moreover, GO enrichment analysis of these target genes showed that 1,319, 301 and 790 GO terms were significantly enriched in biological process, cellular component and molecular function, respectively. Functions of the SDE miRNAs were concentrated in the ubiquitin-dependent protein catabolic process, phospholipid metabolic process, macromolecule catabolic process, chromatin modification, positive regulation of cellular metabolic process, endoplasmic reticulum membrane, bounding membrane of organelle, ubiquitinyl hydrolase activity, and guanyl-nucleotide exchange factor activity (Fig 4A). In particular, many significant enrichment GO terms were related to hematopoietic, immune response and lipid metabolism in the category of biological process, for instance, erythrocyte differentiation, hematopoietic or lymphoid organ development, immune response, immune system processes, lipid metabolic processes, and fatty acid metabolic processes (Fig 4B).

thumbnail
Fig 3. The significantly enriched pathways for potential genes targeted by SDE miRNAs in grass carp spleen.

The abscissa indicates the enriched KEGG pathways and class. The ordinate indicates the enrichment ratio in each KEGG pathway, and the enrichment ratio is calculated as follows: sample number/background number. ***p-value < 0.001; **p-value < 0.01; *p-value < 0.05.

https://doi.org/10.1371/journal.pone.0266189.g003

thumbnail
Fig 4. The significantly enriched GO terms for potential genes targeted by SDE miRNAs in grass carp spleen.

(A) Part of significantly enriched GO terms in the category of biological process, cellular component, and molecular function. (B) Significantly enriched GO terms associated with immunity and lipid metabolism in the category of biological process. The abscissa indicates the enriched GO terms and category. The ordinate indicates the enrichment ratio for each GO term, and the enrichment ratio is calculated as follows: sample number/background number. ***p-value < 0.001; **p-value < 0.01; *p-value < 0.05.

https://doi.org/10.1371/journal.pone.0266189.g004

Interaction networks of mRNAs in grass carp spleen

To further reveal the potential regulatory mechanisms of these SDE miRNAs, negatively correlated SDE miRNA-mRNA pairs were identified based on the KEGG enrichment analysis results for the SDE miRNA prediction target genes, and then miRNA-mRNA interaction networks were constructed by Cytoscape software (Version 3.2). Thereinto, the immune-related interaction network contains 10 miRNA-mRNA pairs (Fig 5), which were screened from 13 immune and hematopoietic pathways including toll-like receptor signaling pathway, RIG-I-like receptor signaling pathway, NOD-like receptor signaling pathway, cytosolic DNA-sensing pathway, T cell receptor signaling pathway, B cell receptor signaling pathway, hematopoietic cell lineage, Fc epsilon RI signaling pathway, complement and coagulation cascades, natural killer cell mediated cytotoxicity, leukocyte transendothelial migration, chemokine signaling pathway, and platelet activation. In this network, cid-miR-144-5p is the key node of this network. The gene NFKBA targeted by cid-miR-144-5p is a common component of seven pathways. Similarly, another gene P38 targeted by cid-miR-144-5p is also a common component of seven pathways. Therefore, cid-miR-144-5p appears to have a wide range of regulatory effects. In contrast, other miRNAs in this network have specific roles.

thumbnail
Fig 5. The miRNA-mRNA interaction networks are associated with immunity and hematopoiesis.

The dot indicates miRNA, and the box indicates a target gene that had a negative correlation with a given miRNA.

https://doi.org/10.1371/journal.pone.0266189.g005

The development of spleen is accompanied by a series of cell proliferation and differentiation events. To this end, a miRNA-mRNA interaction network associated with cell proliferation and differentiation was constructed (Fig 6), which contained 28 miRNA-mRNA pairs involving TNF signaling pathway, MAPK signaling pathway, PI3K-Akt signaling pathway, AMPK signaling pathway, Ras signaling pathway, Jak-STAT signaling pathway, NF-kappa B signaling pathway, mTOR signaling pathway, ECM-receptor interaction, cGMP-PKG signaling pathway, cAMP signaling pathway, FoxO signaling pathway, VEGF signaling pathway, and Rap1 signaling pathway. In this network, genes such as SOCS3, LAMA3_5, DDIT4, ATP1A, TNFB, HGF, NFKBIA, CFLAR, P38, and RASGRF2 are involved in more than two pathways. Likewise, cid-miR-122, cid-miR-18c, cid-miR-20b-5p, cid-miR-205-5p, and cid-miR-363-3p, as well as members of mir-144, mir-203, and mir-9 families, all acted on more than two genes. These miRNAs, together with their predicted target genes, constituted the key nodes of this network and appear to have a broad range of functions. In contrast, some miRNAs, such as cid-miR-182-5p and cid-miR-375 in the MAPK signaling pathway, cid-miR-725-3p in the AMPK signaling pathway, and cid-miR-216a in the FoxO signaling pathway, which target only one gene in this network and have specific roles.

thumbnail
Fig 6. The miRNA-mRNA interaction networks are associated with cell proliferation and differentiation.

The dot indicates miRNA, and the box indicates a target gene that had a negative correlation with a given miRNA.

https://doi.org/10.1371/journal.pone.0266189.g006

In addition, we found that the predicted target genes of SDE miRNAs were significantly enriched in some pathways associated with lipid metabolism. Hence, an interaction network was constructed using 23 miRNA-mRNA pairs selected from five pathways, such as fatty acid metabolism, glycerolipid metabolism, insulin signaling pathway, adipocytokine signaling pathway, and PPAR signaling pathway (Fig 7). In this network, cid-miR-182-5p, cid-miR-18c, cid-miR-725-3p, cid-miR-144-5p, and members of the mir-203 family all target more than two genes. Similarly, SOCS3, SCD, and glpK are common components of two pathways and were targeted by three miRNAs. These genes became the interrelated centers in the miRNA-mRNA network maintaining collectively the stability of this interaction network. It is noteworthy that the gene PPP1R3 in the insulin signaling pathway was targeted by four miRNAs, and the gene ACSL was a common component of fatty acid metabolism, the PPAR signaling pathway, and the adipocytokine signaling pathway targeted by three miRNAs.

thumbnail
Fig 7. The miRNA-mRNA interaction networks are associated with lipid metabolism or deposition.

The dot indicates miRNA, and the box indicates a target gene that had a negative correlation with a given miRNA.

https://doi.org/10.1371/journal.pone.0266189.g007

Discussion

Due to the pivotal role played in modulating immune responses, the spleen is gaining more attention. In the past decade, spleen research has focused on morphology, organogenesis, and function as well as the genetic regulation underlying these processes. In recent years, with the rapid development and application of high-throughput sequencing technologies, considerably more attention has been focused on elucidating the molecular mechanisms regulating spleen function. The primary goal is to identify genes involved in spleen function by transcriptome analysis in different species, such as the mouse [8, 82], goose [83], pig [84], giant panda [85], common carp [86] as well as bighead carp and silver carp [87]. In addition, to elucidate epigenetic regulation of spleen function, non-coding RNAs such as miRNA, lncRNA or circRNA were identified by transcriptome sequencing in different species, such as the mouse [82], pig [83], giant panda [84], Chinese giant salamander [88], green-spotted puffer fish [89], rainbow trout [90], and common carp [19]. These studies improve our understanding of the regulation of spleen function in animals. However, some investigations are limited to spleen function under stress conditions such as pathogen infection. As such, little is known about the molecular biological mechanisms regulating the spleen’s physiological function. Although the transcriptome profiles involved in grass carp spleen have been studied [35], miRNA expression profiles in normal spleen tissues have not been reported. In this study, 324 conserved miRNAs and nine novel miRNAs were identified from the spleen tissues of one-year-old and three-year-old grass carp (S2 Table). The expression profile characteristics of these miRNAs were described, and the miRNA regulatory network related to the spleen’s physiological function was constructed. These results expand the number of known miRNAs associated with spleen function and provide a meaningful framework to understand the molecular mechanisms regulating spleen function in grass carp and closely related species.

Evolutionary conservation is one key characteristic of miRNAs. Many homologous miRNAs are highly conserved in sequence across metazoans [91]. In addition, some miRNA families widely exist in multiple species and exhibit either expansion or contraction in different species, thus playing the corresponding regulatory role by generating different members of miRNAs [92]. Accumulating evidence indicates that conserved miRNAs across species may have similar functions in regulating different biological processes from lower to higher organisms. Thus, the identification of conserved miRNAs or the conservation analysis of miRNAs may help us infer the functions of miRNAs in a specific biological context on their known functions in other species. Due to the lack of adequate and available genomic information and the absence of known grass carp miRNAs in the miRBase database at present, the main approach for miRNA identification of grass carp was to match the sequences obtained by high-throughput screening with the known miRNA sequences of all animals in the miRBase database, thereby identifying known conserved miRNAs. Likewise, this approach was also used in this study to identify 324 known conserved miRNAs from the grass carp spleen. These identified miRNAs belong to 105 miRNA gene families, which are found in 94 species ranging from lower Platyhelminthes to the highest mammals. Among them, 28 families such as mir-1, mir-124, mir-7, and mir-9, were shared by both invertebrates and vertebrates; 63 miRNA families such as mir-103, mir-128, mir-130, mir-132, and mir-138, were present only in vertebrates; 14 families such as mir-459, mir-430, mir-462, and mir-7147, were identified only in fish (Fig 1). Members of 15 miRNA families, including mir-10, let-7, mir-15, mir-17, mir-130, mir-30, mir-8, mir-19, mir-153, mir-181, mir-27, mir-132, mir-25, mir-221, and mir-29, were expressed in grass carp spleen, accounting for 41% of the conserved miRNAs identified in this study. These data demonstrate that miRNAs expressed in grass carp spleen were dominated by conserved miRNAs. This is consistent with previous results that reported highly conserved miRNAs in the spleen of Chinese Giant Salamander [88]. Previously, the identification of miRNAs in the spleens of Chinese giant salamander [88], common carp [19], rainbow trout [90] and green-spotted puffer fish [89] has been reported. We compared miRNAs expressed in the spleen of grass carp with those of the four species and found that 26.5%, 33.3%, 36.1% and 53.7% of the grass carp spleen miRNAs were identical to miRNAs in the spleen of Chinese giant salamander, common carp, rainbow trout, and green-spotted puffer fish, respectively, although the number of miRNAs identified by different studies varied greatly. In particular, 22 miRNAs from the families of mir-140, mir-144, mir-145, mir-181, mir-210, mir-2188, mir-221, and mir-455 were expressed in the spleen of the above five species. These studies demonstrated that miRNAs expressed in the vertebrate spleen are phylogenetically conserved. These conserved miRNAs may regulate the spleen physiology in the vertebrate. Therefore, the known conserved miRNAs identified in this study can provide a valuable resource for further research on the molecular mechanism of the spleen’s physiological function in grass carp and other vertebrates.

It is estimated that miRNA genes represent approximately 1% of the known eukaryotic genomes [92]. Thus, the identification of novel miRNAs has become an essential step for understanding the variety of functions of miRNAs. Currently, next-generation sequencing such as Illumina deep sequencing is an efficient method for the identification of novel miRNAs in animals lacking adequate and available genome information. Unfortunately, only nine novel miRNAs in the grass carp spleen were identified using Illumina deep sequencing in this study. This is similar to previous studies where twenty-four and twelve novel miRNAs were obtained in the spleen of the Chinese giant salamander [88] and the common carp [19], respectively. In contrast, another study identified 926 novel miRNAs in the spleen of rainbow trout [90]. Although these studies all adopted next-generation sequencing, the number of novel miRNAs obtained was significantly different due to the use of different criteria for miRNA identification. In fact, 298 potential novel miRNAs were initially identified using the miRDeep2 software (Version 0.1.3). However, after strict sequence alignment, many of these potential novel miRNAs were found to be the length and/or end-sequence variations of known miRNA sequences identified in this study. Therefore, these potential novel miRNAs that match known miRNA but have a difference of 1–4 nucleotides at the 3’ end were treated as corresponding homologue miRNAs. Finally, 157 conserved miRNAs were identified that match known zebrafish miRNAs in the miRBase database but are different at the end of mature sequences (S1 Table). The above criteria were used to identify miRNA based on the following facts: (i) the length and sequence heterogeneity of miRNA are often detected in high-throughput sequencing [19, 9395]; (ii) the length and sequence heterogeneity in miRNA datasets are often the result of sequencing errors and post-transcriptional modifications of RNA [96]; and (ⅲ) a difference at the last base of the miRNA between two species caused by sequencing errors is generally not treated as a change [97]. Therefore, based on the above results and analysis, it was inferred that the number of novel miRNAs or organ-specific miRNAs in the spleen of grass carp may be relatively small, and are mainly widespread and conserved miRNAs. The majority of the conserved miRNAs showed length and/or end-sequence variations that affect the formation of the miRNA/ target mRNA hybrid duplex, which increases the diversity of miRNAs and their targets [98] and thus perform distinct roles in the grass carp spleen.

The expression profile of miRNAs reflects their different roles in specific tissues or developmental stages as well as corresponding biological mechanisms. In this study, we analyzed the expression profile characteristics of miRNAs in the spleen of grass carp at two developmental stages and found that 308 known miRNAs were shared at the two developmental stages. The TPM density distribution showed that the expression levels of these miRNAs were highly consistent between one-year-old and three-year-old spleens (S2 Fig). This temporal expression feature of miRNAs in the grass carp spleen is completely consistent with previous results in the spleen of giant pandas [88]. In grass carp spleen, known miRNAs are highly expressed, and most of them have more than one family member. However, novel miRNAs are usually weakly expressed as previously reported in the spleen of Chinese giant salamander [88] and common carp [19]. Most known miRNAs are conserved and they are often widespread. Five miRNAs, including miR-21, let-7a, miR-26a-5p, miR-451, and miR-142a-5p, were among the top 10 miRNAs with the highest abundance in the grass carp spleen in this study (Table 1) and in the spleen of the Chinese giant salamander in a previous report [88]. The above results indicate that the miRNA expression profile of spleen may have a high consistency among vertebrates.

Generally, abundant miRNAs play fundamental and broad regulatory functions in maintaining biological processes. Interestingly, many of the known miRNAs identified in this study were highly expressed where miRNAs with TPM > 100 reached more than 39%. In particular, ten, seven, five, four and four members in the family of let-7, mir-10, mir-30, mir-126 and mir-15, respectively, were expressed in high abundance in grass carp spleen (Table 1). Examination of the literature for miRNAs of TPM > 1000, it was found that these miRNAs were reported in other species and closely involved in animal immunity (Table 1). This included innate and adaptive immune responses, T-cell activation and cytotoxicity, innate anti-viral response, phagocytosis of myeloid inflammatory cells, early iNKT cell development, inflammatory factor secretion from macrophages, macrophage polarization, hypoxia-induced immunosuppression, Treg-mediated immunological tolerance, and inflammatory response to bacterial infection (Table 1). These suggest that the abundance miRNAs in the grass carp spleen may function similar to those of their orthologs and play important roles in regulating the spleen’s physiological function.

The expression of miRNAs usually exhibits temporal and spatial specificity. In this study, twenty differentially expressed miRNAs between one-year-old and three-year-old grass carp spleens were identified and closely related to hematopoiesis and immunity as determined by functional enrichment analysis (Figs 3 and 4). Most of the differentially expressed miRNAs were down-regulated compared with the one-year-old spleen. Except for miR-192, miR-144-5p, miR-2188-5p, miR-451 and miR-144-3p, these miRNAs were expressed in relatively low abundance during grass carp spleen development (S3 Table). Interestingly, in these differentially expressed miRNAs, target gene prediction revealed potential interactions between miR-144-3p and miR-205-5p with Tlx1, miR-144-5p with Bapx1, and miR-192 with Wt1. It is known that Tlx1, Bapx1 and Wt1 are major marker genes for spleen development. In particular, transcription factor Tlx1 is a marker for early spleen mesenchymal cells controlling cell fate specification and organ expansion during spleen development [99]. These studies suggest that differentially expressed miRNAs play important roles in regulating spleen immune function in grass carp.

A single miRNA can regulate multiple target genes or a single gene can be targeted by multiple miRNAs simultaneously. These miRNAs construct a complex network and collectively regulate gene expression. The spleen has important functions in haematopoiesis and immunity. The spleen’s development is achieved through complex interactions between spleen mesenchymal cells and invading hematopoietic and endothelial cells, involving a series of cellular events. To reveal the regulatory mechanism of miRNAs underlying the spleen’s physiological function in grass carp, miRNA–mRNA interaction networks related to haematopoiesis and immunity, and cell proliferation and differentiation were constructed (Figs 5 and 6). In these interaction networks, some genes, as a component of multiple pathways, link multiple pathways to each other. For instance, NFKBIA links across seven immune-related pathways such as toll-like receptor signaling pathway, cytosolic DNA-sensing pathway, B cell receptor signaling pathway, chemokine signaling pathway, RIG-I-like receptor signaling pathway, NOD-like receptor signaling pathway and T cell receptor signaling pathway. Likewise, P38 links across five pathways closely related to cell proliferation and differentiation including TNF signaling pathway, MAPK signaling pathway, VEGF signaling pathway, Rap1 signaling pathway and FoxO signaling pathway. Some miRNAs regulate the functions of multiple pathways by targeting multiple genes in these interaction networks. Thus, it is proposed that a complex regulatory network, formed by the interaction between miRNAs and their target gene and between pathways, regulates the spleen’s physiological function in grass carp. Among the miRNAs constituting the above interaction network, some miRNAs including miR-144-3p, miR-192, and miR-100-3p have previously been demonstrated to be involved in immunity [42, 54, 69]. In particular, two members of the mir-144 family have potential target relationships with seven genes in three interaction networks, including TLR5, IL11, NFKBIA, P38, SOCS3, ACSL, and HGF. MiRNAs from the mir-144 family have been shown not only to participate in immunity, but also play an important role in hematopoietic regulation [100, 101]. Therefore, the interactive networks constructed in this study may truly reflect the complexity of the post-transcriptional regulation of the spleen’s physiological function in grass carp.

Supporting information

S1 Fig. Post-sequencing analysis of sRNAs from the spleen of grass carp.

(A) Length distribution of clean reads. (B) Venn diagram of total sRNA between the cis1 and cis3. (C) Venn diagram of unique sRNA between the cis1 and cis3. (D) Annotation of unique sRNA in the RFam 11.0 database. (E) Number and distribution of clean reads mapped to the reference genome sequence.

https://doi.org/10.1371/journal.pone.0266189.s001

(TIF)

S2 Fig. MiRNA expression analysis in grass carp spleen.

(A) Histogram of miRNA expression distribution in cis1. (B) Histogram of miRNA expression distribution in cis3. (C) Boxplot of miRNA expression distribution in cis1 and cis3. (D) Scatterplot comparing the number of miRNA reads for cis1 and cis3. (E) Heatmap of the differentially expressed miRNAs between cis1 and cis3.

https://doi.org/10.1371/journal.pone.0266189.s002

(TIF)

S1 Table. The conserved miRNAs expressed in the spleen of grass carp.

https://doi.org/10.1371/journal.pone.0266189.s003

(DOCX)

S2 Table. List of novel miRNAs identified in grass carp spleen.

https://doi.org/10.1371/journal.pone.0266189.s004

(DOCX)

S3 Table. The significant differentially expressed miRNAs identified in this study.

https://doi.org/10.1371/journal.pone.0266189.s005

(DOCX)

Acknowledgments

We thank Majorbio (Shanghai, China) for help with bioinformatic analysis and Dr. Xiaojun Liu (College of Animal Science and Technology, Henan Agricultural University) for proofreading the manuscript.

References

  1. 1. Zapata A, Diez B, Cejalvo T, Gutiérrez-de Frías C, Cortés A. Ontogeny of the Immune System of Fish. Fish Shellfish Immunol. 2006; 20(2): 126–36. pmid:15939627
  2. 2. Pearson MP, Stevens ED. Size and Hematological Impact of the Splenic Erythrocyte Reservoir in Rainbow trout,Oncorhynchus Mykiss. Fish Physiol Biochem. 1991; 9(1): 39–50. pmid:24214608
  3. 3. Raida MK, Buchmann K. Development of Adaptive Immunity in Rainbow Trout, Oncorhynchus Mykiss (Walbaum) Surviving an Infection With Yersinia Ruckeri. Fish Shellfish Immunol. 2008; 25(5): 533–41. Epub 2008 Jul 23. pmid:18706505
  4. 4. Chaves-Pozo E, Muñoz P, López-Muñoz A, Pelegrín P, Ayala AG, Mulero V, et al. Early Innate Immune Response and Redistribution of Inflammatory Cells in the Bony Fish Gilthead Seabream Experimentally Infected With Vibrio Anguillarum. Cell Tissue Res. 2005; 320(1): 61–8. pmid:15714279
  5. 5. Espenes A, Press CM, Dannevig BH, Landsverk T. Immune-complex Trapping in the Splenic Ellipsoids of Rainbow Trout (Oncorhynchus Mykiss). Cell Tissue Res. 1995; 282(1): 41–8. pmid:8581925
  6. 6. Whyte SK. The Innate Immune Response of Finfish—A Review of Current Knowledge. Fish Shellfish Immunol, 2007; 23(6): 1127–51. pmid:17980622
  7. 7. Brendolan A, Rosado MM, Carsetti R, Selleri L, Dear TN. Development and Function of the Mammalian Spleen. Bioessays. 2007; 29(2): 166–77. pmid:17226804
  8. 8. Hoffman BG, Williams KL, Tien AH, Lu V, Algara TR, Ting JPY, et al. Identification of Novel Genes and Transcription Factors Involved in Spleen, Thymus and Immunological Development and Function. Genes Immun. 2006; 7(2): 101–12. pmid:16355110
  9. 9. Jones JF. Development of the Spleen. Lymphology. 1983; 16(2): 83–9. pmid:6350736
  10. 10. Langenau DM, Palomero T, Kanki JP, Ferrando AA, Zhou Y, Zon LI, et al. Molecular Cloning and Developmental Expression of Tlx (Hox11) Genes in Zebrafish (Danio Rerio). Mech Dev. 2002; 117(1–2): 243–8. pmid:12204264
  11. 11. Petrie-Hanson L, Ainsworth AJ. Differential Cytochemical Staining Characteristics of Channel Catfish Leukocytes Identify Cell Populations in Lymphoid Organs. Vet Immunol Immunopathol. 2009; 73(2): 129–44.
  12. 12. Santos NM dos, Romano N, Sousa M de, Ellis AE, Rombout JH. Ontogeny of B and T Cells in Sea Bass (Dicentrarchus Labrax, L.). Fish Shellfish Immunol. 2000; 10(7): 583–96. pmid:11081436
  13. 13. Erik J Sontheimer, Richard W Carthew . Silence From Within: Endogenous siRNAs and miRNAs. Cell. 2005; 122(1): 9–12. pmid:16009127
  14. 14. Carrington JC, Ambros V. Role of microRNAs in Plant and Animal Development. Science. 2003; 301(5631): 336–8. pmid:12869753
  15. 15. Bizuayehu TT, Babiak I. MicroRNA in Teleost Fish. Genome Biol Evol. 2014; 6(8): 1911–37. pmid:25053657
  16. 16. Williams AE. Functional Aspects of Animal microRNAs. Cell Mol Life Sci. 2008; 65(4): 545–62. pmid:17965831
  17. 17. Gong W, Huang Y, Xie J, Wang G, Yu D, Zhang K, et al. Identification and Comparative Analysis of the miRNA Expression Profiles From Four Tissues of Micropterus Salmoides Using Deep Sequencing. Genomics. 2018; 110(6): 414–22. pmid:30287401
  18. 18. Guo C, Cui H, Ni S, Yan Y, Qin Q. Comprehensive Identification and Profiling of Host miRNAs in Response to Singapore Grouper Irido virus (SGIV) Infection in Grouper (Epinephelus Coioides). Dev Comp Immunol. 2015; 52(2): 226–35. pmid:26027797
  19. 19. Li G, Zhao Y, Wen L, Liu Z, Yan F, Gao C. Identification and Characterization of MicroRNAs in the Spleen of Common Carp Immune Organ. Journal of Cellular Biochemistry. 2014; 115(10): 1768–78. pmid:24819892
  20. 20. Sha Z, Gong G, Wang S, Lu Y, Wang L, Wang Q, et al. Identification and Characterization of Cynoglossus Semilaevis microRNA Response to Vibrio Anguillarum Infection Through High-Throughput Sequencing. Dev Comp Immunol. 2014; 44(1): 59–69. pmid:24296438
  21. 21. He L, Zhao M, Yu X, Zhao Y, Fu L, Qiao X, et al. MicroRNA-182-3p Negatively Regulates Cytokines Expression by Targeting TLR5M in Orange-Spotted Grouper, Epinephelus Coioides. Fish Shellfish Immunol. 2019; 93: 589–96. pmid:31351112
  22. 22. Guan XL, Zhang BC, Sun L. pol-miR-194a of Japanese Flounder (Paralichthys Olivaceus) Suppresses Type I Interferon Response and Facilitates Edwardsiella Tarda Infection. Fish Shellfish Immunol. 2019; 87: 220–5. pmid:30641186
  23. 23. Xu X, Shen Y, Fu J, Lu L, Li J. De Novo Assembly of the Grass Carp Ctenopharyngodon Idella Transcriptome to Identify miRNA Targets Associated With Motile Aeromonad Septicemia. PLoS One. 2014; 9(11): e112722. pmid:25409340
  24. 24. Gong W, Huang Y, Xie J, Wang G, Yu D, Sun X. Genome-wide Identification and Characterization of Conserved and Novel microRNAs in Grass Carp (Ctenopharyngodon Idella) by Deep Sequencing. Comput Biol Chem. 2017; 68: 92–100. pmid:28282565
  25. 25. Chen L, Huang R, Zhu D, Yang C, He L, Li Y, et al. Deep Sequencing of Small RNAs From 11 Tissues of Grass Carp Ctenopharyngodon Idella and Discovery of Sex-Related microRNAs. J Fish Biol. 2019; 94(1): 132–41. pmid:30471229
  26. 26. He L, Zhang A, Chu P, Li Y, Huang R, Liao L, et al. Deep Illumina Sequencing Reveals Conserved and Novel microRNAs in Grass Carp in Response to Grass Carp Reovirus Infection. BMC Genomics. 2017; 18(1): 195. pmid:28219339
  27. 27. Xu X, Shen Y, Fu J, Lu L, Li J. Next-generation Sequencing Identified microRNAs That Associate With Motile Aeromonad Septicemia in Grass Carp. Fish Shellfish Immunol. 2015; 45(1): 94–103. pmid:25698074
  28. 28. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012; 9(4): 357–9. pmid:22388286
  29. 29. Griffiths-Jones S, Saini HK, van Dongen S, Enright AJ. miRBase: tools for microRNA genomics. Nucleic Acids Res. 2008; 36: D154–158. pmid:17991681
  30. 30. Hofacker IL. Vienna RNA secondary structure server. Nucleic Acids Res. 2003; 31(13): 3429–31. pmid:12824340
  31. 31. Friedländer MR, Mackowiak SD, Li N, Chen W, Rajewsky N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. 2012; 40(1): 37–52. pmid:21911355
  32. 32. Ambros V, Bartel B, Bartel DP, Burge CB, Carrington JC, Chen X, et al. A uniform system for microRNA annotation. RNA. 2003; 9(3): 277–9. pmid:12592000
  33. 33. Wang L, Feng Z, Wang X, Wang X, Zhang X. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010; 26(1): 136–8. pmid:19855105
  34. 34. Storey J. D. The positive false discovery rate: A Bayesian interpretation and the q-value. The Annals of Statistics. 2003; 31: 2013–35.
  35. 35. Li G, Zhao Y, Wang J, Liu B, Sun X, Guo S, et al. Transcriptome Profiling of Developing Spleen Tissue and Discovery of Immune-Related Genes in Grass Carp (Ctenopharyngodon Idella). Fish Shellfish Immunol. 2017; 60: 400–10. pmid:27965162
  36. 36. Enright AJ, John B, Gaul U, Tuschl T, Sander C, Marks DS. MicroRNA targets in Drosophila. Genome Biol. 2003; 5(1): R1. pmid:14709173
  37. 37. Friedman RC, Farh KK, Burge CB, Bartel DP. Most mammalian mRNAs are conserved targets of microRNAs. Genome Res. 2009; 19(1): 92–105. pmid:18955434
  38. 38. Klopfenstein DV, Zhang L, Pedersen BS, Ramírez F, Vesztrocy AW, Naldi A, et al. GOATOOLS: A Python library for Gene Ontology analyses. Sci Rep. 2018; 8(1): 10872. pmid:30022098
  39. 39. Xie C, Mao X, Huang J, Ding Y, Wu J, Dong S, et al. KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic Acids Res. 2011; 39: W316–22. pmid:21715386
  40. 40. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods. 2001; 25(4): 402–8. pmid:11846609
  41. 41. Jiang S. Recent findings regarding let-7 in immunity. Cancer Lett. 2018; 434: 130–1. pmid:30053607
  42. 42. Lee HM, Kim TS, Jo EK. MiR-146 and miR-125 in the regulation of innate immunity and inflammation. BMB Rep. 2016; 49(6): 311–8. pmid:26996343
  43. 43. Zhu Y, Zhang S, Li Z, Wang H, Li Z, Hu Y, et al. miR-125b-5p and miR-99a-5p downregulate human γδ T-cell activation and cytotoxicity. Cell Mol Immunol. 2019; 16(2): 112–5. pmid:29429995
  44. 44. Fordham JB, Naqvi AR, Nares S. Regulation of miR-24, miR-30b, and miR-142-3p during macrophage and dendritic cell differentiation potentiates innate immunity. J Leukoc Biol. 2015; 98(2): 195–207. pmid:25990241
  45. 45. Chu F, Hu Y, Zhou Y, Guo M, Lu J, Zheng W, et al. MicroRNA-126 deficiency enhanced the activation and function of CD4+ T cells by elevating IRS-1 pathway. Clin Exp Immunol. 2018; 191(2): 166–79. pmid:28987000
  46. 46. Ferretti C, La Cava A. miR-126, a new modulator of innate immunity. Cell Mol Immunol. 2014; 11(3): 215–17. pmid:24531618
  47. 47. Lindner SE, Lohmüller M, Kotkamp B, Schuler F, Knust Z, Villunger A, et al. The miR-15 family reinforces the transition from proliferation to differentiation in pre-B cells. EMBO Rep. 2017; 18(9): 1604–17. pmid:28705801
  48. 48. Gagnon JD, Kageyama R, Shehata HM, Fassett MS, Mar DJ, Wigton EJ, et al. miR-15/16 Restrain Memory T Cell Differentiation, Cell Cycle, and Survival. Cell Rep. 2019; 28(8): 2169–81.e4. pmid:31433990
  49. 49. Du H, Cui S, Li Y, Yang G, Wang P, Fikrig E, et al. MiR-221 negatively regulates innate anti-viral response. PLoS One. 2018; 13(8): e0200385. pmid:30089112
  50. 50. Naqvi AR, Fordham JB, Nares S. miR-24, miR-30b, and miR-142-3p regulate phagocytosis in myeloid inflammatory cells. J Immunol. 2015; 194(4): 1916–27. pmid:25601927
  51. 51. Zhou Q, Haupt S, Prots I, Thümmler K, Kremmer E, Lipsky PE, et al. miR-142-3p is involved in CD25+ CD4 T cell proliferation by targeting the expression of glycoprotein A repetitions predominant. J Immunol. 2013; 190(12): 6579–88. pmid:23650616
  52. 52. Lim SP, Ioannou N, Ramsay AG, Darling D, Gäken J, Mufti GJ. miR-181c-BRK1 axis plays a key role in actin cytoskeleton-dependent T cell function. J Leukoc Biol. 2018; 103(5): 855–66. pmid:29656550
  53. 53. Ziętara N, Łyszkiewicz M, Witzlau K, Naumann R, Hurwitz R, Langemeier J, et al. Critical role for miR-181a/b-1 in agonist selection of invariant natural killer T cells. Proc Natl Acad Sci USA. 2013; 110(18): 7407–12. pmid:23589855
  54. 54. Hsu AY, Wang D, Liu S, Lu J, Syahirah R, Bennin DA, et al. Phenotypical microRNA screen reveals a noncanonical role of CDK2 in regulating neutrophil migration. Proc Natl Acad Sci USA. 2019; 116(37): 18561–70. pmid:31451657
  55. 55. Liu HY. Down-regulation of miR-144 after Mycobacterium tuberculosis infection promotes inflammatory factor secretion from macrophages through the Tpl2/ERK pathway. Cell Mol Biol (Noisy-le-grand). 2016; 62(2): 87–93. pmid:26950457
  56. 56. Huang Y, Du KL, Guo PY, Zhao RM, Wang B, Zhao XL, et al. IL-16 regulates macrophage polarization as a target gene of mir-145-3p. Mol Immunol. 2019; 107: 1–9. pmid:30634164
  57. 57. Schulte LN, Westermann AJ, Vogel J. Differential activation and functional specialization of miR-146 and miR-155 in innate immune sensing. Nucleic Acids Res. 2013; 41(1): 542–53. pmid:23143100
  58. 58. Labbaye C, Testa U. The emerging role of MIR-146A in the control of hematopoiesis, immune function and cancer. J Hematol Oncol. 2012; 5: 13. pmid:22453030
  59. 59. Wu MZ, Cheng WC, Chen SF, Nieh S, O’Connor C, Liu CL, et al. miR-25/93 mediates hypoxia-induced immunosuppression by repressing cGAS. Nat Cell Biol. 2017; 19(10): 1286–96. pmid:28920955
  60. 60. Cruz LO, Hashemifar SS, Wu CJ, Cho S, Nguyen DT, Lin LL, et al. Excessive expression of miR-27 impairs Treg-mediated immunological tolerance. J Clin Invest. 2017; 127(2): 530–42. pmid:28067667
  61. 61. Cho S, Wu CJ, Yasuda T, Cruz LO, Khan AA, Lin LL, et al. miR-23~27~24 clusters control effector T cell differentiation and function. J Exp Med. 2016; 213(2): 235–49. pmid:26834155
  62. 62. Zhu QY, Liu Q, Chen JX, Lan K, Ge BX. MicroRNA-101 targets MAPK phosphatase-1 to regulate the activation of MAPKs in macrophages. J Immunol. 2010; 185(12): 7435–42. pmid:21068409
  63. 63. Yang Y, Xu J, Chen H, Fei X, Tang Y, Yan Y, et al. MiR-128-2 inhibits common lymphoid progenitors from developing into progenitor B cells. Oncotarget. 2016; 7(14): 17520–31. pmid:27008703
  64. 64. Zhang T, Zhang Z, Li F, Ping Y, Qin G, Zhang C, et al. miR-143 Regulates Memory T Cell Differentiation by Reprogramming T Cell Metabolism. J Immunol. 2018; 201(7): 2165–75. pmid:30150287
  65. 65. Fang Y, Xu XY, Shen Y, Li J. miR-148 targets CiGadd45ba and CiGadd45bb to modulate the inflammatory response to bacterial infection in grass carp. Dev Comp Immunol. 2020; 106: 103611. pmid:31953153
  66. 66. Liu X, Zhan Z, Xu L, Ma F, Li D, Guo Z, et al. MicroRNA-148/152 impair innate response and antigen presentation of TLR-triggered dendritic cells by targeting CaMKIIα. J Immunol. 2010; 185(12): 7244–51. pmid:21068402
  67. 67. Chen Z, Stelekati E, Kurachi M, Yu S, Cai Z, Manne S, et al. miR-150 Regulates Memory CD8 T Cell Differentiation via c-Myb. Cell Rep. 2017; 20(11): 2584–97. pmid:28903040
  68. 68. Sang W, Sun C, Zhang C, Zhang D, Wang Y, Xu L, et al. MicroRNA-150 negatively regulates the function of CD4(+) T cells through AKT3/Bim signaling pathway. Cell Immunol. 2016; 306–307: 35–40. pmid:27329362
  69. 69. Bezman NA, Chakraborty T, Bender T, Lanier LL. miR-150 regulates the development of NK and iNKT cells. J Exp Med. 2011; 208(13): 2717–31. pmid:22124110
  70. 70. Chu Q, Xu T. miR-192 targeting IL-1RI regulates the immune response in miiuy croaker after pathogen infection in vitro and in vivo. Fish Shellfish Immunol. 2016; 54: 537–43. pmid:27164215
  71. 71. Bi D, Cui J, Chu Q, Xu T. MicroRNA-21 contributes to suppress cytokines production by targeting TLR28 in teleost fish. Mol Immunol. 2017; 83: 107–14. pmid:28129531
  72. 72. Wang L, He L, Zhang R, Liu X, Ren Y, Liu Z, et al. Regulation of T lymphocyte activation by microRNA-21. Mol Immunol. 2014; 59(2): 163–71. pmid:24631982
  73. 73. Ruan Q, Wang P, Wang T, Qi J, Wei M, Wang S, et al. MicroRNA-21 regulates T-cell apoptosis by directly targeting the tumor suppressor gene Tipe2. Cell Death Dis. 2014; 5(2): e1095. pmid:24577093
  74. 74. Xu T, Chu Q, Cui J. Rhabdovirus-Inducible MicroRNA-210 Modulates Antiviral Innate Immune Response via Targeting STING/MITA in Fish. J Immunol. 2018; 201(3): 982–94. pmid:29967101
  75. 75. Li C, Zhao M, Zhang C, Zhang W, Zhao X, Duan X, et al. miR210 modulates respiratory burst in Apostichopus japonicus coelomocytes via targeting Toll-like receptor. Dev Comp Immunol. 2016; 65: 377–81. pmid:27545641
  76. 76. Mok Y, Schwierzeck V, Thomas DC, Vigorito E, Rayner TF, Jarvis LB, et al. MiR-210 is induced by Oct-2, regulates B cells, and inhibits autoantibody production. J Immunol. 2013; 191(6): 3037–48. pmid:23960236
  77. 77. Yuan X, Berg N, Lee JW, Le TT, Neudecker V, Jing N, et al. MicroRNA miR-223 as regulator of innate immunity. J Leukoc Biol. 2018; 104(3): 515–24. pmid:29969525
  78. 78. Chen L, Song Y, He L, Wan X, Lai L, Dai F, et al. MicroRNA-223 Promotes Type I Interferon Production in Antiviral Innate Immunity by Targeting Forkhead Box Protein O3 (FOXO3). J Biol Chem. 2016; 291(28): 14706–16. pmid:27226534
  79. 79. Steiner DF, Thomas MF, Hu JK, Yang Z, Babiarz JE, Allen CD, et al. MicroRNA-29 regulates T-box transcription factors and interferon-γ production in helper T cells. Immunity. 2011; 35(2): 169–81. pmid:21820330
  80. 80. Ma F, Xu S, Liu X, Zhang Q, Xu X, Liu M, et al. The microRNA miR-29 controls innate and adaptive immune responses to intracellular bacterial infection by targeting interferon-γ. Nat Immunol. 2011; 12(9): 861–9. pmid:21785411
  81. 81. Chapman LM, Ture SK, Field DJ, Morrell CN. miR-451 limits CD4+ T cell proliferative responses to infection in mice. Immunol Res. 2017; 65(4): 828–40. pmid:28378118
  82. 82. Reza AMMT Cho SK, Choi Y, Hong K, Kim J. Microarray Profiling of miRNA and mRNA Expression in Rag2 Knockout and Wild-Type Mouse Spleens. Sci Data. 2018; 5: 170199. pmid:29313843
  83. 83. Wang A, Liu F, Chen S, Wang M, Jia R, Zhu D, et al. Transcriptome Analysis and Identification of Differentially Expressed Transcripts of Immune-Related Genes in Spleen of Gosling and Adult Goose. Int J Mol Sci. 2015; 16(9): 22904–26. pmid:26402676
  84. 84. Che T, Li D, Jin L, Fu Y, Liu Y, Liu P, et al. Long Non-Coding RNAs and mRNAs Profiling During Spleen Development in Pig. PLoS One. 2018; 13(3): e0193552. pmid:29538394
  85. 85. Peng R, Liu Y, Cai Z, Shen F, Chen J, Hou R, et al. Characterization and Analysis of Whole Transcriptome of Giant Panda Spleens: Implying Critical Roles of Long Non-Coding RNAs in Immunity. Cell Physiol Biochem. 2018; 46(3): 1065–77. pmid:29669315
  86. 86. Li G, Zhao L, Liu Z, Gao C, Yan F, Liu B, et al. De novo assembly and characterization of the spleen transcriptome of common carp (Cyprinus carpio) using Illumina paired-end sequencing. Fish & Shellfish Immunology. 2015; 44(2): 420–9. pmid:25804493
  87. 87. Li G, Zhao Y, Guo S, Liu B, Chen Y, Sun X, et al. Comparative analysis of spleen transcriptome detects differences in evolutionary adaptation of immune defense functions in bighead carp and silver carp. Fish and Shellfish Immunology. 2019; 84: 148–57. pmid:30287346
  88. 88. Huang Y, Ren HT, Xiong JL, Gao XC, Sun XH. Identification and Characterization of Known and Novel microRNAs in Three Tissues of Chinese Giant Salamander Base on Deep Sequencing Approach. Genomics. 2017; 109(3–4): 258–64. pmid:28476431
  89. 89. Yi S, Lu D, Peng W, Wang T, Zhang Y, Lin H. Differential Expression Profiling of Spleen microRNAs in Response to Two Distinct Type II Interferons in Tetraodon Nigroviridis. PLoS One. 2014; 9(5): e96336. pmid:24800866
  90. 90. Cao Y, Wang D, Li S, Zhao J, Xu L, Liu H, et al. A Transcriptome Analysis Focusing on Splenic Immune-Related mciroRNAs of Rainbow Trout Upon Aeromonas Salmonicida Subsp. Salmonicida Infection. Fish Shellfish Immunol. 2019; 91: 350–7. pmid:31128295
  91. 91. Wheeler BM, Heimberg AM, Moy VN, Sperling EA, Holstein TW, Heber S, et al. The deep evolution of metazoan microRNAs. Evol Dev. 2009; 11(1): 50–68. pmid:19196333
  92. 92. Bartel DP. MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004; 116(2): 281–97. pmid:14744438
  93. 93. Li G, Li Y, Li X, Ning X, Li M, Yang G. MicroRNA identity and abundance in developing swine adipose tissue as determined by Solexa sequencing. J Cell Biochem. 2011; 112(5): 1318–28. pmid:21312241
  94. 94. Landgraf P, Rusu M, Sheridan R, Sewer A, Iovino N, Aravin A, et al. A mammalian microRNA expression atlas based on small RNA library sequencing. Cell. 2007; 129(7): 1401–14. pmid:17604727
  95. 95. Nielsen M, Hansen JH, Hedegaard J, Nielsen RO, Panitz F, Bendixen C, et al. MicroRNA identity and abundance in porcine skeletal muscles determined by deep sequencing. Anim Genet. 2010; 41(2): 159–68. pmid:19917043
  96. 96. Ebhardt HA, Tsang HH, Dai DC, Liu Y, Bostan B, Fahlman RP. Meta-analysis of small RNA-sequencing errors reveals ubiquitous post-transcriptional RNA modifications. Nucleic Acids Res. 2009; 37(8): 2461–70. pmid:19255090
  97. 97. Kadri S, Hinman VF, Benos PV. RNA deep sequencing reveals differential microRNA expression during development of sea urchin and sea star. PLoS One. 2011; 6(12): e29217. pmid:22216218
  98. 98. Jazdzewski K, Liyanarachchi S, Swierniak M, Pachucki J, Ringel MD, Jarzab B, et al. Polymorphic mature microRNAs from passenger strand of pre-miR-146a contribute to thyroid cancer. Proc Natl Acad Sci USA. 2009; 106(5): 1502–5. pmid:19164563
  99. 99. Lenti E, Farinello D, Yokoyama KK, Penkov D, Castagnaro L, Lavorgna G, et al. Transcription Factor TLX1 Controls Retinoic Acid Signaling to Ensure Spleen Development. J Clin Invest. 2016; 126(7): 2452–64. pmid:27214556
  100. 100. Su Z, Si W, Li L, Zhou B, Li X, Xu Y, et al. MiR-144 regulates hematopoiesis and vascular development by targeting meis1 during zebrafish development. Int J Biochem Cell Biol. 2014; 49: 53–63 pmid:24448023
  101. 101. Xu P, Palmer LE, Lechauve C, Zhao G, Yao Y, Luan J, et al. Regulation of gene expression by miR-144/451 during mouse erythropoiesis. Blood. 2019; 133(23)::2518–28. pmid:30971389