Next Article in Journal
Discovery of GABA Aminotransferase Inhibitors via Molecular Docking, Molecular Dynamic Simulation, and Biological Evaluation
Next Article in Special Issue
LncRNA and Protein Expression Profiles Reveal Heart Adaptation to High-Altitude Hypoxia in Tibetan Sheep
Previous Article in Journal
YUCCA2 (YUC2)-Mediated 3-Indoleacetic Acid (IAA) Biosynthesis Regulates Chloroplast RNA Editing by Relieving the Auxin Response Factor 1 (ARF1)-Dependent Inhibition of Editing Factors in Arabidopsis thaliana
Previous Article in Special Issue
Novel Transcriptomic Interactomes of Noncoding RNAs in the Heart under Altered Thyroid Hormonal States
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
IJMS
Volume 24
Issue 23
10.3390/ijms242316991
ijms-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Transcriptome Analysis of mRNA and lncRNA Related to Muscle Growth and Development in Gannan Yak and Jeryak

by
Yali Wei
,
Dashan Guo
,
Yanbin Bai
,
Zhanxin Liu
,
Jingsheng Li
,
Zongchang Chen
,
Bingang Shi
,
Zhidong Zhao
*,
Jiang Hu
*,
Xiangmin Han
,
Jiqing Wang
,
Xiu Liu
,
Shaobin Li
and
Fangfang Zhao
Gansu Key Laboratory of Herbivorous Animal Biotechnology, College of Animal Science and Technology, Gansu Agricultural University, Lanzhou 730070, China
*
Authors to whom correspondence should be addressed.
Int. J. Mol. Sci. 2023, 24(23), 16991; https://doi.org/10.3390/ijms242316991
Submission received: 23 October 2023 / Revised: 22 November 2023 / Accepted: 26 November 2023 / Published: 30 November 2023
(This article belongs to the Special Issue New Horizon for Non-coding RNAs)
Browse Figures
Versions Notes

Abstract

:
The production performance of Jeryak, resulting from the F1 generation of the cross between Gannan yak and Jersey cattle, exhibits a significantly superior outcome compared with that of Gannan yak. Therefore, we used an RNA-seq approach to identify differentially expressed mRNAs (DEMs) and differentially expressed lncRNAs (DELs) influencing muscle growth and development in Gannan yaks and Jeryaks. A total of 304 differentially expressed lncRNAs and 1819 differentially expressed mRNAs were identified based on the screening criteria of |log 2 FC| > 1 and FDR < 0.05. Among these, 132 lncRNAs and 1081 mRNAs were found to be down-regulated, while 172 lncRNAs and 738 mRNAs were up-regulated. GO and KEGG analyses showed that the identified DELs and DEMs were enriched in the entries of pathways associated with muscle growth and development. On this basis, we constructed an lncRNA–mRNA interaction network. Interestingly, two candidate DELs (MSTRG.16260.9 and MSTRG.22127.1) had targeting relationships with 16 (MYC, IGFBP5, IGFBP2, MYH4, FGF6, etc.) genes related to muscle growth and development. These results could provide a basis for further studies on the roles of lncRNAs and mRNAs in muscle growth in Gannan yaks and Jeryak breeds.

1. Introduction

Yak is a specialty germplasm resource mainly distributed in Tibet, Sichuan, Yunnan, and Gansu and is an important resource for local herders to survive. However, the slow growth of yak, slow muscle development, low fat deposition efficiency, and poor meat tenderness significantly impede the production performance and economic benefits of yak. Jeryak is a hybrid resulting from the crossbreeding of Jersey cattle and Gannan yak. Phenotypic analysis revealed that under the same feeding conditions, Jeryak exhibited significantly superior production performance in terms of body weight, shoulder height, body length, and chest circumference compared with Gannan yak [1,2]. Hence, it is of great relevance to study the molecular mechanism of muscle growth and development in Jeryak and Gannan yak.
The developmental process of skeletal muscle determines the yield and quality of livestock and poultry meat products [3]. The growth and development of muscle are determined by the number of pre-birth fibers of the myofiber type, while the hypertrophy process of the area of myofibers mainly occurs in the post-natal period [4,5,6]. Muscle fibers are the fundamental building blocks of muscles, and their properties have a direct impact on animal meat production and quality [7]. The diameter, density, area, and type of muscle fibers are crucial indicators that reflect the characteristics of muscle fibers [8,9]. It has been confirmed that an increase in muscle fiber area and diameter leads to a reduction in muscle tenderness and tensile strength; on the contrary, a decrease in muscle fiber area and diameter leads to increased meat tenderness [10]. In addition to myofiber area and diameter, the myofiber type also affects muscle growth and development [11,12]. Furthermore, the expression of specific mRNAs in skeletal muscle tissues, including IGFs [13], PAXs [14], MRFs [15], and MEF2 [16], exerts a significant influence on both muscle development and meat quality.
Long non-coding RNA (lncRNA), as one key member of non-coding RNAs, has been extensively investigated in recent years to substantiate its regulatory role in muscle growth and development. Among the various approaches employed to explore lncRNA targets, cis-regulation and trans-regulation of genes have emerged as commonly utilized methods [17,18]. For instance, in skeletal muscle, lncMAAT exhibits cis-regulation to enhance the expression of the Mbnl1 gene while concurrently exerting trans-regulation control over SOX6 to suppress miR-29b expression through a regulatory module, thereby delaying muscle atrophy [19]. Additionally, lncRNAs serve as molecular sponges for microRNAs, forming a bond with miRNAs to alleviate the suppressive impact of miRNAs on target genes. Through its role of molecular sponge, lncIRS1 absorbs the miR-15 family, thereby activating the IGF1-PI3K/AKT pathway and exerting an impact on muscle atrophy [20]. An additional instance involves the functioning of lncMD, which enhances bovine myoblast differentiation by increasing the expression of IGF2 through the adsorption of miR-125b [21]. It has also been discovered that lncRNAs exert a regulatory function in the process of muscle regeneration following injury. For that instance, lncRNA MAR1 exhibits a potential positive impact on skeletal muscle hypertrophy and may also contribute to mitigating aging-related muscle atrophy [22]. During the pre-natal and early post-neonatal stages in porcine skeletal muscle, lncRNA MEG3 demonstrates elevated levels of expression, and it plays a crucial part in the development and upkeep of early skeletal muscle [23]. Many researchers have identified the crucial genes that influence muscle development through the comprehensive sequencing of lncRNAs across distinct species. For example, Yang et al. sequenced the longest dorsal muscles of Leiqiong and Lufeng cattle and found 12 genes regulated by DE lncRNA cis; they suggested that the DKK1 and SIGMAR1 genes may play a crucial role in muscle development regulation [24]. However, limited research has been conducted regarding the influence of lncRNAs on the process of muscle growth and development in Jeryaks.
Therefore, in this study, adult Gannan yaks and Jeryaks were selected to investigate the regulatory mechanisms of the differences in muscle fiber tissue characteristics between Gannan yaks and Jeryaks by using HE staining, and fast and slow muscle fluorescence staining. Further, RNA-seq was used to study the expression of lncRNA in the dorsal longissimus muscle tissues of the two breeds, and the mRNA–lncRNA interaction network was constructed to identify the key genes for muscle growth and development. These results could provide a basis for further excavation of important lncRNAs in the skeletal muscle development of Jeryaks, and subsequent characterization and functional studies.

2. Results

2.1. Observable Disparities in Muscle Fibers among Various Breeds of Beef Cattle

First, we performed fast and slow myofiber fluorescence staining (Figure 1A) and HE staining (Figure 1B) on Gannan yaks and Jeryaks. The results showed that there were highly significant differences in myofiber characteristics between Gannan yaks and Jeryaks. IPP (Image-Pro Plus 6.0) software analysis showed that the average myofiber area, myofiber diameter, and myofiber length of Gannan yaks were significantly greater than those of Jeryaks (p < 0.01), The ratio of fast-to-slow muscle fibers and the ratio of fast muscle fibers to muscle fiber area were significantly greater than in Jeryaks (p < 0.01), while the ratio of slow muscle fibers to muscle fiber area was significantly smaller than that of Jeryaks (p < 0.01) (Figure 1C).

2.2. Sequencing Data Summary

For this study, six sequencing libraries were constructed, and the cDNA libraries of Gannan yak (M) and Jeryak (P) were sequenced using the Illumina NovaSeq 6000 sequencing platform. On average, each group obtained raw data of 86,609,737 and 91,797,352 reads. On average, each group obtained 86,111,288.67 (99.41%) and 88,349,502.1 (99.50%) high-quality sequences (clean reads) after eliminating the splice sequences and considering N (modular base) reads and reads with low quality at the 3’ end of the reads. Additionally, all six sequenced libraries exhibited Q30 values exceeding 91%, indicating excellent sequencing quality. Furthermore, the GC/AT content was approximately 60% in alignment with the theoretical genome distribution (Table 1). The comparison between pristine reads and the yak’s reference genome revealed a genome comparison rate of over 70% for each sample, demonstrating a favorable comparison effect (Table 2). More than 55% of the reads were uniquely mapped to the genome, indicating a high level of accuracy and reliability in the pristine reads, which are suitable for further analysis. We identified a substantial number of lncRNAs through comprehensive screening and filtering. PCA analysis showed good reproducibility within sample groups and large differences between groups (Figure 2A).
Among these, a total of 3082 lncRNA were recently identified, exhibiting an average size of 1165 base pairs. The largest percentage (30.61%) was observed in the lncRNA category, with a length ranging from 1000 to 1500 base pairs (Figure 2C), while approximately 61.77% of lncRNA consisted of two exons (Figure 2B). Additional categorization of the identified lncRNA indicated that intergenic lncRNA constituted 53.35%, while sense lncRNA made up 13.51%, antisense lncRNA comprised 9.32%, bidirectional lncRNA represented 6.5%, intronic lncRNA constituted 1.44%, and other lncRNAs accounted for 15.87% of the total (Figure 2D).

2.3. Expression of Long Non-Coding RNA (lncRNA) and Messenger RNA (mRNA) at the Genomic Level

With an analysis of gene structure and expression, we aimed to identify discrepancies between lncRNA and mRNA obtained with RNA-seq. Based on the FPKM value, our findings revealed comparatively elevated expression levels for both mRNA and lncRNA (Figure 3A). The average number of exons in mRNA was higher (23.04) compared with lncRNA (average 2.8). In addition, 86.78% of messenger RNA contained four or more exons, whereas 81.35% of long non-coding RNA had three or fewer exons (as shown in Figure 3B). The length of mRNA (averaging 2291 bp) was significantly greater than that of lncRNA (averaging 1363 bp) (Figure 3C). In the longissimus dorsi muscle of Gannan yak and Jeryak, the measurement was conducted on a combined count of 3743 lncRNAs and 16,143 mRNAs. Using Gannan yak as the reference group, a total of 304 differentially expressed lncRNAs (DELs) and 1819 differentially expressed mRNAs (DEMs) were identified. Among these, 132 lncRNAs and 1081 mRNAs were found to be down-regulated, while 172 lncRNAs and 738 mRNAs were up-regulated (Figure 3D,G). Furthermore, a separate cluster analysis was conducted on DELs and DEMs, revealing that the identical sets of DELs and DEMs were grouped together, thus confirming the precision and dependability of the samples (Figure 3E,F). Tables S1 and S2 show the DEMs and DELs, respectively.

2.4. Performing GO and KEGG Analyses on Differentially Expressed Genes (DEGs)

The 783 up-regulated DEMs in Class-up were significantly enriched in 1480 GO entries, including 1264 GO-BP, 167 GO-MF, and 49 GO-CC (Table S3). In Figure 4A, the top 60 GO entries that were significantly enriched are depicted. The GO entry that was enriched to the greatest extent was cell adhesion (GO:0007155), and there were 17 other GO entries that pertained to the growth and development of muscles, including skeletal muscle cell differentiation (GO:0035914). There were 1081 down-regulated DEMs in Class-down that were significantly enriched in 470 GO entries (317 GO-BPs, 62 GO-MFs, and 91 GO-CCs) (Table S4). In Figure 4C, the top 60 GO entries that were significantly enriched are depicted. The intracellular ribonucleoprotein complex (GO:0030529) emerged as the most prominently enriched GO term, while muscle growth and development were represented by four GO entries, including muscle cell migration (GO:0014812). There were 36 KEGG pathways significantly enriched by the 738 up-regulated DEMs in Class-up, with the most significantly enriched pathway being alcoholism, followed by Systemic lupus erythematosus (Table S5). Furthermore, muscle growth and development were linked to five KEGG pathways, specifically ECM–receptor interaction, cAMP signaling pathway, MAPK signaling pathway, focal adhesion, and cell adhesion molecules. There were 25 KEGG pathways significantly enriched by the 1081 down-regulated DEMs in Class-down, with the ribosome pathway being the most significantly enriched (Table S6). Muscle growth and development were linked to two KEGG pathways, specifically the MAPK signaling pathway (fly) and the regulation of actin cytoskeleton. Figure 4B,D exhibit the notable enhancement in KEGG pathways in Class-up, along with the top 20 KEGG pathways that exhibited significant enrichment in Class-down.

2.5. Cis-/Trans-Regulation of the lncRNA Target Genes

Long non-coding RNA can control the expression of adjacent messenger RNA. By examining mRNA located within a 100 kb vicinity of the lncRNA and protein-coding genes, we successfully detected a total of 150 differentially expressed lncRNAs and 125 differentially expressed mRNAs (Table S7). Subsequently, the potential target genes of lncRNA were subjected to GO and KEGG analyses (Figure 5A,B). A total of 1224 GO entries (Table S8) were identified as significantly enriched (p < 0.05). Among them, the three most significantly enriched GO entries were cell (GO:0005623), cell part (GO:0044464), and intracellular (GO:0005622). Furthermore, we identified 90 KEGG pathways that exhibited significant enrichment (Table S9). Among these, the MAPK signaling pathway, the oxytocin signaling pathway, and Human cytomegalovirus infection emerged as the top three pathways with the highest level of enrichment.
Afterward, we detected 304 differentially expressed lncRNAs (DELs) and 1796 differentially expressed mRNAs (DEMs) by applying a threshold of r > ±0.95 to investigate the target genes of trans-lncRNA (Table S10). Figure 5C,D show the GO and KEGG results for DEL target genes. There were 1005 GO terms that showed significant enrichment (Table S11). Among these, ribosome (GO:0005840), the structural component of the ribosome (GO:0003735), and regulation of response to injury (GO:1903034) were the top three most significantly enriched GO entries. Out of 1005, 17 findings (p < 0.05) were associated with the enhancement and formation of muscles, including the migration of muscle cells. Additionally, we discovered 20 KEGG pathways that exhibited significant enrichment in KEGG (Table S12). Among these, the three most notably enriched pathways were ribosome, alcoholism, and fluid shear stress and atherosclerosis. Muscle growth and development were linked to four KEGG pathways, including the MAPK signaling pathway and ECM–receptor interaction.

2.6. Building a Network of Interactions between lncRNA and mRNA

Here, we selected eight signaling pathways associated with muscle growth and development, among which the genes enriched in the MAPK signaling pathway were MEF2C, FGF6, GADD45A, FGFR4, and MYC. By conducting a literature search, we identified DEMs that impact muscle growth and development. We then obtained DELs and their corresponding target genes from cis and trans prediction results. Finally, we integrated these findings with DEMs from sequencing results to create a DEL–DEM interaction network. The network of DEL–DEM co-expression included 8 DELs and 55 trans-targets (Figure 6B), as well as 3 DELs and 3 cis-targets (Figure 6C). Several lncRNAs (MSTRG.16260.9, MSTRG.1959.1, ENSBGRT00000010172, and MSTRG.26608.2) were found to be up-regulated, while others (MSTRG.22127.1, MSTRG.9040.1, MSTRG.6643.2, and MSTRG.18522.2) were down-regulated, potentially linked to muscle growth and development. The results of the co-expression network indicated that certain DELs interacted with multiple DEMs. For instance, a specific DEL (MSTRG.16260.9 and MSTRG.22127.1) exhibited targeting associations with 16 genes associated with muscle growth and development. This suggests that lncRNA and mRNA are part of the essential non-coding/coding RNA that plays a crucial role in regulating the development of muscle growth.

2.7. Verification of Sequencing Results

We randomly chose seven lncRNAs and eight mRNAs from DEMs and DELs, respectively, and confirmed the precision of the sequencing findings with real-time quantitative PCR. The qPCR data in this study (Figure 7A,B) were found to be in agreement with the expression pattern observed in the sequencing data.

3. Discussion

The meat production and meat quality of Jeryak are superior to that of Gannan yak. With the aim of studying the mechanism of rapid growth in Jeryak, we explored the role of genes and lncRNA in the growth mechanism from the perspectives of myofiber characterization and transcriptome sequencing. Myofibers are the most basic units that make up a muscle, and their characteristics have a direct influence on the animal’s meat production performance and meat quality [8,9]. Meat tenderness was found to be greater when the muscle fibers had a smaller diameter and area [10]. In the present study, we found that the average muscle fiber area, average muscle fiber diameter, and average muscle fiber length of Jeryak were significantly smaller than those of Gannan yak. Moreover, the quality of muscle was found to be superior in cases where the proportion of slow-twitch muscle fibers exceeded that of fast-twitch muscle fibers [11,12]. We also observed that the area of fast muscle fiber was significantly increased in Gannan yaks compared with Jeryaks, while the area of slow muscle fibers exhibited a significant decrease relative to Jeryaks. These variations may serve as pivotal factors contributing to the disparities in post-partum meat production performance and meat quality between Gannan yaks and Jeryaks, thereby establishing the foundation for investigating their molecular regulatory mechanisms in this study.
In general, the regulation of muscle growth and development is orchestrated by core genes and signaling pathways. lncRNA, which is characterized by a length exceeding a certain value, has emerged as crucial a regulator of gene function and expression levels in this context [25,26,27,28]. According to recent research, lncRNA plays a crucial role in the regulation of skeletal muscle growth and development [29,30,31]. In the present study, a total of 304 DELs and 1819 DEMs were identified, among which 132 lncRNAs and 1081 mRNAs were found to be down-regulated, while 172 lncRNAs and 738 mRNAs were up-regulated. In particular, DEMs exhibited significant enrichment of GO terms associated with muscle growth and myoblast differentiation, including processes such as muscle cell migration. Based on the gene function analysis, we identified five up-regulated genes related to muscle development: MYH4, FRZB, IGFBP5, MYH3, and FOXO1. Additionally, three down-regulated genes associated with muscle development were also identified: IGFBP3, SIXI, and MEF2C.
IGFBP3 and IGFBP5 are part of the Insulin-like growth factor binding proteins (IGFBPs) family. During the process of myoblast differentiation, the inhibition of muscle cell differentiation occurs due to the presence of IGFBPs. These IGFBPs hinder the expression of IGFs, resulting in a decrease in protein synthesis. Moreover, the overexpression of IGFBP5 can either enhance or hinder the IGF-mediated myoblast differentiation or survival by facilitating or impeding the action of IGFs [32,33]. Phosphorylation of the FOXO family occurs through the PI3K/AKT pathway in response to insulin, resulting in exclusion from the nucleus and subsequent inactivation [34]. Activation of FOXO1 inhibits C2C12 differentiation and prevents AKT-induced myotube specialization [34]. Accili et al. similarly demonstrated that FOXO1 exerts inhibitory effects on myofibroblast differentiation [35]. Furthermore, it has been demonstrated that FOXO1 facilitates the degradation of muscle wasting [36]. Controlling the expression of genes related to the growth of skeletal muscles, such as MYOG, MYHC, MYOD, IGF1R, and INSR, is a vital function carried out by the Six1 gene [37,38]. Six1 has a substantial impact on the growth and development of skeletal muscles from the embryonic stage to post-natal periods, as stated in this regulation [39,40]. Several studies carried out in zebrafish have demonstrated that the deletion of or reduction in Six1 expression exerts a substantial impact on skeletal muscle development [41]. Additionally, it has been observed that knockdown of either Six1a or Six1b in zebrafish leads to a significant increase in apoptosis or cell death within embryonic myoblasts [39,40]. This highlights the critical role of Six1 in promoting myoblast survival and normal function. MEF2C is commonly present in muscle tissues, particularly in skeletal and cardiac muscles [42]. MEF2C was the first DNA-binding factor identified to possess muscle characteristics and possesses the ability to regulate muscle growth by directly binding to promoters of muscle-specific genes such as MCK and MYHC. Additionally, MEF2C plays a crucial role in regulating skeletal myogenesis through direct interaction with MYHC [43].
The analysis results of the KEGG pathways revealed significant enrichment of DEMs in the Class-up category in 36 pathways. Among these pathways, five were found to be associated with growth and development, including ECM–receptor interaction, cAMP signaling pathway, and several others. In the Class-down category, DEMs were notably enhanced in 25 pathways related to muscle growth and development, such as the MAPK signaling pathway (fly) and regulation of actin cytoskeleton.
The metabolic adaptation of skeletal muscle to exercise is facilitated by the MAPK signaling pathway, which plays a significant role in this process. During exercise, the MAPK signaling pathway enhances insulin-dependent glucose uptake, increases oxidative metabolism, and promotes oxidative phosphorylation of mitochondria in muscle [44,45]. Moreover, this pathway is involved in the regulation of cell cycle arrest in mature myoblasts, a crucial step in initiating muscle differentiation [46]. In skeletal muscle cells, MEF2C activates the MAPK P38 pathway to induce differentiation of skeletal muscle [47]. In this study, we demonstrated that MEF2C, GADD45A, FGF6, FGFR4, and MYC may play crucial roles in regulating myogenesis in Gannan yak and Jeryak. These genes could be potential candidates for further enhancing the meat production performance of Gannan yak. The extracellular matrix (ECM), widely distributed in animal tissues, is a complex network synthesized and secreted by cells. It serves as both a physical scaffold for embedded cells and an essential support system for muscle differentiation [48,49]. Moreover, the ECM not only provides a structural framework for cell embedding but also regulates various cellular processes, such as growth, differentiation, migration, and survival. Additionally, it acts as a crucial reservoir of growth factors, facilitating intercellular signaling through attachment and release of numerous growth and signaling factors [50]. Focal adhesion plays a pivotal role in the regulation of skeletal muscle development, being intricately involved in several crucial signaling pathways, such as the WNT, MAPK, and PI3K-AKT intracellular pathways [51]. Moreover, it collaborates with the ECM to finely regulate intracellular pathways associated with cell movement, proliferation, and differentiation [52].
To better understand the potential functions of DELs, we performed GO and KEGG analyses on their cis-target and trans-target genes. In this study, we successfully identified 304 DELs and predicted their target genes, which were mainly enriched in GO entries related with muscle growth and development, particularly muscle structure development. These compelling findings strongly imply that DELs may exert a substantial influence on the intricate processes underlying muscle growth and development. Similarly, some lncRNAs have been recognized for their crucial involvement in the regulation of muscle growth and development. The KEGG analysis conducted on the target genes of DELs also revealed a significant association with key pathways, such as ECM–receptor interaction and the MAPK signaling pathway, which are known to play vital roles in muscle growth and development (p < 0.05). In addition, it has been observed that the Rap1 signaling pathway interacts with the β-adrenergic signaling pathway, thereby making a substantial contribution to skeletal muscle growth and development [53]. Kosuru et al. investigated the expression and localization of the Rap1 protein during different developmental stages of skeletal muscle and discovered that the Rap1 protein accumulates at the neuromuscular and tendon junctions in the early and late stages of myogenesis, respectively [54]. To screen lncRNAs and mRNAs that may be related to muscle growth and development more intuitively and deeply, we constructed a DEL–DEM co-expression network. Here, we selected eight signaling pathways associated with muscle growth and development, among which the genes enriched in the MAPK signaling pathway were MEF2C, FGF6, GADD45A, FGFR4, and MYC. The genes in these pathways that were present in the background were identified by intersecting them with the differentially expressed mRNAs that had been previously sequenced in the transcriptome [55]. These intersecting mRNAs and differential lncRNAs were used to construct an interaction network. In the interaction network, some DELs were found to interact with multiple DEMs. For example, two candidate DELs (MSTRG.16260.9 and MSTRG.22127.1) were found to have a targeting relationship with 16 genes related to muscle growth and development (MYC, IGFBP5, IGFBP2, MYH4, FGF6, FGFR4, etc.). We found that FGF (fibroblast growth factor) promotes proliferation and inhibits differentiation of skeletal muscle satellite cells [56,57]. FGF6, which is expressed at high levels, enhances the survival of C2C12 and primary myoblasts by activating Cyclin D1 through FGFRI. Down-regulation of FGF6 leads to enhanced myoblast differentiation via FGFR4-dependent activation of ERK1/2, resulting in increased expression of MYHC, MYOD, and MYOG [58]. In our study, FGF6 was found to be a down-regulated DEM that targets the up-regulated MSTRG.16260.9.
In addition, a large number of studies have found that lncRNA are important cis-acting regulatory elements that can regulate the expression of protein-coding genes (examples of such lncRNAs include Malat1 [59], Yam1 [60], and lncRNA-NEF [61]) through cis-acting. Cai et al. discovered a protein-coding gene (Six1 gene) situated less than 100 kb away from lncRNA-Six1 [62]. According to functional validation of lncRNA-Six1, it regulates the expression of the Six1 gene and promotes cell proliferation and muscle growth [62]. A protein-coding gene (FGF6) was discovered to be situated less than 100 kb away from MSTRG.26306.4 and exhibited differential expression in Gannan yaks compared to Jeryaks. However, there are some limitations to our study. We focused on identifying differentially expressed genes in the constructed network but did not perform adequate functional analyses. Therefore, it is important to validate our findings with in vitro and in vivo functional studies. Further studies should include validation of the functions of MSTRG.16260.9 and MSTRG.22127.1 in trans-acting, as well as cis-acting MSTRG.26306.4 and the target gene FGF6, to explore their modes of action and specific roles in muscle growth and development.

4. Materials and Methods

4.1. Preparation of Animals and Tissues

The three Gannan yaks (Gannan male yak × Gannan female yak) and three Jeryaks (Jersey male cattle × Gannan female yaks) involved in this study were sourced from Hezuo City, Gannan Tibetan Autonomous Prefecture, Gansu Province. All animals were 4 years old and healthy. They were fed and watered ad libitum under identical conditions. The Experimental Animal Ethics Committee of Gansu Agricultural University approved the selection of adult bulls with comparable body weights, which were subsequently slaughtered in an abattoir. After slaughter, longest dorsal muscle samples were collected and stored in a refrigerator at −80 °C.

4.2. Staining of Muscle Tissue Using Hematoxylin and Eosin, as Well as Fluorescence Staining of Fast/Slow Muscle Fibers

To enhance the observation of the histological morphology of the muscle tissue, staining techniques, including HE staining and fluorescent staining for fast/slow muscle, were employed. The method by Zhang [63] was used to carry out HE staining. Initially, the muscle tissue was extracted from the fixative, followed by the preparation of paraffin sections. The only procedures conducted were embedding, sectioning, staining with hematoxylin and eosin, dehydration, mounting, and image acquisition. Next, the tissue sections were placed inside a container filled with EDTA antigen recovery buffer. Once cooled, the slides were immersed in PBS (pH 7.4) and cleansed by shaking them with a decolorizing shaker. Subsequently, a drop of blocking serum was introduced, followed by overnight application of a primary antibody and subsequently a secondary antibody. Subsequently, the nuclei were re-stained by adding DAPI staining solution drop by drop, followed by incubation of the cells at room temperature in the absence of light for a duration of 10 min. They were then cleaned using PBS. We added autofluorescence quencher, sealed the slides, observed the slices under a fluorescence microscope, and collected images. DAPI-stained cell nuclei appeared blue, while fast muscle fibers were indicated by the color red, and slow muscle fibers were indicated by the color green.
CaseViewer2.2 scanning browser software was utilized to select target regions for data analysis, with eight regions per slice. Afterward, the pictures were brought into image-pro Plus 6.0 software for image analysis to determine their diameter, area, and type. The statistical analyses were conducted utilizing SPSS 25.0 software.

4.3. Extraction of RNA, Construction of Library, and RNA Sequencing

Total RNA was extracted from muscle tissue using the Trizol reagent kit (Invitrogen, Carlsbad, CA, USA), and RNA integrity was precisely assessed using an Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA). The analysis of RNA integrity and the detection of DNA contamination was conducted by employing electrophoresis on 1% agarose gel. Following the isolation of complete RNA, the ribosomal RNA was eliminated, and subsequently, the RNA was fragmented into smaller pieces. Random oligonucleotides were used as primers to synthesize the initial cDNA strand, while the second strand was synthesized using the DNA polymerase I system. Eventually, the process of PCR amplification was carried out, resulting in the acquisition of the library. PCR was used to amplify the libraries, which were then sequenced using Illumina NovaSeq 6000 by Gene Denovo Biotechnology (located in Guangzhou, China). The library-sequencing parameter was PE150 (paired-end 150 bp) (located in Guangzhou, China).

4.4. Data Quality Control and Reference Genome Alignment

For obtaining reads of excellent quality and free from impurities, it is recommended to utilize fastp (version 0.18.0) for filtering the reads. The settings are as follows: (1) eliminate reads that have adapters; (2) discard reads that have over 10% unidentified nucleotides (N); (3) eliminate low-quality reads that have more than 50% bases with low quality (Q value ≤ 20).
Using HiSAT2 software (version 2.1.0) [64], the high-quality clean reads of each sample were compared with the yak reference genome (LU_Bosgru_v3.0) to determine the remaining matches. Genome and gene alignment were used to evaluate all samples. The software Stringtie (version 1.3.4) [65] was utilized for transcript reconstruction.

4.5. lncRNA Identification and lncRNA Target Gene Prediction

Furthermore, for acquiring lncRNA data and conducting transcript screening, transcripts falling under categories “i”, “u”, “x”, and “o” were initially chosen. Subsequently, transcripts with lengths below 200 nt and exon counts less than 2 were eliminated. Following that, the software applications CPC (version 0.9-r2) [66], CNCI (version 2) [67], and Feelnc (version v0.2) [68] were utilized to predict the coding ability. Ultimately, the intersection of the outcomes from the three software applications was considered the identified lncRNA. There are two main approaches to predicting lncRNA targets: cis-acting and trans-acting [69]. To predict cis-acting effects, the co-location threshold was established as 100 kb in both the upstream and downstream regions of lncRNA in order to identify cis-target genes. Afterwards, the trans-target genes were filtered by examining the association between lncRNA and mRNA expression in the samples (Pearson correlation coefficient, r > 0.95).

4.6. Analysis of Genes with Differential Expression and Pathways

The Cuffdiff program [70] was utilized to conduct a comparative analysis of the 6 transcripts collectively, with FPKM values serving as indicators of gene or transcript expression levels. DESeq software (version 1.24.0) [71,72]. was utilized to normalize the counts of each sample. DEMs and DELs were identified by screening for mRNAs and lncRNAs with absolute log2FC values greater than 1 and FDR values less than 0.05.
To comprehend the role of lncRNA target genes, an analysis of GO (Gene Ontology) enrichment was conducted using the top GO. The hypergeometric distribution method was utilized to calculate significantly enriched GO entries, with a criterion of p < 0.05 for determining significant enrichment. Cluster profiler software (version 4.10.0) was used to conduct the pathway enrichment analysis of KEGG (Kyoto Encyclopedia of Genes and Genomes), with a significance criterion of p < 0.05.

4.7. DEL–DEM Co-Expression Network Construction

We anticipated lncRNA cis- and trans-target genes, intersected the forecast outcomes with mRNA sequencing outcomes, and subsequently examined the literature to filter the genes associated with muscle growth and development. Finally, we used Cytoscape (v3.6.0) software to map the target relationship between lncRNA and its target genes.

4.8. Verification of DEMs and DELs Using qRT-PCR

To verify the reliability of the sequencing results, 7 DELs and 8 DEMs were randomly selected for qRT-PCR. Primer 5.0 was utilized to create the primers, with GAPDH (glyceraldehyde-3-phosphate dehydrogenase) serving as the internal reference gene. The primer sequence details can be found in Table 3 and Table 4. First, total RNA was extracted from the samples using TRIzol reagent (Invitrogen, Carlsbad, CA, USA). Second, we used a reverse transcription kit (TransGen Biotech, Beijing, China) and a real-time fluorescence quantification kit (TransGen Biotech, Beijing, China) for cDNA synthesis and real-time fluorescence quantitative PCR. All qRT-PCR reactions were performed with an ABI Prism 7500 real-time PCR system (Thermo Fisher Scientific, Waltham, MA, USA). The total reaction system was 20 μL: 10 μL of 2× PerfectStart Green qPCR SuperMix, 0.4 μL of Passive Reference Dye (50×) (optional), 0.8 μL of primers, 2 μL of cDNA templates, and 6.8 μL of RNase-free water. The relative quantification of lncRNAs was determined utilizing the 2−ΔΔCt approach and compared to the findings from transcriptome sequencing.

5. Conclusions

In the overview, we analyzed DELs and DEMs in the longest dorsal muscle tissues using Gannan yak as a control group. The identified DELs and DEMs were enriched in the pathways related to muscle growth and development. Based on these findings, we constructed a network of interactions between DELs and DEMs and discovered that MSTRG.16260.9 and MSTRG.22127.1 might play a role in regulating myogenesis in Gannan yak and Jeryak. These lncRNAs could potentially be utilized as candidates to further improve the production performance of yak.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ijms242316991/s1.

Author Contributions

Conceptualization, Y.W. and Z.Z.; Methodology, Y.B., Z.C., B.S. and J.L.; Validation, Y.W., D.G., Z.L., Y.B. and J.L.; Formal analysis, Y.W., F.Z. and D.G.; Investigation, Y.W., J.L., D.G., Z.L. and Z.C.; Sources, Z.Z., J.H. and X.H.; Writing—original manuscript preparation, Y.W.; Writing—review and editing, Y.W. and Z.Z.; Supervision, J.W. and Z.Z.; Project management, X.L., S.L. and Z.Z.; Funding acquisition, J.H. and Z.Z. All authors have read and agreed to the published version of the manuscript.

Funding

This research was supported by National Natural Science Foundation (32360821); Discipline Team Project of Gansu Agricultural University (GAU-XKTD-2022-22); Development and demonstration of high-efficiency production technology of yaks and cattle farming and animal husbandry cycle in pastoral areas of Qilian Mountains (2022CYZC-43); Gansu Provincial Department of Education: Young PhD Support Program (2023QB-128).

Institutional Review Board Statement

The animal study was reviewed and approved by the Ethical Commission of Gansu Agricultural University as well as the Ministry of Science and Technology of the People’s Republic of China (ethical approval file No. GSAU-Eth-AST-2023-001).

Informed Consent Statement

Not applicable.

Data Availability Statement

The data presented in the study are deposited in the NCBI repository, https://www.ncbi.nlm.nih.gov/bioproject/PRJNA1023693, URL (accessed on 10 June 2023), accession number: PRJNA1023693.

Acknowledgments

We thank all the members of the Gansu Key Laboratory of Herbivorous Animal Biotechnology of the College of Animal Science and Technology, Gansu Agricultural University, who contributed with their efforts to these experiments.

Conflicts of Interest

The authors declare that this research was conducted in the absence of any commercial or financial relationships that could be construed as potential conflict of interest.

References

  1. Guo, S.Z.; Bao, Y.Q.; Ma, D.L.; Li, B.M.; Wang, W.B.; Xu, G.Q.; Hu, J.; Bao, Z.X.J.C.; Wang, L.B.; Zhang, H.X.; et al. Determination of slaughter performance and meat quality of Jeryak in alpine pastures. Chin. Herbiv. Sci. 2019, 39, 72–74. [Google Scholar]
  2. Guo, S.Z.; Ma, D.L.; Yu, S.J.; Li, B.M.; Bao, Z.X.J.C.; Wang, L.B.; Ma, Z.T.; Zhang, Y.Z.; Niu, X.Y.; Zhou, J.; et al. Observations on the effect of crossbreeding between Jersey cattle and Gannan yaks in alpine pastures. China Cattle Sci. 2018, 44, 32–35. [Google Scholar]
  3. Oksbjerg, N.; Therkildsen, M. Chapter 3—Myogenesis and Muscle Growth and Meat Quality. In New Aspects of Meat Quality; Purslow, P.P., Ed.; Woodhead Publishing: Cambridge, UK, 2017; pp. 33–62. [Google Scholar]
  4. Te, K.G.; Reggiani, C. Skeletal muscle fibre type specification during embryonic development. J. Muscle Res. Cell Motil. 2002, 23, 65–69. [Google Scholar] [PubMed]
  5. Albrecht, E.; Teuscher, F.; Ender, K.; Wegner, J. Growth- and breed-related changes of muscle bundle structure in cattle. J. Anim. Sci. 2006, 84, 2959–2964. [Google Scholar] [CrossRef]
  6. Agarwal, M.; Sharma, A.; Kumar, P.; Kumar, A.; Bharadwaj, A.; Saini, M.; Kardon, G.; Mathew, S.J. Myosin heavy chain-embryonic regulates skeletal muscle differentiation during mammalian development. Development 2020, 147, 184507. [Google Scholar] [CrossRef] [PubMed]
  7. Buonaiuto, G.; Lopez-Villalobos, N.; Niero, G.; Degano, L.; Dadati, E.; Formigoni, A.; Visentin, G. The application of Legendre Polynomials to model muscularity and body condition score in primiparous Italian Simmental cattle. Ital. J. Anim. Sci. 2022, 21, 350–360. [Google Scholar] [CrossRef]
  8. Wang, T.; Xu, Y.Q.; Yuan, Y.X.; Xu, P.W.; Zhang, C.; Li, F.; Wang, L.N.; Yin, C.; Zhang, L.; Cai, X.C.; et al. Succinate induces skeletal muscle fiber remodeling via SUNCR1 signaling. EMBO Rep. 2019, 20, e47892. [Google Scholar] [CrossRef]
  9. Xu, M.; Chen, X.; Huang, Z.; Chen, D.; Li, M.; He, J.; Chen, H.; Zheng, P.; Yu, J.; Luo, Y.; et al. Effects of dietary grape seed proanthocyanidin extract supplementation on meat quality, muscle fiber characteristics and antioxidant capacity of finishing pigs. Food Chem. 2022, 367, 130781. [Google Scholar] [CrossRef]
  10. Xu, X.; Mishra, B.; Qin, N.; Sun, X.; Zhang, S.; Yang, J.; Xu, R. Differential Transcriptome Analysis of Early Postnatal Developing Longissimus Dorsi Muscle from Two Pig Breeds Characterized in Divergent Myofiber Traits and Fatness. Anim. Biotechnol. 2019, 30, 63–74. [Google Scholar] [CrossRef]
  11. Kim, G.D.; Jeong, J.Y.; Yang, H.S.; Hur, S.J. Differential abundance of proteome associated with intramuscular variation of meat quality in porcine longissimus thoracis et lumborum muscle. Meat Sci. 2019, 149, 85–95. [Google Scholar] [CrossRef]
  12. Kim, G.D.; Jeong, J.Y.; Jung, E.Y.; Yang, H.S.; Lim, H.T.; Joo, S.T. The influence of fiber size distribution of type IIB on carcass traits and meat quality in pigs. Meat Sci. 2013, 94, 267–273. [Google Scholar] [CrossRef]
  13. Heijmans, B.T.; Kremer, D.; Tobi, E.W.; Boomsma, D.I.; Slagboom, P.E. Heritable rather than age-related environmental and stochastic factors dominate variation in DNA methylation of the human IGF2/H19 locus. Hum. Mol. Genet. 2007, 16, 547–554. [Google Scholar] [CrossRef] [PubMed]
  14. Gunther, S.; Kim, J.; Kostin, S.; Lepper, C.; Fan, C.M.; Braun, T. Myf5-positive satellite cells contribute to Pax7-dependent long-term maintenance of adult muscle stem cells. Cell Stem Cell 2013, 13, 590–601. [Google Scholar] [CrossRef] [PubMed]
  15. Zammit, P.S. Function of the myogenic regulatory factors Myf5, MyoD, Myogenin and MRF4 in skeletal muscle, satellite cells and regenerative myogenesis. Semin. Cell Dev. Biol. 2017, 72, 19–32. [Google Scholar] [CrossRef] [PubMed]
  16. Desjardins, C.A.; Naya, F.J. The Function of the MEF2 Family of Transcription Factors in Cardiac Development, Cardiogenomics, and Direct Reprogramming. J. Cardiovasc. Dev. Dis. 2016, 3, 26. [Google Scholar] [CrossRef] [PubMed]
  17. Bridges, M.C.; Daulagala, A.C.; Kourtidis, A. LNCcation: lncRNA localization and function. J. Cell Biol. 2021, 220, e202009045. [Google Scholar] [CrossRef] [PubMed]
  18. Qian, X.; Zhao, J.; Yeung, P.Y.; Zhang, Q.C.; Kwok, C.K. Revealing lncRNA Structures and Interactions by Sequencing-Based Approaches. Trends Biochem. Sci. 2019, 44, 33–52. [Google Scholar] [CrossRef]
  19. Li, J.; Yang, T.; Tang, H.; Sha, Z.; Chen, R.; Chen, L.; Yu, Y.; Rowe, G.C.; Das, S.; Xiao, J. Inhibition of lncRNA MAAT Controls Multiple Types of Muscle Atrophy by cis- and trans-Regulatory Actions. Mol. Ther. 2021, 29, 1102–1119. [Google Scholar] [CrossRef]
  20. Li, Z.; Cai, B.; Abdalla, B.A.; Zhu, X.; Zheng, M.; Han, P.; Nie, Q.; Zhang, X. LncIRS1 controls muscle atrophy via sponging miR-15 family to activate IGF1-PI3K/AKT pathway. J. Cachexia Sarcopenia 2019, 10, 391–410. [Google Scholar] [CrossRef]
  21. 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. Biochim. Biophys. Acta 2016, 1863, 2835–2845. [Google Scholar] [CrossRef]
  22. Zhang, Z.K.; Li, J.; Guan, D.; Liang, C.; Zhuo, Z.; Liu, J.; Lu, A.; Zhang, G.; Zhang, B.T. A newly identified lncRNA MAR1 acts as a miR-487b sponge to promote skeletal muscle differentiation and regeneration. J. Cachexia Sarcopenia 2018, 9, 613–626. [Google Scholar] [CrossRef]
  23. Yu, X.; Wang, Z.; Sun, H.; Yang, Y.; Li, K.; Tang, Z. Long non-coding MEG3 is a marker for skeletal muscle development and meat production traits in pigs. Anim. Genet. 2018, 49, 571–578. [Google Scholar] [CrossRef]
  24. Yang, C.; Wu, L.F.; Liu, G.B.; Li, Y.K.; Liu, D.W.; Sun, B.L. Characterization of lncRNA expression in the longest dorsal muscle of Leiqiong and Lufeng cattle and analysis of its associated ceRNA network. Acta Vet. Zootech. Sin. 2023, 54, 1951–1963. [Google Scholar]
  25. Guttman, M.; Donaghey, J.; Carey, B.W.; Garber, M.; Grenier, J.K.; Munson, G.; Young, G.; Lucas, A.B.; Ach, R.; Bruhn, L.; et al. lincRNAs act in the circuitry controlling pluripotency and differentiation. Nature 2011, 477, 295–300. [Google Scholar] [CrossRef] [PubMed]
  26. Fatica, A.; Bozzoni, I. Long non-coding RNAs: New players in cell differentiation and development. Nat. Rev. Genet. 2014, 15, 7–21. [Google Scholar] [CrossRef] [PubMed]
  27. Li, Y.; Chen, X.; Sun, H.; Wang, H. Long non-coding RNAs in the regulation of skeletal myogenesis and muscle diseases. Cancer Lett. 2018, 417, 58–64. [Google Scholar] [CrossRef] [PubMed]
  28. Song, C.; Wang, J.; Ma, Y.; Yang, Z.; Dong, D.; Li, H.; Yang, J.; Huang, Y.; Plath, M.; Ma, Y.; et al. Linc-smad7 promotes myoblast differentiation and muscle regeneration via sponging miR-125b. Epigenetics 2018, 13, 591–604. [Google Scholar] [CrossRef] [PubMed]
  29. Archacka, K.; Ciemerych, M.A.; Florkowska, A.; Romanczuk, K. Non-Coding RNAs as Regulators of Myogenesis and Postexercise Muscle Regeneration. Int. J. Mol. Sci. 2021, 22, 11568. [Google Scholar] [CrossRef]
  30. Chen, R.; Lei, S.; She, Y.; Zhou, S.; Shi, H.; Li, C.; Jiang, T. Lnc-GD2H Promotes Proliferation by Forming a Feedback Loop With c-Myc and Enhances Differentiation Through Interacting With NACA to Upregulate Myog in C2C12 Myoblasts. Front. Cell Dev. Biol. 2021, 9, 671857. [Google Scholar] [CrossRef]
  31. Alessio, E.; Buson, L.; Chemello, F.; Peggion, C.; Grespi, F.; Martini, P.; Massimino, M.L.; Pacchioni, B.; Millino, C.; Romualdi, C.; et al. Singlecellanalysisrevealstheinvolvementofthelongnon-codingRNAPvt1inthemodulationofmuscleatrophy and mitochondrialnetwork. Nucleic Acids Res. 2019, 47, 18. [Google Scholar]
  32. Jackman, R.W.; Kandarian, S.C. The molecular basis of skeletal muscle atrophy. Am. J. Physiol.-Cell Physiol. 2004, 287, C834–C843. [Google Scholar] [CrossRef] [PubMed]
  33. Mukherjee, A.; Wilson, E.M.; Rotwein, P. Insulin-like growth factor (IGF) binding protein-5 blocks skeletal muscle differentiation by inhibiting IGF actions. Mol. Endocrinol. 2008, 22, 206–215. [Google Scholar] [CrossRef]
  34. Hribal, M.L.; Nakae, J.; Kitamura, T.; Shutter, J.R.; Accili, D. Regulation of insulin-like growth factor-dependent myoblast differentiation by Foxo forkhead transcription factors. J. Cell Biol. 2003, 162, 535–541. [Google Scholar] [CrossRef] [PubMed]
  35. Accili, D.; Arden, K.C. FoxOs at the crossroads of cellular metabolism, differentiation, and transformation. Cell 2004, 117, 421–426. [Google Scholar] [CrossRef] [PubMed]
  36. Xu, M.; Chen, X.; Chen, D.; Yu, B.; Huang, Z. FoxO1: A novel insight into its molecular mechanisms in the regulation of skeletal muscle differentiation and fiber type specification. Oncotarget 2017, 8, 10662–10674. [Google Scholar] [CrossRef]
  37. Wu, W.; Huang, R.; Wu, Q.; Li, P.; Chen, J.; Li, B.; Liu, H. The role of Six1 in the genesis of muscle cell and skeletal muscle development. Int. J. Biol. Sci. 2014, 10, 983–989. [Google Scholar] [CrossRef] [PubMed]
  38. Liu, Y.; Chakroun, I.; Yang, D.; Horner, E.; Liang, J.; Aziz, A.; Chu, A.; De Repentigny, Y.; Dilworth, F.J.; Kothary, R.; et al. Six1 regulates MyoD expression in adult muscle progenitor cells. PLoS ONE 2013, 8, e67762. [Google Scholar] [CrossRef] [PubMed]
  39. Hetzler, K.L.; Collins, B.C.; Shanely, R.A.; Sue, H.; Kostek, M.C. The homoeobox gene SIX1 alters myosin heavy chain isoform expression in mouse skeletal muscle. Acta Physiol. 2014, 210, 415–428. [Google Scholar] [CrossRef] [PubMed]
  40. O’Brien, J.H.; Hernandez-Lagunas, L.; Artinger, K.B.; Ford, H.L. MicroRNA-30a regulates zebrafish myogenesis through targeting the transcription factor Six1. J. Cell Sci. 2014, 127, 2291–2301. [Google Scholar]
  41. Le Grand, F.; Grifone, R.; Mourikis, P.; Houbron, C.; Gigaud, C.; Pujol, J.; Maillet, M.; Pages, G.; Rudnicki, M.; Tajbakhsh, S.; et al. Six1 regulates stem cell repair potential and self-renewal during skeletal muscle regeneration. J. Cell Biol. 2012, 198, 815–832. [Google Scholar] [CrossRef]
  42. Loumaye, A.; Lause, P.; Zhong, X.; Zimmers, T.A.; Bindels, L.B.; Thissen, J.P. Activin A Causes Muscle Atrophy through MEF2C-Dependent Impaired Myogenesis. Cells-Basel 2022, 11, 1119. [Google Scholar] [CrossRef]
  43. Shen, L.; Chen, L.; Zhang, S.; Zhang, Y.; Wang, J.; Zhu, L. MicroRNA-23a reduces slow myosin heavy chain isoforms composition through myocyte enhancer factor 2C (MEF2C) and potentially influences meat quality. Meat Sci. 2016, 116, 201–206. [Google Scholar] [CrossRef] [PubMed]
  44. Bengal, E.; Aviram, S.; Hayek, T. p38 MAPK in Glucose Metabolism of Skeletal Muscle: Beneficial or Harmful? Int. J. Mol. Sci. 2020, 21, 6480. [Google Scholar] [CrossRef] [PubMed]
  45. Somwar, R.; Koterski, S.; Sweeney, G.; Sciotti, R.; Djuric, S.; Berg, C.; Trevillyan, J.; Scherer, P.E.; Rondinone, C.M.; Klip, A. A dominant-negative p38 MAPK mutant and novel selective inhibitors of p38 MAPK reduce insulin-stimulated glucose uptake in 3T3-L1 adipocytes without affecting GLUT4 translocation. J. Biol. Chem. 2002, 277, 50386–50395. [Google Scholar] [CrossRef] [PubMed]
  46. Lee, J.; Hong, F.; Kwon, S.; Kim, S.S.; Kim, D.O.; Kang, H.S.; Lee, S.J.; Ha, J.; Kim, S.S. Activation of p38 MAPK induces cell cycle arrest via inhibition of Raf/ERK pathway during muscle differentiation. Biochem. Biophys. Res. Commun. 2002, 298, 765–771. [Google Scholar] [CrossRef]
  47. Al, M.A.; Mehta, V.; Li, G.; Figeys, D.; Wiper-Bergeron, N.; Skerjanc, I.S. Skeletal myosin light chain kinase regulates skeletal myogenesis by phosphorylation of MEF2C. EMBO J. 2011, 30, 2477–2489. [Google Scholar]
  48. Muncie, J.M.; Weaver, V.M. The Physical and Biochemical Properties of the Extracellular Matrix Regulate Cell Fate. Curr. Top. Dev. Biol. 2018, 130, 1–37. [Google Scholar]
  49. Frantz, C.; Stewart, K.M.; Weaver, V.M. The extracellular matrix at a glance. J. Cell Sci. 2010, 123, 4195–4200. [Google Scholar] [CrossRef]
  50. Birch, H.L. Extracellular Matrix and Ageing. Subcell. Biochem. 2018, 90, 169–190. [Google Scholar]
  51. Romer, L.H.; Birukov, K.G.; Garcia, J.G. Focal adhesions: Paradigm for a signaling nexus. Circ. Res. 2006, 98, 606–616. [Google Scholar] [CrossRef]
  52. Lassiter, D.G.; Nylen, C.; Sjogren, R.; Chibalin, A.V.; Wallberg-Henriksson, H.; Naslund, E.; Krook, A.; Zierath, J.R. FAK tyrosine phosphorylation is regulated by AMPK and controls metabolism in human skeletal muscle. Diabetologia 2018, 61, 424–443. [Google Scholar] [CrossRef]
  53. Lynch, G.S.; Ryall, J.G. Role of beta-adrenoceptor signaling in skeletal muscle: Implications for muscle wasting and disease. Physiol. Rev. 2008, 88, 729–767. [Google Scholar] [CrossRef]
  54. Kosuru, R.; Chrzanowska, M. Integration of Rap1 and Calcium Signaling. Int. J. Mol. Sci. 2020, 21, 1616. [Google Scholar] [CrossRef] [PubMed]
  55. Li, X.; Bai, Y.; Li, J.; Chen, Z.; Ma, Y.; Shi, B.; Han, X.; Luo, Y.; Hu, J.; Wang, J.; et al. Transcriptional analysis of microRNAs related to unsaturated fatty acid synthesis by interfering bovine adipocyte ACSL1 gene. Front. Genet. 2022, 13, 994806. [Google Scholar] [CrossRef] [PubMed]
  56. Pawlikowski, B.; Vogler, T.O.; Gadek, K.; Olwin, B.B. Regulation of skeletal muscle stem cells by fibroblast growth factors. Dev. Dyn. 2017, 246, 359–367. [Google Scholar] [CrossRef] [PubMed]
  57. Zofkie, W.; Southard, S.M.; Braun, T.; Lepper, C. Fibroblast growth factor 6 regulates sizing of the muscle stem cell pool. Stem Cell Rep. 2021, 16, 2913–2927. [Google Scholar] [CrossRef] [PubMed]
  58. Cai, Q.; Wu, G.; Zhu, M.; Ge, H.; Xue, C.; Zhang, Q.; Cheng, B.; Xu, S.; Wu, P. FGF6 enhances muscle regeneration after nerve injury by relying on ERK1/2 mechanism. Life Sci. 2020, 248, 117465. [Google Scholar] [CrossRef]
  59. Zhang, B.; Arun, G.; Mao, Y.S.; Lazar, Z.; Hung, G.; Bhattacharjee, G.; Xiao, X.; Booth, C.J.; Wu, J.; Zhang, C.; et al. The lncRNA Malat1 Is Dispensable for Mouse Development but Its Transcription Plays a cis-Regulatory Role in the Adult. Cell Rep. 2012, 2, 111–123. [Google Scholar] [CrossRef]
  60. Lu, L.; Sun, K.; Chen, X.; Zhao, Y.; Wang, L.; Zhou, L.; Sun, H.; Wang, H. Genome-wide survey by ChIP-seq reveals YY1 regulation of lincRNAs in skeletal myogenesis. EMBO J. 2013, 32, 2575–2588. [Google Scholar] [CrossRef]
  61. Liang, W.C.; Ren, J.L.; Wong, C.W.; Chan, S.O.; Waye, M.M.; Fu, W.M.; Zhang, J.F. LncRNA-NEF antagonized epithelial to mesenchymal transition and cancer metastasis via cis-regulating FOXA2 and inactivating Wnt/beta-catenin signaling. Oncogene 2018, 37, 1445–1456. [Google Scholar] [CrossRef]
  62. Cai, B.; Li, Z.; Ma, M.; Wang, Z.; Han, P.; Abdalla, B.A.; Nie, Q.; Zhang, X. LncRNA-Six1 Encodes a Micropeptide to Activate Six1 in Cis and Is Involved in Cell Proliferation and Muscle Growth. Front. Physiol. 2017, 8, 230. [Google Scholar] [CrossRef]
  63. Zhang, J.; Ma, G.; Guo, Z.; Yu, Q.; Han, L.; Han, M.; Zhu, Y. Study on the apoptosis mediated by apoptosis-inducing-factor and influencing factors of bovine muscle during postmortem aging. Food Chem. 2018, 266, 359–367. [Google Scholar]
  64. Kim, D.; Langmead, B.; Salzberg, S.L. HISAT: A fast spliced aligner with low memory requirements. Nat. Methods 2015, 12, 357–360. [Google Scholar] [CrossRef]
  65. Pertea, M.; Pertea, G.M.; Antonescu, C.M.; Chang, T.C.; Mendell, J.T.; Salzberg, S.L. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol. 2015, 33, 290–295. [Google Scholar] [CrossRef]
  66. Kang, Y.J.; Yang, D.C.; Kong, L.; Hou, M.; Meng, Y.Q.; Wei, L.; Gao, G. CPC2: A fast and accurate coding potential calculator based on sequence intrinsic features. Nucleic Acids Res. 2017, 45, W12–W16. [Google Scholar] [CrossRef]
  67. Sun, L.; Luo, H.; Bu, D.; Zhao, G.; Yu, K.; Zhang, C.; Liu, Y.; Chen, R.; Zhao, Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013, 41, e166. [Google Scholar] [CrossRef] [PubMed]
  68. Wucher, V.; Legeai, F.; Hédan, B.; Rizk, G.; Lagoutte, L.; Leeb, T.; Jagannathan, V.; Cadieu, E.; David, A.; Lohi, H.; et al. FEELnc: A tool for long non-coding RNA annotation and its application to the dog transcriptome. Nucleic Acids Res. 2017, 45, e57. [Google Scholar] [CrossRef] [PubMed]
  69. Mercer, T.R.; Dinger, M.E.; Mattick, J.S. Long non-coding RNAs: Insights into functions. Nat. Rev. Genet. 2009, 10, 155–159. [Google Scholar] [PubMed]
  70. Trapnell, C.; Williams, B.A.; Pertea, G.; Mortazavi, A.; Kwan, G.; van Baren, M.J.; Salzberg, S.L.; Wold, B.J.; Pachter, L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat. Biotechnol. 2010, 28, 511–515. [Google Scholar] [CrossRef] [PubMed]
  71. Anders, S.; Huber, W. Differential expression analysis for sequence count data. Genome Biol. 2010, 11, R106. [Google Scholar] [CrossRef]
  72. Love, M.I.; Huber, W.; Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014, 15, 550. [Google Scholar] [CrossRef] [PubMed]
Ijms 24 16991 g001
Figure 1. Fast/slow muscle fluorescence staining and paraffin section HE staining of muscle tissue. (A) Red represents fast muscle fiber; green represents slow muscle fiber; and blue represents the nucleus of the muscle cell. (B) Red represents a single muscle fiber. (C) Comparison of myofiber characteristics of different breeds of cattle. Gannan Yak in blue and Jeryak in red.
Figure 1. Fast/slow muscle fluorescence staining and paraffin section HE staining of muscle tissue. (A) Red represents fast muscle fiber; green represents slow muscle fiber; and blue represents the nucleus of the muscle cell. (B) Red represents a single muscle fiber. (C) Comparison of myofiber characteristics of different breeds of cattle. Gannan Yak in blue and Jeryak in red.
Ijms 24 16991 g001
Ijms 24 16991 g002
Figure 2. Analysis of lncRNA sequencing data. (A) Principal component analysis (PCA) of sample.The figure shows Gannan yaks (M1, M2, M3) in orange and Jeryaks (M1, M2, M3) in blue. (B) Exon number distribution of novel lncRNAs. (C) Length distribution of the novel lncRNAs. (D) Classification of novel lncRNAs.
Figure 2. Analysis of lncRNA sequencing data. (A) Principal component analysis (PCA) of sample.The figure shows Gannan yaks (M1, M2, M3) in orange and Jeryaks (M1, M2, M3) in blue. (B) Exon number distribution of novel lncRNAs. (C) Length distribution of the novel lncRNAs. (D) Classification of novel lncRNAs.
Ijms 24 16991 g002
Ijms 24 16991 g003
Figure 3. Genomic expression of lncRNA and mRNA. (A) The expression of lncRNA and mRNA. (B,C) Comparison of length and exon count between lncRNA and mRNA. (E,F) Cluster heatmap of DEMs and DELs. The FPKM array was centered and scaled in the row direction using R package pheatmap (v1.0.8). Red indicates higher expression, and blue represents lower expression. The log10 (FPKM + 1) value was converted (scale number) and clustered. (D,G) Volcano plot of DEMs and DELs.
Figure 3. Genomic expression of lncRNA and mRNA. (A) The expression of lncRNA and mRNA. (B,C) Comparison of length and exon count between lncRNA and mRNA. (E,F) Cluster heatmap of DEMs and DELs. The FPKM array was centered and scaled in the row direction using R package pheatmap (v1.0.8). Red indicates higher expression, and blue represents lower expression. The log10 (FPKM + 1) value was converted (scale number) and clustered. (D,G) Volcano plot of DEMs and DELs.
Ijms 24 16991 g003
Ijms 24 16991 g004
Figure 4. GO and KEGG analyses of DEMs. (A) GO terms significantly enriched by DEMs in Class-up. (B) Top 20 KEGG pathways significantly enriched by DEMs in Class-up. (C) GO terms significantly enriched by DEMs in Class-down. (D) Top 20 KEGG pathways significantly enriched by DEMs in Class-down.
Figure 4. GO and KEGG analyses of DEMs. (A) GO terms significantly enriched by DEMs in Class-up. (B) Top 20 KEGG pathways significantly enriched by DEMs in Class-up. (C) GO terms significantly enriched by DEMs in Class-down. (D) Top 20 KEGG pathways significantly enriched by DEMs in Class-down.
Ijms 24 16991 g004
Ijms 24 16991 g005
Figure 5. GO and KEGG analyses of cis- and trans-target genes of DELs. (A) The GO analysis of the top 60 cis-target genes with the corresponding categories. (B) KEGG pathways analysis of top 20 cis-target genes. (C) GO analysis of top 60 trans-target genes along with the categories of BP, CC, and MF. (D) The KEGG pathways of the top 20 trans-target genes.
Figure 5. GO and KEGG analyses of cis- and trans-target genes of DELs. (A) The GO analysis of the top 60 cis-target genes with the corresponding categories. (B) KEGG pathways analysis of top 20 cis-target genes. (C) GO analysis of top 60 trans-target genes along with the categories of BP, CC, and MF. (D) The KEGG pathways of the top 20 trans-target genes.
Ijms 24 16991 g005
Ijms 24 16991 g006
Figure 6. The DEL–DEM interaction network. (A) Network analysis of DE-lncRNA target genes and their enrichment pathways. (B) DEL and trans-target interaction network. The trans-target was also a DEM. (C) DEL and cis-target interaction network. The cis-target was also a DEM. The circle type and inverted-triangle type represent mRNA and lncRNA, respectively. Yellow represents the enriched signaling pathways; red represents up-regulation; and green represents down-regulation.
Figure 6. The DEL–DEM interaction network. (A) Network analysis of DE-lncRNA target genes and their enrichment pathways. (B) DEL and trans-target interaction network. The trans-target was also a DEM. (C) DEL and cis-target interaction network. The cis-target was also a DEM. The circle type and inverted-triangle type represent mRNA and lncRNA, respectively. Yellow represents the enriched signaling pathways; red represents up-regulation; and green represents down-regulation.
Ijms 24 16991 g006
Ijms 24 16991 g007
Figure 7. qRT-PCR results of DEMs and DELs. (A) Changes in mRNA expression in M and P treatment groups. (B) qRT-PCR validation of changes in lncRNA expression in M and P treatment groups. Three biological replicates were used for each group. The qRT-PCR data were determined with the 2−ΔΔCt method using GAPDH as an internal reference. Validation data from mRNA and lncRNA sequencing results were further normalized to log2 (fold change). Data represent means ± standard errors.
Figure 7. qRT-PCR results of DEMs and DELs. (A) Changes in mRNA expression in M and P treatment groups. (B) qRT-PCR validation of changes in lncRNA expression in M and P treatment groups. Three biological replicates were used for each group. The qRT-PCR data were determined with the 2−ΔΔCt method using GAPDH as an internal reference. Validation data from mRNA and lncRNA sequencing results were further normalized to log2 (fold change). Data represent means ± standard errors.
Ijms 24 16991 g007
Table 1. Output statistics and quality assessment of cDNA library sequencing.
Table 1. Output statistics and quality assessment of cDNA library sequencing.
SampleRaw ReadsClean ReadsLow Quality (%)Q20 (%)Q30 (%)GC Content (%)
M191,597,21891,113,0840.3396.8492.1465.18
M290,970,93490,561,7900.2796.9392.4065.54
M377,261,06076,658,9920.5996.8892.2064.48
P191,772,58891,351,0060.2596.9592.2472.24
P293,352,14492,856,1460.3196.8992.1871.10
P390,267,32489,794,2080.3096.6891.7772.18
Table 2. Clean reads and reference genome alignment results.
Table 2. Clean reads and reference genome alignment results.
SampleTotalTotal Mapped (%)Unique Mapped (%)Multiple Mapped (%)
M168,296,27055,971,679 (81.95%)44,550,776 (65.23%)11,420,903 (16.72%)
M268,471,61655,259,551 (80.70%)44,632,547 (65.18%)10,627,004 (15.52%)
M359,586,39049,279,578 (82.70%)39,697,584 (66.62%)9,581,994 (16.08%)
P150,960,52436,586,723 (71.79%)28,892,599 (56.70%)7,694,124 (15.10%)
P256,700,50641,707,068 (73.56%)32,919,039 (58.06%)8,788,029 (15.50%)
P351,396,22836,426,543 (70.87%)28,481,656 (55.42%)7,944,887 (15.46%)
Table 3. lncRNA primers used in qRT-PCR.
Table 3. lncRNA primers used in qRT-PCR.
lncRNAForward (5′→3′)Reverse (5′→3′)
MSTRG.1959.1GGAGGCTGAGACTGGAGGATGGGGTGTCGCTATGTTGC
MSTRG.16260.9GTCCATCCATCCGCATCTCCCTCCACTCTGACCATCCCT
MSTRG.4090.7CTTGAAAAGTGGCTGTGGTTACGGATGGTGTCTGGAGGT
MSTRG.12211.4CAGTAAAGTGTCTGCCTGTAATGCTATCTCCCGAGTGGGTTGC
MSTRG.18038.1ATTCGGGAAGGAGGTACAGGTGGGGTTTACTCAAGGCACT
MSTRG.31101.2GGCGGTCTACGGGTTATTCCGCTTTCGTTCACAGGCTAA
MSTRG.22127.1TGGCTGCTCTACTGTCCTCTGCTGGTTGCCCCTGAATACG
GAPDHAGTTCAACGGCACAGTCAAGGACCACATACTCAGCACCAGCA
Table 4. mRNA primers used in qRT-PCR.
Table 4. mRNA primers used in qRT-PCR.
mRNAForward (5′→3′)Reverse (5′→3′)
ADARB1AGCTGAACGAGATCAAGCCCCTCGAACACCTGTCCGTTGA
FGFR4CCTTGCTTCTGCACAACGTCCCTGTCCATCCTTGAGCCAG
PFKFB3AACCGTGTGCAGGATCACATGGTCCTTCAGGTTCTGCTCC
HEYLACTTCCGGAGCATCGGTTTTTGCATAGCTGTTGAGGTGGG
MEF2CCACTGGGAAACCCCAACCTTAGCAGACCTGGTGAGTTTCG
IGFBP3CGCTACAAGCGTTGTTGGACTGCTGTGGTCTTCTTCCGAC
RPL14TTCAGGCGCTTCGTAGAGGGCAGGTGCTTTCTGGGATG
RPS3AGCGGTCGGCAAGAATAAGTGGTGGCTCGTATCCATC
GAPDHAGTTCAACGGCACAGTCAAGGACCACATACTCAGCACCAGCA
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Wei, Y.; Guo, D.; Bai, Y.; Liu, Z.; Li, J.; Chen, Z.; Shi, B.; Zhao, Z.; Hu, J.; Han, X.; et al. Transcriptome Analysis of mRNA and lncRNA Related to Muscle Growth and Development in Gannan Yak and Jeryak. Int. J. Mol. Sci. 2023, 24, 16991. https://doi.org/10.3390/ijms242316991

AMA Style

Wei Y, Guo D, Bai Y, Liu Z, Li J, Chen Z, Shi B, Zhao Z, Hu J, Han X, et al. Transcriptome Analysis of mRNA and lncRNA Related to Muscle Growth and Development in Gannan Yak and Jeryak. International Journal of Molecular Sciences. 2023; 24(23):16991. https://doi.org/10.3390/ijms242316991

Chicago/Turabian Style

Wei, Yali, Dashan Guo, Yanbin Bai, Zhanxin Liu, Jingsheng Li, Zongchang Chen, Bingang Shi, Zhidong Zhao, Jiang Hu, Xiangmin Han, and et al. 2023. "Transcriptome Analysis of mRNA and lncRNA Related to Muscle Growth and Development in Gannan Yak and Jeryak" International Journal of Molecular Sciences 24, no. 23: 16991. https://doi.org/10.3390/ijms242316991

APA Style

Wei, Y., Guo, D., Bai, Y., Liu, Z., Li, J., Chen, Z., Shi, B., Zhao, Z., Hu, J., Han, X., Wang, J., Liu, X., Li, S., & Zhao, F. (2023). Transcriptome Analysis of mRNA and lncRNA Related to Muscle Growth and Development in Gannan Yak and Jeryak. International Journal of Molecular Sciences, 24(23), 16991. https://doi.org/10.3390/ijms242316991

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop