Reference gene selection for quantitative PCR in liver, skeletal muscle, and jejunum of Bos indicus cattle

- The objective of the present study was to evaluate the stability of candidate reference genes and select the genes that can be used for normalizing real-time polymerase chain reaction (PCR) in the liver, skeletal muscle, and jejunum tissues of Nellore or Nellore × Angus steers fed different diets. Fourteen purebred and 14 crossbred steers were used, in which half of the animals of each genetic group received a diet containing whole shelled corn (WSC) and the other half whole shelled corn and sugarcane bagasse (WSCB). Stability was analyzed by the RefFinder program. To validate the selection of candidate reference genes, the expression of target genes was evaluated using the different groups of reference genes. The most stable genes were 18S , ACTB , and CASC3 for skeletal muscle; HMBS , ACTB , and 18S for the liver; and GAPDH , ACTB , and CASC3 for the jejunum, regardless of breed and diet provided. Possible errors caused in data analyses were clarified comparing the more and less stable genes as reference for normalization of the target genes FASN , ACOX , SCD1 , MGAM , and SLC2A1 . The use of the more stable and less stable sets of reference genes may lead to different conclusions in respect to the expression profile of the target studied gene. The selection of more suitable reference genes for each experiment is of utmost importance to ensure the reliability of gene expression studies so that they can be applied in practice.


Introduction
Elucidation of various molecular mechanisms, previously unknown, occurred after discovery of the polymerase chain reaction (PCR) method in the scientific community, which has become prominent in numerous research studies involving animal production (Oliveira et al., 2014;Teixeira et al., 2017;Lopes et al., 2020).
The relative quantification of reverse transcription PCR (RT-qPCR) is an accurate and sensitive method to quantify mRNA levels of target genes. However, data processing can affect interpretation of results (Lisowski et al., 2008), being necessary a normalization to correct inevitable experimental variations. The most common method of data normalization is the reference genes, which generally have a relatively stable expression pattern in tissues, since they are required for basic cellular processes (Chapman and Waldenström, 2015). However, it is becoming increasingly apparent that the stability of some usual reference genes can change according to species, tissue types, and experimental conditions (Lisowski et al., 2008;Kaur et al., 2018;Lozano-Villegas et al., 2021). Therefore, selecting reference genes with less variation has become a requirement for RT-qPCR analysis, and an important factor for obtaining reliable results. If gene expression is normalized with a reference gene that fluctuate under experimental conditions, relevant changes in expression of the target genes can easily be missed or overemphasized (Chapman and Waldenström, 2015).
Various studies have already described selection of reference genes associated with different experimental conditions in bovines and have selected reference genes related to the mammary gland, fat, liver (Saremi et al., 2012;Bonnet et al., 2013), skeletal muscle (Mberema and Sparagano, 2017), somatic cells (Muhaghegh-Dolatabady et al., 2017), blood (Lozano-Villegas et al., 2021), and kidney, pituitary gland, and thyroid (Lisowski et al., 2008). However, there are no studies regarding the reference genes more suitable to muscle, liver, and jejunum in Zebu cattle under nutritional stress.
The hypothesis of this study is that the stability of reference genes differs among genetic groups and diets, and that this may result in errors in analysis of gene expression. Thus, this study aimed to evaluate the stability of candidate reference genes and select the genes that can be used for normalizing RT-PCR in the liver, skeletal muscle, and jejunum tissues of Nellore or Nellore × Angus steers fed different diets.

Animals, diets, and sample collection
Research on animals was conducted according to ethical standards and approved by the local Ethics and Animal Welfare Committee (protocol 056/15).
Fourteen Nellore purebred steers and 14 Angus × Nellore crossbred steers with mean age of 25 months and body weight of 353±25.3 kg were housed in individual pens into a completely randomized design in a 2 × 2 factorial arrangement (two diets and two breeds). Half of the animals of each genetic group received a diet containing whole shelled corn (WSC) and the other half received whole shelled corn and sugarcane bagasse (WSCB) ( Table 1).
The animals were slaughtered after the experimental period of 96 days preceded by a 20-day adaptation period. Muscle, liver, and jejunum samples were collected soon after slaughter, labeled, and placed in 5-mL cryogenic tubes, frozen, and transported in liquid nitrogen until being stored in an ultra-freezer (−80 °C).

Design and efficiency of the primers
The primers of the reference genes were designed by means of sequences registered and published in the GenBank public database, at the National Center for Biotechnology Information (NCBI) platform. Characterization of the genes and the Open Reading Frames (ORF) of the selected sequences were obtained by the ORFinder tool of the NCBI. The primers were designed using the Oligo Perfect Designer software based on the sequence accessed in the GenBank and analyzed using OligoAnalizer 3.1.
2.3. RNA extraction and cDNA synthesis RNA extraction was performed using two methods, in accordance with the tissue studied to avoid degradation and also to get better quality and quantity. RNA was extracted from muscle samples (80-100 mg) using the reagent Trizol (Invitrogen, Gaithersburg, MD, USA) according to manufacturer's instructions and was homogenized (IKA ® T18 ULTRA-TURRAX ® Basic Homogenizer, Wilmington, NC, USA) together with 1 mL of Trizol reagent. Adaptations were made in relation to the amounts Reference gene selection for quantitative PCR in liver, skeletal muscle, and jejunum of Bos indicus cattle Coelho et al. 4 of chloroform (300 µL) and isopropanol (400 µL), and the step of chloroform addition was performed twice. Samples were treated with DNase DNA free (Ambion, Austin, TX, USA) to eliminate possible contaminations.
The liver and jejunum RNA were extracted using the Promega SV RNA Isolation kit (Promega, Madison, WI, USA) according to manufacturers. Tissue samples were weighed (130 mg) and homogenized (IKA ® T18 ULTRA-TURRAX ® Basic Homogenizer, Wilmington, NC, USA with 1 mL of RNA Lysis Buffer (RLA).
Samples were then quantified in a nanospectrophotometer (DeNovix DS-11, Wilmington, DE, USA) at A260 nm for the purpose of determining the amount (ng/μL) and purity (260/280 and 260/230) of the extracted RNA. After these procedures, the RNA samples were stored at −80 ℃.
Two methods were used for cDNA synthesis: the High-Capacity cDNA Reverse Transcription kit (Applied Biosystems, Foster City, CA, USA) for muscle samples, and the GoScript™ reverse transcription kit (Promega, Madison, WI, USA) for liver and jejunum samples, according to manufacturers' specifications.

Quantitative gene expression analysis
For quantitative gene expression analysis by RT-qPCR, the Mastercycler ® ep Realplex (Eppendorf) model was used with the SYBR Green detection system (Applied Biosystems, Foster City, CA, USA) and the cDNA obtained from RNA extracted from the liver, muscle, and jejunum of bovines. The thermal reaction conditions were 2 min at 50 ℃ and 10 min at 95 ℃, followed by 40 cycles of 15 s at 95 ℃ and 1 min at 60 ℃, ending with the melting curve. For each reaction, 1 μL of cDNA at the concentration of 1:5, 0.3 μL of each primer (1.2 µM), and 5.0 μL of Master Mix SYBR Green were used, for a final volume of 10.0 μL/sample on a MicroAmp Optical 96-well reaction plate (Applied Biosystems, Foster City, CA, USA).
The entire RT-qPCR experiment for each gene under study was conducted from cDNAs obtained from seven biological replications, with three technical replicates for each one. The Cq was determined by the number of cycles in which the fluorescence generated within a reaction crosses the baseline (threshold).

Evaluation of reference genes
The data were calibrated by the GenEX program (MultiD, Gothenburg, Sweden) for correcting the values according to the efficiency found for each gene evaluated. The SigmaPlot 12.0 software was employed for construction of boxplot diagrams to illustrate the levels and variations in expression of candidate reference genes. Data were then analyzed using the RefFinder online program, and thus the more stable genes for each condition were determined. Pairwise evaluation was carried out by the GeNorm algorithm to determine the minimum number of genes necessary; we needed to add one gene as a reference when the standard deviation (SD) was greater than 0.15.

Validation of expression of reference genes
For validation of the reference genes, a completely randomized design in a 2 × 2 factorial arrangement was used, with four treatments [Nellore steers fed a whole shelled corn diet (NWSC); Nellore steers fed a whole shelled corn and sugarcane bagasse diet (NWSCB); Nellore × Angus steers fed a whole shelled corn diet (NAWSC); and Nellore × Angus steers fed a whole shelled corn and sugarcane bagasse diet (NAWSCB)] and seven replications per treatment. Relative expression of target genes were calculated according to the method described by Pfaffl (2001), which is based on Cq values that are corrected for the amplification efficiency of each primer pair.
Reference gene selection for quantitative PCR in liver, skeletal muscle, and jejunum of Bos indicus cattle Coelho et al. 5 Gene expressions were analyzed using the PROC MIXED procedure of SAS (Statistical Analysis System, version 9.2), through the model: in which γ ij is the response variable measured, μ is the overall mean, α i is the i-th level of the fixed effect of diet (two levels), β j is the j -th level of the fixed effect of breed (two levels), αβ ij is the interaction of breed (β j ) and diet (α i ), and e ij is the random error associated with γ ij .
A Shapiro-Wilk test was used to analyze data normality, when gene expression data did not exhibit normal distribution, the data were transformed using the RANK procedure of SAS. After that, analysis of variance (ANOVA) and Tukey's test were used on the validation data of reference genes to determine statistical differences among the means of the different treatments.

Results
Levels of gene expression were determined from values of the Quantification cycle (Cq) of the seven reference genes evaluated in the different tissues of the animals, regardless of their genetic group and diets. Simultaneous evaluation of the effect of breeds and diets on muscle, liver, and jejunum tissue ( Figure 1) showed that the lowest mean Cq values were 11.9, 11.9, and 13.7 and the highest were 25.7, 27, and 27.4, respectively.
The ranking of candidate reference genes was obtained by the RefFinder program from the geometric mean of the results from the software programs GeNorm, BestKeeper, NormFinder, and  Delta Cq. The minimum number of genes for tissues, breeds, and the combination between variables were determined by pairwise analysis by the GeNorm algorithm (Vandesompele et al., 2002). Vandesompele et al. (2002) proposed 0.15 as the cutoff value for the paired variation. After inclusion of the second reference gene, a substantial reduction in the V value (0.15) was observed in GeNorm, showing the need for at least two and at most four reference genes for normalization of the study for the different conditions of this study (Tables 3-5). For the muscle of pure and crossbred animals fed WSC four reference genes were needed (Table 3); on the other hand, for liver (Table 4) and jejunum (Table 5), the minimum number of genes was three.
The results obtained showed that the ranking (more stable to the less stable gene) for the muscle according to the RefFinder program were Ribosomal Protein 18S (18S), Beta-actin (ACTB), Cancer susceptibility candidate 3 (CASC3), Eukaryotic translation elongation factor 1 alpha 2 (EEF1A2), Glyceraldehyde-3-phosphate (GAPDH), Hydroxymethylbilane synthase mRNA (HMBS), and Ubiquitin C (UBC) ( Table 6), being 18S, ACTB, CASC3, and EEF1A2 (Table 3) the most stable genes. The use of CASC3 in purebred animals and of EEF1A2 was expendable in comparing the different genetic groups and diets simultaneously (Table 3). For the liver, the ranking of the reference genes (RefFinder method) was HMBS, ACTB, 18S, CASC3, GAPDH, UBC, and EEF1A2 (Table 7) and the most stable genes were HMBS,  ACTB, 18S, and GAPDH (Table 4). In the crossbred animals receiving the two diets or only WSC, GAPDH was more stable than the 18S gene that was selected in purebred animals. Finally, in the jejunum, the most stable reference genes were GAPDH, ACTB, CASC3, and 18S (Table 5), considering that the ranking of the reference genes were GAPDH, ACTB, CASC3, UBC, HMBS, 18S, and EEF1A2 (Table 8). In general, there was a small effect of diet on the reference genes in the tissues.

Discussion
A reference gene should have expression near the target gene; however, that is not always possible. Thus, the general directives recommend that the Cq values be from 15 to 30 (Wan et al., 2010). In spite of the wide variation found in the Cq of the genes tested (Figure 1), the values found in    Some studies have already described the selection of reference genes associated with different tissues and conditions in bovines ( Table 9) and showed that the choice of genes is highly dependent on experimental conditions and that most of the studies were carried out with milk cattle and taurine cattle. Saremi et al. (2012) found different reference genes recommended for the liver and muscle tissues of bovines fed conjugated linoleic acid (CLA). Lisowski et al. (2008) evaluated different cattle breeds in different developmental stages and likewise found that two or more reference genes are necessary and that they differ according to the tissue/organ evaluated (liver, kidney, pituitary gland, or thyroid gland).
Studies carried out prior to publication of the Minimum Information for Publication of Quantitative Real-Time PCR Experiments (MIQE guidelines) used only one reference gene (Hamalainen et al., 2001;Deindl et al., 2002;Glare et al., 2002) to normalize RT-qPCR data. However, after its publication, the use of two or more reference genes has been recommended to ensure stable expression under all experimental conditions (Bustin et al., 2009). In this study, as a consequence of the data analysis by the GeNorm software, the need for at least two, and up to four, reference genes was shown for normalization of the study for the different breeds and diets studied and for the skeletal muscle, liver, and jejunum tissues. Likewise, Pérez et al. (2008) and Saremi et al. (2012) showed the need for three genes for muscle; and Lisowski et al. (2008) and Saremi et al. (2012) showed the need for two to three reference genes for the liver. However, we did not find studies that investigated the reference genes more suitable for the jejunum of bovines, which is important for studies of gene expression of membrane transport proteins.
A review of literature regarding reference genes in gene expression studies showed that since the publication of the MIQE, the most common choice of reference genes were ACTB (used in 38% of the studies) and GAPDH (37%); the 18S gene (12%) was a little less used in vertebrates (Chapman and Waldenström, 2015). For the tissues evaluated, at least one of the three genes cited above was stable, making them strong candidates for reference genes in studies with RT-qPCR in zebu cattle and their crosses. The ACTB was among the more stable reference genes under all the conditions in this study (Tables 3-5). One of the reasons for this is that ACTIN is the most abundant protein in eukaryotic cells (Nakajima-Iijima et al., 1985). There are at least six ACTIN isoforms in vertebrates: four types in muscle (skeletal, cardiac, vascular smooth, and stomach smooth actins) and two non-muscular types (β-cytoplasmic and γ-cytoplasmic actins) (Vandekerckhove and Weber, 1979), and they are expressed in all types of tissues, as shown in the results of this study, making this gene an optimal normalizer (Tables 3-5). Reference gene selection for quantitative PCR in liver, skeletal muscle, and jejunum of Bos indicus cattle Coelho et al. 9 The EEF1A2 gene was considered the least stable reference gene under most of the conditions studied in the liver and jejunum, and UBC was the least stable in muscle (Tables 6-8). Both are genes expressed under diverse conditions; however, in this study, there was greater variation in relation to the tissues studied. The EEF1A gene has two isoforms, EEF1A1 and EEF1A2, that codify the second most abundant protein (after ACTIN) and catalyzes bonding of tRNA to the ribosomal site. The isoform EEF1A2 is present especially in the brain, heart, and muscle (Lee and Surh, 2009). The expression of this gene in muscle was greater and more stable, with a mean Cq value of 17.6 (Figure 1a), compared with liver ( Figure 1b) and jejunum (Figure 1c) tissues, with mean values of 27 and 26.2, respectively. In the latter two tissues, this gene was selected as one of the three least stable (Tables 7 and 8), whereas for muscle, the gene was stable in the purebred and crossbred animals, regardless of the diet fed (Table 6).
In muscle, the UBC gene was highly variable in the breeds and diets studied, and thus, it is not considered a good normalizer for this tissue (Figure 1a). A possible explanation for this is that the gene codifies a protein that participates in selective degradation of short-lived proteins in eukaryotic cells and is involved in immune response and cell development and death (Hershko and Ciechanover, 1998). UXT, RPS9, and RPS15 ACTB and GAPDH Bionaz and Loor (2007) Based on the results that show differences in the stability of candidate reference genes, statistical analyses were performed using the more stable and less stable genes found in this study to reveal possible errors in analysis of gene expression of target genes.
Fatty Acid Synthase (FASN) and Acyl-CoA Oxidase 1 (ACOX) are genes involved in synthesis and oxidation of saturated fatty acids, and the protein coded by ACOX is the first enzyme of beta-oxidation of fatty acids in the peroxisomes (Schulz, 1996). The gene FASN, in its turn, codifies the enzyme responsible for the final steps of de novo synthesis of long-chain fatty acids from acetyl-CoA, malonyl-CoA, and NADPH (Ladeira et al., 2016). When the more stable genes were used for normalization of the study (Figure 2), both FASN and ACOX target genes in the muscle tissue were more expressed in WSC diets, and FASN was also more expressed in crossbred animals. However, upon normalizing the data with the less stable genes, statistical differences were not found for the diets studied, and only tendencies were found for the breeds in both genes. Therefore, these results show that the wrong choice of reference genes would lead to misinterpretation and distortion of results.
The target gene Stearoyl-CoA Desaturase (SCD1) was evaluated in the liver (Figure 3), and Maltase-Glucoamylase 2 (MGAM) and Solute Carrier Family 2 Member 1 (SLC2A1) in the jejunum (Figure 4). The SCD1 gene is present throughout the organism and is highly expressed in the liver (Ntambi and Miyazaki, 2004). It is responsible for codifying an enzyme responsible for desaturation of long-chain fatty acids, thus acting in regulation of membrane fluidity (Ntambi et al., 1988;Ntambi and Miyazaki, 2004). Evaluation of expression of the SCD1 gene in the liver showed that upon using the more stable genes, statistical effects of the diets and breeds would not be found. In contrast, with the use of the less stable genes, the inclusion of sugarcane bagasse in the diet would increase expression of the gene in the N -Nellore; NA -Nellore × Angus; WSC -whole shelled corn; WSCB -whole shelled corn and sugarcane bagasse; B×D -breed and diet interaction; SEM -standard error of the mean.   (Figure 3). In the jejunum, the MGAM gene is expressed in the brush border of the small intestine membrane and codifies an enzyme responsible for hydrolyzing maltose into glucose (Mochizuki et al., 2010). Statistical differences were not found in expression of the MGAM gene in the breeds and diets evaluated when the less stable reference genes were used. However, upon using the more stable ones, there was a tendency that the diets containing sugarcane bagasse would result in greater expression N -Nellore; NA -Nellore × Angus; WSC -whole shelled corn; WSCB -whole shelled corn and sugarcane bagasse; B×D -breed and diet interaction; SEM -standard error of the mean.   (Figure 4). Finally, the SLC2A1, which is expressed in the enterocytes of the small intestine and codifies the sodium-glucose linked transporter (SGLT1), has its expression similar for diets and breeds, regardless of selection of more stable or less stable genes (Figure 4). In this case, there was statistical difference for diets in both groups of reference genes used as normalizers in this study. However, the relative expression with the use of the less stable genes was much greater, which can also lead to poor interpretation of data by authors.
It can be inferred that the choice of reference genes highly affects the results of analyses of gene expression that can occur due to differences brought about under determined metabolic conditions of determined tissues. Thus, there is the possibility of results exhibiting wide variation in gene expression. In the present study, there were few variations in the expression and stability of reference genes among the genetic groups and the diets studied within a specific tissue. However, analysis of the expression of the reference genes among the tissues showed that expressions of GAPDH and EEF1A2 were greater and less stable in muscle compared with the liver and jejunum. Therefore, these genes would likely not be considered ideal for normalization of studies comparing these three tissues.

Conclusions
The more stable genes for zebu animals and their crosses were 18S, ACTB, and CASC3 for muscle tissue; HMBS, ACTB, and 18S for liver tissue; and GAPDH, ACTB, and CASC3 for jejunum tissue, regardless of the breed and diet. Comparison of the more stable and less stable genes as a reference for normalization of target genes shows that the wrong choice of these genes will lead to misleading interpretation of the results.

Conflict of Interest
The authors declare no conflict of interest.

Author Contributions
Conceptualization