Genome-wide identification of new reference genes for RT-qPCR normalization in CGMMV-infected Lagenaria siceraria

Plant preparation, virus inoculation

The cultivation and management of bottle gourd (L. siceraria, accession “Hangzhou Gourd”) were performed as follow: after soaking and germination, the seeds were first transplanted into 10 cm nutrient preparations with soil rich in organic matter, and when the seedlings grew to two and a half leaf stage, transplanted them into 20 L PVC drums (1 per barrel). The mixed substrate used was peat: vermiculite: perlite: organic fertilizer = 4:4:1:1, The pH of the culture substrate was about 7.0 and the water content was maintained at about 70% relative humidity. A “flower-free” nutrient solution was used once a week (N:P:K = 20:20:20) (Shanghai Yongtong Chemical Co., Ltd., Shanghai, China). The greenhouse conditions were daily temperature 25–28 °C, night temperature 18–20 °C; photoperiod 14 h/d (light intensity is greater than 87.5 µmol m⁻² s⁻¹). Scaffolding, topping, pruning and pollination were carried out according to routine management. The fruits were harvested 10 days after pollination.

CGMMV inoculum (CGMMV-ZJ) was sap from L. siceraria plants with typical symptoms that had been infected with a CGMMV infectious clone 14 d earlier (Zheng et al., 2015). At least six plants were inoculated with CGMMV-infected sap at the two and a half leaf stage on the two expanding leaves. Approximately 1 g of plant tissue was homogenized in 20 volumes of inoculation buffer (0.1 M phosphate buffer, pH 7.5, 0.2% sodium sulfite and 0.01 M 2-mercaptoethanol), while the mock plants were only inoculated with buffer.

RNA sequencing

According to the protocol of TruSeq Small RNA Sample Prep Kits (Illumina, San Diego, CA, USA), the total RNA of about 5 µg was extracted for the preparation of small RNA library. Sequencing of the RNA-Seq libraries was carried out on an Illumina Hiseq2500 at LC-BIO (Hangzhou, China) following the manufacturer’s protocol.

RNA and first strand cDNA preparation

Three replicate samples of flesh tissue of the ripe fruits and newly expanding leaves from both inoculated and control plants were collected for RNA extraction. Total RNA was extracted from each sample using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) according to the manufacturer’s instructions. The RNA quantity and quality from each sample was evaluated by denaturing agarose gel electrophoresis and microfluidic capillary electrophoresis with the Agilent 2100 bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). Only RNA samples with a complete band and A260/A280 ratio in the range 1.8–2.0 were used for the next step. All RNA samples were stored at −70 °C. For virus detection, the first strand cDNA was synthesized using ReverTra Ace-α-^® kit (TOYOBO, Japan) following the product’s protocol. The infection of CGMMV in the tissues was confirmed by CGMMV specific primers, and the primers of ZYMV and WMV were also used to monitor the presence of these two common viruses occurred in cucurbit crops (Heeju et al., 2015). For RT-qPCR, first strand cDNA was synthesized from 1 µg total RNA using PrimeScript™ RT reagent Kit with gDNA Eraser (Perfect Real Time) kit (TaKaRa, Dalian, China) according to the manufacturer’s instructions. The negative controls without PrimeScript RT Enzyme Mix I were analyzed in parallel to detect the presence of genomic DNA contamination in the RNA samples.

Selection of candidate RGs

Partial candidate RGs in leaves and fruits were selected from publicly available references, but most were from our RNA sequencing data. To obtain RGs that are stably and moderately or highly expressed in CGMMV-infected leaves, we kept the Reads Per Kilobases per Million reads (RPKM), ratio of the maximum to the minimum RPKM (RPKMmax/min), and coefficient of variation (CV) to >40, <2.0, and <0.3 at p < 0.05, respectively. In fruit, the RPKM, RPKMmax/min, and CV were maintained at >40, <2.0, and <0.2 at p < 0.05, respectively. All selected internal RGs have only one transcript and were ranked from small to large according to their RPKMmax/min values.

To select better RGs for both leaves and fruits, the RPKM and RPKMmax/min were kept at >40 and <2.0 for the RGs commonly used in cucurbit plants in keeping with the RNA-seq data. The RGs for leaves and fruits were screened and analyzed simultaneously with 14 common RGs of cucurbit crops in previous studies. A total of 16 RGs from bottle gourd leaves and fruits were screened and analyzed.

Primer design and verification of selected gene amplicons

The fourteen common RGs were amplified according to the references or based on primers designed by Primer-BLAST of the RNA sequence data of leaves and fruits (Table 1). Specific primers for the candidate RGs from our RNA-sequencing data were designed using Primer 3 (http://primer3.ut.ee/) (Table 2). All PCR amplicon lengths were between 80–200 bp. All primers were synthesized by a commercial supplier (Biosune, Hangzhou, China).

To check the specificity of all primers, the cDNA of each sample was amplified by PCR, and the amplified products were separated by electrophoresis on 3% agarose gel and purified using a QIAquick Gel Extraction Kit (Qiagen, Hilden, Germany) according to the manufacturer’s instructions, and cloned into pEASY-Blunt zero (Transgen, Beijing, China) followed by sequencing.

The quantification cycle (Cq) values obtained by qRT-PCR on a standard curve generated from a fourfold dilution series of one sample at six dilution points for three technical replicates were used to draw the standard curve to get R2 and slope values. The PCR amplification efficiency of each primer was calculated using the equation: E(Efficiency)% = (10 ^{[−1∕slope]} −1) ×100%.

Quantitative real-time PCR

qRT–PCR was carried out in 384-well plates using the QuantStudio 6 Flex real-time PCR detection system (ABI, USA). Each reaction mixture consisted of 5 µL SYBR Green Realtime PCR Master Mix (TaKaRa, Dalian, China), 0.5 µL cDNA diluted fivefold, 0.5 µL (10 mM) each of forward and reverse primers, and 3.5 µL RNA-free H₂O, equating to a final volume of 10 µL in each well. The qPCR reaction was as follows: initial denaturation at 95 °C for 5 min and 40 cycles of amplification (95 °C 15 s, 58 °C 20 s and 72 °C 20 s). Subsequently, fluorescence acquisition was performed after each cycle. A melting curve was generated after 40 cycles of amplification by heating at 65–95 °C. Cq values and baseline were set automatically by the QuantStudio Real-Time PCR Software v1.2 (ABI, USA) using default parameters.

Gene expression stability analysis

The programs geNorm (Vandesompele et al., 2002), NormFinder (Andersen, Jensen & Orntoft, 2004), BestKeeper (Pfaffl et al., 2004) and RefFinder (http://150.216.56.64/referencegene.php?type=reference) were used to analyze the stability of the candidate RGs under CGMMV infection conditions. All software packages were used according to the manufacturer’s instructions.

Validation of the selected RGs

LsIPT and LsDdRP genes in CGMMV-infected bottle gourd leaf and fruit tissue were selected to detect the effectiveness of these identified RGs. Primers for the two genes were designed as described above and listed in Table 1. The best RGs identified by the algorithms above were used for normalization.

Transcriptome analysis of Lagenaria siceraria under CGMMV infection based on RNA-seq

Bottle gourd leaves and fruits infected by CGMMV were collected from three replicate virus-inoculated plants, and the presence of CGMMV in each sample was further confirmed by RT-PCR and western blot, and the contamination of ZYMV and WMV was excluded by RT-PCR (Fig. S1). The control leaves and fruits samples were in parallel collected from three mock bottle gourds. The analysis of bottle gourd transcripts before and after CGMMV infection showed 639 and 3,930 non-differentially expressed genes (— log2 fold_change — <1, P ≤ 0.05) in the leaves and fruits of bottle gourd, respectively (Tables S1 and S2). And these non-differentially expressed genes were used as the source of candidate RGs from the RNA-Seq dataset.

Selection of candidate RGs

In the present study, 11 and 86 candidate RGs respectively from bottle gourd leaves and fruits were screened from our RNA-Seq dataset by setting up a series of conditions (Tables S3 and S4). Only two genes, LsH3 and LsWD (Table S5), met the criteria to be candidate RGs for both leaves and fruits. The primers were designed based on the gene sequences in the database. To select candidate RGs that could be used in both bottle gourd leaves and fruits, these two novel genes with other 14 traditional candidate RGs were used to compare their expression stability. To select candidate RGs that could be used in bottle gourd leaves or fruits, separately, the commonly used reference gene sequences were then compared with the bottle gourd transcriptome data (Table S6). All the 11 candidate RGs of bottle gourd leaves screened from our RNA-Seq data, and seven traditional candidate RGs, LsPP2A, LsADP, LsEF1α, LsCYP, LsRPS15, LsTBP, and LsRPL23 screened from transcriptome comparison data were selected as candidate RGs of bottle gourd leaves. Of the 86 genes screened from the bottle gourd fruit transcriptome data, the first 11 (based on the ratio of RPKMmax/min) were selected for further analysis in addition to LsH3 and LsWD, and five traditional candidate RGs, LsPP2A, LsADP, LsTBP, LsTUA, and LsGAPDH screened from transcriptome comparison data were selected as candidate RGs of bottle gourd fruits. Therefore, a total of 16 common RGs of both bottle gourd leaves and fruits (Table 1) and 18 RGs of bottle gourd leaves and fruits separately were screened and analyzed (Tables 1 and 2).

Evaluation of target specificity and amplification efficiency in RT-qPCR reactions

Preliminary evaluation of candidate reference gene primers was performed by evaluating primer specificity and efficiency. The single peak in melting curve analyses following RT-qPCR confirmed the specific amplification of each gene (Fig. S2). Each amplicon was detected by agarose gel electrophoresis, only a single fragment of the expected size (80–200 bp) was observed (Fig. S3). Further sequencing results showed all genes sequences were exactly 100% identical to those of the corresponding genes in bottle gourd transcriptome databases. Amplification efficiencies of bottle gourd leaves ranged from 90.1% to 110.1%, whereas those of the fruits ranged from 89.9% to 116.2% (Tables 1 and 2). Furthermore, the standard curves showed good linear relationships (>0.981) between the Cq values and the log-transformed copy numbers of all candidated RGs (Tables 1 and 2). There was no band detected in the negative controls, indicating that the genomic DNA contamination does not exist (Fig. S3).

Expression intensity of candidate RGs

In order to fully understand the relative expression intensity of all candidate RGs in bottle gourd, three biological and three technical replicates (n = 9 for each gene) were used to determine the Cq values for all RGs. From the graph for bottle gourd leaves, the mean Cq values of these candidate genes we selected ranged from 16.51 (LsCYP) to 25.63 (LsTBP) (Fig. 1A), which represented the highest and lowest accumulation levels, respectively. The minimal variation in gene expression was LsUBC (<0.81 cycles) in bottle gourd leaves. The lowest and highest median Cq value of the mRNA accumulation levels in the fruits of bottle gourd was 20.94 (LsH3) and 28.88 (LsTBP), respectively (Fig. 1B). LsPLA expression exhibited the least amount of variation (<1.14 cycles) in bottle gourd fruits. The Cq value 15.18 of LsCYP in leaves was the lowest and 28.04 of LsPP2A in fruits was the highest (Fig. 1C). The minimal variation in gene expression observed in both leaves and fruits were 2.95 cycles (LsACT), 3.06 cycles (LsADP) and 3.78 cycles (LsWD). However, only the comparison of raw Cq values is not sufficient to evaluate the expression stability of candidate RGs, a further intensive statistical analysis was required for more accurate assessment.

Comparison of the expression stability of the universal traditional candidate RGs and novel candidate RGs in both bottle gourd leaves and fruits

To select universal candidate RGs expressed stably both in bottle gourd leaves and fruits infected by CGMMV, we screened from our RNA-Seq dataset by setting up a series of conditions (RPKMmax/min < 2, and CV < 0.3, p < 0.05, RPKM min > 40). Only two genes, LsH3 and LsWD, met the criteria to be candidate RGs and we then compared the expression stability of these two novel genes with 14 other commonly-used cucurbit RGs. The expression profiles (Fig. 2A) and the variations of the 16 genes (Fig. 2B; Tables S5 and S6) in bottle gourd suggested that LsWD, LsH3 and LsTBP and LsGAPDH were expressed more stably (CV < 0.3), and the variation in respective expression levels was the lowest among all the genes during CGMMV infection, indicating that these genes, especially the two novel genes LsWD and LsH3, may be more suitable for normalization than other traditional candidate RGs. In addition, LsWD, LsH3 and LsGAPDH each had only one transcript, which could facilitate primer design and ensure the accuracy and reliability of the RT-qPCR compared to other genes.

The stability of these 16 potential RGs was further evaluated with two statistical methods. By geNorm analysis, the average gene expression stability (M) of all of the universal candidate RGs were compared, among them, LsH3 and LsGAPDH showed the lowest M value (M = 0.284) in both leaves and fruits, followed by LsWD (M = 0.340) (Fig. 3A), indicating that these genes displayed the most stable profiles. Pairwise variation V_n∕n+1 was less than 0.15 in all leaf and fruit samples (Fig. 3B), indicating that adding other RGs was not necessary, and demonstrating that at least two reference genes were required for more reliable normalization, the top two gene were LsH3 and LsGAPDH. The raw Cq values were also transformed into Q values for NormFinder analysis. The lowest stability value of NormFinder analysis indicates the most stably expressed gene. By NormFinder analysis, the best three universal RGs in both leaf and fruit were LsTBP, LsWD and LsH3 (Fig. 4A). So, both geNorm and NormFinder analysis suggested that the two novel candidate RGs LsWD and LsH3 were suitable to evaluate the gene expression stability of bottle gourd leaves and fruits infected by CGMMV.

Expression stability of candidate RGs in CGMMV infected leaves and fruits of bottle gourd separately from transcriptome analysis

To further select more suitable candidate RGs expressed stably in bottle gourd leaves or fruits infected by CGMMV separately, 18 RGs of bottle gourd leaves or fruits separately obtained from the RNA-Seq dataset were compared using different algorithms. By geNorm analysis, the average gene expression stability of the 18 candidate RGs of bottle gourd leaves and fruits screened from the RNA-Seq were all less than 1.5, respectively (Figs. 3C and 3E; Table 2). For all the tested leaf samples, LsTBP and LsCYP showed the lowest M value (M = 0.213) (Fig. 3C) of all the candidate RGs while in fruit samples, LsP4H and LsVAMP had the lowest M value (M = 0.212) (Fig. 3E). The pairwise variation V_n∕n+1 of each sample was also less than 0.15 (Figs. 3D and 3F), and the top two reference genes those had the lowest M value were needed for more reliable normalization at least. Similarly, by NormFinder analysis, the best three genes screened from the RNA-Seq dataset were LsCYP, LsH3 and LsPP2A in leaves (Fig. 4B), and LsTBP, LsP4H and LsXRN1 in fruits (Fig. 4C).

For further clarification, the expression stability of these candidate RGs was examined by two more algorithms. BestKeeper software can only compare the expression levels of up to 10 internal control genes in 100 samples, so only the top ten genes identified by geNorm and NormFinder were selected for subsequent assessment. Of these, the top three candidate internal RGs were LsCYP (r = 0.995, p-value =0.001), LsH3 (r = 0.984, p-value =0.001) and LsTBP (r = 0.964, p-value =0.002) in leaf samples, and LsP4H (r = 0.978, p-value =0.001), Ls VAMP (r = 0.975, p-value =0.001) and LsTBP (r = 0.959, p-value =0.002) in fruit samples (Table S7). The results of BestKeeper were therefore broadly consistent with geNorm and NormFinder. We also compared and ranked the tested candidate RGs based on a web-based comprehensive analysis tool, RefFinder, which suggested that the top three candidate RGs screened from the RNA-Seq dataset in bottle gourd leaves were LsCYP, LsH3 and LsTBP, while those in bottle gourd fruits were LsP4H, LsADP, and LsTBP (Table S7). These should therefore be the best RGs to use in RT-qPCR.

Validation of the candidate RGs

According to the transcriptional data, expression of the IPT and DdRP genes of L. siceraria changed significantly in response to CGMMV. LsIPT increased 1.78 fold in leaves and decreased 1.2 fold in fruits compared with their mock-inoculated tissues, while LsDdRP increased 1.4 fold in leaves and increased 1.63 fold in fruits (Table S8). These genes were therefore chosen to evaluate the reliability of the top candidate RGs as indicated by the previous analysis. The top rank RGs LsH3, LsGAPDH, LsWD, LsCYP, LsTBP and LsP4H selected by geNorm and RefFinder were used as candidate RGs. Among them, LsTBP, LsWD, LsH3, LsGAPDH and LsCYP were selected to use as RGs in leaves, and all these genes together with LsP4H for fruit.

LsIPT increased 2.75 fold when LsWD was used as the reference gene in leaves, with 2.95, 3.82, 3.93 and 4.43 fold increases using LsGAPDH, LsH3, LsCYP and LsTBP respectively. The values of LsIPT normalized fold expression in fruits were 0.54 (LsTBP), 0.52 (LsWD), 0.60 (LsH3), 0.57 (LsGAPDH), 0.59 (LsCYP) and 0.56 (LsP4H) (Fig. 5). LsDdRP increased 1.78, 1.95, 2.51, 2.57 and 2.88 fold in leaves when LsWD, LsGAPDH, LsH3, LsCYP and LsTBP were the internal RGs respectively, while the corresponding values in fruits were 1.42 (LsWD), 1.53 (LsTBP), 1.60 (LsH3), 1.62 (LsGAPDH), 1.80 (LsP4H) and 2.31 (LsCYP) (Fig. 5). The RT-qPCR results showed that all candidate RGs gave results consistent with the data from the transcriptome database. Overall, the most suitable internal RGs chosen for use in leaves were LsWD, LsGAPDH and LsH3 and those for fruits were LsH3, LsGAPDH, LsP4H and LsCYP.