Oncotarget

Research Papers:

Genome-wide analysis of circular RNAs in prenatal and postnatal muscle of sheep

PDF |  HTML  |  Supplementary Files  |  How to cite  |  Order a Reprint

Oncotarget. 2017; 8:97165-97177. https://doi.org/10.18632/oncotarget.21835

Metrics: PDF 883 views  |   HTML 1751 views  |   ?  

Cunyuan Li, Xiaoyue Li, Yang Yao, Qiman Ma, Wei Ni _, Xiangyu Zhang, Yang Cao, Wureli Hazi, Dawei Wang, Renzhe Quan, Xiaoxu Hou, Zhijin Liu, Qianqian Zhan, Li Liu, Mengdan Zhang, Shuting Yu and Shengwei Hu

Abstract

Cunyuan Li1,*, Xiaoyue Li1,*,Yang Yao1,*, Qiman Ma1, Wei Ni1, Xiangyu Zhang1, Yang Cao1, Wureli Hazi2, Dawei Wang1, Renzhe Quan1, Xiaoxu Hou1, Zhijin Liu1, Qianqian Zhan1, Li Liu1, Mengdan Zhang1, Shuting Yu1 and Shengwei Hu1

1College of Life Sciences, Shihezi University, Shihezi, Xinjiang, 832003, China

2College of Animal Science and Technology, Shihezi University, Shihezi, Xinjiang, 832003, China

*These authors contributed equally to this work

Correspondence to:

Wei Ni, email: niweiwonderful@sina.com

Shengwei Hu, email: hushengwei@163.com

Keywords: sheep, circRNAs, muscle, expression profiles, RNA-seq

Received: June 05, 2017     Accepted: September 23, 2017     Published: October 12, 2017

ABSTRACT

Circular RNAs (circRNAs), a type of non-coding RNA with circular structure, were generated by back splicing and widely expressed in animals and plants. Recent studies have shown that circRNAs extensively participate in cell proliferation, cell differentiation, cell autophagy and other biological processes. However, the role and expression of circRNAs in the development and growth of muscle have not been studied in sheep. In our study, we first used RNA-seq to study the circRNAs in prenatal and postnatal longissimus dorsi muscle of sheep. A total of 6113 circRNAs were detected from the RNA-seq data. Several circRNAs were identified using reverse transcription PCR, DNA sequencing and RNase R digestion experiments. The expression levels of several circRNAs in prenatal and postnatal muscle were confirmed by Real-Time RT-PCR. The gene ontology (GO) and KEGG enrichment analysis of the host gene of the circRNAs showed that these circRNAs were mainly involved in the growth and development of muscle related signaling pathways. These circRNAs might sponge microRNAs (miRNAs) in predicted circRNA-miRNA-mRNA networks. The circRNAs expression profiles in muscle provided an important reference for the study of circRNAs in sheep.


Genome-wide analysis of circular RNAs in prenatal and postnatal muscle of sheep | Li | Oncotarget

INTRODUCTION

CircRNAs were firstly discovered in plant viruses and hepatitis δ virus [1, 2]. Although circRNA has been observed in eukaryotes a few decades ago, it was mostly considered assplicing errors at best and lacking biological function [3, 4]. With breakthroughs in high-throughput sequencing techniques [5], several recent studies have shown that there are many single-stranded exonic circRNAs present in mammals [68], and a large number of studies have shown that circRNAs are generated by head-to-tail splicing [9, 10]. There are several proposed models involved in the biogenesis of circRNAs, including spliceosome-dependent circulation path, intron-pairing-driven circularization path, lariat-driven circularization path, protein factors associated circulation path [1114]. CircRNAs were stable, conserved, tissue- and stage-specific expression [7]. The biological functions of circRNAs are still being indicated, but they are shown to sponge miRNAs [15], regulate gene transcription [16], translate into proteins [14], interact with RNA binding proteins [17] and modulate the stability of mRNAs [16]. A growing number of circRNAs are aberrantly expressed in colorectal cancer [18], hepatocellular carcinoma [19], gastric cancer [20], and laryngeal squamous cell carcinoma [21] and have the potential to become new diagnostic or prognostic biomarkers.

The skeletal muscle has become an important material for studying the mechanism of specific cell differentiation and proliferation in mammal [22]. Classic developmental studies have understood the origin and development of skeletal muscle cells in embryo. It has been known for many years that vertebrate skeletal muscle derived from embryonic mesodermal precursor cells [23]. The medial and lateral domains of somites can produce different muscle groups [2226]. Longissimus dorsi muscle as the largest part of the ridge of the spine, it is a prime cut of high motion relevance for fresh and cured meet production [27]. Recent studies on skeletal muscle of monkeys in a range of ages [28] and chicken embryo skeletal muscle [29] have shown that circRNAs may affect muscle growth. Previous studies have shown that miRNAs such as miR-133, miR-206 and miR-1 expressed specifically in muscle [30], and long non-coding RNA (LINCMD1) can act as a sponge of miR-133 to regulate muscle development [31]. Recent studies have shown that circRNAs can also act as sponges for adsorbing miRNAs [32].

Sheep has been used as a powerful model organism for human disease studies and is an important animal in agriculture because of its utility in meat production. However, until now the expression and function of circRNAs in the muscle of sheep is unclear. In our study, we first used RNA-seq to study the circRNAs of prenatal and postnatal muscle in sheep. The circRNAs expression profiles of the longissimus dorsi muscle was given by our study, which will facilitate better understanding roles of circRNAs in development and growth of muscle.

RESULTS

Deep sequencing of circRNAs in sheep muscle

In order to understand the identity, abundance and difference of circRNAs in the muscle of embryonic and adult sheep, we performed paired-end ribominus RNA sequencing (RNA-seq) according to the flow chart (Figure 1), and used the FindCirc calculation pipeline to detect the circRNAs [7]. Through the RNA-seq of the embryo longissimus dorsi muscle (LDM_E) and adult longissimus dorsi muscle (LDM_A), and we obtained 15.91GB and 13.94GB clear data from the results respectively. After processing, a total of 6113 circRNAs were identified from these data, and the length of the circRNAs were mainly focused below 13000 nt (Figure 2A). These circRNAs were distributed over 26 autosomes and X chromosomes, and these circRNAs consisted mainly of introns and exons, and only a small fraction are composed of intergenic sequences (Figure 2B). The annotation, chromosomal location, hosting mRNA and so on of all circRNAs are shown in Supplementary File 1.

Pipeline for circRNAs identification.

Figure 1: Pipeline for circRNAs identification.

The information of circRNAs in sheep muscle by deep sequencing.

Figure 2: The information of circRNAs in sheep muscle by deep sequencing. (A) The length of circRNAs in sheep muscle. Light red on behalf of LDM_A and light blue on behalf of LDM_E. (B) The number of introns, exons and intergenic of circRNAs in sheep muscle.

Validation of sheep circRNAs

We have further done some experimental validation in order to confirm the RNA-seq and expression of circRNAs in sheep. We randomly selected 10 circRNAs and designed divergent primers across the circRNAs junctions (Figure 3A). RT-PCR amplification and DNA sequencing techniques were used to confirm the circRNAs data of RNA-seq. The results of the RT-PCR amplification show a single band of the expected size (Figure 3B). DNA sequencing results confirm the head-to-tail circularization splicing of these circRNAs (Figure 3C). For circRNA 0003451 and circRNA 0005256, there is a single base difference between results of normal DNA and deep sequencing, which may be single nucleotide polymorphism (SNP) similar to some SNP in miRNAs [33]. Real-Time RT-PCR was further used for verifying circRNAs resistance to RNase R digestion experiments. All selected 10 circRNAs showed resistance to RNase R digestion. As expected, linear control of β-Actin was sensitive to the RNase R digestion, and it was completely digested (Figure 3D).

Verification of circRNAs data from RNA-sequencing.

Figure 3: Verification of circRNAs data from RNA-sequencing. (A) Divergent primers to amplify the circular junctions. Red arrows represent divergent primers. (B) RT-PCR amplification of circRNAs with divergent primers. PCR products of circRNA 0002456, circRNA 0005179, circRNA 0000666, circRNA 0000552, circRNA 0005250, circRNA 0005256, circRNA 0003541, circRNA 0005243, circRNA 0004690 and circRNA 0004676 were analysed by gel electrophoresis. (C) Head-to-tail junctions were confirmed by DNA sequencing. The unmatched bases in circRNA 0005256 and circRNA 0003541 are indicated in green. (D) Resistance test of circRNAs to RNase R digestion by Real-Time RT-PCR. β-Actin was used as a linear control. Error bars indicate ± SD.

Analysis and validation of differentially expressed circRNAs

A total of 5086 differentially expressed circRNAs were found between the embryo longissimus dorsi muscle (LDM_E) and adult longissimus dorsi muscle (LDM_A), including 2940 up-regulated and 2146 down-regulated circRNAs (Figure 4A). For verify the accuracy of these differentially expressed, we randomly selected 10 significant differentially expressed circRNAs using Real-Time RT-PCR (Figure 4B). The expression levels of 9 circRNAs detected by RT-PCR analysis (except for circRNA 0005179) were consistent with the RNA-seq data. CircRNA 0002456, circRNA 0000552, circRNA 0004676, circRNA 0004690 and circRNA 0000666 were higher in the embryonic group (LDM_E) than in adult group (LDM_A), and in contrast, circRNA 0005256, circRNA 0003451, circRNA 0005243 and circRNA 0005250 were highly expressed in the adult group (LDM_A) relative to the embryonic group (LDM_E). Real-Time RT-PCR results showed that the expression level of circRNA 0005179 was higher in the embryonic group (LDM_E) than in adult group (LDM_A), which was contrary with the RNA-seq data. The results show that there is a strong agreement between real-time RT-PCR and RNA-seq data, and it is shown that the identified circRNAs largely reflect the true differential expression in vivo.

Analysis and validation of differentially expressed circRNAsin LDM_E and LDM_A group.

Figure 4: Analysis and validation of differentially expressed circRNAsin LDM_E and LDM_A group. (A) The volcano plot analysis of all circRNAs between LDM_E and LDM_A group. The logarithm of the significant difference between the two samples was analyzed by log2 (fold change) as the abscissa, and the negative logarithm-log10 (P-value) of the P value was calculated as the ordinate (P < 0.05).Red dots indicate up-regulated genes; green dots indicate down-regulated genes. (B) Heapmap analysis of ten differentially expressed genes was performed with data from RNA-seq and Real-Time RT-PCR. Green represents low expression, and red represents high expression, P < 0.05.

Enrichment of differentially expressed circRNAs

For the hosting genes of differentially expressed circRNAs, we performed enrichment analysis of them by GO and KEGG pathways. A total of 897 hosting genes of differentially expressed circRNAs were analyzed by GO analysis (P < 0.05, which were enriched separately in biological process, molecular function and cellular components category (Figure 5A). The top 20 genes from each GO category were related in metabolic process, position in cell or organelle, and function of nucleic binding, which suggests that some circRNAs were involved to the basic process of muscle growth and development. A total of 270 terms were enriched by KEGG pathway in which AMPK signaling pathway, ECM-receptor interactions, ErbB signaling pathway, focal adhesion, valine, leucine and isoleucine degradation, ubiquitin mediated proteolysis, mTOR signaling pathway and glutamatergic synapse were related to muscle differentiation and proliferation (Figure 5B) (Supplementary File 2). Analysis of these results indicate that circRNAs may play a vital role in muscle development and growth.

Annotations and enrichment of differentially expressed circRNAs.

Figure 5: Annotations and enrichment of differentially expressed circRNAs. (A) The GO analysis showed 897 significantly (LDM_E VS LDM_A) enriched terms (P < 0.05) in biological process, molecular function and cellular components category. (B) A total of 270 terms were enriched of differentially expressed circRNAs in LDM_E and LDM_A group by KEGG pathway analysis. The AMPK signaling pathway, ECM-receptor interactions, ErbB signaling pathway, focal adhesion, valine, leucine and isoleucine degradation, ubiquitin mediated proteolysis, mTOR signaling pathway and glutamatergic synapse, which were involved muscle growth and development, were labeled with red lines.

Predicted functions of sheep circRNAs

Functional study of circRNAs are mainly focused on that circRNAs as an adsorbed miRNAs sponge that interacts with miRNAs [3437]. Many studies have indicated that cicrRNAs affect gene expression at post-transcriptional levels through miRNA response elements (MRE) [38, 39]. In our study, a total of 286956 interactions were observed with the various miRNAs in the identified 6113 circRNAs, and the interaction between the circRNAs and the miRNAs were analyzed by using miRanda and psRobot (Supplementary File 3). It is noteworthy that several well-known miRNAs are closely related to growth and development of skeletal muscle, which are considered to be promising candidates for future research. For example, circRNAs (circRNA 0000385, circRNA 0000582 and circRNA 0001099 etc) have multiple conservative target sites for muscle development-related miRNAs (miR-143, miR-133 and miR-23etc, respectively). We further analyzed the number of miRNAs and MREs that interacted with these circRNAs (Table 1). This result suggests that it is likely that the process of the development and growth of the muscles were affected by the circRNAs.

Table 1: Potential miRNA (muscle development-related miRNAs) targets of circRNAs

CircRNA ID

Regulation (LDM_E VS LDM_A)

Gene symbol

Potential miRNA targets (No. MREs)

oar_circ_0000385

Up

RAD51C

let-7a (12); let-7b (6); let-7c (9); let-7d (3); let-7f (15); let-7g (8); let-7i (1); miR-103 (15); miR-107 (7); miR-10a (1); miR-10b (4); miR-1185-3p (2); miR-1185-5p (1); miR-1193-5p (6); miR-1197-3p (2); miR-1197-5p (2); miR-133 (3); miR-134-3p (4); miR-134-5p (1); miR-136 (4); miR-148a (4); miR-150 (4); miR-152 (5); miR-154b-3p (4); miR-16b (6); miR-181a (9); miR-191 (5); miR-194 (4); miR-199a-3p (5); miR-200a (5); miR-200b (4); miR-200c (7); miR-218a (6); miR-21 (1); miR-23a (5); miR-26a (6); miR-26b (3); miR-299-5p (4); miR-29a (5); miR-29b (5); miR-30a-3p (5); miR-30a-5p (8); miR-30b (7); miR-30c (4); miR-30d (10); miR-323a-3p (4); miR-323a-5p (1); miR-323b (5); miR-323c (1); miR-329a-3p (2); miR-329a-5p (1); miR-329b-3p (1); miR-329b-5p (1); miR-369-3p (3); miR-370-3p (3); miR-376a-3p (2); miR-376a-5p (4); miR-376b-3p (3); miR-376b-5p (1); miR-376d (1); miR-376e-5p (1); miR-380-5p (2); miR-381-3p (1); miR-382-5p (6); miR-3955-3p (1); miR-3955-5p (1); miR-3956-3p (1); miR-3958-5p (1); miR-3959-3p (1); miR-409-3p (4); miR-410-3p (2); miR-411b-3p (2); miR-412-5p (3); miR-485-3p (1); miR-485-5p (3); miR-487a-3p (3); miR-487b-5p (3); miR-493-3p (1); miR-493-5p (6); miR-494-3p (5); miR-495-3p (4); miR-495-5p (1); miR-496-3p (5); miR-496-5p (1); miR-539-3p (5); miR-539-5p (1); miR-541-5p (4); miR-543-3p (3); miR-544-3p (8); miR-654-5p (2); miR-665-3p (3); miR-668-3p (1); miR-758-3p (1);

oar_circ_0000582

Down

YME1L1

miR-103 (2); miR-107 (1); miR-10b (2); miR-1185-3p (1); miR-1197-3p (1); miR-125b (4); miR-134-5p (4); miR-143 (3); miR-150 (1); miR-154a-3p (1); miR-194 (2); miR-200b (2); miR-200c (1); miR-221 (5); miR-23a (4); miR-23b (2); miR-25 (1); miR-26a (8); miR-26b (4); miR-29a (1); miR-29b (2); miR-30a-3p (2); miR-30a-5p (1); miR-30b (6); miR-30c (5); miR-30d (1); miR-323a-3p (2); miR-323b (2); miR-323c (1); miR-329a-5p (3); miR-329b-5p (1); miR-362 (2); miR-374a (4); miR-376a-3p (1); miR-376a-5p (2); miR-376b-3p (1); miR-376b-5p (1); miR-376c-3p (4); miR-376c-5p (6); miR-376d (1); miR-376e-3p (1); miR-376e-5p (1); miR-377-3p (3); miR-382-5p (2); miR-3955-3p (1); miR-3959-5p (1); miR-409-3p (1); miR-412-5p (1); miR-487a-5p (2); miR-494-3p (3); miR-494-5p (2); miR-495-5p (3); miR-496-3p (2); miR-496-5p (2); miR-539-3p (4); miR-539-5p (1); miR-544-3p (3); miR-668-5p (1); miR-758-3p (2)

oar_circ_0001099

Down

C7

miR-134-3p (2); miR-23b (2); miR-377-3p (2); miR-323a-3p (1); miR-329a-3p (1); miR-409-3p (3); miR-493-3p (1); miR-106b (1); miR-106a (1); miR-22-3p (1); miR-23a (3); miR-148a (2); miR-26a (1); miR-380-5p (2); miR-200b (6); miR-376b-5p (1); miR-362 (1); miR-376e-3p (1); miR-221 (1); miR-27a (1); miR-409-5p (1); miR-26b (1); miR-654-3p (1); miR-376e-5p (1); miR-381-3p (2); miR-495-3p (3); miR-30b (1); miR-125b (1); miR-382-3p (2); miR-3956-3p (2); miR-379-5p (1); miR-1193-5p (1); miR-665-3p (2); miR-152 (2); miR-3955-3p (4); miR-376a-5p (1); miR-30c (1); miR-3957-3p (1); miR-30d (1); miR-17-5p (1); miR-30a-5p (1); miR-665-5p (1); miR-1185-3p (1); miR-411b-3p (1); miR-1185-5p (1); miR-376c-3p (1); miR-181a (2); miR-30a-3p (2); miR-376c-5p (1); miR-323c (2); miR-494-3p (2); miR-543-3p (1); miR-487a-5p (1); miR-329b-3p (3); miR-544-5p (1); miR-3959-3p (1); miR-323b (2); miR-544-3p (3); miR-323a-5p (2); miR-199a-3p (3); miR-370-3p (1); miR-410-5p (2); miR-134-5p (1); miR-200c (5); miR-377-5p (1); miR-99a (1); miR-143 (1)

oar_circ_0001413

Up

SEMA7A

miR-181a (1); miR-329a-3p (1); miR-133 (1); miR-329b-3p (1); miR-493-3p (2); miR-107 (1); miR-1185-3p (1); miR-125b (1); miR-103 (1); miR-362 (1); miR-668-5p (1); miR-134-5p (1); miR-665-3p (1); miR-218a (1)

oar_circ_0003451

Down

TTN

miR-410-5p (1); miR-154b-5p (1); miR-758-5p (1); miR-544-3p (1)

oar_circ_0005250

Down

MYH7

miR-134-3p (10); miR-493-3p (6); miR-539-3p (1); miR-376b-5p (1); miR-654-5p (1); miR-409-3p (2); miR-376e-5p (1); let-7d (2); miR-379-3p (1); let-7f (2); let-7g (2); let-7a (2); let-7b (2); let-7c (2); miR-1193-5p (4); miR-323c (1); let-7i (2); miR-3955-3p (3); miR-376a-5p (3); miR-22-3p (3); miR-412-5p (7); miR-3956-3p (3); miR-1197-3p (2); miR-758-5p (1); miR-541-3p (5); miR-3958-5p (5); miR-29b (1); miR-3959-5p (2); miR-370-3p (8); miR-154b-3p (1); miR-377-5p (2); miR-493-5p (2); miR-380-5p (3); miR-495-5p (2); miR-27a (3); miR-433-5p (2); miR-1193-3p (4); miR-152 (3); miR-150 (2); miR-382-5p (1); miR-1197-5p (5); miR-148a (3); miR-411a-3p (1); miR-376c-5p (1); miR-3959-3p (2); miR-323a-5p (1); miR-154b-5p (2); miR-3958-3p (1); miR-143 (2); miR-29a (1); miR-154a-5p (1); miR-23a (1); miR-668-3p (2); miR-543-5p (8); miR-362 (2); miR-433-3p (1); miR-16b (2); miR-380-3p (1); miR-125b (9); miR-495-3p (1); miR-382-3p (1); miR-106b (2); miR-487a-5p (1); miR-106a (2); miR-3957-3p (1); miR-411b-3p (1); miR-1185-5p (2); miR-181a (3); miR-200a (2); miR-200b (3); miR-200c (4); miR-299-3p (3); miR-544-3p (7); miR-30a-3p (2); miR-17-5p (2); miR-23b (1); miR-103 (2); miR-107 (2); miR-431 (1); miR-432 (10); miR-221 (1); miR-154a-3p (1); miR-654-3p (5); miR-3955-5p (2); miR-379-5p (1); miR-487b-5p (1); miR-412-3p (1); miR-136 (3); miR-133 (1); miR-411b-5p (1); miR-665-5p (1); miR-1185-3p (1); miR-10a (2); miR-10b (2); miR-758-3p (1); miR-3957-5p (2); miR-485-5p (5); miR-494-3p (2); miR-370-5p (4); miR-377-3p (2)

MREs: miRNA response elements.

Bioinformatics analysis of circRNA-miRNA-mRNA networks

To elucidate the function of circRNAs, we used TargetScan and miRanda to predict the differentially expressed circRNA target and downstream regulated mRNA, and establish the basic circRNA-miRNA connectivity. The circRNA-miRNA-mRNA interaction network was established using Cytoscape version 3.5.1 (Figure 6A, Supplementary File 4). In order to further understand the function of circRNA in skeletal muscle growth and development, we had further established a circRNA-miRNA-mRNA network with potentially effective circRNAs (circRNA 0000385, circRNA 0000582 and circRNA 0001099 etc), their predicted miRNA targets and downstream regulated mRNA (Figure 6B6G). We used KEGG pathway analysis to functionally annotate the predicted mRNA targets within the above networks (As shown in Figure 7). These networks had several significantly enriched pathways related to muscle growth and development. This information may help us to explore the basic mechanism of circRNA in the growth and development of sheep skeletal muscle.

The circRNA-miRNA co-expression network.

Figure 6: The circRNA-miRNA co-expression network. (A) The circRNA-miRNA-mRNA interaction network between circRNA and its target gene. Ellipse nodes represent miRNAs, rectangle nodes represent genes and triangle nodes represent circRNAs. Red color and green color represents up and down regulation respectively. (B-G) The circRNA-miRNA-mRNA network with potentially effective circRNAs (circRNA 0000385 (B), circRNA 0000582 (C), circRNA 0001099 (D), circRNA 0001413 (E), circRNA 0003451 (F), circRNA 0005250 (G). Ellipse nodes represent miRNAs (blue), rectangle nodes represent genes (green) and triangle nodes represent circRNAs (red).

Bioinformatics analysis of circRNA-miRNA-mRNA network with potentially effective circRNAs.

Figure 7: Bioinformatics analysis of circRNA-miRNA-mRNA network with potentially effective circRNAs. Target mRNAs of the 6 networks (circRNA 0000385 (A), circRNA 0000582 (B), circRNA 0001099 (C), circRNA 0001413 (D), circRNA 0003451 (E), circRNA 0005250 (F)) are functionally annotated by KOBAS and KEGG pathway analysis. Top 10 significant (P < 0.05) enriched pathway terms of the 6 networks. Involving muscle growth and development of the signal pathway, marked with a red line.

DISCUSSION

With the development of high-throughput techniques, the deep sequencing of non-polyadenylated RNA populations indicates cumulative signals from certain exons (called excised exons) [47]. Further studies have shown that it is circRNAs [48]. In recent years, the genome of some breeds of sheep has been assembled and annotated, but sheep transcriptome analysis still needs to obtain more data from different tissue samples. Although thousands of unique circRNAs have been identified in cell types of different species, especially in humans and mice [4953], and some of them have been demonstrated play an important role in animal development and growth. However, the circRNAs of sheep almost do not know. In this study, we used a comprehensive sequencing and analysis of the longest muscle ribosomal deletion of RNA from the embryonic and adult sheep, and 6113 circRNAs were successfully identified in numerous circRNAs.

As an important regulators of gene expression, circRNAs have an important role in growth and development of animal. A total of 5086 differentially expressed circRNAs between sheep embryos and adult muscle were identified by our study. These circRNAs may have specific biological effects on the development of muscle. The development and growth of muscle contains a number of changes in the expression of many genes and non-coding RNAs [54]. Sun et al. reported that lncRNAs played an important role in the development of bovine longissimus dorsi muscle [55]. Recently, Sun et al. reported that 5,566 lncRNA and 4,360 circRNA were differentially expressed in longissimus dorsi muscle of Landrace and Lantang pigs, indicating that there is a potentially post-transcriptional regulation of non-coding RNA involved in the development of pig muscle [56]. Therefore, we predicted that circRNAs may play a new post-transcriptional regulation during growth and development of muscle in sheep.

The current study of the function of circRNAs is mainly focused on competitive endogenous RNA (ceRNA). CircRNAs can function as a sponge for adsorbing miRNAs to influence post-transcriptional regulation [5760]. In this study, miRNAs target sites were found in sheep circRNAs by using miRanda and psRobot. There are many circRNAs that interact with a lot of muscle-related miRNAs (miR-143, miR-133 and miR-23 etc). Interestingly, we found that some circRNAs contain simultaneously multiple target sites of different miRNAs. For example, oar_circ_0001413 contains target sites of miRNA-133, miRNA-125, miRNA-134, miRNA-103, miRNA-107, miRNA-1185, miRNA-181, miRNA-218, miRNA-329, simultaneously. One study has shown that circRNA (circHIPK3) regulates cell growth by sponging to 9 miRNAs with 18 potential binding sites [61]. Thus, as potential ceRNA, these circRNAs which can interact with multiple miRNAs may play a broad endogenous regulatory role in the growth and development of muscle. Oar_circ_0001413, as potential ceRNA, may be center regulator for growth and development of skeletal muscle, which will be analyzed in future study.

Current studies had shown that some pathways were involved in muscle growth, development and degradation, including AMPK signaling pathway, ECM-receptor interactions, ErbB signaling pathway, focal adhesion, valine, leucine and isoleucine degradation, ubiquitin mediated proteolysis, mTOR signaling pathway and Glutamatergic Synapse [6269]. However, the report on the role of muscle circRNAs in sheep is limited. Our study clearly showed that circRNAs played an important role in regulating growth and development of muscle in sheep through enriched KEGG pathways and GO pathways. In particular, AMPK signaling pathway, valine, leucine and isoleucine degradation and focal adhesion are significantly enriched.

In conclusion, a number of circRNAs in the muscle of sheep were identified in our study. The differentially expressed circRNAs were screened in prenatal and postnatal muscle of sheep. At the same time, we found several circRNAs involved in the regulation of growth and development of muscle by KEGG pathway analysis. Our study provides valuable resources for circRNAs biology, especially in the muscle of the sheep, and helps for understanding the function of circRNA in sheep.

MATERIALS AND METHODS

Ethics statement

All procedures involving animals were approved by the Animal Care Committee of Shihezi University. The study was performed in accordance with the ethical standards laid down in the 1964 declaration of Helsinki and its later amendments.

Animals

We collected longissimus dorsi muscle samples from three adult kazakh sheep (female) that had been slaughtered at the slaughterhouse. The longissimus dorsi muscle samples of three sheep embryos were carefully collected by surgery from estrus of three pregnant ewes after mating 100 days. Muscle samples were quickly dispensed into the cryopreserved tubes without RNase and immediately placed in liquid nitrogen until RNA isolation. All experiments involving animals were carried out under the protocol approved by the Animal Care Committee of Shihezi University.

Library construction and circRNAs sequencing analysis

According to the manufacturer’s protocol, total RNAs were extracted from frozen tissues after grinding with liquid nitrogen using TRIzol (Invitrogen, CA, USA). Using Bioanalyzer 2100 and RNA6000 Nano Kit (Agilent, CA, USA) detected quantity and purity of total RNAs. Take the same amount of RNAs from three adult sheep and three embryos sheep were respectively pooled to construct library. In order to establish a cDNA library of circRNAs, we take approximately 10 μg of total RNA and use Rnase R (RNR-07250, epicentre) to digest linear RNA. We used IlluminaHiseq 2000/2500 to perform paired-end sequencing on the cDNA library of circRNAs according to the recommended protocol. CircRNAs analysis was performed following the steps in the pipeline (Figure 1). Low quality data was deleted, and reference genome and gene annotation was downloaded from the genome website (http://genome.ucsc.edu/). Index of the reference genome was built using Bowtie v2.0.6, and paired-end clean reads were aligned to the reference genome by using TopHat v2.0.9. CircRNAs candidates were predicted as described previously [40]. Only those who contain more than two independent junction-spanning reads and correspond to the GU / AG intron rules can be determined as a circRNA.

Target site prediction and enrichment analysis

The interaction of the circRNAs and miRNAs were analyzed using miRanda and psRobot [41, 42]. Gene Ontology (GO) analysis was used to host genes of the circRNAs using DAVID software [43], and A KEGG enrichment analysis was used to host genes of circRNAs using KOBAS software [44, 45]. P < 0.05 was used as a criterion for the determination of whether the enrichment analysis was significant

RT-PCR analysis and DNA sequencing

Total RNAs were extracted from sheep muscle using TRIzol (Invitrogen, CA, USA), and 10 μg of purified RNAs were used to synthesize cDNA using RT-PCR kit (Takara, Dalian, China). Ten specific primers (circRNA 0002456, circRNA 0000552, circRNA 0005179, circRNA 0005256, circRNA 0003451, circRNA 0005243, circRNA 0005250, circRNA 0004676, circRNA 0004690 and circRNA 0000666) were used for circRNAs PCR and the primer sequences for these circRNAs were listed in Supplementary Table 1. The PCR products were analyzed by gel electrophoresis and DNA sequencing. We used DNAMAN software to compare the sequence data of PCR products obtained from DNA sequencing with the sheep reference genome and RNA-seq data.

Real-Time RT-PCR analysis

We used Real-Time RT-PCR to detect the expression levels of 10 circRNAs (circRNA 0002456, circRNA 0000552, circRNA 0005179, circRNA 0005256, circRNA 0003451, circRNA 0005243, circRNA 0005250, circRNA 0004676, circRNA 0004690 and circRNA 0000666). Total RNA was treated with RNase R (RNR-07250, epicenter) prior to cDNA synthesis to detect the resistance of circRNAs to RNase R digestion. Total RNA was synthesized by RT-PCR kit (Takara, Dalian, China) and real-time RT-PCR was used to verify the differentially expressed. All Real-Time RT-PCR analyzes were performed using SYBR Green (TaKaRa Biotech, Dalian) according to the manufacturer’s protocol. We selected linear β-actin as an internal reference to normalize the expression of circRNA [46]. Three independent experiments were performed on triplicate samples.

Data availability

The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

Author contributions

CYL, XYL, YY, XYZ, WN and SWH conceived and designed the experiment, analyzed and interpreted the data. CYL, XYL, DWW, WN and SWH wrote the manuscript. QMM, YC, XYZ, RZQ, YY, SY, XXH, ZJL collected sheep skeletal muscle sample. CYL and XYL participated in RNA extraction and RT-PCR analysis. All authors read and approved the final manuscript.

ACKNOWLEDGMENTS

We thank all contributors of the present study.

CONFLICTS OF INTEREST

The authors have declared that no competing interests exist.

FUNDING

This work was supported by the National Natural Science Foundation of China (NSFC) (31660644, 31460037, 31660718, 31360615), Outstanding youth (2015ZRKXJQ01 and 2016BC001), and the Recruitment Program of Global Young Experts (1000Plan).

REFERENCES

1. Sanger HL, Klotz G, Riesner D, Gross HJ, Kleinschmidt AK. Viroids are single-stranded covalently closed circular RNA molecules existing as highly base-paired rod-like structures. Proceedings of the National Academy of Sciences. 1976; 73:3852–3856.

2. Kos A, Dijkema R, Arnberg AC, Vander Meide PH, Schellekens H. The hepatitis delta (δ) virus possesses a circular RNA Nature. 1986; 323:558–560.

3. Zhang XO, Wang HB, Zhang Y, Lu X, Chen LL, Yang L. Complementary sequence-mediated exon circularization. Cell. 2014; 159:134–147.

4. Salzman J, Gawad C, Wang PL, Lacayo N, Brown PO. Circular RNAs are the predominant transcript isoform from hundreds of human genes in diverse cell types. PloS one. 2012; 7:e30733.

5. Jeck WR, Sharpless NE. Detecting and characterizing circular RNAs. Nature biotechnology. 2014; 32:453–461.

6. Jeck WR, Sorrentino JA, Wang K, Slevin MK, Burd CE, Liu J, William F, Marzluff Sharpless NE. Circular RNAs are abundant, conserved, and associated with ALU repeats. RNA. 2013; 19:141–157.

7. Memczak S, Jens M, Elefsinioti A, Torti F, Krueger J, Rybak A, Maier L, Mackowiak SD, Gregersen LH, Munschauer M, Loewer A, Ziebold U, Landthaler M, et al. Circular RNAs are a large class of animal RNAs with regulatory potency. Nature. 2013; 495:333–338.

8. Wang PL, Bao Y, Yee MC, Barrett SP, Hogan GJ, Olsen MN, Dinneny JR, Brown PO, Salzman J. Circular RNA is expressed across the eukaryotic tree of life. PloS one. 2014; 9:e90859.

9. Braun S, Domdey H, Wiebauer K. Inverse splicing of a discontinuous pre-mRNA intron generates a circular exon in a HeLa cell nuclear extract. Nucleic acids research. 1996; 24:4152–4157.

10. Ashwal-Fluss R, Meyer M, Pamudurti NR, Ivanov A, Bartok O, Hanan M, Evantal N, Memczak S, Rajewsky N, Kadener S. circRNA biogenesis competes with pre-mRNA splicing. Molecular cell. 2014; 56:55–66.

11. Chen L, Huang C, Wang X, Shan G. Circular RNAs in eukaryotic cells. Current genomics. 2015; 16:312–318.

12. Huang C, Shan G. What happens at or after transcription: insights into circRNA biogenesis and function. Transcription. 2015; 6:61–64.

13. Hsiao KY, Sun HS, Tsai SJ. Circular RNA–New member of noncoding RNA with novel functions. Experimental Biology and Medicine. 2017; 242:1136–1141.

14. Liu J, Liu T, Wang X, He A. Circles reshaping the RNA world: from waste to treasure. Molecular cancer. 2017; 16:58.

15. Hansen TB, Jensen TI, Clausen BH, Bramsen JB, Finsen B, Damgaard CK, Kjems J. Natural RNA circles function as efficient microRNA sponges. Nature. 2013; 495:384–388.

16. Li Y, Zheng Q, Bao C, Li S, Guo W, Zhao J, Chen D, Gu J, He X, Huang S. Circular RNA is enriched and stable in exosomes: a promising biomarker for cancer diagnosis. Cell research. 2015; 25:981–984.

17. Qu S, Yang X, Li X, Wang J, Gao Y, Shang R, Sun W, Dou KF, Li H. Circular RNA: a new star of noncoding RNAs. Cancer letters. 2015; 365:141–148.

18. Huang G, Zhu H, Shi Y, Wu W, Cai H, Chen X. cir-ITCH plays an inhibitory role in colorectal cancer by regulating the Wnt/β-catenin pathway. PloS one. 2015; 10:e0131225.

19. Yu L, Gong X, Sun L, Zhou Q, Lu B, Zhu L. The circular RNA Cdr1as act as an oncogene in hepatocellular carcinoma through targeting miR-7 expression. PloS one. 2016; 11:e0158347.

20. Li P, Chen S, Chen H, Mo X, Li T, Shao Y, Xiao BX, Guo J. Using circular RNA as a novel type of biomarker in the screening of gastric cancer. Clinica chimica acta. 2015; 444:132–136.

21. Xuan L, Qu L, Zhou H, Wang P, Yu H, Wu T, Wang X, Li Q, Tian L, Liu M, Sun Y. Circular RNA: a novel biomarker for progressive laryngeal cancer.American journal of translational research. 2016; 8:932–939.

22. Molkentin JD, Olson EN. Defining the regulatory networks for muscle development. Current opinion in genetics & development. 1996; 6:445–453.

23. Christ B, Ordahl CP. Early stages of chick somite development. Anatomy and embryology. 1995; 191:381–396.

24. Buffinger N, Stockdale FE. Myogenic specification in somites: induction by axial structures. Development. 1994; 120:1443–1452.

25. Stern HM, Hauschka SD. Neural tube and notochord promote in vitro myogenesis in single somite explants. Developmental biology. 1995; 167:87–103.

26. Stern HM, Brown AM, Hauschka SD. Myogenesis in paraxial mesoderm: preferential induction by dorsal neural tube and by cells expressing Wnt-1. Development. 1995; 121:3675–3686.

27. Óvilo C, Benítez R, Fernández A, Núñez Y, Ayuso M, Fernández AI, Rodríguez C, Isabel B, Rey AI, López-Bote C, Silió L. Longissimus dorsi transcriptome analysis of purebred and crossbred Iberian pigs differing in muscle characteristics. BMC genomics. 2014; 15:413.

28. Abdelmohsen K, Panda AC, De S, Grammatikakis I, Kim J, Ding J, Noh JH, Kim KM, Mattison JA, Cabo R, Gorospe M. Circular RNAs in monkey muscle: age-dependent changes. Aging (Albany NY). 2015; 7:903–910. http://doi.org/10.18632/aging.100834.

29. Ouyang H, Nie Q, Zhang X. P3052 Characterization of Circular RNAs in relation to embryonic muscle development in chicken.Journal of Animal Science. 2016; 94:79.

30. Horak M, Novak J, Bienertova-Vasku J. Muscle-specific microRNAs in skeletal muscle development. Developmental biology. 2016; 410:1–13.

31. Cesana M, Cacchiarelli D, Legnini I, Santini T, Sthandier O, Chinappi M, Tramontano A, Bozzoni I. A long noncoding RNA controls muscle differentiation by functioning as a competing endogenous RNA. Cell. 2011; 147:358–369.

32. Wang K, Long B, Liu F, Wang JX, Liu CY, Zhao B, Zhou LY, Sun T, Wang M, Ying Gong TY, Liu J, Dong YH, Li N, et al. A circular RNA protects the heart from pathological hypertrophy and heart failure by targeting miR-223. European heart journal. 2016; 37:2602–2611.

33. Gong J, Tong Y, Zhang HM, Wang K, Hu T, Shan G, Sun J, Guo AY. Genome-wide identification of SNPs in microRNA genes and the SNP effects on microRNA target binding and biogenesis.Human mutation. 2012; 33:254–263.

34. Rybak-Wolf A, Stottmeister C, Glažar P, Jens M, Pino N, Giusti S, Hanan M, Behm M, Bartok O, Ashwal-Fluss R, Herzog M, Schreyer L, Papavasileiou P, et al. Circular RNAs in the mammalian brain are highly abundant, conserved, and dynamically expressed.Molecular cell. 2015; 58:870–885.

35. Brinkmeier ML, Davis SW, Carninci P, MacDonald JW, Kawai J, Ghosh D, Hayashizakib Y, Lyons RH, Camper SA. Discovery of transcriptional regulators and signaling pathways in the developing pituitary gland by bioinformatic and genomic approaches. Genomics. 2009; 93:449–460.

36. Shan L, Wu Q, Li Y, Shang H, Guo K, Wu J, Wei H, Zhao J, Yu J, Li MH. Transcriptome profiling identifies differentially expressed genes in postnatal developing pituitary gland of miniature pig. DNA research. 2013; 21:207–216.

37. Hansen TB, Jensen TI, Clausen BH, Bramsen JB, Finsen B, Damgaard CK, Kjems J. Natural RNA circles function as efficient microRNA sponges. Nature. 2013; 495:384–388.

38. Cao S, Wei D, Li X, Zhou J, Li W, Qian Y, Wang Z, Li G, Pan X, Lei D. Novel circular RNA expression profiles reflect progression of patients with hypopharyngeal squamous cell carcinoma. Oncotarget. 2017; 8:45367–45379. http://doi.org/10.18632/oncotarget.17488.

39. Zou M, Huang C, Li X, He X, Chen Y, Liao W, Liao Y, Sun J, Liu Z, Zhong L, Bin J. Circular RNA expression profile and potential function of hsa_circRNA_101238 in human thoracic aortic dissection. Oncotarget. 2017; 8:81825–81837. http://doi.org/10.18632/oncotarget.18998.

40. Hansen TB, Venø MT, Damgaard CK, Kjems J. Comparison of circular RNA prediction tools. Nucleic acids research. 2016; 44:e58.

41. Li JH, Liu S, Zhou H, Qu LH, Yang JH. starBase v2.0: decoding miRNA-ceRNA, miRNA-ncRNA and protein–RNA interaction networks from large-scale CLIP-Seq data. Nucleic acids research. 2013; 42:D92-D97.

42. Betel D, Wilson M, Gabow A, Marks DS, Sander C. The microRNA. org resource: targets and expression. Nucleic acids research. 2008; 36:D149–D153.

43. Huang DW, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nature protocols. 2009; 4:44–57.

44. Ghosal S, Das S, Sen R, Chakrabarti J. HumanViCe: host ceRNA network in virus infected cells in human. Frontiers in genetics. 2014; 5:249.

45. Xie C, Mao X, Huang J, Ding Y, Wu J, Dong S, Kong L, Gao G, Li CY, Wei L. KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic acids research. 2011; 39:W316-W322.

46. Peletto S, Bertuzzi S, Campanella C, Modesto P, Maniaci MG, Bellino C, Ariello D, Quasso A, Caramelli M, Acutis PL. Evaluation of internal reference genes for quantitative expression analysis by real-time PCR in ovine whole blood. International journal of molecular sciences. 2011; 12:7732–7747.

47. Zhang XO, Wang HB, Zhang Y, Lu X, Chen LL, Yang L. Complementary sequence-mediated exon circularization. Cell. 2014; 159:134–147.

48. Yang L, Duff MO, Graveley BR, Carmichael GG, Chen LL. Genomewide characterization of non-polyadenylated RNAs. Genome biology. 2011; 12:R16.

49. Liu J, Damon M, Guitton N, Guisle I, Ecolan P, Vincent A, Cherel P, Gondret F. Differentially-expressed genes in pig Longissimus muscles with contrasting levels of fat, as identified by combined transcriptomic, reverse transcription PCR, and proteomic analyses. Journal of agricultural and food chemistry. 2009; 57:3808–3817.

50. Li L, Guo J, Chen Y, Chang C, Xu C. Comprehensive CircRNA expression profile and selection of key CircRNAs during priming phase of rat liver regeneration. BMC genomics. 2017; 18:80.

51. Gruner H, Cortés-López M, Cooper DA, Bauer M, Miura P. CircRNA accumulation in the aging mouse brain. Scientific reports. 2016; 6:38907.

52. Bonizzato A, Gaffo E, Te Kronnie G, Bortoluzzi S. CircRNAs in hematopoiesis and hematological malignancies. Blood cancer journal. 2016; 6:e483.

53. Hanan M, Soreq H, Kadener S. CircRNAs in the brain. RNA biology. 2016; 1–7.

54. Chen J, Yang XJ, Xia D, Wegner J, Jiang Z, Zhao RQ. Sterol regulatory element binding transcription factor 1 expression and genetic polymorphism significantly affect intramuscular fat deposition in the longissimus muscle of Erhualian and Sutai pigs1. Journal of animal science. 2008; 86:57–63.

55. Sun X, Li M, Sun Y, Cai H, Lan X, Huang Y, Bai Y, Qi X, Chen H. The developmental transcriptome sequencing of bovine skeletal muscle reveals a long noncoding RNA, lncMD, promotes muscle differentiation by sponging miR-125b. Biochimica et Biophysica Acta. 2016; 1863:2835–2845.

56. Sun J, Xie M, Huang Z, Li H, Chen T, Sun R, Wang J, Xi Q, Wu T, Zhang Y. Integrated analysis of non-coding RNA and mRNA expression profiles of 2 pig breeds differing in muscle traits. Journal of Animal Science. 2017; 95:1092–1103.

57. Huang M. Zhong Z, Lv M, Shu J, Tian Q, Chen J. Comprehensive analysis of differentially expressed profiles of lncRNAs and circRNAs with associated co-expression and ceRNA networks in bladder carcinoma. Oncotarget. 2016; 7:47186–47200. http://doi.org/10.18632/oncotarget.9706.

58. Kartha RV, Subramanian S. Competing endogenous RNAs (ceRNAs): new entrants to the intricacies of gene regulation. Frontiers in genetics. 2014; 5:8.

59. Xu H, Guo S, Li W, Yu P. The circular RNA Cdr1as, via miR-7 and its targets, regulates insulin transcription and secretion in islet cells. Scientific reports. 2015; 5:12453.

60. Liu Q, Zhang X, Hu X, Dai L, Fu X, Zhang J, Ao Y. Circular RNA related to the chondrocyte ECM regulates MMP13 Expression by functioning as a MiR-136 ‘Sponge’ in human cartilage degradation. Scientific reports. 2016; 6:22572.

61. Zheng Q, Bao C, Guo W, Li S, Chen J, Chen B, Luo Y, Lyu D, Li Y, Shi G, Liang L, Gu J, He X, et al. Circular RNA profiling reveals an abundant circHIPK3 that regulates cell growth by sponging multiple miRNAs. Nature communications. 2016; 7:11215.

62. Chen ZP, Stephens TJ, Murthy S, Canny BJ, Hargreaves M, Witters LA, Kemp BE, McConell GK. Effect of exercise intensity on skeletal muscle AMPK signaling in humans. Diabetes. 2003; 52:2205–2212.

63. Flück M, Carson JA, Gordon SE, Ziemiecki A, Booth FW. Focal adhesion proteins FAK and paxillin increase in hypertrophied skeletal muscle. American Journal of Physiology-Cell Physiology. 1999; 277:C152-C162.

64. Llovera M, Garcia-Martinez C, Agell N, Lopez-Soriano FJ, Argiles JM. Muscle wasting associated with cancer cachexia is linked to an important activation of the ATP-dependent ubiquitin-mediated proteolysis. International journal of cancer. 1995; 61:138–141.

65. Nair KS, Schwartz RG, Welle S. Leucine as a regulator of whole body and skeletal muscle protein metabolism in humans. American Journal of Physiology-Endocrinology And Metabolism. 1992; 263:E928-E934.

66. Bodine SC, Stitt TN, Gonzalez M, Kline WO, Stover GL, Bauerlein R, Zlotchenko E, Scrimgeour A, Lawrence JC, Glass DJ, Yancopoulos GD. Akt/mTOR pathway is a crucial regulator of skeletal muscle hypertrophy and can prevent muscle atrophy in vivo. Nature cell biology. 2001; 3:1014–1019.

67. Brunelli G, Spano P, Barlati S, Guarneri B, Barbon A, Bresciani R, Pizzi M. Glutamatergic reinnervation through peripheral nerve graft dictates assembly of glutamatergic synapses at rat skeletal muscle.Proceedings of the National Academy of Sciences of the United States of America. 2005; 102:8752–8757.

68. Lebrasseur NK, Coté GM, Miller TA, Fielding RA, Sawyer DB. Regulation of neuregulin/ErbB signaling by contractile activity in skeletal muscle. American Journal of Physiology-Cell Physiology. 2003; 284:C1149–C1155.

69. Li Y, Xu Z, Li H, Xiong Y, Zuo B. Differential transcriptional analysis between red and white skeletal muscle of Chinese Meishan pigs. International journal of biological sciences. 2010; 6:350–360.


Creative Commons License All site content, except where otherwise noted, is licensed under a Creative Commons Attribution 3.0 License.
PII: 21835