http://orcid.org/0000-0002-4354-6236
Figures
Abstract
Domestic yak (Bos grunniens) is the most crucial livestock in the Qinghai-Tibetan Plateau, providing meat and other necessities for local people. The skeletal muscle of adult livestock is composed of muscle fibers, and fiber composition in muscle has influence on meat qualities, such as tenderness, pH, and color. Real-time quantitative polymerase chain reaction (RT-qPCR) is a powerful tool to evaluate the gene expression of muscle fiber, but the normalization of the data depends on the stability of expressed reference genes. Unfortunately, there is no consensus for an ideal reference gene for data normalization in muscle tissues of yak. In this study, we aimed to assess the stability of 14 commonly used candidate reference genes by using five algorithms (GeNorm, NormFinder, BestKeeper, Delat Ct and Refinder). Our results suggested UXT and PRL13A were the most stable reference genes, while the most commonly used reference gene, GAPDH, was most variably expressed across different muscle tissues. We also found that the extensor digitorum lateralis (EDL), trapezius pars thoracica (TPT), and psoas major (PM) muscle had the higher content of type I muscle fibers and the lowest content of type IIB muscle fibers, while gluteobiceps (GB) muscle had the highest content of type IIB muscle fibers. Our study provides the suitable reference genes for accurate analysis of yak muscle fiber composition.
Citation: Wu X, Zhou X, Ding X, Chu M, Liang C, Pei J, et al. (2020) Reference gene selection and myosin heavy chain (MyHC) isoform expression in muscle tissues of domestic yak (Bos grunniens). PLoS ONE 15(2): e0228493. https://doi.org/10.1371/journal.pone.0228493
Editor: Marcio de Souza Duarte, Universidade Federal de Viçosa, BRAZIL
Received: August 11, 2019; Accepted: January 16, 2020; Published: February 6, 2020
Copyright: © 2020 Wu et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: All relevant data are within the manuscript and its Supporting Information files.
Funding: This study was supported by the Agricultural Science and Technology Innovation Program (CAAS-ASTIP-2014-LIHPS-01) and National Beef Cattle Industry Technology & System (CARS-37). Ping Yan received both of the awards. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The yak (Bos grunniens) is one of the livestock that can adapt to the harsh environment of the Qinghai-Tibetan Plateau and adjacent Alpine regions, which provides meat, milk, hair, transport and fuel for residents [1]. There are around 14 million yaks in the world, 95% of them are distributed in China [2]. As the consumer demand for high quality meat products continues to increase, there is a growing demand for high quality yak meat in most areas of China. Therefore, understanding the characteristics of meat quality is a major priority in yak production.
Skeletal muscle is a heterogeneous tissue that consists of a large variety of physiologically and biochemically diverse fiber types [3]. According to the different structure, contractile properties, metabolic and morphological traits, skeletal muscle fibers of adult livestock are classified into four types: type I (slow-twitch, oxidative metabolism), type IIA (fast-twitch, oxidative metabolism), type IIB (fast twitch, glycolytic metabolism) and type IIX (fast-twitch, intermediate metabolism) [4]. Many studies indicate that muscle fiber type composition has profound influence on the biochemical properties of the muscle, and the meat quality traits of livestock [5,6]. Muscle fiber type composition is commonly used to predict the meat quality traits of livestock. Myosin heavy chain (MyHC) is the most abundant muscle structural protein, comprising about 35% of the protein pool [7]. Four adult MyHC isoforms (I, IIA, IIX and IIB) have been identified in the skeletal muscle of representative mammalian [8]. The slow-twitch muscle fibers of skeletal muscle primarily express the MyHC I isoform, while fast-twitch muscle fibers express other faster isoforms (MyHC 2A, MyHC 2X and MyHC 2B) [9]. Therefore, muscle fiber types are characterized by isoforms of the MyHC genes which are expressed. The composition of each muscle fiber type can generally be determined by the protein and mRNA expression of four MyHC isoforms using immunohistochemistry [10], SDS-PAGE [11], in situ hybridization [12], and Reverse transcription PCR (RT-PCR) approaches [13], but these approaches are semiquantitative. Currently, the real-time quantitative polymerase chain reaction (RT-qPCR) has become a powerful method used to analyze the muscle fiber composition due to its high sensitivity, specificity, simplicity, and reproducibility [14].
RT-qPCR is considered as the gold standard for defining accurate expression profiles of selected genes [15,16]. The reliability of RT-qPCR results can be easily influenced by several variables, including the initial sample amount, RNA integrity and quantity, accurate reverse transcription, and primer efficiency. To minimize the influence of these factors, reference genes are usually used to calibrate and normalize the expression level of target genes [17]. Ideal reference genes should have a ubiquitous expression and very low variance across various development stages, types of tissues, environmental conditions, and health conditions [18]. Recently, some commonly used reference genes for RT-qPCR have been analyzed in different tissues of yak [19–22]. Unfortunately, as of now, there are no available endogenous control genes for RT-qPCR data normalization that have been identified for various muscle of yak.
The aim of this study was to select and evaluate the stable reference genes that can be used to normalize mRNA expression by RT-qPCR in skeletal muscle tissue of yak. To attain this aim, we selected 14 commonly used candidate reference genes from published literature for gene expression studies. The expression stabilities of these reference genes were analyzed using the GeNorm [23], NormFinder [24], BestKeeper [25], Delat Ct [26] and Refinder programs [27]. Additionally, the validated reference genes were used to identify the muscle fiber composition in ten skeletal muscles of yak.
Materials and methods
Muscle samples
Three male Ashidan yaks were raised at the Qinghai Datong yak farm of Province with the same feeding conditions. Yaks were slaughtered at the body weight of 176.00 ± 8.54 kg (30 months of age). After slaughter, ten different skeletal muscles were immediately frozen in liquid nitrogen and stored at -80 °C until further processing. These ten skeletal muscles were extensor digitorum lateralis (EDL), psoas major (PM), latissimus dorsi (LD), gluteobiceps (GB), gastrocnemius (GC), fibularis longus (FL), semitendinosus (ST), trapezius pars thoracica (TPT), supraspinatus (SP) and longissimus dorsi muscle (LDM). This study was carried out in strict accordance with the recommendations in “Guidelines for Experimental Animals” of the Ministry of Science and Technology (Beijing, China), and all of the experimental protocols and procedures were approved by the Animal Administration and Ethics Committee of Lanzhou Institute of Husbandry and Pharmaceutical Sciences of CAAS (Permit No. SYXK-2014-0002).
RNA extraction and cDNA synthesis
Total RNA was isolated using the animal tissue RNA isolation kit (ZDGSY, Beijing, China) following the manufacturer’s protocol. RNA concentration and purity (A260/A280) was determined using NanoDrop2000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). Integrity of RNA was evaluated using 1% agarose gel. Total RNA (1000 ng) was reverse transcribed according to the instructions of PrimeScript™ RT reagent kit with gDNA Eraser (TaKaRa, Dalian, China). The cDNA was stored at -20 °C until required.
Primers design
Based on literature review, the following genes were chosen to evaluate using the RT-qPCR method for expression of stability: glyceraldehyde-3-phosphate dehydrogenase (GAPDH), actin-beta (ACTB), ubiquitously expressed prefoldin-like chaperone (UXT), TATA box binding protein (TBP), tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein zeta (YWHAZ), ribosomal protein L13a (RPL13A), dehydrogenase complex subunit A (SDHA), ribosomal protein S15 (RPS15), hypoxanthine phosphoribosyl transferase 1 (HPRT1), peptidylprolyl isomerase A (PPIA), hydroxymethylbilane synthase (HMBS), mitochondrial ribosomal protein L39 (MRPL39), beta-2-microglobulin (B2M) and protein phosphatase 1 regulatory inhibitor subunit 11 (PPP1R11) (S1 Table). Among these candidate reference genes, PPP1R11, UXT, MRPL39, RPS15, HMBS, YWHAZ, TBP and ACTB have been reported to be optimal reference genes for normalization in yak [19–22]. In addition, HPRT1, PPIA, SDHA, B2M and RPL13A were identified as reference gene for expression studies in muscle tissue of pig, chicken and human [28–31]. GAPDH have been utilized as reference gene in numerous studies of yak. Primers for TBP, ACTB, PPIA, HPRT1, GAPDH and SDHA were used from Li et al [21]. The sequences of other genes were obtained from NCBI. All the primer pairs were designed by Primer 5.0 software with the length of 20 ± 3 bases, and size of amplicon ranging from 79–198 bp (Table 1). Primer specificity of each reference gene was confirmed by melting curve analysis and 1.5% agarose gel electrophoresis.
RT-qPCR analysis
All RT-qPCR analyses were conducted in 96-well plates using the BioRad CFX96™ Real-Time PCR system (Bio-Rad Laboratories, CA, USA). The PCR program for amplification was as follows: an initial denaturation step of 3 min at 95 °C, and 45 cycles of 95 °C for 10 s, 60 °C for 10 s, and 72 °C for 10 s. A final melting program ranging from 65 °C to 95 °C at increments of 0.5 °C and acquiring the fluorescence after each step. The reaction mixture contained 1 μl of diluted cDNA, 10 μl of SYBR TB Green mix (TaKaRa, Dalian, China), 1 μl of each sense and anti-sense, and 7 μL of ddH2O in 20 μl of total reaction volume. To ensure repeatability of the experiments, each reaction was run in triplicate. A ten-fold dilution series of cDNA samples was used to generate a relative standard curve. The correlation coefficient (R2) and slope were obtained from the linear regression model created by the BioRad CFX96™ Real-Time PCR system. PCR amplification efficiency (E) was calculated according to the formula: E = (10(−1/slope) − 1) × 100.
Reference genes expression stability analysis
The expression stability of the selected reference genes was evaluated using four Microsoft Excel-based programs: geNorm [23], Normfinder [24], Bestkeeper [25], Delta Ct [26] and a comprehensive online tool RefFinder [27] (http://150.216.56.64/referencegene.php) according to the developer’s instructions. GeNorm was also used to calculate the pairwise variation (Vn/Vn+1 value) to determine the suitable number of reference genes required for accurate normalization.
Validation of reference genes
The muscle fiber differentiation-related gene Myogenin (MyoG) was selected as the target for validating the two most stable reference genes, UXT and RPL13A. The reference genes that varied the most, particularly, GAPDH, was also used to normalize the expression of target gene. The primers sequences for MyoG are shown in Table 1. RT-qPCR was performed as described above. The 2-ΔΔCt method was used to calculate the relative normalized expression in EDL and GB muscle [32]. Data were expressed as mean ± SD of three biological replicates for each muscle. The differences on the MyoG expression levels were evaluated by Student’s t test using IBM SPSS Statistics 23.0 (SPSS Inc., Armonk, NY, USA).
Measurement of muscle fiber characteristics
MyHC expression.
Gene-specific primers for the quantification of mRNA expression of MyHC I and MyHC IIB gene were designed using Primer 5.0 software (S2 Table). RT-qPCR was performed as described above. The relative expression of MyHC I and MyHC IIB mRNA were calculated via the 2-ΔΔCt method [32].
Mitochondrial DNA (mtDNA) copy number.
Genomic DNA was isolated using the animal tissue Genomic DNA Kit (ZDGSY, Beijing, China) following the manufacturer’s instructions. The copy number of mtDNA per cell was determined using RT-qPCR according to the previous method [33,34]. ATP synthase subunit a (ATP6), cytochrome C oxidase I (COX1) and NADH dehydrogenase, subunit 1 (ND1) and glucagon (GCG) gene were selected to calculate the number of mtDNA copies per cell. All the primers used are described in S2 Table. RT-qPCR was performed as described above. The following modified equation was used to determine the mtDNA copy numbers of per cell: mtDNA copy numbers = [(No. of copies of the mtDNA gene)/(No. of copies of GCG)/2].
Results
Selection of candidate reference genes
To select suitable reference genes for the ten skeletal muscles of yak, we chose 14 candidate reference genes (GAPDH, ACTB, UXT, TPB, YWHAZ, RPL13A, SDHA, RPS15, HPRT1, PPIA, HMBS, MRPL39, B2M and PPP1R11) based on previous literatures. Candidate reference genes were distributed among different chromosomes that belong to different function classes to minimize the risk of coregulation.
Primer specificity and amplification efficiency of candidate genes
Gene name, primer pair sequences, amplicon size, PCR efficiency (E), and correlation coefficients (R2) of the reference genes are provided in Table 1. The sizes of all PCR products were between 79–198 bp, and melting temperatures ranged from 84 °C to 88 °C. All PCR assays produced a specific amplification of the expected size, which exhibited a single sharp peak in the melting curve (S1 Fig). The amplification efficiencies of all the primers for the reference genes ranged from 93–109%, with R2 varied between 0.9921 and 1.0000 (Table 1).
Expression levels of candidate reference genes
The Ct values of 14 candidate reference genes in all of the tested samples were shown in Fig 1. The mean Ct values of the 14 candidate reference genes exhibited relatively different variation in the ten different muscle tissues (S3 Table). Mean Ct values of candidate reference genes were between 20.86 to 32.64. The smaller the Ct value, the higher gene expression, and vice versa. The RPL13A gene had the lowest mean Ct value (21.73), corresponding to its highest level, while HMBS was the lowest expressed gene with a mean Ct value of 30.24.
Each box indicates the 25th and 75th percentiles. Whiskers caps correspond to the maximum and minimum values. A line within the boxes depicts the median.
Stability of the reference genes
In this study, the stability of expression of the reference gene was analyzed by GeNorm, NormFinder, BestKeeper and Delta-Ct program. Then, the genes were re-ranked based on geometric mean by using the RefFinder online tool. The results of each program are provided in Table 2.
GeNorm analysis.
The GeNorm program evaluated the stability of reference genes based on expression stability value (M). Genes with the lower the M value have the higher the expression stability. All of the analyzed reference genes had an M-value below the recommended cut-off value of 1.5, indicating that each reference gene has a relative stability across the tested samples (Table 2). According to the GeNorm result, the UXT and RPL13A had the highest stability, with the lowest M-value (0.242). In contrast, GAPDH exhibited the least stable (M value of 0.730) (Table 2).
GeNorm was also used to calculate the pairwise variation (Vn/Vn+1) between sequential normalization factors (NF) to determine the minimum number of genes necessary for normalization. A cutoff value of 0.15 is the recommended threshold to identify the benefit of additional reference gene for the normalization [23]. As evident in Fig 2, pairwise variations V2/3 in ten skeletal muscles were less than 0.15, which indicates that two reference genes were suitable for gene expression analysis (Table 2).
The V cutoff value was set to 0.15, below which the inclusion of an additional reference gene would not have a beneficial effect on the normalization.
NormFinder analysis.
According to the NormFinder program, genes with lower average expression stability values present with more stable expression. NormFinder ranking indicated that PPP1R11 and UXT had the most stable expression with stability values of 0.255 and 0.293, while B2M (0.739) and GAPDH (1.159) were ranked the least stable genes with the higher stability values (Table 2).
BestKeeper analysis.
The BestKeeper calculates the stability of the candidate reference genes on the basis of standard deviation (SD) of Ct values and the correlation coefficient (r) of expression among the reference genes. Determining the gene stability based on the SD value is the more conservative approach [35]. A gene with the lowest SD value is ranked as the most stable reference gene, while a gene with a SD above 1 is not a stable reference gene. BestKeeper indicated that SD values of all reference genes were less than 1. UXT and YWHAZ had the lowest SD value (0.351) and thus were the most stable genes in the tissue, while HMBS (0.925) and GAPDH (0.988) were the least expressed genes (Table 2).
The Delta Ct analysis.
The Delta Ct program ranks the most stable reference genes by comparing relative expression value of tested reference genes in pairs in all of the samples. Genes with lower SD values represent the higher expression stability. The ranking of reference genes was similar to what was obtained using the GeNorm algorithm. UXT was the most stable reference gene in ten skeletal muscles with the lowest SD value (0.579), whereas GAPDH (1.247) had the lowest stability (Table 2).
RefFinder analysis.
To avoid the divergent stability rankings of different software, the online RefFinder tool was also used to construct a comprehensive ranking of the most stable candidate reference genes. Reference genes with the lowest geometric mean are recommended by RefFinder. According to the recommended comprehensive ranking, UXT and RPL13A were the most suitable genes for endogenous controls because they had the lowest geometric mean of ranking value, indicating they were the most stable (Table 2).
Reference gene validation
The relative expression of MyoG was used to validate the stability of the best-ranked candidate reference genes. When normalized using the stable reference genes UXT and RPL13A either alone or in combination, there was no significant difference in MyoG expression between EDL and GB muscle. On the contrary, when the least stable reference gene GAPDH was used for normalization, the expression of MyoG in EDL muscle was significantly higher than that in GB muscle (Fig 3).
The expression levels were normalized using UXT, RPL13A or GAPDH as reference genes individually and with the mean of UXT+ RPL13A. **P < 0.01.
The characteristics of skeletal muscle fiber types
The myofiber characteristics of ten different muscles from yak are shown in Fig 4. The muscle had higher percentages of type I fibers and lower percentages of type IIb fibers, and vice versa. Among these muscles, EDL, TPT and PM muscle had the higher expression of MyHC I mRNA and copy number of mtDNA per cell, indicating that their oxidative capacity was higher compared with those of other skeletal muscles. In contrast, GB muscle had the highest expression of MyHC IIB mRNA and lowest copy number of mtDNA per cell, suggesting it to be more proficient in anaerobic glycolytic metabolism.
(A and B) Normalized expression of MyHC I and MyHC IIB in ten different muscle tissues, (C) MtDNA copies per cell in ten different muscle tissues.
Discussion
Improving meat quality is part of the main breeding objectives for domestic animal production. Previous studies have shown that muscle fiber composition directly influences the eating quality of meat, as they are the basic component of skeletal muscles [36]. One of the important aspects in the study of muscle fiber composition is the expression analysis of MyHC genes by RT-qPCR. However, RT-qPCR data can be easily influenced by the stability of the reference genes chosen for RT-qPCR data normalization. Using unsuitable reference genes can lead to deviations in expression normalization and thus generate incorrect results [37]. The optimal reference gene should be constantly expressed in all types of organisms throughout the different development stages [18]. In reality, there is not an ideal set of genes that are consistent in all organisms. Therefore, it is necessary to select and validate suitable reference genes in different cases for expression analysis prior to RT-qPCR. In this study, 14 reference genes of yak were selected, and their expression was evaluated using a multi-algorithm analysis in ten different muscles of yak. To the best of our knowledge, this is the first study to select and validate the reference gene for different muscle tissues of yak.
To evaluate the gene expression stability of the selected genes among different skeletal muscles, we used four programs (GeNorm, NormFinder, BestKeeper and Delta CT) to obtain rankings of 14 candidate reference genes. GeNorm uses average pairwise variation of each reference gene to calculate its gene-stability value (M) [23]. BestKeeper estimates gene expression stability on the basis of standard deviation (SD) The NormFinder calculates expression stability according to the estimated intragroup and intergroup variation [24]. BestKeeper estimates gene expression stability on the basis of standard deviation (SD) [25]. The Delta-Ct method compares the relative expression of reference gene in pairs to rank the stability based on reproducibility of gene expression variation. [26]. In previous studies, as in ours, some discrepancies of stability rankings for reference genes were presented using different programs. The Best Keeper and Delta Ct methods identified UXT as the most stable reference gene in different skeletal muscles, but NormFinder recommended PPP1R11 as the most stable reference gene. GeNorm recommended RPL13A and UXT as the most stable genes, respectively. These divergent stability rankings for the reference genes may be attributed to the distinct mathematical models of each algorithm [38]. To address this problem, the web-based tool RefFinder was used to select the most stable reference genes for accurate normalization. It is necessary to unify the inconsistent results in previous scientific papers [18].
According to the literature, GAPDH, which is involved in primary metabolism and other cellular processes, has been mostly used to normalize the quantitative expression in yak [39,40]. It was assumed that the expression GAPDH was quite stable in many other mammalian species. However, several studies have shown that GAPDH is not an ideal reference gene in developing skeletal muscle. For example, GAPDH was considered as the optimal reference gene in human muscle [41], however GAPDH exhibited poor performance in the muscle tissues in several animals including mouse, cattle, goat, and pig [38,42–44]. Our data showed that GAPDH was identified as the least stable reference gene with any of the algorithms. These results highlight the fact that the most commonly used reference genes may vary considerably under some experiments and selection of the best reference genes is necessary under specific experimental conditions and species. The ACTB gene usually presents as a high stability gene, and it was considered as stable reference gene in skeletal muscle tissue of cattle and pig [43,45]. Li et al. investigated stabilities of ten reference genes, including GAPDH and ACTB, in six tissues of fetal yak. They found that TBP and ACTB were the two most stable, while GAPDH performed poorly [21]. However, our results showed that ACTB could be regarded as the third most stable gene based on the BestKeeper program while its ranking position indicated moderate stability in other three programs and RefFinder, this could be more applicable in various other tissue types other than the skeletal muscle of yak.
Previous studies have demonstrated that using at least two reference genes is sufficient to improve the accuracy and reliability of data normalization [17]. The GeNorm V values showed that all reference genes with a mean pairwise variation value <0.15, demonstrating that two reference genes were enough for gene expression analysis, and the addition of one more reference gene would not significantly improve accuracy [23]. In the present study, the comprehensive ranking of RefFinder analysis showed that UXT and RPL13A were the best reference genes that were expressed with the highest stability. UXT serves as a cofactor which regulates gene transcription, and it also plays an essential role in concert with the corepressor URI1 to regulate androgen receptor AR-mediated transcription [46]. UXT has previously been found suitable as a reference gene for cells and tissues of some ruminants, including taurine cattle, zebu, riverine buffalo, and goat [47–49]. Recent studies showed that UXT was also selected as a suitable reference gene in yak. During pregnancy and lactation stage, MRPS15, RPS23, and UXT were the most stable gene in yak mammary tissue [19]. Bai et al. investigated the stability of thirteen reference genes and found that RPS9, PPP1R11, UXT, and MRPL39 were the most stable reference genes to normalize RT-qPCR data from milk somatic cells of yak [20]. RPL13A encodes one of the protein components of the large 60S ribosomal subunit, and it has an important role in regulation of translation of specific mRNAs as part of a non-ribosomal complex [50,51]. Notably, some ribosomal protein genes (RPL13A, RPS5, RPL32 and RPS19) have been highly regarded as suitable reference gene in canine metencephalon and body tissues [52,53]. RPL13A was the reference gene in 11 distinct brain regions of normal humans [54] and in mouse lung tissues [55]. Therefore, to obtain the reliable gene expression results from various skeletal muscle of yak, we selected UXT and RPL13A as the most appropriate pair of reference genes.
To validate the best selected reference genes, the relative expression level of the MyoG gene was measured in muscle samples from different muscle fiber types. MyoG is a muscle-specific transcription factor that belongs to myogenic regulatory factor family. It plays a role in muscle differentiation and cell cycle exit [56]. Significant differences were revealed in the expression patterns of the MyoG gene when it was normalized with the two most stable genes (UXT and RPL13A) compared to when it was normalized with one of the least stable genes (GAPDH). These results showed that selection of reliable reference genes is necessary to improve the accuracy of normalization results.
Many previous reports showed that MyHC I, MyHC IIA, MyHC IIX and MyHC IIB isforms are expressed respectively in muscle fibers, and that the MyHC IIB isform is not present in skeletal muscles of cattle and buffalo [57,58]. We used the RT-qPCR to perform MyHC isoform mRNA expression in ten different skeletal muscles of yak, and found that MyHC IIB is expressed in yak muscle. Wang also found that the MyHC IIB isoform protein was expressed in yak muscle [59]. These differences may be attributable to species variation. In this study, there were significant differences in muscle fiber characteristics among ten different muscles. Different metabolic capacity based on the copy number of mtDNA was consistent with the muscle fiber composition. These differences may be attributed to the different location of muscle and variable requirements for oxygen diffusion in muscle fiber [7]. PM muscle is a typical oxidative muscle fiber type. Our results showed that PM muscle had the higher type I fibers and mitochondrial content, which is consistent with measurements made by Lang et al [36]. Hwang et al. reported that LD muscle of Hanwoo (Korean native cattle) steers had a higher percentage of type IIB fibers [60]. In contrast, a moderate fat content was observed in LD muscle of yak in the present study. It is therefore reasonable to assume that the muscle types and species-specific molecular regulation network is due to species variation (taurine cattle vs. yak).
Conclusions
To validate suitable reference genes for gene expression normalization in skeletal muscle tissue of yak, we selected 14 candidate reference genes using four systematic statistical algorithms (geNorm, NormFinder, Delat Ct, and BestKeeper). The obtained results were compared and ranked using RefFinder. Based on our comprehensive analysis, we identified UXT and PRL13A as the most stable reference genes for normalization of gene expression across different muscle tissues. We also examined the gene expression of muscle fiber in ten different skeletal muscles, providing a basis for possible selection towards improved meat quality of yak.
Supporting information
S1 Fig. Melting curves of ten genes of RT-qPCR.
https://doi.org/10.1371/journal.pone.0228493.s001
(TIF)
S1 Table. Ct values of candidate reference genes in all muscle samples.
https://doi.org/10.1371/journal.pone.0228493.s002
(XLS)
S2 Table. The nucleotide sequences of these 14 reference genes.
https://doi.org/10.1371/journal.pone.0228493.s003
(DOC)
S3 Table. Primer sequences for muscle fiber composition related genes.
https://doi.org/10.1371/journal.pone.0228493.s004
(DOCX)
References
- 1.
Wiener G, Jianlin H, Ruijun L. The Yak. 2nd ed: FAO Regional Office for Asia and the Pacific Food and Agriculture Organization of the United Nations, Bangkok, Thailand; 2003.
- 2. Shi Q; Guo Y; Engelhardt SC; Weladji RB; Zhou Y; Long M, et al. Endangered wild yak (Bos grunniens) in the Tibetan plateau and adjacent regions: Population size, distribution, conservation perspectives and its relation to the domestic subspecies. J Nat Conserv. 2016;32: 35–43.
- 3. Pette D; Staron RS. Transitions of muscle fiber phenotypic profiles. Histochem Cell Biol. 2001;115: 359–372. pmid:11449884
- 4. Schiaffino S; Reggiani C. Fiber types in mammalian skeletal muscles. Physiol Rev. 2011;91: 1447–531. pmid:22013216
- 5. Joo ST; Kim GD; Hwang YH; Ryu YC. Control of fresh meat quality through manipulation of muscle fiber characteristics. Meat Sci. 2013;95: 828–836. pmid:23702339
- 6. Ryu YC; Choi YM; Lee SH; Shin HG; Choe JH; Kim JM, et al. Comparing the histochemical characteristics and meat quality traits of different pig breeds. Meat Sci. 2008;80: 363–369. pmid:22063341
- 7. Lefaucheur L. A second look into fibre typing–Relation to meat quality. Meat Sci. 2010;84: 257–270. pmid:20374784
- 8. Kim GD; Jeong JY; Yang HS; Joo ST. Identification of myosin heavy chain isoforms in porcine longissimus dorsi muscle by electrophoresis and mass spectrometry. Electrophoresis. 2013;34: 1255–61. pmid:23463416
- 9. Choi YM; Kim BC. Muscle fiber characteristics, myofibrillar protein isoforms, and meat quality. Livest Sci. 2009;122: 105–118.
- 10. Sazili AQ; Parr T; Sensky PL; Jones SW; Bardsley RG, et al. The relationship between slow and fast myosin heavy chain content, calpastatin and meat tenderness in different ovine skeletal muscles. Meat Sci 2005;69: 17–25. pmid:22062635
- 11. Picard B; Barboiron C; Duris MP; Gagniére H; Jurie C, et al.Geay Y. Electrophoretic separation of bovine muscle myosin heavy chain isoforms. Meat Sci. 1999;53: 1–7. pmid:22062926
- 12. Serrano AL; Murgia M; Pallafacchina G; Calabria E; Coniglio P; Lomo T, et al. Calcineurin controls nerve activity-dependent specification of slow skeletal muscle fibers but not muscle growth. Proc Natl Acad Sci U S A. 2001;98: 13108–13113. pmid:11606756
- 13. Ennion S; Sant’Ana Pereira J; Sargeant AJ; Young A; Goldspink G. Characterization of human skeletal muscle fibers according to the myosin heavy chains they express. J Muscle Res Cell Motil. 1995;16: 35–43. pmid:7751403
- 14. Zurmanová J; Soukup T. Comparison of myosin heavy chain mRNAs, protein isoforms and fiber type proportions in the rat slow and fast muscles. Physiol Res. 2013;62: 445. pmid:23590611
- 15. Derveaux S; Vandesompele J; Hellemans J. How to do successful gene expression analysis using real-time PCR. Methods. 2010;50: 227–230. pmid:19969088
- 16. Chen G; Kronenberger P; Teugels E; De Greve J. Influence of RT-qPCR primer position on EGFR interference efficacy in lung cancer cells. Biol Proced Online. 2010;13: 1. pmid:21369532
- 17. Gomes AÉI; Stuchi LP; Siqueira NMG; Henrique JB; Vicentini R; Ribeiro ML, et al. Selection and validation of reference genes for gene expression studies in Klebsiella pneumoniae using Reverse Transcription Quantitative real-time PCR. Sci Rep. 2018;8: 9001. pmid:29899556
- 18. Sarker N; Fabijan J; Emes RD; Hemmatzadeh F; Meers J; Moreton J, et al. Identification of stable reference genes for quantitative PCR in koalas. Sci Rep. 2018;8: 3364. pmid:29463845
- 19. Jiang M; Lee JN; Bionaz M; Deng XY; Wang Y. Evaluation of Suitable Internal Control Genes for RT-qPCR in Yak Mammary Tissue during the Lactation Cycle. PLoS One 2016;11: e0147705. pmid:26808329
- 20. Bai WL; Yin RH; Zhao SJ; Jiang WQ; Yin RL; Ma ZJ, et al. Technical note: Selection of suitable reference genes for studying gene expression in milk somatic cell of yak (Bos grunniens) during the lactation cycle. J Dairy Sci. 2014;97: 902–910. pmid:24342693
- 21. Li M, Wu X, Guo X, Bao P, Ding X, Chu M, et al. Identification of optimal reference genes for examination of gene expression in different tissues of fetal yaks. Czech J Anim Sci. 2017;62: 426–434.
- 22. Wu X; Zhou X; Ding X; Chu M; Liang C; Pei J, et al. The selection of reference genes for quantitative real-time PCR in the Ashidan yak mammary gland during lactation and dry period. Animals. 2019;9: 943.
- 23. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, et al. Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biol. 2002;3: research0034.1–0034.11.
- 24. Andersen CL, Jensen JL, Ørntoft TF. Normalization of real-time quantitative reverse transcription-PCR data: a model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets. Cancer Res. 2004;64: 5245–5250. pmid:15289330
- 25. Pfaffl MW, Tichopad A, Prgomet C, Neuvians TP. Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper—Excel-based tool using pair-wise correlations. Biotechnol Lett. 2004;26: 509–515. pmid:15127793
- 26. Silver N; Best S; Jiang J; Thein SL. Selection of housekeeping genes for gene expression studies in human reticulocytes using real-time PCR. BMC Mol Biol. 2006;7: 33. pmid:17026756
- 27. Xie F; Xiao P; Chen D; Xu L; Zhang B. miRDeepFinder: a miRNA analysis tool for deep sequencing of plant small RNAs. Plant Mol Biol. 2012;80: 75–84.
- 28. Feng X; Xiong Y; Qian H; Lei M; Xu D; Ren Z. Selection of reference genes for gene expression studies in porcine skeletal muscle using SYBR green qPCR. J Biotechnol. 2010;150: 288–293. pmid:20887758
- 29. McBryan J; Hamill RM; Davey G; Lawlor P; Mullen AM. Identification of suitable reference genes for gene expression analysis of pork meat quality and analysis of candidate genes associated with the trait drip loss. Meat Sci. 2010;86: 436–439. pmid:20579813
- 30. Nascimento CS; Barbosa LT; Brito C; Fernandes RPM; Mann RS; Pinto APG, et al. Identification of suitable reference genes for real time quantitative polymerase chain reaction assays on pectoralis major muscle in chicken (Gallus gallus). Plos One. 2015;10: e127935.
- 31. Hildyard JCW; Taylor-Brown F; Massey C; Wells DJ; Piercy RJ. Determination of qPCR Reference genes suitable for normalizing gene expression in a canine model of duchenne muscular dystrophy. Journal of Neuromuscular Diseases. 2018;5: 177–191. pmid:29614692
- 32. Livak KJ; Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods. 2001;25: 402–408. pmid:11846609
- 33. Ma J; Wang H; Liu R; Jin L; Tang Q; Wang X, et al. The miRNA Transcriptome Directly Reflects the Physiological and Biochemical Differences between Red, White, and Intermediate Muscle Fiber Types. Int J Mol Sci. 2015;16: 9635–9653. pmid:25938964
- 34. Facucho-Oliveira JM; Alderson J; Spikings EC; Egginton S; St John JC. Mitochondrial DNA replication during differentiation of murine embryonic stem cells. J Cell Sci. 2007;120: 4025–4034. pmid:17971411
- 35. Zhao H; Ma T; Lin J; Liu L; Sun W; Guo L, et al. Identification of valid reference genes for mRNA and microRNA normalisation in prostate cancer cell lines. Sci Rep. 2018;8: 1949. pmid:29386530
- 36. Lang Y; Wang Y; Sun B; Zhang S; Xie P; Sha K, et al. Myofiber characteristics and eating quality of three major muscles from Chinese Simmental cattle. Can J Anim Sci. 2016;97: 101–108.
- 37. Liu X; Guan H; Song M; Fu Y; Han X; Lei M, et al. Reference gene selection for qRT-PCR assays in Stellera chamaejasme subjected to abiotic stresses and hormone treatments based on transcriptome datasets. PeerJ. 2018;6: e4535. pmid:29632740
- 38. Zhu W; Lin Y; Liao H; Wang Y. Selection of reference genes for gene expression studies related to intramuscular fat deposition in Capra hircus skeletal muscle. Plos One. 2015;10: e0121280. pmid:25794179
- 39. Xin J; Chai Z; Zhang C; Zhang Q; Zhu Y; Cao H, et al. Transcriptome profiles revealed the mechanisms underlying the adaptation of yak to high-altitude environments. Sci Rep. 2019;9: 7558. pmid:31101838
- 40. Ding Y; Zhang N; Li J; Jin Y; Shao B. Molecular cloning and expression of ghrelin in the hypothalamus-pituitary-gastrointestinal tract axis of the Yak (Bos grunniens) in the Qinghai-Tibetan Plateau. Anatomia, Histologia, Embryologia. 2018;47: 583–590. pmid:30178622
- 41. Mahoney DJ; Carey K; Fu M; Snow R; Cameron-Smith D; Parise G, et al. Real-time RT-PCR analysis of housekeeping genes in human skeletal muscle following acute exercise. Physiol Genomics. 2004;18: 226–231. pmid:15161965
- 42. Uddin MJ; Cinar MU; Tesfaye D; Looft C; Tholen E; Schellander K. Age-related changes in relative expression stability of commonly used housekeeping genes in selected porcine tissues. BMC research notes. 2011;4: 441. pmid:22023805
- 43. Pérez R; Tupac-Yupanqui I; Dunner S. Evaluation of suitable reference genes for gene expression studies in bovine muscular tissue. BMC Mol Biol. 2008;9: 79. pmid:18786244
- 44. Thomas KC; Zheng XF; Garces Suarez F; Raftery JM; Quinlan KGR; Yang N, et al. Evidence based selection of commonly used RT-qPCR reference genes for the analysis of mouse skeletal muscle. Plos One. 2014;9: e88653. pmid:24523926
- 45. Erkens T; Van Poucke M; Vandesompele J; Goossens K; Van Zeveren A; Peelman LJ. Development of a new set of reference genes for normalization of real-time RT-PCR data of porcine backfat and longissimus dorsi muscle, and evaluation with PPARGC1A. BMC Biotechnol. 2006;6: 41. pmid:17026777
- 46. Mita P; Savas JN; Djouder N; Yates JR; Ha S; Ruoff R, et al. Regulation of Androgen Receptor-Mediated Transcription by RPB5 Binding Protein URI/RMP. Mol Cell Biol. 2011;31: 3639–3652. pmid:21730289
- 47. Bonnet M; Bernard L; Bes S; Leroux C. Selection of reference genes for quantitative real-time PCR normalisation in adipose tissue, muscle, liver and mammary gland from ruminants. Animal. 2013;7: 1344–1353. pmid:23552195
- 48. Choudhary R; Kumar S; Singh SV; Sharma AK; Goud TS; Srivastava AK, et al. Validation of putative reference genes for gene expression studies in heat stressed and α-MSH treated melanocyte cells of Bos indicus using real-time quantitative PCR. Mol Cell Probe. 2016;30: 161–167.
- 49. Kaur R; Sodhi M; Sharma A; Sharma VL; Verma P; Swami SK, et al. Selection of suitable reference genes for normalization of quantitative RT-PCR (RT-qPCR) expression data across twelve tissues of riverine buffaloes (Bubalus bubalis). Plos One. 2018;13: e0191558. pmid:29509770
- 50. Xue S; Barna M. Specialized ribosomes: a new frontier in gene regulation and organismal biology. Nat Rev Mol Cell Bio. 2012;13: 355–369.
- 51. Mazumder B; Sampath P; Seshadri V; Maitra RK; Dicorleto PE; Fox PL. Regulated release of L13a from the 60S ribosomal subunit as a mechanism of transcript-specific translational control. Cell. 2003;115: 187–198. pmid:14567916
- 52. Park S; Huh J; Kim Y; Lee S; Kim S; Kim S, et al. Selection of Internal Reference Genes for Normalization of Quantitative Reverse Transcription Polymerase Chain Reaction (qRT-PCR) Analysis in the Canine Brain and Other Organs. Mol Biotechnol. 2013;54: 47–57. pmid:22531949
- 53. Peters IR; Peeters D; Helps CR; Day MJ. Development and application of multiple internal reference (housekeeper) gene assays for accurate normalization of canine gene expression studies. Vet Immunol Immunop. 2007;117: 55–66.
- 54. Grube S; Göttig T; Freitag D; Ewald C; Kalff R; Walter J. Selection of suitable reference genes for expression analysis in human glioma using RT-qPCR. J Neuro-Oncol. 2015;123: 35–42.
- 55. Yin R; Tian F; Frankenberger B; de Angelis MH; Stoeger T. Selection and evaluation of stable housekeeping genes for gene expression normalization in carbon nanoparticle-induced acute pulmonary inflammation in mice. Biochem Bioph Res Co. 2010;399: 531–536.
- 56. Hernández-Hernández JM; García-González EG; Brun CE; Rudnicki MA. The myogenic regulatory factors, determinants of muscle development, cell identity and regeneration. Semin Cell Dev Biol. 2017;72: 10–18. pmid:29127045
- 57. Chikuni K; Muroya S; Nakajima I. Absence of the functional myosin heavy chain 2b isoform in equine skeletal muscles. Zoolog Sci. 2004;21: 589–596. pmid:15170063
- 58. Francisco CL; Jorge AM; Dal-Pai-Silva M; Carani FR; Cabeço LC; Silva SR. Muscle fiber type characterization and myosin heavy chain (MyHC) isoform expression in Mediterranean buffaloes. Meat Sci. 2011;88: 535–541. pmid:21371827
- 59.
Wang, L. The effect of muscle fiber types and their differences in metabolic enzymes activity on postmortem tenderness of yak meat. M.Sc. Thesis, Gansu Agricultural University.2016.
- 60. Hwang Y; Kim G; Jeong J; Hur S; Joo S. The relationship between muscle fiber characteristics and meat quality traits of highly marbled Hanwoo (Korean native cattle) steers. Meat Sci. 2010;86: 456–461. pmid:20598446