Identification and characterization of microRNAs in the intestinal tissues of sheep (Ovis aries)

Lei Hou; Zhibin Ji; Guizhi Wang; Jin Wang; Tianle Chao; Jianmin Wang

doi:10.1371/journal.pone.0193371

Abstract

Sheep are small ruminants, and their long intestines exhibit high digestive and absorptive capacity in many different rearing conditions; however, the genetic bases of this characteristic remains unclear. MicroRNAs (miRNAs) play a major role in maintaining both intestinal morphological structure as well as in regulating the physiological functions of this organ. However, no study has reported on the miRNA expression profile in the intestinal tissues of sheep. Here, we analyzed and identified the miRNA expression profile of three different intestinal tissues (i.e., duodenum, cecum, and colon) of sheep (Ovis aries) using high-throughput sequencing and bioinformatic methods. In total, 106 known miRNAs were identified, 458 conserved miRNAs were detected, 192 unannotated novel miRNAs were predicted, and 195 differentially expressed miRNAs were found between the different tissues. Additionally, 3,437 candidate target genes were predicted, and 17 non-redundant significantly enriched GO terms were identified using enrichment analysis. A total of 99 candidate target genes were found to significantly enriched in 4 KEGG biological pathways. A combined regulatory network was constructed based on 92 metabolism-related candidate target genes and 65 differentially expressed miRNAs, among which 7 miRNAs were identified as hub miRNAs. Via these mechanisms, miRNAs may play a role in maintaining intestinal homeostasis and metabolism. This study helps to further explain the mechanisms that underlie differences in tissue morphology and function in three intestinal segments of sheep.

Citation: Hou L, Ji Z, Wang G, Wang J, Chao T, Wang J (2018) Identification and characterization of microRNAs in the intestinal tissues of sheep (Ovis aries). PLoS ONE 13(2): e0193371. https://doi.org/10.1371/journal.pone.0193371

Editor: Michael Schubert, Laboratoire de Biologie du Développement de Villefranche-sur-Mer, FRANCE

Received: March 29, 2017; Accepted: February 9, 2018; Published: February 28, 2018

Copyright: © 2018 Hou et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: All relevant data are within the paper and its Supporting Information files.

Funding: This work was financially supported by the Funds of Shandong “Double Tops” Program (No. SYL2017YSTD12), Shandong Provincial Modern Agriculture Industry Technology System (no. SDAIT-0901101), National Key Technology Support Program (2015BAD03B05) and Study on Super-high Reproductive New Breed and Supporting technological system of Small-Tailed Han Sheep (no. SDLZ2015131x). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Competing interests: The authors have declared that no competing interests exist.

Introduction

Ruminants have a unique digestive system compared with monogastric animals. Specifically, the gastrointestinal tract of ruminants harbors a large number of microflora that can transform and utilize dietary crude fiber [1]. Conventional research has primarily focused on the complex forestomach of ruminants, while there are few studies on the function of the small and large intestines [2]. Sheep are small ruminants whose small intestine is long and thin, being 27-fold longer than the body of the animal. The volume of the small intestine constitutes 22% of the total volume of the digestive tract and contains the most digestive enzymes in sheep. The duodenum is the most proximal segment and is characterized by the densest villi, the largest diameter, and the deepest position in the small intestine. The duodenum receives pancreatic juice from the pancreas and bile from the liver, aiding the intestines in the digestion of starch, fat, and protein. Following the rumen, the cecum and proximal colon of ruminants are the second largest sites of fermentation. Large amounts of cell wall carbohydrates, cellulose, and hemicellulose enter the cecum and colon, wherein they are fermented and utilized. The produced volatile fatty acids are rapidly taken up by the cecum and colon. The energy conversion efficiency of these two organs has been reported to be higher compared with that of the rumen [3, 4, 5]. The intestines of ruminants perform unique physiological functions, indicating that molecular regulation in the intestines is far more complicated than previously thought.

In recent years, the functions of miRNAs in the intestines have been extensively explored, with many novel miRNAs being identified. MiRNAs have been found to participate widely in the maintenance of intestinal homeostasis and the regulation of metabolic homeostasis [6, 7]. Moreover, miRNAs play an important regulatory role in intestinal mucosal immunity, stress resistance, the development and maintenance of intestinal morphological structure, as well as the regulation of various intestinal functions [7–10]. Studies have found that miRNAs are widely expressed in the intestinal tissues of pigs [11–13], mice [14, 15], and cattle [16, 17]. Some miRNAs are only found in particular intestinal segments, suggesting that miRNAs may play important roles in particular intestinal segments. The miRBase21.0 database lists a total of 793 bovine miRNAs, while only 153 miRNAs are identified in sheep. There have been previous studies of miRNAs in different sheep tissues, such as muscle [18], ovary [19], fat [20,21], wool follicle [22], and heart [23]. Currently, there is a lack of research regarding the miRNAs expression profile in the duodenum, cecum, and colon of sheep. Thorough and in-depth studies of the regulation of miRNAs expression in the intestines of sheep have vital scientific implications for revealing the regulatory mechanisms that affect the intestinal digestion and metabolism in sheep.

In the present study, we compared the miRNAs expression profile between the duodenum, cecum, and colon of sheep (Ovis aries) using Illumina-Solexa high-throughput sequencing technology. These results provide a reference for understanding the regulatory mechanism of digestion, absorption, and metabolism in different intestinal regions of sheep and for the generation of high-yield, high-efficiency breeds.

Materials and methods

Ethics statement of experimental animals

All of the animal experiments in this study were approved by the Ethics Committee of Animal Protection and Use at Shandong Agricultural University (License No.: 2004006). The experiments were performed in accordance with the “Experimental Animal Guidelines” put forth by the Ministry of Science and Technology (Beijing, China). All of the surgical procedures followed the recommendations proposed by the European Commission in 1997. Animals used in this study were the same as those used in our previous study [24], and were slaughtered in the same way.

Sample collection and total RNA extraction

The experimental animals were three 11-month-old purebred ewes (Ovis aries) of similar weight. The animals were obtained from the same farm (Linqv Huanong Sheep Farm, Weifang, Shandong, China) and reared in consistent environments, i.e., the sheep were fed in stables, in natural light and without temperature control, free access to water and food.

The animals were slaughtered quickly to collect approximately intestinal tissues. All animals used in this study were handled in strict accordance with good clinical practices and all efforts were made to minimize suffering. Intestinal tissues were taken from the intestinal epithelium of the duodenum, the cecum, and the proximal colon immediately after slaughter. The samples were preserved in liquid nitrogen until use. Total RNA was extracted from the intestinal tissues using TRIzol. The integrity and concentration of the extracted RNA were determined using an Agilent 2100 Bioanalyzer. Samples with a RIN score >7.5 were used for library construction.

Library construction and sequencing

The libraries were constructed following the mRNA preparation procedures for Illumina sequencer samples, with appropriate modifications. Specifically, total RNA from the duodenum, cecum and colon of three sheep was separately pooled, Illumina TruSeq samll RNA Sample Prep Kit (Cat#RS-200-0012) was used with 1 ug of total RNA for the construction of sequencing libraries. three libraries were constructed (Duodenum, DU; Cecum, CE; and Colon, CO), which were prepared for sequencing using standard Illumina protocols. All of the sequencing was performed in the Beijing Genomics Institute (BGI, Shenzhen, China). Sequencing data had been submitted to GEO database (GSE100649).

Sequence data analysis and miRNAs identification

Raw sequencing image data were converted to raw reads via base calling. The following reads were removed: low-quality reads, reads with contaminated 5'-adapters, reads lacking 3'-adapters, reads lacking an insertion sequence, reads containing polyA sequences, and small fragments (<18 nt). A set of clean reads was ultimately obtained for analysis. Generally, the whole workflow is shown in S1 Fig.

To analyze sequence distribution and expression, the clean reads were mapped to the genome of Ovis aries (assembly Oar_v4.0, ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/298/735/GCF_000298735.2_Oar_v4.0/GCF_000298735.2_Oar_v4.0_genomic.fna.gz) using SOAP aligner v2.21 [25]. To identify the known miRNAs, the clean reads were compared with known mammalian precursor and mature miRNAs in the miRBase21.0 database [26] using BLAST. The criteria were as follows: Clean reads can be completely aligned to precursor (no mismatch), then compared with the mature miRNA, which allow dislocation, but at least there were 16 base coverage without mismatch. Those miRNAs whose reads number less than 10 in each sample were removed as false positive. The clean reads were compared with the Genbank (ftp://ftp.ncbi.nlm.nih.gov/genbank/) and Rfam v11.0 (http://rfam.janelia.org/) [27] databases to remove rRNA, scRNA, snoRNA, snRNA, and tRNA sequences. The clean reads that were not matched in the databases but were matched to the antisense strand of exons, introns, and intergenic regions in the genome were used to predict novel miRNAs with miReap (http://mireap.sourceforge.net/). Using miReap, we predicted novel miRNA precursor hairpins with a minimum free energy ≤ 18 (kcal/mol), while all other parameters were settled as default. Then, retained clean reads were compared with predicted novel miRNA precursors using a same method and parameters of our known miRNA identification, and all successfully mapped reads were automatic constructed into novel miRNAs.

For the comparison analysis of novel miRNAs, we collected 2,381 sheep unannotated novel miRNAs from other 8 sheep miRNA sequencing studies [20, 23, 28–33]. With at least 18 base coverage and not more than one mismatch, and e-value < 1e-5, the comparison analysis of novel miRNAs was performed with BLAST.

Differential expression of miRNAs and prediction of target genes

The expression levels of known and novel miRNAs were normalized to the same order of magnitude to get the expression of transcript per million (TPM).

Normalization formula: Normalized expression = Actual miRNA count/Total count of clean reads*1000000.

The standardized results were then used to calculate fold-change and P-value to analyze the fold differences between different samples.

Fold-change formula: Fold_change = log2(treatment/control)

P-value formula:

The P-value was calculated based on Poisson distribution to determine the significance of any difference in miRNA expression between different samples. And the n, the Bonferroni corrected P-value was calculated. The significantly differentially expressed miRNAs were declared at a fold change ≥ 2 and a Bonferroni corrected P-value <0.01.

Differentially expressed miRNAs were identified and their predicted target genes were predicted using miRanda [34], and TargetScan v7.0 [35] software. Gene sequences used for target prediction analysis were downloaded from (ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/298/735/GCF_000298735.2_Oar_v4.0/GCF_000298735.2_Oar_v4.0_rna.fna.gz). For the TargetScan v7.0, a seed match ≥ 7 and context+ score percentile ≥ 99 was settled as the threshold. And for the miRanda, a free energy ≤ -20 (kcal/mol) was chosen as a cutoff. The intersection of these prediction results was taken as the set of candidate target genes.

GO annotation and KEGG pathway analysis

To generate a comprehensive description of the biological characteristics of candidate target genes, we mapped all of the candidate target genes to the Gene Ontology (GO) [36] database using DAVID v6.8 online software with the FDR <0.1, and background of genes with FPKM>1 from previous report [37], other parameters are in accordance with the DAVID default settings [38]. The GO terms were subjected to classification and enrichment analysis regarding molecular function, cellular component, and biological process. Significantly enriched GO terms were calculated using the hypergeometric test with all genes in the reference genome as the background. In addition, we mapped all of the candidate target genes to the Kyoto Encyclopedia of Genes and Genomes (KEGG, http://www.genome.jp/kegg/) [39] database to obtain significantly enriched metabolism and signal transduction pathways. Significant enrichment was identified at FDR <0.1.

Regulatory network construction

The interactions between filtered candidate target genes (filter criteria: combined score ≥0.4) were explored using the STRING 10.0 databas [40]. An interactive regulatory network between miRNAs and candidate target genes was generated using Cytoscape v3.4 [41].

Real-time quantitative RT-PCR validation

We randomly selected 14 differentially expressed miRNAs (4 known miRNAs and 10 novel miRNAs) for quantitative validation by real-time quantitative PCR. The primers were designed with mature miRNA sequences as the templates. Base modification was performed at the 5'-end if necessary (S1 Table). U6 was chosen as the reference gene.

Each miRNA was tested using three biological replicates, and each biological replicate had three technical replicates. The expression of each miRNA was calculated based on a standard curve using the 2^-ΔΔCT method, and the expression is presented as 2^-ΔΔCT ± SEM. The significance of difference in miRNA expression between three intestinal segments was tested using one-way ANOVA (P ≤0.05).

Results

Sequence data from three segments of intestinal tissues

In the DU, CE, and CO libraries, we obtained 28,493,897; 27,836,933, and 26,057,470 high-quality reads, respectively. After removal of reads lacking 3'-adapters, reads with 5'-contamination, reads <18 nt, reads lacking insertion sequences, and reads containing polyA sequences, we obtained 27,498,891, 26,490,229, and 25,351,117 clean reads in the three libraries (96.51%, 95.16%, and 97.29% of the high quality reads, respectively, Table 1). The length of sRNAs obtained by sequencing generally ranged from 19 to 24 nt, accounting for 91.62%, 91.76%, and 93% of the total reads in the DU, CE, and CO libraries, respectively. Among these reads, those of 22 nt were the most abundant, constituting 37.26%, 33.62%, and 33.48% of the DU, CE, and CO libraries, respectively. Reads of 21 and 23 nt in length were the next most abundant. The distribution of read length was in agreement with the length distribution pattern of miRNA sequences in mammals (Fig 1).

Download:

Fig 1. Length distribution of small RNAs in three libraries.

Sequence length distribution of clean reads based on the abundance and distinct sequences; the most abundant size class was 22 nt, followed by 21 nt, 23 nt and 20nt.

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

Download:

Table 1. Distribution of total sequence for small RNAs in small-tailed han sheep.

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

In the DU, CE, and CO libraries, the total numbers of clean reads mapped to the genome were 16,733,785, 17,005,544, and 16,567,397, respectively (60.85%, 64.20%, and 65.35% of the total reads, respectively, Table 1). Chromosomes Chr5, Chr2, and Chr11 harbored the largest number of reads, whereas Chr25 and Chr26 harbored the smallest number of reads (S2 Table).

To obtain the components contained in clean reads, we classified and annotated the sRNAs. In the DU, CE, and CO libraries, a total of 13,432,974, 14,651,426, and 12,438,925 reads, respectively, were annotated as various types of RNAs (rRNAs, tRNAs, snRNAs, snoRNAs, and miRNAs). Annotated repeat reads represented 123,851, 133,191, and 171,894 reads in the DU, CE, and CO libraries, respectively. Moreover, there were 10,310,790, 10,565,745, and 12,137,476 unannotated reads in the three DU, CE, and CO libraries. These unannotated reads were used in the subsequent prediction of novel miRNAs (Fig 2).

Download:

Fig 2. Classification of small RNAs in three libraries.

The clean reads were annotated and classified as miRNA, rRNA, tRNA, snRNA and snoRNA in GenBank and Rfam databases, and partial reads were not annotated.

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

The identification of known miRNAs and novel miRNAs

To identify known miRNAs, we aligned the clean reads against 106 precursor miRNAs and 153 mature miRNAs of sheep in the miRBase 21.0 database. A total of 106 mature miRNAs were detected in the three libraries (S5 Table). These mature miRNAs were derived from 92 precursor miRNAs and distributed in 51 miRNA families. Of the clean reads, 40%, 50.4%, and 49.3% were matched to precursor miRNAs. But these miRNAs only accounted for 0.18%, 0.27%, and 0.18% of total unique sRNAs in the DU, CE, and CO libraries, respectively. In the DU, CE, and CO libraries, miR-143 had the largest number of reads, while the miR-21, miR-148a, miR-26, and let-7 families were the next most abundant.

The sRNAs that were unannotated to sequences of known miRNAs in the miRBase21.0 database but mapped to the antisense strand of exons, introns, and intergenic regions were selected for predicting novel miRNAs using miReap. The reads with predicted hairpin structures were mapped to known mammalian precursor and mature miRNAs in the miRBase21.0 database, and the remaining unmapped but hairpin-structured reads were considered potential new miRNAs without known miRNA annotations. In total, 458 miRNAs were conserved with other species, 192 miRNAs were identified as unannotated novel miRNAs in the three libraries (S5 Table). Moreover, to get further information on new miRNAs, we collected sheep unannotated novel miRNAs from other studies [20, 23, 28–33], and a comparison analysis has been performed. As a result, 71 novel miRNAs obtained acceptable results with other studies (S6 Table).

Analysis of differentially expressed miRNAs and prediction of target genes

A total of 195 known, conserved and novel differentially expressed miRNAs were identified in the DU, CE, and CO libraries (S7 Table). There were 131 known and novel differentially expressed miRNAs between the DU and CE libraries, with 81 being expressed at higher levels in the DU library and 50 being expressed more highly in the CE library. 123 known and novel differentially expressed miRNAs between the CO and DU libraries, with 49 being expressed more highly in the CO library and 74 being expressed at higher levels in the DU library. Lastly, there were 61 known and novel miRNAs that were differentially expressed between the CO and CE libraries, with 29 being expressed at higher levels in the CO library and 32 being expressed more highly in the CE library.

The 195 known, conserved and novel differentially expressed miRNAs were used to predict target genes using miRanda, and Targetscan. A total of 3,437 non-redundant high-confidence candidate target genes were obtained from the intersection of the obtained gene sets (S8 Table). Newly predicted novel_miR_258 had the largest number of predicted target genes, at 178. There were 9 miRNAs with over 100 predicted target genes, and 97 miRNAs with 10–100 predicted target genes, while 74 miRNAs with less than 10 predicted target genes. LOC101109384 and USP54 were the candidate target genes with the largest number of miRNAs (8 miRNAs), followed by DPF3, GFM2, MAP2, and DOCK5 (6 miRNAs).

The results of qRT-PCR validation for differentially expressed miRNAs showed that the expression trends of 12 miRNAs were consistent with the sequencing results. Only two miRNAs (oar-miR-200b and bta-miR-2478) showed different trends between the qPCR and sequencing data (Fig 3). Specifically, the sequencing results showed high-level expression of oar-miR-200b in the colon, whereas the quantitative PCR results showed relatively high expression levels of this miRNA in the duodenum. Meanwhile, the sequencing results showed the highest expression level of bta-miR-2478 in the cecum, whereas the qPCR results showed relatively high expression levels of this miRNA in the duodenum.

Download:

Fig 3. Validation of differentially expressed miRNAs by qRT-PCR.

A. Result of RNAseq. B. Result of qRT-PCR. oar-miR-10a, oar-miR-148a, oar-miR-200b oar-miR-133 and oar-miR-194 were known sheep miRNAs, cgr-miR-142-5p, chi-miR-215-5p, bta-miR-1388-5p, age-miR-15b, bta-miR-2478, and hsa-miR-141-3p were conserved miRNAs, novel_miR_17, novel_miR_147 and novel_miR_217 were novel miRNA candidates. The relative quantification of expression was calculated using the 2^-△△CT method after the threshold cycle (Ct) and was normalized with the Ct of U6. The relative expression levels were presented as the 2^-△△CT means ± SE. The error bars indicate the standard error of the 2^-△△CT mean values. * represents p <0.05, ** represents p <0.01.

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

GO enrichment for candidate target genes of differentially expressed miRNAs

To better illustrate the function of differentially expressed miRNAs and reveal their associated biological functions in different intestinal segments, we conducted GO enrichment analysis for 3,437 candidate target genes.

The 2,709 candidate target genes that corresponded to the differentially expressed miRNAs between the DU and CE samples were significantly enriched in 16 GO terms (FDR < 0.1). Specifically, 126 candidate target genes were significantly enriched in 2 biological processes, primarily “negative regulation of transcription from RNA polymerase II promoter” and “positive regulation of transcription from RNA polymerase II promoter”; 1,020 candidate target genes were significantly enriched in 8 cellular component terms, with most in “exosome” and “cytoplasm”; and 469 candidate target genes were significantly enriched in 6 molecular function terms, primarily “ATP binding” and “poly(A) RNA binding”.

The 2,443 candidate target genes corresponding to miRNAs that were differentially expressed between the DU and CO libraries were significantly enriched in 15 GO terms (FDR < 0.1). Specifically, 113 candidate target genes were significantly enriched in 2 biological processes, primarily “negative regulation of transcription from RNA polymerase II promoter” and “positive regulation of transcription from RNA polymerase II promoter”, 932 genes candidate target genes were significantly enriched in 8 cellular component terms, mostly “cytoplasm” and “extracellular exosome”; 408 candidate target genes were significantly annotated to 5 molecular function terms, with most in “ATP binding” and “poly(A) RNA binding”.

The 1,076 candidate target genes that were targeted by differentially expressed miRNAs between the CE and CO libraries were significantly enriched in 6 GO terms (FDR < 0.1). Specifically, 257 candidate target genes were significantly enriched in 4 cell components, primarily “cytoplasm” and “nucleoplasm”; 31 candidate target genes were significantly enriched in 1 biological processes, “negative regulation of transcription from RNA polymerase II promoter”; and 90 candidate target genes were significantly enriched in 1 molecular function term, “ATP binding”. A subset of the GO annotation results are provided in Fig 4 and S9 Table.

Download:

Fig 4. GO classification of 3,437 candidate target genes for 195 differentially expressed miRNAs.

The figure shows GO enrichment for the candidate target genes in molecular function, cellular component and biological processes.

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

KEGG analysis and regulatory network construction

To further understand the biological pathways that may be regulated by miRNAs, we performed KEGG pathway enrichment analysis using all 3,437 candidate target genes. KEGG pathway analysis showed that a total of 99 candidate target genes (DU vs CE, 44 genes; DU vs CO, 33 genes; CE vs CO, 29 genes) were significantly annotated to 4 pathways (DU vs CE, “MAPK signaling pathway”; DU vs CO, “Bacterial invasion of epithelial cells” and “Pancreatic cancer”; CE vs CO, “PI3K-Akt signaling pathway”). (S10 Table).

In this study, we used the significantly enriched pathways as the base point to reveal the network regulatory relationship between miRNAs and their target genes. Using protein–protein interaction analysis in the STRING database, we identified 92 candidate target genes with interactions, which were regulated by 65 differentially expressed miRNAs. We cross-referenced these targets with the protein–protein interaction network using Cytoscape 3.4 and constructed a regulatory network of 157 nodes and 615 interactions (Fig 5). In this network, five genes were identified as hub genes (top 5% interaction degree), PIK3CB, VEGFA, MAPK14, PIK3R1, and CASP3, the 5 hub genes may perform important functions in the constructed regulatory network. Therefore, their corresponding predicted regulatory miRNAs were identified as hub miRNAs (PIK3CB: bta-miR-125a and bta-miR-206; VEGFA: novel-miR-70 and novel-miR-378; MAPK14: oar-miR-133; PIK3R1: bta-miR-1343-3p; CASP3: novel_mir_174).

Download:

Fig 5. Regulatory network linking the differential expressed miRNAs and their predicted target genes.

The node size was decided on the basis of the interaction degree value. Blue triangle node: Protein coding gene. Red circular node: miRNA.

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

Discussion

miRNAs are a major class of post-transcriptional regulatory factors that have been heavily studied in recent years. miRNAs are extensively involved in various biological processes, such as development and metabolism. Despite extensive research on miRNAs in intestinal diseases and homeostasis, little work has been dedicated to the intestines of sheep, which exhibit a unique digestive system. Only one recent study examined miRNAs in sheep and investigated the miRNA expression profile of sheep intestines that were infected with parasites [42]. In the present study, based on sequencing of miRNA libraries from the duodenum, cecum, and colon of Ovis aries, we found that miRNAs were primarily concentrated on Chr18, while Chr5 had the highest distribution frequency of reads (S5 Table). Different chromosomes or chromosomal regions can be enriched in genes or other regulatory regions that govern certain functions. Analyzing the distribution of miRNAs on different chromosomes has great implications for investigating the regulatory functions of miRNAs. In the present study, only three known sheep miRNAs were detected on Chr5, which exhibited the highest distribution frequency of reads and accounted for more than 41.9% of the total reads. Our result indicates that there is indeed a bias in terms of the chromosomal distribution of differentially expressed miRNAs.

A typical feature of miRNAs is their strict spatiotemporal expression pattern. In this study, we identified 998,667, 731,745, and 1,035,225 unique sRNAs in the DU, CE, and CO libraries, respectively. There were large differences in the types of miRNAs between the different libraries. This result suggests that different miRNAs participate in the regulation of physiological functions of different intestinal segments by regulating their target genes, thereby mediating the unique physiological functions of these tissues. In the three libraries, miR-143 showed the highest expression level, with 4,606,500, 7,297,109, and 6,614,344 reads in the DU, CE, and CO libraries. According to previous research, miR-143 primarily performs important functions as a tumor suppressor. In the intestines, miR-143 primarily functions in colorectal cancer and post-cancer repair of intestinal epithelial injury [43] and is likely involved in nutrient digestion and absorption. Other highly expressed miRNAs, e.g., miR-10a and miR-10b, were upregulated in the CE and CO libraries. A previous study reported that miR-10a and miR-10b target Hox transcripts in several species and may play a major role in development [44]. Additionally, miR-10b is involved in cell cycle regulation and cancer recognition by the immune system [45, 46]. Taken together, the above information indicates that miR-10a and miR-10b may be related to the rapid regeneration of the intestinal mucosal epithelium in sheep. Moreover, bta-miR-215, hsa-miR-194-5p and bta-miR-192 were differentially expressed between the three types of intestinal tissues. The variation in miRNA expression levels between different intestinal segments may indicate the important regulatory roles of miRNAs in the development and normal maintenance of intestinal morphology and structure, as well as the performance of various tissue-specific functions.

Studies on miRNAs and their target genes are important for understanding the mechanisms of miRNAs in various biological processes [47, 48]. A total of 4,422 candidate target genes were obtained in the present study. Based on the observed regulatory relationship between miRNAs and their target genes, one miRNA might regulate multiple target genes; in turn, a given target gene could be regulated by numerous miRNAs. This result is consistent with those of previous studies [49].

The GO annotation and KEGG pathway enrichment analysis for target genes can provide a comprehensive description of the attributes of genes and gene products in organisms, helping to understand the biological functions and associated signaling pathways of candidate target genes. In the GO annotation, DU vs CE and DU vs CO showed almost exactly the same enrichment results, while the GO analysis result of CE vs CO is much more different. Such results are consistent with functional differences between intestinal segments. In the KEGG analysis, only 4 pathways achieved significant enrichment level. Although the three differentially comparison groups obtained a completely different KEGG enrichment result, (DU vs CE: MAPK signaling pathway; DU vs CO: Bacterial invasion of epithelial cells, Pancreatic cancer; CE vs CO: PI3K-Akt signaling pathway), their related functions are highly similar. All of the four pathways are functionally related to cell proliferation and apoptosis, which indicates the differentially expressed miRNAs might affect cell proliferation and apoptosis of different intestinal segments via targeted regulation of related genes, resulting in different intestinal functions.

The regulation of gene expression is extremely complex. To explore the intricate regulatory relationships between differentially expressed miRNAs, candidate target genes, and gene products, we constructed an interactive regulatory network comprising 65 miRNAs and 92 candidate target genes by exploring gene interactions in the STRING database. In this combined regulatory network, gene PIK3CB, VEGFA, MAPK14, PIK3R1, and CASP3 were identified as hub genes, and their related miRNAs were detected as hub miRNAs. As a potential related miRNA of PIK3CB, miR-206 was reported as a regulatory miRNA on cell proliferation and apoptosis [50, 51], and was identified with PIK3CB showing sensitive to initial injury intensity induced by freeze damage [52]. Another potential related miRNA of PIK3CB, miR-125a, was reported as important regulation miRNA on metabolism [53] and cell development [54]. However, the targeting relationships between miR-125a, miR-206 and PIK3CB are still not reported and verified. Furthermore, in the present study, miR-125a and miR-206 showed opposite expression patterns, miR-125a showed higher expression in cecum and colon, while miR-206 showed highest expression in duodenum. For VEGFA, potential related miRNAs novel-miR-70 and novel-miR-378 were both novel identified miRNAs. Comparative analysis reveals that novel-miR-70 and novel-miR-378 showed homology with unannotated novel miRNA “chr5_15893_mature” [20] and “gi-417531957-ref-NC_019462.1-:15688451:15688536:+” [33], while their detailed functions kept unknown. For gene MAPK14, oar-miR-133 has been identified as the only related miRNA. miR-133 has been reported as an important regulatory factor on cell apoptosis [55] and conversion [56]. The predicted relationship between miR-133 and MAPK14 may indicates that mir-133 may affect the function of different intestine segments by regulating the expression of MAPK14, while more experiments will be required to confirm this assumption. For gene PIK3R1, the only related miRNA is miR-1343-3p, which showed highest expression in colon. Recent studies have shown that miR-1343-3p has tumor suppressor effect on multiple cancers [57, 58]. And for gene CASP3, a cecum specifically expressed novel-miR-174 has been identified as the only related miRNA. To our knowledge, there is currently no report about novel-miR-174, which indicates novel-miR-174 could be a tissue-specific miRNA. In summary, the hub miRNAs detected in this regulatory network may participate in intestinal epithelial morphogenesis, development, and growth, playing a role in maintaining intestinal homeostasis and metabolism. However, the targeting relationships between hub miRNAs and hub genes were not verified, and their specific mechanisms of functions are still unknown. Further researches will be performed to confirm the specific functions of these miRNAs.

Conclusions

In this study, we successfully constructed three sRNA libraries (DU, CE, and CO) in Ovis aries. In total, 106 known miRNAs, 458 conserved miRNAs, 192 unannotated novel miRNAs, and 195 differentially expressed miRNAs were identified by high-throughput sequencing and bioinformatics analysis. Additionally, 3,437 candidate target genes were significantly annotated to 17 GO terms and enriched in 4 KEGG biological pathways, such as “MAPK signaling pathway”, “Bacterial invasion of epithelial cells”, “Pancreatic cancer” and “PI3K-Akt signaling pathway”. Lastly, a combined regulatory network containing 92 candidate target genes and 65 differentially expressed miRNAs was constructed, among which 7 miRNAs were identified as hub miRNAs. The results provide a reference for elucidating, at the post-transcriptional level, the regulatory mechanisms that differentiate the morphological structures and physiological functions of the duodenum, cecum, and colon in sheep.

Supporting information

S1 Fig. Workflow for the differential expression analysis and functional annotation of miRNAs.

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

(TIF)

S2 Fig. Overview of miRNA annotation, identification, differential expression analysis and functional annotation results.

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

(TIF)

S3 Fig. Log2 transformed read distribution box plots for the three libraries.

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

(TIF)

S1 Table. Primer sequences for qRT-PCR validation.

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

(XLSX)

S2 Table. Distribution of reads and sheep known miRNAs on chromosome.

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

(XLSX)

S3 Table. Mireap novel miRNA precursor hairpins prediction results.

https://doi.org/10.1371/journal.pone.0193371.s006

(XLSX)

S4 Table. Mireap novel miRNA construction results.

https://doi.org/10.1371/journal.pone.0193371.s007

(XLSX)

S5 Table. Known, conserved and novel miRNAs identified in three libraries.

https://doi.org/10.1371/journal.pone.0193371.s008

(XLSX)

S6 Table. Comparison results of unannotated novel miRNAs.

https://doi.org/10.1371/journal.pone.0193371.s009

(XLSX)

S7 Table. Differentially expressed miRNAs between different libraries.

https://doi.org/10.1371/journal.pone.0193371.s010

(XLSX)

S8 Table. Targetscan and miRanda differentially expressed miRNAs target prediction results.

https://doi.org/10.1371/journal.pone.0193371.s011

(XLSX)

S9 Table. GO functional enrichment for candidate target genes.

https://doi.org/10.1371/journal.pone.0193371.s012

(XLSX)

S10 Table. KEGG pathway annotations for candidate target genes.

https://doi.org/10.1371/journal.pone.0193371.s013

(XLSX)

Acknowledgments

We thank Prof. Jianmin Wang for revising the language and valuable comments on the manuscript. We thank the Shandong Huanong Mutton Sheep Farm for providing the experimental sheep.

References

1. Wang J, Fan H, Han Y, Zhao J, Zhou Z. Characterization of the microbial communities along the gastrointestinal tract of sheep by 454 pyrosequencing analysis. Asian-Australasian Journal of Animal Sciences. 2017; 30(1):100–110. pmid:27383798
- View Article
- PubMed/NCBI
- Google Scholar
2. Hofmann R R. Evolutionary steps of ecophysiological adaptation and diversification of ruminants: a comparative view of their digestive system. Oecologia.1989; 78(4):443–457. pmid:28312172
- View Article
- PubMed/NCBI
- Google Scholar
3. Lewis S M, Dehority B A. Microbiology and ration digestibility in the hindgut of the ovine. Appl Environ Microbiol. 1985; 50: 356–363. pmid:4051484
- View Article
- PubMed/NCBI
- Google Scholar
4. Degregorio RM, Tucker R E, Mitchell G E Jr, Gill W W. Carbohydrate fermentation in the large intestine of lambs.J Anim Sci. 1982; 54:855–862. pmid:6282799
- View Article
- PubMed/NCBI
- Google Scholar
5. Faichney G J. Production of volatile fatty acids in the sheep caecum. Aust J Agric Res.1969; 20:491–498.
- View Article
- Google Scholar
6. Vienberg S, Geiger J, Madsen S, Dalgaard L T. MicroRNAs in metabolism. Acta Physiologica.2016.
- View Article
- Google Scholar
7. Runtsch M C, Round J L, O’Connell R M. MicroRNAs and the regulation of intestinal homeostasis. Front Genet. 2014; 5(5):347.
- View Article
- Google Scholar
8. Mckenna L B, Schug J, Vourekas A, et al. MicroRNAs control intestinal epithelial differentiation, architecture, and barrier function. Gastroenterology. 2010; 139(5):1654–1664. pmid:20659473
- View Article
- PubMed/NCBI
- Google Scholar
9. Liao Y, Zhang M, Lonnerdal B. Growth factor TGF-beta induces intestinal epithelial cell (IEC-6) differentiation: miR-146b as a regulatory component in the negative feedback loop. Genes & nutrition. 2013; 8(1):69–78.
- View Article
- Google Scholar
10. Ye D, Guo S, Alsadi R, Ma TY. MicroRNA regulation of intestinal epithelial tight junction permeability. Gastroenterology. 2011; 141(4):1323–1333. pmid:21763238
- View Article
- PubMed/NCBI
- Google Scholar
11. Sharbati S, Friedländer M R, Sharbati J, Hoeke L, Chen W, Keller A et al. Deciphering the porcine intestinal microRNA transcriptome. BMC Genomics. 2010; 11(1):275.
- View Article
- Google Scholar
12. Sharbati J, Hanisch C, Pieper R, Einspanier R, Sharbati S. Small molecule and RNAi induced phenotype transition of expanded and primary colonic epithelial cells. Scientific Reports. 2015; 5:12681. pmid:26223582
- View Article
- PubMed/NCBI
- Google Scholar
13. Tao X, Xu Z. MicroRNA transcriptome in swine small intestine during weaning stress. Plos One. 2013; 8(11):e79343. pmid:24260202
- View Article
- PubMed/NCBI
- Google Scholar
14. Yu J, Liu F, Yin P, Zhu X, Cheng G, Wang N, et al. Integrating miRNA and mRNA expression profiles in response to heat stress-induced injury in rat small intestine. Funct Integr Genomics. 2011; 11(2):203–213 pmid:21057845
- View Article
- PubMed/NCBI
- Google Scholar
15. Gao Y, Schug J, McKenna LB, Le Lay J, Kaestner KH, Greenbaum L E, et al. Tissue specific regulation of mouse microRNA genes in endoderm-derived tissues. Nucleic Acids Res. 2011; 39: 454–463. pmid:20843784
- View Article
- PubMed/NCBI
- Google Scholar
16. Liang G, Malmuthuge N, Mcfadden T B, Bao H, Griebel P J, Stothard P, et al. Potential Regulatory Role of MicroRNAs in the Development of Bovine Gastrointestinal Tract during Early Life. Plos One. 2014; 9(3):e92592. pmid:24682221
- View Article
- PubMed/NCBI
- Google Scholar
17. Coutinho LL, Matukumalli LK, Sonstegard TS, Van Tassell CP, Gasbarre LC, Capuco AV, et al. Discovery and profiling of bovine microRNAs from immune-related and embryonic tissues. Physiol Genomics. 2007; 29: 35–43. pmid:17105755
- View Article
- PubMed/NCBI
- Google Scholar
18. Zhang S, Zhao F, Wei C, Sheng X, Ren H, Xu L, et al. Identification and Characterization of the miRNA Transcriptome of Ovis aries. Plos One. 2013; 8(3):e58905. pmid:23516575
- View Article
- PubMed/NCBI
- Google Scholar
19. Chang W, Wang J, Tao D, Yong Z, Jianzhong H E, Shi C, et al. Identification of a novel miRNA from the ovine ovary by a combinatorial approach of bioinformatics and experiments. Journal of Veterinary Medical Science. 2015; 77(12):1617–1624. pmid:26268666
- View Article
- PubMed/NCBI
- Google Scholar
20. Miao X, Luo Q, Qin X, Guo Y. Genome-wide analysis of microRNAs identifies the lipid metabolism pathway to be a defining factor in adipose tissue from different sheep. Scientific Reports. 2015; 5:18470. pmid:26690086
- View Article
- PubMed/NCBI
- Google Scholar
21. Meale S J, Romao J M, He M L, Chaves A V, Mcallister T A, Guan L L. Effect of diet on microRNA expression in ovine subcutaneous and visceral adipose tissues. Journal of Animal Science. 2014; 92(8):3328–37. pmid:24893997
- View Article
- PubMed/NCBI
- Google Scholar
22. Xiaohui Tang, Guangbin Liu, Xiaoyong Du, Jianhua Cao, Xinyun Li, Zhang Luo, et al. Expression Variation of microRNA-31 in Wool Follicle of Tibetan Sheep(Ovis aries). Journal of Agricultural Biotechnology. 2012; 20(10):1178–1183.
- View Article
- Google Scholar
23. Laganà A, Veneziano D, Spata T, et al. Identification of General and Heart-Specific miRNAs in Sheep (Ovis aries). Plos One. 2015; 10(11):e0143313. pmid:26599010
- View Article
- PubMed/NCBI
- Google Scholar
24. Zhang Chunlan, Wang Guizhi, Wang Jianmin, Ji Zhibin, Dong Fei, Chao Tianle. Analysis of Differential Gene Expression and Novel Transcript Units of Ovine Muscle Transcriptomes. Plos One, 2014, 9(2):e89817. pmid:24587058
- View Article
- PubMed/NCBI
- Google Scholar
25. Li R, Yu C, Li Y, Lam T W, Yiu S M, Kristiansen K, et al. Soap2: an improved ultrafast tool for short read alignment. Bioinformatics, 2009, 25(15), 1966–1967. pmid:19497933
- View Article
- PubMed/NCBI
- Google Scholar
26. Kozomara A, Griffithsjones S. miRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Research, 2014, 42(Databaseissue): 68–73.
- View Article
- Google Scholar
27. Nawrocki EP; Burge SW; Bateman A; Daub J; Eberhardt RY; Eddy SR; Floden EW; Gardner PP; Jones TA; Tate J; Finn RD. Rfam 12.0: updates to the RNA families database. Nucleic Acids Research, 2015, 43(Database issue):D130. pmid:25392425
- View Article
- PubMed/NCBI
- Google Scholar
28. Miao X, Luo Q, Qin X. Genome-wide analysis reveals the differential regulations of mRNAs and miRNAs in dorset and small tail han sheep muscles. Gene. 2015;562: 188–196. pmid:25732516
- View Article
- PubMed/NCBI
- Google Scholar
29. Miao X, Qin QLX. Genome-wide transcriptome analysis of mRNAs and microRNAs in Dorset and Small Tail Han sheep to explore the regulation of fecundity. Mol Cell Endocrinol. 2015;402: 32–42. pmid:25573241
- View Article
- PubMed/NCBI
- Google Scholar
30. Yang H, Lin S, Lei X, Yuan C, Tian Z, Yu Y, et al. Identification and profiling of microRNAs from ovary of estrous Kazakh sheep induced by nutritional status in the anestrous season. Anim Reprod Sci. 2016;175: 18–26. pmid:27773477
- View Article
- PubMed/NCBI
- Google Scholar
31. Liu G, Liu R, Li Q, Tang X, Yu M, Li X, et al. Identification of microRNAs in Wool Follicles during Anagen, Catagen, and Telogen Phases in Tibetan Sheep. PLoS One. 2013;8. pmid:24204975
- View Article
- PubMed/NCBI
- Google Scholar
32. Wong LL, Rademaker MT, Saw EL, Lew KS, Ellmers LJ, Charles CJ, et al. Identification of novel microRNAs in the sheep heart and their regulation in heart failure. Sci Rep. 2017;7. pmid:28811555
- View Article
- PubMed/NCBI
- Google Scholar
33. Zhou G, Wang X, Yuan C, Kang D, Xu X, Zhou J, et al. Integrating miRNA and mRNA Expression Profiling Uncovers miRNAs Underlying Fat Deposition in Sheep. Biomed Res Int. 2017;2017. pmid:28293627
- View Article
- PubMed/NCBI
- Google Scholar
34. John B, Enright AJ, Aravin A, Tuschl T, Sander C, Marks DS. Human MicroRNA targets. PLoS Biol, 2005, Jul;3(7):e264.
- View Article
- Google Scholar
35. Agarwal V, Bell G W, Nam J W, Bartel D P. Predicting effective microRNA target sites in mammalian mRNAs. Elife Sciences, 2015, 4(e05005).
- View Article
- Google Scholar
36. Consortium T G O. The Gene Ontology Consortium: The Gene Ontology project in 2008. Nucleic Acids Research, 2007, 36(suppl 1):391.
- View Article
- Google Scholar
37. Wang J, Chao T, Wang G, Ji Z, Liu Z, Hou L, et al. Transcriptome Analysis of Three Sheep Intestinal Regions reveals Key Pathways and Hub Regulatory Genes of Large Intestinal Lipid Metabolism. Sci Rep. 2017;7.
- View Article
- Google Scholar
38. Huang DW, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID Bioinformatics Resources. Nature Protoc, 2009,4(1):44–57.
- View Article
- Google Scholar
39. Minoru K, Yoko S, Masayuki K, Miho F, Mao T. KEGG as a reference resource for gene and protein annotation. Nucleic Acids Research, 2016, 44(D1):D457. pmid:26476454
- View Article
- PubMed/NCBI
- Google Scholar
40. Szklarczyk D, Morris JH, Cook H, Kuhn M, Wyder S, Simonovic M, et al. The STRING database in 2017: quality-controlled protein-protein association networks, made broadly accessible. Nucleic Acids Res. 2017 Jan; 45:D362–68. pmid:27924014
- View Article
- PubMed/NCBI
- Google Scholar
41. Shannon P. Cytoscape: A Software Environment for Integrated Models of Biomolecular Interaction Networks. Genome Research, 2003, 13, 2498–2504. pmid:14597658
- View Article
- PubMed/NCBI
- Google Scholar
42. Jiang S, Li X, Wang X, Ban Q, Hui W, Jia B. MicroRNA profiling of the intestinal tissue of Kazakh sheep after experimental Echinococcus granulosus infection, using a high-throughput approach. Parasite-journal De La Societe Francaise De Parasitologie. 2016; 23:23.
- View Article
- Google Scholar
43. Bai J W, Xue H Z, Zhang C. Down-regulation of microRNA-143 is associated with colorectal cancer progression. European Review for Medical & Pharmacological Sciences. 2016; 20(22).
- View Article
- Google Scholar
44. Lund A H. miR-10 in development and cancer. Cell Death & Differentiation. 2010; 17(2):209–14.
- View Article
- Google Scholar
45. Ibrahim S A, Yip G W, Stock C, Pan J, Neubauer C, Poeter M, et al. Targeting of syndecan‐1 by microRNA miR‐10b promotes breast cancer cell motility and invasiveness via a Rho‐GTPase‐ and E‐cadherin‐dependent mechanism. International Journal of Cancer. 2012; 131(6):E884. pmid:22573479
- View Article
- PubMed/NCBI
- Google Scholar
46. Bauman Y, Yamin R, Vitenshtein A, Mandelboim O, Tsukerman P, Stern-Ginossar N, et al. MiR-10b downregulates the stress-induced cell surface molecule MICB, a critical ligand for cancer cell recognition by natural killer cells. Cancer Research. 2012; 72(21):5463–72. pmid:22915757
- View Article
- PubMed/NCBI
- Google Scholar
47. Lewis B P, Shih I H, Jonesrhoades M W, Bartel D P, Burge C B, et al. Prediction of mammalian microRNA targets. Cell. 2003; 115, 787–798. pmid:14697198
- View Article
- PubMed/NCBI
- Google Scholar
48. Lewis B P, Burge C B, Bartel D P. Conserved seed pairing, often flanked by adenosines, indicates that thousands of human genes are microRNA targets. Cell. 2005; 120(1):15–20. pmid:15652477
- View Article
- PubMed/NCBI
- Google Scholar
49. Bartel D P. MicroRNAs: target recognition and regulatory functions. Cell. 2009; 136(2):215–33. pmid:19167326
- View Article
- PubMed/NCBI
- Google Scholar
50. Sun H, Cai S, Zhang M, Zhao J, Wei S, Luo Y, et al. MicroRNA-206 regulates vascular smooth muscle cell phenotypic switch and vascular neointimal formation. Cell Biol Int. 2017;41: 739–748. pmid:28328152
- View Article
- PubMed/NCBI
- Google Scholar
51. Wang C-N, Wang Y-J, Wang H, Song L, Chen Y, Wang J-L, et al. The Anti-dementia Effects of Donepezil Involve miR-206-3p in the Hippocampus and Cortex. Biol Pharm Bull. 2017;40: 465–472. pmid:28123152
- View Article
- PubMed/NCBI
- Google Scholar
52. R M Jr, Carrigan C T, Abdalla M N, Geddis A V, Leandry L A, Aguilar C A, et al. Rna transcript expression of igf-i/pi3k pathway components in regenerating skeletal muscle is sensitive to initial injury intensity. Growth Hormone & Igf Research, 2016; 32, 14–21.
- View Article
- Google Scholar
53. Sun Z, Zhang W, Li Q. Mir-125a suppresses viability and glycolysis and induces apoptosis by targeting hexokinase 2 in laryngeal squamous cell carcinoma. Cell & Bioscience, 2017; 7(1), 51.
- View Article
- Google Scholar
54. Qin Y, Wu L, Ouyang Y, Zhou P, Zhou H, Wang Y, et al. MiR-125a Is a critical modulator for neutrophil development. PLoS Genet. 2017; 13. pmid:28976973
- View Article
- PubMed/NCBI
- Google Scholar
55. Wang L, Li X, Zhou Y, Shi H, Xu C, He H, et al. Downregulation of mir-133 via mapk/erk signaling pathway involved in nicotine-induced cardiomyocyte apoptosis. Naunyn-Schmiedeberg's Archives of Pharmacology, 2014; 387(2), 197–206. pmid:24190542
- View Article
- PubMed/NCBI
- Google Scholar
56. Morozzi G, Beccafico S, Bianchi R, Riuzzi F, Bellezza I, Giambanco I, et al. Oxidative stress-induced s100b accumulation converts myoblasts into brown adipocytes via an nf-κb/yy1/mir-133 axis and nf-κb/yy1/bmp-7 axis. Cell Death & Differentiation. 2017.
- View Article
- Google Scholar
57. Kim H, Yang JM, Jin Y, Jheon S, Kim K, Lee CT, et al. MicroRNA expression profiles and clinicopathological implications in lung adenocarcinoma according to EGFR, KRAS, and ALK status. Oncotarget. 2017;8: 8484–8498. pmid:28035073
- View Article
- PubMed/NCBI
- Google Scholar
58. Yuan T, Huang X, Woodcock M, Du M, Dittmar R, Wang Y, et al. Plasma extracellular RNA profiles in healthy and cancer patients. Sci Rep. 2016;6. pmid:26786760
- View Article
- PubMed/NCBI
- Google Scholar

[ref1] 1. Wang J, Fan H, Han Y, Zhao J, Zhou Z. Characterization of the microbial communities along the gastrointestinal tract of sheep by 454 pyrosequencing analysis. Asian-Australasian Journal of Animal Sciences. 2017; 30(1):100–110. pmid:27383798
View Article
PubMed/NCBI
Google Scholar

[2] View Article

[3] PubMed/NCBI

[4] Google Scholar

[ref2] 2. Hofmann R R. Evolutionary steps of ecophysiological adaptation and diversification of ruminants: a comparative view of their digestive system. Oecologia.1989; 78(4):443–457. pmid:28312172
View Article
PubMed/NCBI
Google Scholar

[6] View Article

[7] PubMed/NCBI

[8] Google Scholar

[ref3] 3. Lewis S M, Dehority B A. Microbiology and ration digestibility in the hindgut of the ovine. Appl Environ Microbiol. 1985; 50: 356–363. pmid:4051484
View Article
PubMed/NCBI
Google Scholar

[10] View Article

[11] PubMed/NCBI

[12] Google Scholar

[ref4] 4. Degregorio RM, Tucker R E, Mitchell G E Jr, Gill W W. Carbohydrate fermentation in the large intestine of lambs.J Anim Sci. 1982; 54:855–862. pmid:6282799
View Article
PubMed/NCBI
Google Scholar

[14] View Article

[15] PubMed/NCBI

[16] Google Scholar

[ref5] 5. Faichney G J. Production of volatile fatty acids in the sheep caecum. Aust J Agric Res.1969; 20:491–498.
View Article
Google Scholar

[18] View Article

[19] Google Scholar

[ref6] 6. Vienberg S, Geiger J, Madsen S, Dalgaard L T. MicroRNAs in metabolism. Acta Physiologica.2016.
View Article
Google Scholar

[21] View Article

[22] Google Scholar

[ref7] 7. Runtsch M C, Round J L, O’Connell R M. MicroRNAs and the regulation of intestinal homeostasis. Front Genet. 2014; 5(5):347.
View Article
Google Scholar

[24] View Article

[25] Google Scholar

[ref8] 8. Mckenna L B, Schug J, Vourekas A, et al. MicroRNAs control intestinal epithelial differentiation, architecture, and barrier function. Gastroenterology. 2010; 139(5):1654–1664. pmid:20659473
View Article
PubMed/NCBI
Google Scholar

[27] View Article

[28] PubMed/NCBI

[29] Google Scholar

[ref9] 9. Liao Y, Zhang M, Lonnerdal B. Growth factor TGF-beta induces intestinal epithelial cell (IEC-6) differentiation: miR-146b as a regulatory component in the negative feedback loop. Genes & nutrition. 2013; 8(1):69–78.
View Article
Google Scholar

[31] View Article

[32] Google Scholar

[ref10] 10. Ye D, Guo S, Alsadi R, Ma TY. MicroRNA regulation of intestinal epithelial tight junction permeability. Gastroenterology. 2011; 141(4):1323–1333. pmid:21763238
View Article
PubMed/NCBI
Google Scholar

[34] View Article

[35] PubMed/NCBI

[36] Google Scholar

[ref11] 11. Sharbati S, Friedländer M R, Sharbati J, Hoeke L, Chen W, Keller A et al. Deciphering the porcine intestinal microRNA transcriptome. BMC Genomics. 2010; 11(1):275.
View Article
Google Scholar

[38] View Article

[39] Google Scholar

[ref12] 12. Sharbati J, Hanisch C, Pieper R, Einspanier R, Sharbati S. Small molecule and RNAi induced phenotype transition of expanded and primary colonic epithelial cells. Scientific Reports. 2015; 5:12681. pmid:26223582
View Article
PubMed/NCBI
Google Scholar

[41] View Article

[42] PubMed/NCBI

[43] Google Scholar

[ref13] 13. Tao X, Xu Z. MicroRNA transcriptome in swine small intestine during weaning stress. Plos One. 2013; 8(11):e79343. pmid:24260202
View Article
PubMed/NCBI
Google Scholar

[45] View Article

[46] PubMed/NCBI

[47] Google Scholar

[ref14] 14. Yu J, Liu F, Yin P, Zhu X, Cheng G, Wang N, et al. Integrating miRNA and mRNA expression profiles in response to heat stress-induced injury in rat small intestine. Funct Integr Genomics. 2011; 11(2):203–213 pmid:21057845
View Article
PubMed/NCBI
Google Scholar

[49] View Article

[50] PubMed/NCBI

[51] Google Scholar

[ref15] 15. Gao Y, Schug J, McKenna LB, Le Lay J, Kaestner KH, Greenbaum L E, et al. Tissue specific regulation of mouse microRNA genes in endoderm-derived tissues. Nucleic Acids Res. 2011; 39: 454–463. pmid:20843784
View Article
PubMed/NCBI
Google Scholar

[53] View Article

[54] PubMed/NCBI

[55] Google Scholar

[ref16] 16. Liang G, Malmuthuge N, Mcfadden T B, Bao H, Griebel P J, Stothard P, et al. Potential Regulatory Role of MicroRNAs in the Development of Bovine Gastrointestinal Tract during Early Life. Plos One. 2014; 9(3):e92592. pmid:24682221
View Article
PubMed/NCBI
Google Scholar

[57] View Article

[58] PubMed/NCBI

[59] Google Scholar

[ref17] 17. Coutinho LL, Matukumalli LK, Sonstegard TS, Van Tassell CP, Gasbarre LC, Capuco AV, et al. Discovery and profiling of bovine microRNAs from immune-related and embryonic tissues. Physiol Genomics. 2007; 29: 35–43. pmid:17105755
View Article
PubMed/NCBI
Google Scholar

[61] View Article

[62] PubMed/NCBI

[63] Google Scholar

[ref18] 18. Zhang S, Zhao F, Wei C, Sheng X, Ren H, Xu L, et al. Identification and Characterization of the miRNA Transcriptome of Ovis aries. Plos One. 2013; 8(3):e58905. pmid:23516575
View Article
PubMed/NCBI
Google Scholar

[65] View Article

[66] PubMed/NCBI

[67] Google Scholar

[ref19] 19. Chang W, Wang J, Tao D, Yong Z, Jianzhong H E, Shi C, et al. Identification of a novel miRNA from the ovine ovary by a combinatorial approach of bioinformatics and experiments. Journal of Veterinary Medical Science. 2015; 77(12):1617–1624. pmid:26268666
View Article
PubMed/NCBI
Google Scholar

[69] View Article

[70] PubMed/NCBI

[71] Google Scholar

[ref20] 20. Miao X, Luo Q, Qin X, Guo Y. Genome-wide analysis of microRNAs identifies the lipid metabolism pathway to be a defining factor in adipose tissue from different sheep. Scientific Reports. 2015; 5:18470. pmid:26690086
View Article
PubMed/NCBI
Google Scholar

[73] View Article

[74] PubMed/NCBI

[75] Google Scholar

[ref21] 21. Meale S J, Romao J M, He M L, Chaves A V, Mcallister T A, Guan L L. Effect of diet on microRNA expression in ovine subcutaneous and visceral adipose tissues. Journal of Animal Science. 2014; 92(8):3328–37. pmid:24893997
View Article
PubMed/NCBI
Google Scholar

[77] View Article

[78] PubMed/NCBI

[79] Google Scholar

[ref22] 22. Xiaohui Tang, Guangbin Liu, Xiaoyong Du, Jianhua Cao, Xinyun Li, Zhang Luo, et al. Expression Variation of microRNA-31 in Wool Follicle of Tibetan Sheep(Ovis aries). Journal of Agricultural Biotechnology. 2012; 20(10):1178–1183.
View Article
Google Scholar

[81] View Article

[82] Google Scholar

[ref23] 23. Laganà A, Veneziano D, Spata T, et al. Identification of General and Heart-Specific miRNAs in Sheep (Ovis aries). Plos One. 2015; 10(11):e0143313. pmid:26599010
View Article
PubMed/NCBI
Google Scholar

[84] View Article

[85] PubMed/NCBI

[86] Google Scholar

[ref24] 24. Zhang Chunlan, Wang Guizhi, Wang Jianmin, Ji Zhibin, Dong Fei, Chao Tianle. Analysis of Differential Gene Expression and Novel Transcript Units of Ovine Muscle Transcriptomes. Plos One, 2014, 9(2):e89817. pmid:24587058
View Article
PubMed/NCBI
Google Scholar

[88] View Article

[89] PubMed/NCBI

[90] Google Scholar

[ref25] 25. Li R, Yu C, Li Y, Lam T W, Yiu S M, Kristiansen K, et al. Soap2: an improved ultrafast tool for short read alignment. Bioinformatics, 2009, 25(15), 1966–1967. pmid:19497933
View Article
PubMed/NCBI
Google Scholar

[92] View Article

[93] PubMed/NCBI

[94] Google Scholar

[ref26] 26. Kozomara A, Griffithsjones S. miRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Research, 2014, 42(Databaseissue): 68–73.
View Article
Google Scholar

[96] View Article

[97] Google Scholar

[ref27] 27. Nawrocki EP; Burge SW; Bateman A; Daub J; Eberhardt RY; Eddy SR; Floden EW; Gardner PP; Jones TA; Tate J; Finn RD. Rfam 12.0: updates to the RNA families database. Nucleic Acids Research, 2015, 43(Database issue):D130. pmid:25392425
View Article
PubMed/NCBI
Google Scholar

[99] View Article

[100] PubMed/NCBI

[101] Google Scholar

[ref28] 28. Miao X, Luo Q, Qin X. Genome-wide analysis reveals the differential regulations of mRNAs and miRNAs in dorset and small tail han sheep muscles. Gene. 2015;562: 188–196. pmid:25732516
View Article
PubMed/NCBI
Google Scholar

[103] View Article

[104] PubMed/NCBI

[105] Google Scholar

[ref29] 29. Miao X, Qin QLX. Genome-wide transcriptome analysis of mRNAs and microRNAs in Dorset and Small Tail Han sheep to explore the regulation of fecundity. Mol Cell Endocrinol. 2015;402: 32–42. pmid:25573241
View Article
PubMed/NCBI
Google Scholar

[107] View Article

[108] PubMed/NCBI

[109] Google Scholar

[ref30] 30. Yang H, Lin S, Lei X, Yuan C, Tian Z, Yu Y, et al. Identification and profiling of microRNAs from ovary of estrous Kazakh sheep induced by nutritional status in the anestrous season. Anim Reprod Sci. 2016;175: 18–26. pmid:27773477
View Article
PubMed/NCBI
Google Scholar

[111] View Article

[112] PubMed/NCBI

[113] Google Scholar

[ref31] 31. Liu G, Liu R, Li Q, Tang X, Yu M, Li X, et al. Identification of microRNAs in Wool Follicles during Anagen, Catagen, and Telogen Phases in Tibetan Sheep. PLoS One. 2013;8. pmid:24204975
View Article
PubMed/NCBI
Google Scholar

[115] View Article

[116] PubMed/NCBI

[117] Google Scholar

[ref32] 32. Wong LL, Rademaker MT, Saw EL, Lew KS, Ellmers LJ, Charles CJ, et al. Identification of novel microRNAs in the sheep heart and their regulation in heart failure. Sci Rep. 2017;7. pmid:28811555
View Article
PubMed/NCBI
Google Scholar

[119] View Article

[120] PubMed/NCBI

[121] Google Scholar

[ref33] 33. Zhou G, Wang X, Yuan C, Kang D, Xu X, Zhou J, et al. Integrating miRNA and mRNA Expression Profiling Uncovers miRNAs Underlying Fat Deposition in Sheep. Biomed Res Int. 2017;2017. pmid:28293627
View Article
PubMed/NCBI
Google Scholar

[123] View Article

[124] PubMed/NCBI

[125] Google Scholar

[ref34] 34. John B, Enright AJ, Aravin A, Tuschl T, Sander C, Marks DS. Human MicroRNA targets. PLoS Biol, 2005, Jul;3(7):e264.
View Article
Google Scholar

[127] View Article

[128] Google Scholar

[ref35] 35. Agarwal V, Bell G W, Nam J W, Bartel D P. Predicting effective microRNA target sites in mammalian mRNAs. Elife Sciences, 2015, 4(e05005).
View Article
Google Scholar

[130] View Article

[131] Google Scholar

[ref36] 36. Consortium T G O. The Gene Ontology Consortium: The Gene Ontology project in 2008. Nucleic Acids Research, 2007, 36(suppl 1):391.
View Article
Google Scholar

[133] View Article

[134] Google Scholar

[ref37] 37. Wang J, Chao T, Wang G, Ji Z, Liu Z, Hou L, et al. Transcriptome Analysis of Three Sheep Intestinal Regions reveals Key Pathways and Hub Regulatory Genes of Large Intestinal Lipid Metabolism. Sci Rep. 2017;7.
View Article
Google Scholar

[136] View Article

[137] Google Scholar

[ref38] 38. Huang DW, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID Bioinformatics Resources. Nature Protoc, 2009,4(1):44–57.
View Article
Google Scholar

[139] View Article

[140] Google Scholar

[ref39] 39. Minoru K, Yoko S, Masayuki K, Miho F, Mao T. KEGG as a reference resource for gene and protein annotation. Nucleic Acids Research, 2016, 44(D1):D457. pmid:26476454
View Article
PubMed/NCBI
Google Scholar

[142] View Article

[143] PubMed/NCBI

[144] Google Scholar

[ref40] 40. Szklarczyk D, Morris JH, Cook H, Kuhn M, Wyder S, Simonovic M, et al. The STRING database in 2017: quality-controlled protein-protein association networks, made broadly accessible. Nucleic Acids Res. 2017 Jan; 45:D362–68. pmid:27924014
View Article
PubMed/NCBI
Google Scholar

[146] View Article

[147] PubMed/NCBI

[148] Google Scholar

[ref41] 41. Shannon P. Cytoscape: A Software Environment for Integrated Models of Biomolecular Interaction Networks. Genome Research, 2003, 13, 2498–2504. pmid:14597658
View Article
PubMed/NCBI
Google Scholar

[150] View Article

[151] PubMed/NCBI

[152] Google Scholar

[ref42] 42. Jiang S, Li X, Wang X, Ban Q, Hui W, Jia B. MicroRNA profiling of the intestinal tissue of Kazakh sheep after experimental Echinococcus granulosus infection, using a high-throughput approach. Parasite-journal De La Societe Francaise De Parasitologie. 2016; 23:23.
View Article
Google Scholar

[154] View Article

[155] Google Scholar

[ref43] 43. Bai J W, Xue H Z, Zhang C. Down-regulation of microRNA-143 is associated with colorectal cancer progression. European Review for Medical & Pharmacological Sciences. 2016; 20(22).
View Article
Google Scholar

[157] View Article

[158] Google Scholar

[ref44] 44. Lund A H. miR-10 in development and cancer. Cell Death & Differentiation. 2010; 17(2):209–14.
View Article
Google Scholar

[160] View Article

[161] Google Scholar

[ref45] 45. Ibrahim S A, Yip G W, Stock C, Pan J, Neubauer C, Poeter M, et al. Targeting of syndecan‐1 by microRNA miR‐10b promotes breast cancer cell motility and invasiveness via a Rho‐GTPase‐ and E‐cadherin‐dependent mechanism. International Journal of Cancer. 2012; 131(6):E884. pmid:22573479
View Article
PubMed/NCBI
Google Scholar

[163] View Article

[164] PubMed/NCBI

[165] Google Scholar

[ref46] 46. Bauman Y, Yamin R, Vitenshtein A, Mandelboim O, Tsukerman P, Stern-Ginossar N, et al. MiR-10b downregulates the stress-induced cell surface molecule MICB, a critical ligand for cancer cell recognition by natural killer cells. Cancer Research. 2012; 72(21):5463–72. pmid:22915757
View Article
PubMed/NCBI
Google Scholar

[167] View Article

[168] PubMed/NCBI

[169] Google Scholar

[ref47] 47. Lewis B P, Shih I H, Jonesrhoades M W, Bartel D P, Burge C B, et al. Prediction of mammalian microRNA targets. Cell. 2003; 115, 787–798. pmid:14697198
View Article
PubMed/NCBI
Google Scholar

[171] View Article

[172] PubMed/NCBI

[173] Google Scholar

[ref48] 48. Lewis B P, Burge C B, Bartel D P. Conserved seed pairing, often flanked by adenosines, indicates that thousands of human genes are microRNA targets. Cell. 2005; 120(1):15–20. pmid:15652477
View Article
PubMed/NCBI
Google Scholar

[175] View Article

[176] PubMed/NCBI

[177] Google Scholar

[ref49] 49. Bartel D P. MicroRNAs: target recognition and regulatory functions. Cell. 2009; 136(2):215–33. pmid:19167326
View Article
PubMed/NCBI
Google Scholar

[179] View Article

[180] PubMed/NCBI

[181] Google Scholar

[ref50] 50. Sun H, Cai S, Zhang M, Zhao J, Wei S, Luo Y, et al. MicroRNA-206 regulates vascular smooth muscle cell phenotypic switch and vascular neointimal formation. Cell Biol Int. 2017;41: 739–748. pmid:28328152
View Article
PubMed/NCBI
Google Scholar

[183] View Article

[184] PubMed/NCBI

[185] Google Scholar

[ref51] 51. Wang C-N, Wang Y-J, Wang H, Song L, Chen Y, Wang J-L, et al. The Anti-dementia Effects of Donepezil Involve miR-206-3p in the Hippocampus and Cortex. Biol Pharm Bull. 2017;40: 465–472. pmid:28123152
View Article
PubMed/NCBI
Google Scholar

[187] View Article

[188] PubMed/NCBI

[189] Google Scholar

[ref52] 52. R M Jr, Carrigan C T, Abdalla M N, Geddis A V, Leandry L A, Aguilar C A, et al. Rna transcript expression of igf-i/pi3k pathway components in regenerating skeletal muscle is sensitive to initial injury intensity. Growth Hormone & Igf Research, 2016; 32, 14–21.
View Article
Google Scholar

[191] View Article

[192] Google Scholar

[ref53] 53. Sun Z, Zhang W, Li Q. Mir-125a suppresses viability and glycolysis and induces apoptosis by targeting hexokinase 2 in laryngeal squamous cell carcinoma. Cell & Bioscience, 2017; 7(1), 51.
View Article
Google Scholar

[194] View Article

[195] Google Scholar

[ref54] 54. Qin Y, Wu L, Ouyang Y, Zhou P, Zhou H, Wang Y, et al. MiR-125a Is a critical modulator for neutrophil development. PLoS Genet. 2017; 13. pmid:28976973
View Article
PubMed/NCBI
Google Scholar

[197] View Article

[198] PubMed/NCBI

[199] Google Scholar

[ref55] 55. Wang L, Li X, Zhou Y, Shi H, Xu C, He H, et al. Downregulation of mir-133 via mapk/erk signaling pathway involved in nicotine-induced cardiomyocyte apoptosis. Naunyn-Schmiedeberg's Archives of Pharmacology, 2014; 387(2), 197–206. pmid:24190542
View Article
PubMed/NCBI
Google Scholar

[201] View Article

[202] PubMed/NCBI

[203] Google Scholar

[ref56] 56. Morozzi G, Beccafico S, Bianchi R, Riuzzi F, Bellezza I, Giambanco I, et al. Oxidative stress-induced s100b accumulation converts myoblasts into brown adipocytes via an nf-κb/yy1/mir-133 axis and nf-κb/yy1/bmp-7 axis. Cell Death & Differentiation. 2017.
View Article
Google Scholar

[205] View Article

[206] Google Scholar

[ref57] 57. Kim H, Yang JM, Jin Y, Jheon S, Kim K, Lee CT, et al. MicroRNA expression profiles and clinicopathological implications in lung adenocarcinoma according to EGFR, KRAS, and ALK status. Oncotarget. 2017;8: 8484–8498. pmid:28035073
View Article
PubMed/NCBI
Google Scholar

[208] View Article

[209] PubMed/NCBI

[210] Google Scholar

[ref58] 58. Yuan T, Huang X, Woodcock M, Du M, Dittmar R, Wang Y, et al. Plasma extracellular RNA profiles in healthy and cancer patients. Sci Rep. 2016;6. pmid:26786760
View Article
PubMed/NCBI
Google Scholar

[212] View Article

[213] PubMed/NCBI

[214] Google Scholar

Figures

Abstract

Introduction

Materials and methods

Ethics statement of experimental animals

Sample collection and total RNA extraction

Library construction and sequencing

Sequence data analysis and miRNAs identification

Differential expression of miRNAs and prediction of target genes

GO annotation and KEGG pathway analysis

Regulatory network construction

Real-time quantitative RT-PCR validation

Results

Sequence data from three segments of intestinal tissues

The identification of known miRNAs and novel miRNAs

Analysis of differentially expressed miRNAs and prediction of target genes

GO enrichment for candidate target genes of differentially expressed miRNAs

KEGG analysis and regulatory network construction

Discussion

Conclusions

Supporting information

S1 Fig. Workflow for the differential expression analysis and functional annotation of miRNAs.

S2 Fig. Overview of miRNA annotation, identification, differential expression analysis and functional annotation results.

S3 Fig. Log2 transformed read distribution box plots for the three libraries.

S1 Table. Primer sequences for qRT-PCR validation.

S2 Table. Distribution of reads and sheep known miRNAs on chromosome.

S3 Table. Mireap novel miRNA precursor hairpins prediction results.

S4 Table. Mireap novel miRNA construction results.

S5 Table. Known, conserved and novel miRNAs identified in three libraries.

S6 Table. Comparison results of unannotated novel miRNAs.

S7 Table. Differentially expressed miRNAs between different libraries.

S8 Table. Targetscan and miRanda differentially expressed miRNAs target prediction results.

S9 Table. GO functional enrichment for candidate target genes.

S10 Table. KEGG pathway annotations for candidate target genes.

Acknowledgments

References