Identification of reliable reference genes for quantitative real-time PCR normalization in pitaya

Background A suitable reference gene is an important prerequisite for guarantying accurate and reliable results in quantitative real-time PCR (qRT-PCR) analyses. However, there is no absolute universality in reference genes among different species. It’s hard to find an ideal reference gene to fit for different tissues and growth periods. Pitaya (Hylocereus) is commercially produced as a new fruit crop at a large scale in tropical and subtropical regions. To date, there is no report on the identification of the most reliable reference genes for qRT-PCR normalization in pitaya. Results In this study, six candidate reference genes i.e. Actin(1), GAPDH, UBC(1), UBC(2) EF1-α(1) and histone(1) were selected from thirty-nine typical candidate reference genes to determine the most stable reference genes for qRT-PCR normalization in different tissues, temperature stresses and fruit developmental stages of pitaya. Among the six candidate reference genes, Actin(1) and EF1-α(1) were the most stable gene according to calculations of three statistical methods (GeNorm, NormFinder and BestKeeper) while UBC(1) and UBC(2) showed the lowest expression stability. The six candidate reference genes were further validated by comparing expression profiles of key genes related to betalain biosynthesis at flesh coloration stages of Guanhuahong (Hylocereus monacanthus) and Guanhuabai (H. undatus) pitayas. Actin(1) was recommended the best reference gene for accurate normalization of qRT-PCR data. Conclusions In this study, the stability of the selected reference genes for normalizing the qRT-PCR data were identified from pitaya. Actin(1) was the most stably expressed genes in different tissues and fruit developmental stages in pitaya. The present work provides the first data of reference gene identification for pitaya and will facilitate further studies in molecular biology and gene function on Hylocereus and other closely related species. Electronic supplementary material The online version of this article (10.1186/s13007-019-0455-3) contains supplementary material, which is available to authorized users.


Background
Quantitative real-time PCR (qRT-PCR) is a nucleic acid quantitative technology to study gene expression level in molecular biology research [1]. It combines the sensitivity of conventional PCR with a cost-effective assay using a specific fluorescent signal. Compared with northern blot and conventional RT-PCR, qRT-PCR is proven to be highly sensitive and specific and has been increasingly used in various research fields of biological sciences [2][3][4]. A stable reference gene plays an important role in analyzing the relative expression levels of target genes [5]. Reference gene, i.e. internal reference gene, whose expression is relatively constant in various tissues and cells, can be used as the reference in the detection of gene relative expression analyses. However, transcription levels of reference gene may change in different species, treatments and/or developmental stages. Therefore, it is necessary to validate the expression stability of reference The accuracy of gene expression levels using qRT-PCR depends on the use of stable reference genes to normalize the difference between samples [6]. GAPDH, Actin, 18S rRNA, 28S rRNA, tubulin, EF1-α and ubiquitin are commonly used as reference genes for comparison of samples [7][8][9]. Reference genes are always expressed because they produce proteins that cells must have to function, however, the stability of reference genes is limited in different tissues, growth and development stages of the same or different plants [10][11][12][13][14][15]. Therefore, it is necessary to select the most suitable reference gene for accurate qRT-PCR evaluation according to different species and experimental conditions [16,17]. In addition, calculations of reference genes by geNorm [6], NormFinder [18] and BestKeeper [19] is essential for normalization of qRT-PCR analyses.

Plant Methods
Pitaya or pitahaya, also known as dragon fruit, originated from rainforest in tropical and subtropical areas of Central America and Mexico. Pitaya belongs to the genus Hylocereus in the Cactaceae family. Currently, redflesh and white-flesh pitayas have been commercially produced as a new fruit crop at a large scale in Central America, Southeast Asia and China. Pitaya is a fast growing, perennial, hemi-epiphytic, vine-like, crassulacean acid metabolism (CAM) plant. It is a spiny succulent plant which can adapt to a wide ecological range such environmental cues as drought, heat and poor soil. Pitaya has gained great attention of the world due to its high nutritional value as well as its antioxidant capacity and antiproliferative activities [20][21][22]. Therefore, the knowledge about the molecular basis of pitaya may contribute to the discovery of new and promising genes related to important agricultural traits that could be further transferred to a target crop for the development of new crop cultivars.
Gene expression analysis is an important tool to elucidate the complex regulatory networks of the genetic, signalling and metabolic pathway mechanisms during plant life cycle. qRT-PCR is an ideal tool to verify the results of differential gene expression of interest on a smaller scale. To study mRNA expression levels, a set of appropriate reference genes are fundamental to get reliable results from qRT-PCR. In this study, thirty-nine genes annotated as 18S rRNA (18S ribosomal RNA), CYP (cytochrome), Actin, GAPDH (glyceraldehyde-3-phosphate dehydrogenase), EF1-α (eukaryotic elongation factor 1-alpha), histone, TUA (alpha-tubulin), TUB (beta-tubulin), UBQ (ubiquitin) and UBC (ubiquitin-conjugating enzyme) were selected as candidate reference genes to evaluate their expression stability in different pitaya tissues, temperature stresses and fruit developmental stages using three available statistical algorithms i.e. geNorm [6], Nor-mFinder [18] and BestKeeper [19]. The aim of the present study is to identify appropriate reference genes for qRT-PCR normalization in pitaya and facilitate future molecular studies in Hylocereus. In addition, we validated the expression levels of the key genes related to betalain biosynthesis using the selected reference genes.

Plant materials
Two pitaya cultivars i.e. Guanhuahong (H. monacanthus) and Guanhuabai (H. undatus) pitayas from the same orchard of Dalingshan Forest Park were used as plant materials. Guanhuahong pitaya is a red flesh fruit with red peel while Guanhuabai pitaya is a white flesh fruit with red peel (Additional file 1: Fig. S1). Root, stem, receptacle, calyx, petal, filament, anther, ovary, style, stigma, peel and flesh were collected separately and used for expression analyses of candidate reference genes. Pitaya flesh were collected on the 13rd, 16th, 19th, 23rd, 25th, 27th and 29th days after flowering (DAF) in July and August, 2015 (Additional file 1: Fig. S1). Cutting seedlings were collected for expression analyses of candidate reference genes on the 0, 9th, 18th, 36th, 72th and 144th h after cultivating in different temperature stresses (4 °C, 10 °C and 25 °C). All samples were immediately frozen in liquid nitrogen and stored at − 80 °C until use.

Candidate gene selection and primer design
Thirty-nine candidate reference genes with relatively stable expression levels (based on reported RPKM and fold change values) i.e.  Table S1) [23]. The specific primers (Table 1) were designed using Primer Premier 5.0 with the following parameters: melting temperature (Tm) values ranging from 47 to 63 °C, GC percent of 25-65%, primer lengths of 18-22 bp and product length of 100-300 bp. Tenfold serial dilutions of cDNA were used to determine slope of the standard curve to calculate amplification efficiency of primer (E = 10 (−1/slope of the standard curve) ).

RNA extraction and cDNA synthesis
Total RNA was extracted according to the protocol of Wu et al. [24] and was purified by RNase-free columns (Huayueyang, Beijing, China). Genomic DNA was removed by DNase I (TaKaRa, Dalian, China). RNA quality was analyzed by 1.0% agarose gel and RNA integrity was assessed using an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). Only RNA samples with an A260/A280 ratio of 1.8-2.2 and an A260/A230 ratio > 1.8 were used for further analysis. cDNA was synthesized from 2.0 µg of total RNA samples using M-MLV first strand cDNA synthesis kit according to the manufacturer's instructions (TaKaRa).

Analysis of gene expression stability
NormFinder, geNorm and BestKeeper were used to evaluate expression stability of reference genes in various pitaya tissues, temperature stresses and fruit developmental stages. The geNorm program determines stability and the optimal number of genes required to calculate the M value and pairwise variation Vn/Vn + 1 between two sequential normalization factors [6]. The Nor-mFinder program was performed to obtain the best candidate reference genes by the average expression stability between samples based on the results of variance analysis [18]. The BestKeeper program is used for determination of stable reference genes, which calculates pairwise correlations based on standard deviation (SD) and percentage covariance (Cov) values [19].

Validation of reference genes
To validate the reliability of the selected reference genes, the relative expression levels of a key gene Cytochrome P450 (Cyt P450-like1) involved in betalain biosynthesis were analyzed using qRT-PCR at different flesh coloration stages of Guanhuahong and Guanhuabai pitayas. The specific primers were designed according to Hua et al. [23] (Table 1).

Statistical analyses
Statistical analyses were performed using SPSS 17.0 software.

Selection of candidate reference genes and sequence analyses
Thirty-nine reference genes from the transcriptome of pitaya were selected as putative candidate reference genes (Additional file 2: Table S1).  (3) showed 92% sequence identity with Sghistone (KC415067.1) (Additional file 3: Table S2).

Verification of primer specificity and stability
Melting curve plays an important role in assessing the effects of primer specificity. As shown in Additional file 4: Fig. S2, there were only one distinctive peak for Actin (1) (6), TUB (1) and TUB (2). Therefore, those thirty-four candidate reference genes were selected from the thirty-nine candidate reference genes according to the melting curves (Additional file 4: Fig. S2). The efficiency for each primer pair varied from 90.1% for 18S rRNA to 110.4% for EF1-α(5) ( Table 1). The sequence length ranged from 102 to 300 bp (Table 1). Cycle threshold (Ct) values ranged from 17.45 for GAPDH to 37 for histone (2). The boxplot provided expression levels of the 39 candidate reference genes (Additional file 5: Fig. S3).  (1), UBC (2) and histone(1) were selected as candidate reference genes for further analyses.
The coefficient of variance (CV) of Ct values can determine reference gene stability levels. UBC(2) presented a high degree of transcript level variation with a CV value of 15.61% while Actin(1) presented the lowest degree of variation with a CV value of 4.96% (Additional file 7: Table S3). Among the six candidate reference genes, Actin(1) (4.96%) was the most stable gene while UBC(1) (13.16%) and UBC(2) (15.61%) showed the lowest expression stability according to CV value analyses of six reference genes.

GeNorm analysis
GeNorm was used to select the most stable reference gene by calculating the gene expression stability measure (M) based on the average pairwise expression ratio. The gene with the lowest M value is considered the most stable. A gene with M value below 1.5 can be considered as a reference gene. As shown in Fig. 2, the M values of Actin(1), EF1-α(1), GAPDH, UBC(1), UBC (2) and histone(1) were lower than 1.5, suggesting that all of them conformed to basic requirements for the reference gene. Among the six genes, UBC (2) was the least stable reference gene with the highest M value (> 0.214) detected in the two different colored pitayas, different tissues and fruit developmental stages of Guanhuahong and Guanhuabai pitayas. EF1-α (1) and Actin (1)  A single reference gene often can not meet the requirements of accurate quantification for transcription analyses of gene expression. It is necessary to use two or more than two reference genes to produce accurate and reliable normalization. The optimal number of reference genes can be determined using the pairwise variation value (Vn/Vn + 1, V-value) as normalization factor for reference genes. Additional (n + 1) reference genes are necessary to normalize the genes if V-value is higher than the threshold of 0.15. As shown in Fig. 3, the pairwise variation values of reference genes were examined when all samples were considered. The ratio of V2/V3 was lower than the threshold of 0.15, suggesting that the optimum number of reference gene was two. There is no need to add a third reference gene as an internal control for the normalization of gene expression in pitaya. Therefore, a combination of Actin(1) and EF1-α(1) was the optimum Fig. 2 Average expression stability values (M) of the six candidate reference genes using geNorm algorithm. Total represents a comparison between all tissues/fruits between the two cultivars reference genes for qRT-PCR normalization of various samples in pitaya.

NormFinder analysis
Different from geNorm, the NormFinder can determinate optimal number of reference genes for normalization by intra-and inter-group variations. Average gene expression stability values of the six reference genes in different samples were calculated by NormFinder software and ranking was made based on relative stability values. As shown in Fig. 4, UBC(1) was the least stable reference gene in the total dataset, which is different from the result obtained by NormFinder. In different fruit developmental stages, GAPDH and histone(1) emerged as the most stably expressed gene for Guanhuabai pitaya compared with UBC(1) and GAPDH for Guanhuahong pitayas. In various tissues of Guanhuabai and Guanhuahong pitayas, EF1-α(1) was the most stable reference gene. UBC (1) and Actin(1) were the best reference gene for different temperature stresses. Combined with expression analyses of all samples, histone (1) and GAPDH were the most stable reference genes for normalization in pitaya.

BestKeeper analysis
BestKeeper program was used to grade candidate reference gene stability by calculating SD of the Ct values, CV and pair correlation coefficient. A gene with SD value below 1.0 can be considered as stable expression. In different temperature stresses, SD values were 0.38 for UBC (2), 0.42 for UBC(1), 0.44 for EF1-α(1), 0.47 for Actin(1), 0.59 for GAPDH and 1.08 for histone(1), respectively. As for the different fruit developmental stages of Guanhuahong pitaya, SD values of Actin(1) and EF1-α(1) were 0.76 and 0.97, respectively, suggesting Actin(1) and EF1-α(1) could be used as reference genes. Actin(1) was the most stably expressed gene, followed by EF1-α(1) (Fig. 5). Compared with Guanhuabai pitaya, Actin (1) was the most stable gene for normalization of the different fruit developmental stages of Guanhuahong pitaya. In different tissues of Guanhuabai and Guanhuahong Fig. 3 Pairwise variation (V) analyses of the six candidate reference genes using geNorm algorithm. Total represents a comparison between all tissues/fruits between the two cultivars pitayas, EF1-α(1) (0.97) and Actin(1) (0.9) was the most stable gene, respectively; UBC(2) (2.58) was the least stable reference gene with the highest SD value. In the different colored pitayas, Actin(1) (0.97) was considered to be the most stable gene. Taking into account the total dataset, Actin(1) was the most stable gene in all samples tested. The results of BestKeeper were consistent with the analyses of geNorm.

Reference genes validation
To validate the suitability of the two selected reference genes (Actin(1) and EF1-α(1)), and the least suitably ranked candidates UBC(2) according to BestKeeper and geNorm analyses. Relative expression levels of a key gene Cyt P450-like1 involved in betalain biosynthesis were evaluated using qRT-PCR at different flesh coloration stages of Guanhuahong (H. polyrhizus) and Guanhuabai (H. undatus) pitayas (Fig. 6). The expression levels of Cyt P450-like1 increased gradually during flesh color transition from white (23 days) to red stages (25 days) and decreased at full maturation stage (29 days) when single or a combination of reference genes (Actin(1) and EF1-α(1)) were used for normalization. These results were in accordance with betalain accumulation pattern and expression characteristics of Cyt P450-like1 gene related to betalain biosynthesis at all flesh coloration stages of H. monacanthus [23,26]. However, the expression pattern of Cyt P450-like1 normalized by UBC(2) increased gradually at all flesh coloration stages of Guanhuahong pitaya. Furthermore, when UBC(2) was used to normalize Cyt P450-like1 expression in Guanhuabai pitaya, the expression values of 23 days were overestimated. In general, the RT-qPCR profiling of Cyt P450-like1 gene expression supported the analyses of BestKeeper and geNorm .

Discussion
qRT-PCR is the most sensitive technology of measuring gene transcript levels depending on stable reference gene(s) for normalization [27]. An appropriate reference gene can effectively correct the errors of RNA quantity, reverse transcription efficiency, stability of different reaction channels of PCR instrument and standard operation mode, which can help to obtain the real differential expression of the target genes [28]. The reference genes Fig. 4 Expression stability values and ranking of the six candidate reference genes using NormFinder algorithm. Total represents a comparison between all tissues/fruits between the two cultivars such as Actin, 18S rRNA and GAPDH are commonly used as reference genes. However, these reference genes show instability in various plant species or genotypes [13,14,29]. Their expression levels are different in different patterns of cells and tissues, physiological status and growth periods [30,31]. Therefore, it is essential to select the most suitable reference gene for the specific sample types and experimental conditions. Pitaya has become a popular fruit due to its excellent nutritional, commercial and medical values [20][21][22]. To date, there is no report on the identification of the most suitable reference gene in pitaya. In this study, 39 typical reference genes from pitaya transcriptome database [23] were selected to determine the most stable reference gene for qRT-PCR normalization in pitaya. Based on the melting curve and Ct values of reference genes (Additional file 4: Fig. S2 and Additional file 5: S3), six genes, i.e. Actin(1), EF1-α(1), GAPDH, UBC(1), UBC (2) and histone(1) were selected as candidate reference genes for further analyses.
NormFinder [18], geNorm [6] and BestKeeper [19] were commonly used to evaluate the stability of the reference genes. The geNorm was applied to estimate the stability of a candidate gene by pairwise comparison, while the NormFinder and BestKeeper were used to prevent co-regulation and to further assess the results obtained by the geNorm program. In the present study, geNorm, NormFinder and BestKeeper were applied to define the stability of Actin(1), EF1-α(1), GAPDH, UBC(1), UBC (2) and histone(1) reference genes in different tissues, temperature stresses and fruit developmental stages of pitaya. Different rankings for the six candidate reference genes were detected after comparison to the ranking of the candidates generated by the three algorithms (Fig. 2, 3, 4, 5). From the geNorm evaluation, we found that Actin(1) and EF1-α(1) were the optimum reference genes for qRT-PCR normalization of various samples of pitaya. The results of BestKeeper were consistent with the analyses of geNorm, showed Expression stability values and ranking of the six candidate reference genes using BestKeeper algorithm. Total represents a comparison between all tissues/fruits between the two cultivars that Actin(1) was the most stable reference gene in all of the samples examined. However, results from Nor-mFinder analyses showed that histone (1) and GAPDH were the most stable reference genes in all of the samples tested. The reason for the differences in the topranked reference genes may be due to the discrepancies resulted from different approaches of calculations by geNorm, NormFinder and BestKeeper programs. Similar results were also obtained in Vitis vinifera [32], Citrus [33], Euscaphis konishii [34], Cowpea [35] and Actinidia deliciosa [36]. It has been reported that the most discrepant results in gene stability ranking were obtained by BestKeeper [37]. In our study, although there are some differences in the stability of each reference gene expression of samples from the three different algorithms, Actin(1) and EF1-α(1) are the optimum reference genes. Combining the three algorithms with expression analyses, we suggested Actin(1) was the most adequate reference gene for expression studies in different tissues, temperature stresses and fruit developmental stages of pitaya.
Cytochrome P450 is a key gene involved in betalain biosynthesis [23,[38][39][40][41][42]. To validate the candidate reference genes, the relative expression levels of Cyt P450-like1 gene were analyzed by qRT-PCR during flesh coloration of Guanhuabai and Guanhuahong pitayas. Expression levels of Cyt P450-like1 increased gradually during flesh coloration and decreased at full maturation stage when single or a combination of reference genes (Actin(1) and EF1-α(1)) were used as reference genes. Those results were in consistent with betalain accumulation pattern and expression characteristics of Cyt P450-like1 at all flesh coloration stages of H. monacanthus. Moreover, the expression pattern of Cyt P450-like1 increased with fruit maturation of H. undatus [23,24]. Those results suggested that Actin(1) was the most stable reference gene for expression studies of key genes involved in betalain biosynthesis of pitaya. And the most stable reference gene combination for expression studies in pitaya was Actin(1) and EF1-α(1).

Conclusions
To the best of our knowledge, this study is the first report on systematically evaluating the expression stability of different potential reference genes for qRT-PCR in the Cactaceae family. Thirty-nine typical reference genes were selected to determine the most stable reference genes for qRT-PCR normalization in pitaya. Based on the melting curve and Ct values, Actin(1), EF1-α(1), GAPDH, UBC(1), UBC(2) and histone(1) genes were selected for further analyses in different tissues, temperature stresses and fruit developmental stages of pitayas. Actin(1) and EF1-α(1) were identified as the optimum internal control genes according to calculations made with geNorm, NormFinder and BestKeeper programs. Validation of suitable reference genes was carried out to profile the expression of Cyt P450-like1 gene during flesh coloration stages of Guanhuahong and Guanhuabai. Actin(1) was the best reference gene for qRT-PCR normalization. The present study provides the appropriate reference gene for normalization of reliable qRT-PCR data in different pitaya tissues, temperature stresses and fruit developmental stages which will be useful for the expression profiles of target genes in Hylocereus plant and related species.  (1) and Actin(1) + EF1-α(1) were used as one or two most stable reference genes, UBC(2) was used as the least stable reference gene. Different letters indicate significant difference of the expression of the target gene based on three biological replications (P < 0.05, t test; n = 3)