Reference gene selection for qRT-PCR analyses of luffa (Luffa cylindrica) plants under abiotic stress conditions

Selecting suitable internal reference genes is an important prerequisite for the application of quantitative real-time PCR (qRT-PCR). However, no systematic studies have been conducted on reference genes in luffa. In this study, seven reference genes were selected, and their expression levels in luffa plants exposed to various simulated abiotic stresses [i.e., cold, drought, heat, salt, H2O2, and abscisic acid (ABA) treatments] were analyzed by qRT-PCR. The stability of the reference gene expression levels was validated using the geNorm, NormFinder, BestKeeper, and RefFinder algorithms. The results indicated that EF-1α was the most stably expressed and suitable reference gene overall and for the heat, cold, and ABA treatments. Additionally, UBQ expression was stable following the salt treatment, whereas TUB was identified as a suitable reference gene for H2O2 and drought treatments. The reliability of the selected reference genes was verified by analyzing the expression of copper/zinc superoxide dismutase (Cu/Zn-SOD) gene in luffa. When the most unstable reference genes were used for data normalizations, the resulting expression patterns had obvious biases when compared with the expression patterns for the most ideal reference genes used alone or combined. These results will be conducive to more accurate quantification of gene expression levels in luffa.

Primer specificity and amplification efficiency analysis. Seven genes (ACT , TUA , TUB, EF-1α, GAPDH, UBQ, and 18S) were selected as candidate reference genes. Information regarding these genes and their qRT-PCR primer pairs are summarized in Table 1. Gel electrophoresis and melting curve analyses were used to determine primer specificity. All primers amplified a single amplicon of the expected size (Fig. 2). The seven candidate reference genes under various abiotic stress conditions produced a single peak during the melting curve analysis (Fig. 3). These results indicated the primers for these genes were highly specific. The amplification efficiencies for the seven candidate reference genes ranged from 97.6% (TUB) to 104.2% (18S), and the correlation coefficients varied from 0.988 (ACT ) to 0.997 (18S) ( Table 1). Therefore, all primers were appropriate for qRT-PCR analyses.   Polyubiquitin  KR349345   5′-TGC TTC GTC TCA GGG  GTG G-3′  5′-GTC CTG AAT TTT AGC  TTT CAC-3     www.nature.com/scientificreports/ For optimal data normalization, two or more reference genes are necessary for qRT-PCR analyses. The optimal number of reference genes can be determined with the geNorm algorithm, which calculates the pairwise variation V n /V n+1 . If V n /V n+1 is less than 0.15, n is the most suitable number of internal reference genes. However, if V n /V n+1 is greater than 0.15, n + 1 is the ideal number of internal reference genes. For the cold treatment, V 2 / V 3 was less than 0.15 ( Fig. 6), implying that two reference genes are sufficient for normalizing gene expression data. Regarding the ABA treatment, V 4 /V 5 was less than 0.15, indicating that four reference genes should be used for the data normalization. For the other four treatments, V n /V n+1 exceeded 0.15, and the optimal number of reference genes was not determined.
NormFinder analysis. NormFinder calculates the SV for each reference gene based on the intragroup and intergroup variations. Lower SVs correspond to higher gene expression stability. For the heat and ABA treatments as well as overall (i.e., all treatments), EF-1α and GAPDH were the most and least stably expressed reference genes, respectively. For the cold treatment, EF-1α was the most stably expressed reference gene, whereas TUA was the least stably expressed reference gene. Regarding the H 2 O 2 and drought treatments, TUB and GAPDH were the most and least stably expressed reference genes, respectively. In response to the salt treatment, the most and least stably expressed genes were UBQ and 18S, respectively (Table 2).
BestKeeper analysis. BestKeeper calculates the stability of candidate reference genes based on the SD, r, and CV of the Ct data for all reference genes. Low SD and CV values reflect stable gene expression. The rankings based on the BestKeeper analysis revealed that four reference genes were stably expressed following the heat treatment, of which ACT (CV ± SD = 0.59 ± 0.149) was the most stably expressed. For the cold treatment, all seven ref-  www.nature.com/scientificreports/ erence genes were stably expressed, but ACT (CV ± SD = 1.16 ± 0.298) and EF-1α (CV ± SD = 1.63 ± 0.342) were the most stably expressed. Regarding the salt treatment, the expression of four reference genes was stable, with UBQ (CV ± SD = 2.98 ± 0.56) identified as the most stably expressed gene. For the H 2 O 2 treatment, four reference genes were suitable for normalizing gene expression data, among which TUB (CV ± SD = 2.13 ± 0.494) and EF-1α (CV ± SD = 2.7 ± 0.513) were the most stably expressed. In response to the ABA treatment, five reference genes were stably expressed, with TUA (CV ± SD = 1.9 ± 0.44) and ACT (CV ± SD = 1.87 ± 0.484) detected as the most stably expressed genes. Following the drought treatment, the expression levels of five reference genes were considerably stable, with TUB (CV ± SD = 2.25 ± 0.553) and TUA (CV ± SD = 2.46 ± 0.554) revealed as the most stably expressed genes. An analysis of all treatments indicated that four reference genes were stably expressed, with ACT (CV ± SD = 3.33 ± 0.875) identified as the gene with the most stable expression level (Table 3).
RefFinder analysis. Finally, RefFinder was used for the comprehensive ranking of the reference genes under each stress condition. The results indicated that EF-1α was the most stably expressed gene overall and for the heat, cold, and ABA treatments. The UBQ gene had the most stable expression level following the salt treatment, whereas TUB was the most stably expressed gene in response to the H 2 O 2 and drought treatments. In contrast, GAPDH was the least stably expressed gene overall and for the heat, salt, H 2 O 2 , ABA, and drought treatments. Regarding the cold treatment, TUA was the least stably expressed gene (Table 4).

Reference gene validation
To validate the accuracy of the analysis of reference gene expression stability, the relative expression levels of luffa LcCu/Zn-SOD1 gene were analyzed under various stresses conditions. The two most stable reference genes and one unstable reference gene according to RefFinder were selected for normalizing gene expression data. The results indicated that when the most ideal reference genes were used alone or combined as the internal reference control, the expression of the LcCu/Zn-SOD1 exhibited similar trends with minor changes. Moreover, the  www.nature.com/scientificreports/ use of two reference genes generally improved the quantification of LcCu/Zn-SOD1 gene expression. However, the LcCu/Zn-SOD1 gene expression patterns were considerably different when the unstable reference gene was used to normalize data. Specifically, the expression levels were significantly higher at 9 or 12 h under heat, salt, H 2 O 2 , ABA and drought treatment conditions, reflecting the overestimation of gene expression levels at these two time-points (Fig. 7). Therefore, the use of inappropriate reference genes may lead to biased and incorrect expression analyses.

Discussion
The qRT-PCR assay, which has been widely used for quantifying target gene transcription levels, is an important research tool applicable for characterizing gene functions 17,18 . However, the accuracy of gene expression analyses mainly depends on the selection of an appropriate internal control, which is often referred to as a reference or housekeeping gene. Reference genes were initially selected primarily based on housekeeping gene functions. For example, ACT and TUB genes encode the basic components of the cytoskeleton, whereas the proteins encoded by GAPDH, EF-1α, and UBQ genes contribute to the basic metabolic processes of organisms. Accordingly, these genes were expected to be stably expressed in all cells and physiological states 2 . Subsequent research proved that reference gene expression levels may vary in different tissues or growth and development stages. The expression might also be influenced by biotic and abiotic stresses and hormones. Thus, there are no internal reference genes that are stably expressed under all experimental conditions [19][20][21][22] . For example, GAPDH expression is highly stable in grape berries, but is relatively unstable in wheat 23,24 . Both ACT and UBQ are stably expressed in wheat, but not in tomato 24,25 . Gong et al. 26 analyzed the expression stability of 18 candidate reference genes of Goji at   www.nature.com/scientificreports/ different developmental stages and under drought stress condition. They observed that two novel reference genes LbCML38 and LbRH52 were more stably expressed than the commonly used reference genes. Liu et al. 27 analyzed perennial ryegrass under high-temperature stress and identified HIS3 and eIF4A as the most suitable reference genes. In another study, EF-1α and UBQ were detected as the most stable internal reference genes in potato and Arabis alpina under drought conditions 28,29 . Therefore, verifying the expression stability of reference genes under different experimental conditions is critical when selecting reference genes to standardize gene expression levels.
However, there has been relatively little research regarding the systematic validation of reference genes in luffa.
Only 18S has been reported as a reference gene in luffa, but it has often been used without any proper verification, possibly resulting in inaccurate analyses of luffa gene expression patterns.
In this study, the expression stability of seven candidate reference genes (ACT , TUA , TUB, EF-1α, GAPDH, UBQ, and 18S) after heat, cold, salt, H 2 O 2 , ABA, and drought treatments was analyzed with BestKeeper, geNorm, NormFinder, and RefFinder. The objective of this study was to identify the most stable reference genes for investigating gene expression in luffa plants exposed to abiotic stress factors. The data revealed EF-1α as the most stably expressed gene following the heat, cold, and ABA treatments. This gene, which encodes a polyribosomal protein, is widely expressed in various cells. In many studies on Cucurbitaceae species, EF-1α has been used as an appropriate reference gene. In cucumber, EF-1α, Fbox, CAC , and TIP41 are reportedly stably expressed regardless of the abiotic stress, growth regulator treatment, and nitrogen status 30,31 . In another study on cucumber, EF-1α was confirmed as the most stably expressed reference gene 32 . In an earlier study on zucchini, UFP, EF-1α, RPL36aA, PP2A, and CAC were identified as the ideal reference genes for normalizing expression data 33 . Ye et al. 34 screened and evaluated reference genes for fluorescence-based real-time quantitative PCR analyses of wax gourd and confirmed that EF-1α is the most stably expressed reference gene under cold conditions and in different tissues. In the current study, UBQ was the most stably expressed reference gene following the salt treatment, which is consistent with the results of an earlier study by Jin et al. 35 . The UBQ gene belongs to the ubiquitin gene family, which affects many important biological processes, including cell cycle regulation, growth and apoptosis, signal transduction, and immune responses. Because of its high sequence homology and the fact it is highly conserved, UBQ has been used as a reference gene for rice 36 , black fungus 37 , and cucumber 32 . Our data also indicated that TUB was the most stably expressed reference gene for the H 2 O 2 and drought treatments. This gene is involved in plant growth and development, with an important role in maintaining cell morphology, promoting intracellular transport, and mediating cell movement and cell division, and has often been used as a reference gene. Kiarash et al. 38 determined that TUB is the most appropriate reference gene for wheat plants under drought conditions. Although ACT is the most commonly used reference gene for cucurbit crops [39][40][41] , it was not identified as the best reference gene for any of the conditions analyzed in this study. Furthermore, ACT expression reportedly varies substantially in many other crops 42 . Kong et al. 20 suggested that ACT should be used together with other reference genes for qRT-PCR analyses of melons. In the current study, except for the cold treatment, GAPDH was the most unstably expressed reference gene, implying that GAPDH is not a suitable reference gene for investigations of the effects of abiotic stresses on luffa gene expression. This is consistent with the results of studies by Liu et al. 27 and Yu et al. 43 . Although GAPDH is the most unsuitable reference gene for annual ryegrass and ramie under different stress conditions, it is ideal for studies of Juglans regia L 44 , L. speciosa 45 and Stipa grandis 46 . Therefore, potential internal reference genes must be evaluated and confirmed as appropriate for different experimental systems and materials before they are used to analyze gene expression.
To verify the reliability of the identified reference genes, we analyzed the expression of luffa LcCu/Zn-SOD1 gene. Because of its protective effects on cell membranes, Cu/Zn-SOD can protect plants from oxidative damage due to reactive oxygen species, while also contributing to plant responses to various abiotic stresses 47 . In the present study, the expression of LcCu/Zn-SOD1 gene initially increased in response to the abiotic stress, after which it decreased. This expression pattern was consistent with previously reported stress-induced expression trends [48][49][50] . Earlier studies tended to emphasize the use of a single reference gene for calibrating the expression of target genes. However, a single reference gene may be unstably expressed under different experimental conditions. Thus, applying multiple genes may increase the accuracy and reliability of data normalizations to some extent 51 . In this study, the combined use of two stably expressed reference genes generally enhanced the quantification of Cu/Zn-SOD gene expression when compared with the data obtained with a single reference gene. Deng et al. 22 suggested that potential interactions between reference genes should be examined when using multiple reference genes. A positive or negative relationship between the expression of the reference genes may affect the analysis of the target gene expression level. Reference genes with linearly additive expression may be used together as internal reference controls to improve the normalization of expression data under various experimental conditions. Therefore, the method used to select multiple reference genes should be carefully considered to ensure appropriate combinations are applied.

Methods
Plant materials and treatments. 'Minyan No. 1' , a commercial F 1 luffa hybrid, was used in the present study. The seeds were first wrapped with wet gauze and then germinated in an incubator at 30 °C. The seedlings were sown in 10 × 8 cm pots containing soil and cultivated in an artificial climate incubator with a 16-h artificial light (25 °C):8-h dark (18 °C) cycle and 65-75% relative humidity. Seedlings at the third true leaf stage were used for the following treatments. For the heat and cold treatments, the seedlings were incubated at 40 °C and 4 °C, respectively, for 24 h. For the ABA and H 2 O 2 treatments, the seedlings were sprayed with 200 μM ABA and 100 μM H 2 O 2 , respectively. Regarding the salt and drought stress treatments, the seedlings were transplanted into full-strength Hoagland's nutrient solution containing 100 mM NaCl or 15% PEG 6000, respectively, for 24 h. Leaves were then collected at 0, 2, 6, 9, 12, and 24 h after the treatments. Each experiment was completed Reference gene selection and primer design. The seven genes initially selected as candidate reference genes included the ACT , TUA , TUB, EF-1α, and GAPDH genes as well as two previously cloned reference genes, UBQ 16 and 18S 15 . Details regarding these genes are provided in Table 1. All qRT-PCR primers were designed with the Primer 5.0 software. The specificity of each primer pair was evaluated by 1% agarose gel electrophoresis and a melting curve analysis. The PCR amplification efficiency (E) and correlation coefficient (R 2 ) for each qRT-PCR primer pair were calculated with a fivefold cDNA dilution series 52 .
qRT-PCR analysis. A qRT-PCR assay was completed with SYBR Premix Ex Taq (TaKaRa) and the QuantS- Box plots were drawn to visualize the reference gene expression levels and variations. The geNorm, NormFinder, BestKeeper, and RefFinder algorithms were used to assess the expression stability of selected genes. Specifically, geNorm calculates the expression stability value (M) and pairwise variation (V), and the most stably expressed gene is the one with the lowest M value. Moreover, geNorm determines the optimal number of reference genes according to the relative value V n /V n+1 . NormFinder calculates a stability value (SV) for each gene based on the variance analysis, and the gene with the lowest SV is identified as the most stably expressed gene. For both algorithms, the Ct values should first be transformed by 2 −ΔCt . BestKeeper mainly determines the stability of reference gene expression based on the standard deviation (SD), Pearson correlation coefficient (r), and the coefficient of variation (CV) of the Ct data for all reference genes. The most stably expressed gene is the one with the lowest SD and CV values. Finally, the reference genes were ranked based on the geometric mean (GM) values calculated with RefFinder (http://www.leonx ie.com/refer enceg ene.php).

Validation of reference genes.
To confirm the reliability of the reference genes, the expression profiles of Cu/Zn-SOD gene were determined and normalized against the data for the two most stable and the least stable reference genes 22 . The qRT-PCR amplification conditions were as described above and the primer pairs are listed in Table 1. The relative expression data were calculated with the 2 −ΔΔCt method 53 . Statistical analyses were performed with the SPSS 18 program.  TUB  5′-ACC GTG AGA AGA TGA GGG AA-3' 5′-GTT ACT AAT TGT CGA GGT CC-3'   EF-1α  5′-AGG CTG CTT ATT TCT ACG GA-3'  5′-AAA TGC AGG CTT GGC TGG TT-3'   GAPDH  5′-CCC TGA ATC CCC ACT TCT TT-3'  5′-AAA CTA ACC CAT GGT GTA GG-3'   ACT  5′-CTT CGA GCC AAA TCG CTT TC-3'  5′-TCT CAG GTA AAA CCG TAC CG-3'   TUA  5′-TAT TCT CAG AGG CAC ACT CG-3'  5′-AAC GGC AGG CTC TTG AAC TA-3'