Selection of Reliable Reference Genes for Real-time qRT-PCR Analysis of Zi Geese (Anser anser domestica) Gene Expression
Article information
Abstract
Zi geese (Anser anser domestica) belong to the white geese and are excellent layers with a superior feed-to-egg conversion ratio. Quantitative gene expression analysis, such as Real-time qRT-PCR, will provide a good understanding of ovarian function during egg-laying and consequently improve egg production. However, we still don’t know what reference genes in geese, which show stable expression, should be used for such quantitative analysis. In order to reveal such reference genes, the stability of seven genes were tested in five tissues of Zi geese. Methodology/Principal Findings: The relative transcription levels of genes encoding hypoxanthine guanine phosphoribosyl transferase 1 (HPRT1), β-actin (ACTB), β-tubulin (TUB), glyceraldehyde-3-phosphate-dehydrogenase (GADPH), succinate dehydrogenase flavoprotein (SDH), 28S rRNA (28S) and 18S rRNA (18S) have been quantified in heart, liver, kidney, muscle and ovary in Zi geese respectively at different developmental stages (1 d, 2, 4, 6 and 8 months). The expression stability of these genes was analyzed using geNorm, NormFinder and BestKeeper software. Conclusions: The expression of 28S in heart, GAPDH in liver and ovary, ACTB in kidney and HPRT1 in muscle are the most stable genes as identified by the three different analysis methods. Thus, these genes are recommended for use as candidate reference genes to compare mRNA transcription in various developmental stages of geese.
INTRODUCTION
The research on improving egg production is always the focus on poultry breeding and management, especially for geese because of their poor laying performance and reproductive seasonality (Shi et al., 2008; Kang et al., 2009). Zi geese (Anser anser domestica) are excellent layers with a superior feed-to-egg conversion ratio. This species breeds only in northeast area of China including Heilongjiang and Jilin provinces. In domestic fowl breeding programs, the traditional method to improve egg production is to select those with the biggest egg-laying amount or the rate of lay (Kuhnlein et al., 1997). With the development of molecular biotechnologies such as genomics and proteomics analysis, alternative methods including selection of breeders with egg-laying or other important trait marker genes or proteins are incorporated. Recent progress in molecular breeding technologies has provided tools to study complex biological traits under different physiological conditions using quantitative genes expression analysis (Chen et al., 2007). Among them, selection of breeders by comparing mRNA transcription between different domestic fowl samples has become a very important molecular biological protocol for improving egg production (Kuhnlein et al., 1997; Yen et al., 2006; Chen et al., 2007; Kang et al., 2009). In order to improve production performance of geese, researches on breeding and genetics should be focused on molecular genetic markers mapping, genome analysis and identification of candidate genes for specific performance. Real-time quantitative reverse transcription polymerase chain reaction (Real-time qRT-PCR) is a major development of PCR technology that enables reliable detection and measurement of DNA (cDNA) generated during each cycle of PCR process (Arya et al., 2005).
Real-time qRT-PCR, with its capacity to detect and measure very small amount of nucleic acids in a wide range of samples from numerous sources, has been used extensively in molecular biology, e.g., molecular diagnostics, life sciences, agriculture, and medicine. Real-Time qRT-PCR represents a rapid and reliable method for the detection and quantification of mRNA transcription levels of a selected gene in various biological specimens, or at different developmental stages or different physiological status (McCurley and Callard, 2008; Beekman et al., 2011). However, recently, a growing body of research have demonstrated that these genes expression can change in different tissues, during growth and differentiation, in response to biochemical stimuli, and in disease states (Janovick-Guretzky et al., 2007; Wen and Mao, 2007). The expression levels of the ideal endogenous reference genes should be constant in different experimental conditions. So the limitation for the application of qRT-PCR is the need for suitable internal reference genes which reduce the specimen differences and allow the quantification of this gene expression to be comparable (Huggett et al., 2005). It is well-known that normalization is critically important to reduce sample-to-sample variations, which includes the RT efficiency and RNA integrity, cDNA sample loading, instrumental errors, and the presence of PCR inhibitors etc. (Stahlberg et al., 2004; Bustin et al., 2005). Normalization of target gene expression levels must be performed before doing relative compaisions (Pfaffl et al., 2001). Several normalization strategies have been proposed, and the use of endogenous reference genes is currently the preferred one. An ideal endogenous reference gene used for normalization of data in quantitative real-time PCR should have the following features: constant expression levels among all individuals, organs and cells, during different developmental stages, and various experimental treatments (Jin et al., 2004). Usually, the well-known housekeeping genes are chosen such as the glyceraldehydes-3-phosphate dehydrogenase (GAPDH), β-actins (ACTB), hypoxantine phosphoribosyltransferase (HRPT), β-tubulin (TUB), elongation factor 1 alpha (EF1A) and 18S, 28S rRNA (Shu et al., 2004; Fernandes et al., 2008; Øvergard et al., 2010). Several commonly used reference genes in studies of domestic fowl gene expression include the ACTB, GAPDH and 18S rRNA. Shu et al. (2004) used ACTB as a reference gene to study the expression level of three novel expressed sequence tags (ESTs) in hypothalamus, pituitary, muscle, liver and fat tissues of Shaoxing ducks (Shu et al., 2004). Scholz chose 18S rRNA as the endogenous reference to analyze the sex-dependent gene expression level in early brain development of chicken embryos by Real-time PCR (Scholz1 et al., 2006). Yen used ACTB as a reference gene for normalization of data in transcript analysis of pituitary gland genes in laying geese (Ding et al., 2007). Chen used ACTB as a reference gene to analyze relative mRNA expression levels in the hypothalamus/pituitary glands in the Red-feather Taiwan country chicken which show significantly different reproductive performance (Chen et al., 2007). Ding chose 18S rRNA as a reference gene for normalization of data in transcriptional analysis of Vitellogenin I, apoVLDL-II, ethanolamine kinase, G-protein gamma-5 subunit, and leucyl-tRNA synthase expression level in the livers of the laying and pre-laying geese (Ding et al., 2007).
Thus far, the genes encoding GAPDH, ACTB and 18S rRNA have been used as endogenous reference genes for qRT-PCR in geese (Chen et al., 2007; Ding et al., 2007; Kang et al., 2009), but the stability analysis of these candidate genes in geese has not yet been reported. Based on earlier gene expression studies in domestic fowl, we have tested the stability of expression of seven housekeeping genes, including GAPDH, ACTB, HPRT1, TUB, succinate dehydrogenase flavoprotein A (SDHA), 18S rRNA and 28S rRNA in this study. The geNorm (Vandesompe et al., 2002), Normfinder (Andersen et al., 2004), BestKeeper (Pfaffl et al., 2004), and the comparative ΔCt method are popular algorithms to determine the most stable endogenous reference genes from a set of tested candidate reference genes in a given experimental condition. Hence three different software tools were used to validate the stability of the selected housekeeping genes in different developmental stages (1 d, 2, 4, 6 and 8 months old) of Zi geese tissues (heart, liver, kidney, muscle, ovary) using real-time qRT-PCR.
MATERIALS AND METHODS
Geese and tissue collection
The study protocol was approved by the Animal Care Committee of Jilin University. Thirty female Zi geese were randomly selected from one hundred geese in a local breeding farm and raised according to the standard program used at the farm (Daqing, China). Six geese were sacrificed at the age of 1 d, 2, 4, 6 and 8 months respectively. Geese were sacrificed by electrical stunning followed by exsanguination. Heart, liver, kidney, muscle and ovary samples were rapidly removed, wrapped in foil, frozen in liquid nitrogen, and then stored at −70°C until analysis.
Total RNA isolation
Total RNA was isolated from the Zi geese tissues (hearts, livers, kidneys, muscles and ovaries) at 1 d, 2, 4, 6 and 8 months respectively according to the Trizol® reagent (Invitrogen, USA) manufacturer’s instructions. Total RNAs were treated with DNase I (RNase Free, Takara, Japan) according to the manufacturer’s instructions to remove contaminations of genomic DNA. Total RNA concentration and purity was determined using a SmartSpec™ Plus spectrophotometer (Bio-Rad, USA). The optical density (OD) ratio A260/A280 nm was measured with the spectrophotometer was 1.95±0.12 (OD A260/A280 ratio± SD).
Reverse transcription
Total RNA (1.5 μg) from Zi geese, 500 ng/μl of random hexamers primer (Promega, USA), 10 mM deoxynucleoside triphosphate (dNTP) Mix (Takara) and sterile MilliQ water (to a total volume of 12 μl) were heated to 65°C for 5 min in order to disrupt possible secondary structures and then quickly chilled on ice. Thereafter, 5×First-Strand Buffer was combined with 0.1 M dithiothreitol (DTT) and 40 units/μl of RNaseOUT™ Recombinant Ribonu-clease Inhibitor (Invitrogen). The mixture was mixed gently and incubated at 37°C for 2 min. Then a total of 200 units of M-MLV reverse transcriptase was added and incubated at 25°C for 10 min. Reverse transcription was performed at 37°C for 50 min, and the reaction mixture was heated to 70°C for 15 min. The final cDNA products were diluted 10-fold prior to use in real-time PCR.
Primer design
All primers designed for all reference genes were based on sequences published in Genbank (http://www.ncbi.nlm.nih.gov/Table 1). Primer pairs for qPCR amplification were designed using Primer premier 5.0 (http://www.premierbiosoft.com), BLAST searches were performed to confirm the total gene specificity of the primer sequences (http://www.ncbi.nlm.nih.gov/BLAST/). The specificity of all primers was checked by electrophoresis of RT-PCR products on the 1% agarose gel (Figure 1a).

Confirmation of amplicon size and primer specificity of the selected genes. (a) Agarose gel electrophoresis showing the 28S and 18S rRNA bands. (b) Agarose gel electrophoresis showing specific RT-PCR products of the expected size for each gene. M represents DNA size marker. (c) Melting curves generated for all genes.
Real-time qPCR
The qRT-PCR was performed on the first strand cDNA using the Line-Gene K Real-time PCR Detection System and software (Bioer, China) with SYBR® Premix Ex Taq™ (Takara). Briefly, each reaction (50 μl) consisted of 1 μl 10-fold diluted cDNA template, 25 μl of SYBR® Premix Ex Taq™ (2×Concentration), 0.5 μl of 20 μM of PCR Forward Primer and PCR Reverse Primer, and 23 μl of nuclease-free water. Thermal cycling was performed with an initial denaturation step of 10 s at 94°C, followed by 45 cycles of 5 s at 94°C, and 56°C for 30 s, and then a final extension at 72°C for 20 s. Finally, a dissociation curve was generated by increasing temperature starting from 65 to 95°C to determine the specificity of the reactions. The crossing cycle number (Cp) was automatically determined for each reaction by the Line-Gene K Real-Time PCR Detection System and software (Bioer) with default parameters using the second derivative method. As a control for genomic DNA contamination, an equivalent amount of total RNA without reverse transcription was tested for each reaction. A no-template control (NTC) was also included in each reaction. Relative quantitation of gene expression was performed in three replicates for each sample. The quality of standard curves was judged by the slope of the standard curve and the square of the Pearson correlation coefficient (R2). The PCR amplification efficiency of each primer pair is calculated from the slope of a standard curve using the following equation: Efficiency % = (10(−1/slope)−1)×100% (Bustin et al., 2009).
Statistical analysis
Differences in expression levels of GAPDH, ACTB, HPRT1, TUB, SDHA, 18S rRNA and 28S rRNA with developmental stage were examined by one-way ANOVA. The IBM SPSS statistics 17.0 package (IBM, USA) was used for all analyses. Significance levels were set at p<0.05. Determination of reference gene expression stability were calculated using geNorm, NormFinder and BestKeeper.
For evaluation of expression stability of the candidate reference genes, Ct values for all samples were calculated and the stability of the genes was determined utilizing three different software tools: geNorm, Normfinder and BestKeeper. The gene expression stability (M) and the optimal number of endogenous reference genes for normalization were determined by using the geNorm algorithm as previously described (Vandesompe et al., 2002). The second algorithm utilized was NormFinder as previously described by Andersen (Andersen et al., 2004). It is an algorithm that attempts to find the optimum reference genes out of a group of candidate genes. It can also, in contrast to geNorm, take information of groupings of samples into account, such as treatment/control, sick/healthy, or different developmental stages. The BestKeeper algorithm creates a pairwise correlation coefficient between each gene and the BestKeeper index (BI). This index was then compared to each individual candidate housekeeping gene by pair-wise correlation analyses, with each combination assigned a value for the Pearson correlation coefficient (r) and the probability (p) (Pfaffl et al., 2004). The gene with the highest coefficient of correlation with the BI indicates the highest stability.
RESULTS
Selection of candidate reference genes and primer design
We have investigated seven housekeeping genes commonly used as internal controls in expression studies, including the GAPDH, ACTB, HPRT1, TUB, SDHA, 18S rRNA and 28S rRNA. The primers were designed according to the Zi geese mRNA sequences which are available in GenBank (Table 1).
Quality control of the nucleic acids and qPCR
The optical density (OD) ratio A260/A280 of the RNA was 1.8 to 2.0. Agarose gel electrophoresis (Figure 1b) that 28S:18S ratio was approximately 2:1 indicated that the RNA was intact. Agarose gel electrophoresis (Figure 1a) and melting curve analysis (Figure 1c and Table 1) revealed that all primer pairs amplified a single PCR product with the expected size. Furthermore, sequence analysis of cloned amplicons revealed that all sequenced amplified fragments were identical or nearly identical to the sequences generated by the designed primers. In order to quantitatively determine the transcriptional level of each candidate gene, the average cycle threshold (Ct) value of each gene was calculated (Table 1). As expected, the average Ct value of different gene varied. A standard curve using a dilution series of the cloned amplicons was made to calculate the gene-specific PCR efficiency. The correlation coefficient (R2) of the slope of the standard curve, the PCR amplification efficiency (E) and the PCR efficiency standard deviation (SD) of each gene were listed in Table 1. All primer pairs utilized in this study presented amplification efficiency between 91 to 104% (Table 1).

Confirmation of amplicon size and primer specificity of the selected genes. (a) Agarose gel electrophoresis showing the 28S and 18S rRNA bands. (b) Agarose gel electrophoresis showing specific RT-PCR products of the expected size for each gene. M represents DNA size marker. (c) Melting curves generated for all genes.
Expression stability of candidate reference genes
Three different software tools were used to calculate the expression stability of the candidate reference genes: geNorm, NormFinder and BestKeeper.
geNorm analysis
In present study, average expression stability (M value) of all genes was calculated by geNorm (version 3.5). The M values of the candidate reference genes across Zi geese tissues (heart, liver, kidney, muscle,oaries) are shown in Table 2. GAPDH and 28S had the highest expression stability in Zi geese heart tissues, (the lowest M values). GAPDH and SDH had the highest expression stability in liver tissues, ACTB and SDH had the highest expression stability in kidney tissues, HPRT1 and 18S had the highest expression stability in Zi geese muscle tissues, GAPDH and HPRT1 had the highest expression stability in ovaries tissues. TUB (heart, kidney and muscle), 18S (ovarie) and ACTB (liver) had the highest M values, indicating less stable expression across Zi geese tissues (Table 2). In different developmental stages of Zi geese heart tissue, low V4/5 value = 0.210, the 4 member set GAPDH, 28S, SDH and ACTB is an excellent choice for the calculation of the NF. In liver tissue, low V3/4 value = 0.237, the 3 member set GAPDH, HPRT1 and SDH is an excellent choice for the calculation of the NF. In kidney tissues, the inclusion of a 4th gene has no significant effect (low V3/4 value = 0.148<1.5) on the NF. The 3 member set ATC, GAPDH and SDH is an excellent choice for the calculation of the NF. In muscle tissue, low V2/3 value = 0.124<1.5, the 2 member set HPRT1 and 18S is an excellent choice for the calculation of the NF. In oaries tissue, low V2/3 value = 0.152, the 2 member set GAPDH and HPRT1 is an excellent choice for the calculation of the NF.
NormFinder analysis
Table 2 shows the ranking order of the seven candidate reference genes mentioned above, using the NormFinder program to calculate their expression stability. Genes that are more stably expressed are indicated by lower average expression stability values. The analysis ranks 28S, SDH, GAPDH and HPRT1 were the four most stable genes in different developmental stages of Zi geese heart tissues (Table 2). SDH, GAPDH and HPRT1 were the three most stable genes in liver tissues. SDH, ACTB and GAPDH were the three most stable genes in kidney tissues. GAPDH and 28S were the two most stable genes in muscle tissues. GAPDH and HPTR1 were the two most stable genes in oaries tissues. Thus, both geNorm and NormFinder rank the same genes as the most stable and the entire order is identical.
BestKeeper analysis
The results of reference gene evaluation by the BestKeeper tool are shown in Table 2. The BestKeeper revealed that in different developmental stages of Zi-geese heart tissues the best correlations were obtained for SDH (r = 0.903), GAPDH (r = 0.881), 28S (r = 0.857) and HPRT1 (r = 0.846) with p value of 0.001 (Table 2). TUB are ranked as the least stable genes. In liver tissues the best correlations were obtained for GAPDH (r = 0.951), 28S (r = 0.916) and 18S (r = 0.909) with p value of 0.001. TUB are ranked as the least stable genes. In kidney tissues the best correlations were obtained for ACTB (r = 0.975), SDH (r = 0.943) and GAPDH (r = 0.936) with p value of 0.001. TUB are ranked as the least stable genes. In muscle tissues the best correlations were obtained for HPRT1 (r = 0.980) and 18S (r = 0.979) with p value of 0.001. TUB are ranked as the least stable genes. In ovaries tissues the best correlations were obtained for GAPDH (r = 0.914) and SDH (r = 0.876) with p value of 0.001. 28S are ranked as the least stable genes.
DISCUSSION
Real-time qRT-PCR, with its capacity to detect and measure very small amount of nucleic acids in a wide range of samples from numerous sources, has been used extensively in molecular biology, e.g., molecular diagnostics, life sciences, agriculture, and medicine. Normalization is required to reduce the tube-to-tube variations caused by variable RNA quality or reverse transcription efficiency, inaccurate quantification, and pipetting etc (O’Connell, 2002). Endogenous reference genes are thus commonly used to normalize the expression levels of analyzed genes. ACTB, together with GAPDH, TUB, EF1A and 18S rRNA, are expressed constitutively and are involved in basic housekeeping functions required for cell maintenance. Because of this, they are commonly used as reference genes to normalize gene expression studies (Sturzenbaum and Kille, 2001; Jin et al., 2004; Wen and Mao, 2007; Øvergard et al., 2010). Recently, a growing body of research has demonstrated that these genes expression can change in different tissues, during growth and differentiation, in response to biochemical stimuli, and in disease states (Wen and Mao, 2007; Janovick-Guretzky et al., 2007). The expression levels of the ideal endogenous reference genes should be constant in different experimental conditions. The geNorm, Normfinder, BestKeeper, and the comparative ΔCt method are popular algorithms to determine the most stable endogenous reference genes from a set of tested candidate reference genes in a given experimental condition.
In this study, we have selected seven candidate reference genes from Genbank to analyze their candidacy to be used as reference genes. Also, we have developed a qRT-PCR method for GAPDH, HPRT1, ACTB, 18S, 28S, SDHA and TUB as the target gene. The specificity of the qRT-PCR primer pairs was confirmed by agarose gel electrophoresis, Tm analysis and sequencing of the amplicons. The PCR amplification efficiency was estimated, and the reference genes were ranked according to their expression level stability across various developmental stages in Zi geese tissues using geNorm, Normfinder and BestKeeper algorithms. When gene expression stability in Zi geese was analyzed by geNorm, which had been recently noted as one of the best methods to determine the most stably expressed genes for qRT-PCR analysis (10, 8), the most stable genes in the seven series were different as shown in Table 2 and Figure 2. In different developmental stages of Zi geese heart tissues, the stability rank was GAPDH, 28S>SDH>ACTB >18S>HPRT1>TUB, and the optimal number of reference genes was four. In order to avoid co-regulation, we have also determined the stability of the selected genes using Normfinder and BestKeeper. The 28S, SDH, GAPDH and HPRT1, as identified by both the NormFinder and the BestKeeper tools, were the four most stable genes, which supported the geNorm analysis in this experiment. The stability of the four genes, from the highest to the lowest, was 28S, SDH, GAPDH and HPRT1. In liver tissues, the stability ranking was GAPDH, HPRT1>SDH>TUB>28S>18S>ACTB, and optimal number of reference genes was three. In kidney tissues, the stability ranking was ACTB, SDH>GAPDH>HPRT1>28S>18S>TUB, and the optimal number of reference genes was three. In conclusion, the three algorithms did not rank the candidate reference genes in the same order, but all indicated that ACTB, SDH and GAPDH should be the most stably expressed genes in the experimental conditions applied in this study. In muscle, the stability ranking was HPRT1, 18S>GAPDH>28S>SDH>ACTB>TUB, and the optimal number of reference genes was two. In ovary, the stability ranking was GAPDH, HPRT1>SDH>TUB>18S>ACTB>28S, and the optimal number of reference genes was two. The NormFinder identified that GAPDH and HPRT1 to be the two most stable genes whereas, in contrast, the BestKeeper identified GAPDH and SDH to be the two most stable genes.

Determination of the optimal number of reference genes calculated by geNorm. The reference genes were ranked according to their expression levels stability across various developmental stages in Zi geese tissues using geNorm. When gene expression stability in Zi geese was analyzed by geNorm, the stability were different in heart, liver, kideny, muscle and ovary.
IMPLICATIONS
This research is the first attempt to validate a set of commonly used candidate reference genes in various developmental stages in Zi geese tissues for the normalization of gene expression using qRT-PCR. Analysis of stability using geNorm, NormFinder and BestKeeper reveals that the geometric mean of GAPDH, 28S, SDH and HPRT1 (in heart), GAPDH, HPRT1 and SDH (in liver), ACTB, SDH and GAPDH (in kidney), HPRT1 and 18S (in muscle), and GAPDH and HPRT1 (in ovary) are recommended to be used as reference genes in Zi geese. These methods may further be employed to identify the most stable reference genes in other tissues or under other experimental conditions in the future studies on geese. Also, this study may serve as a good foundation for further studies on how to improve the economic traits of geese, both in egg, fatty livers and meat production.
Acknowledgements
The work was supported by the Natural Science Foundation of Heilongjiang Province of China (ZD201116). The authors would like to thank Mrs. Yibing Han, for her valuable time in providing language help. We also thank the anonymous reviewers for their valuable comments and suggestions.