MicroRNA and circular RNA profiling in the deposited fat tissue of Sunite sheep

Abstract

As the most typical deposited fat, tail fat is an important energy reservoir for sheep adapted to harsh environments and plays an important role as a raw material in daily life. However, the regulatory mechanisms of microRNA (miRNA) and circular RNA (circRNA) in tail fat development remain unclear. In this study, we characterized the miRNA and circRNA expression profiles in the tail fat of sheep at the ages of 6, 18, and 30 months. We identified 219 differentially expressed (DE) miRNAs (including 12 novel miRNAs), which exhibited a major tendency to be downregulated, and 198 DE circRNAs, which exhibited a tendency to be upregulated. Target gene prediction analysis was performed for the DE miRNAs. Functional analysis revealed that their target genes were mainly involved in cellular interactions, while the host genes of DE circRNAs were implicated in lipid and fatty acid metabolism. Subsequently, we established a competing endogenous RNA (ceRNA) network based on the negative regulatory relationship between miRNAs and target genes. The network revealed that upregulated miRNAs play a leading role in the development of tail fat. Finally, the ceRNA relationship network with oar-miR-27a_R-1 and oar-miR-29a as the core was validated, suggesting possible involvement of these interactions in tail fat development. In summary, DE miRNAs were negatively correlated with DE circRNAs during sheep tail fat development. The multiple ceRNA regulatory network dominated by upregulated DE miRNAs may play a key role in this developmental process.

Keywords: Sunite sheep; circular RNA; competing endogenous RNA; microRNA; tail fat.

Figure 1

(A) The length distribution of sequenced miRNAs. 6M: A1, A2, A3; 18M: B1, B2, B3; 30M: C1, C2, C3. (B) Percentage of circRNA host genes. For example, the blue pie chart shows that 39.75% of the genes are transcribed to form 1 circRNA, meanwhile the red pie chart shows that 17.68% of the genes can be transcribed to make 2 circRNAs, and so on.

Figure 2

miRNA differential expression analysis. (A–C) The volcano plot of DE miRNA. (A) 30M vs. 6M; (B) 30M vs. 18M; (C) 18M vs. 6M. Annotated as the top five DE miRNA based on P-value. (D) Venn diagram analysis of DE miRNA. (E) Expression trend analysis of DE miRNA. The upper number indicates the ordinal number of each trend, while the lower number is the number of genes enriched, and those with color are the significantly enriched trends.

Figure 3

circRNA differential expression analysis. (A–C) Volcano plot of DE circRNAs. (A) 30M vs. 6M; (B) 30M vs. 18M; (C) 18M vs. 6M. Annotated as the top five DE circRNAs based on P-value. (D) Venn diagram analysis of DE circRNA. (E) DE circRNA expression trend analysis. The upper number indicates the ordinal number of each trend, while the lower number is the number of genes enriched, and those with color are the significantly enriched trends.

Figure 4

Validation of DE miRNAs and DE circRNAs via qRT-PCR. The blue bar represents qRT-PCR data, and the red line represents RNA-seq data.

Figure 5

The top 15 GO and KEGG terms for DE miRNAs from the three comparison groups. (A) GO analysis; (B) KEGG analysis.

Figure 6

The top 15 enriched GO and KEGG terms for DE circRNAs. (A) GO analysis; (B) KEGG analysis.

Figure 7

ceRNA regulatory network analysis in sheep tail fat. (A) Down-up-down mode. (B) Up-down-up mode. The shapes represent different RNAs and the colors represent different regulations.

Figure 8

ceRNA dual-luciferase reporter gene analysis. (A) Binding site validation for oar-miR-27a_R-1 and ACSL4 as well as oar-27a_R-1 for circRNA1985; (B) Binding site validation for oar-miR-29a and GPAM as well as oar-miR-29a and circRNA3539. Bars without slashes show the result of inserting the mutant sequence into the plasmid, while solid bars show the result of inserting the original sequence. The group with extra miRNA sequences is shown in red, and the group without any miRNA sequences (control group) is shown in blue. ***P < 0.001.