Prediction of candidate reference genes in VBT transcriptomics and primer design

The fruit transcriptome includes green fruits, red fruits and blue fruits, which represent the progressive stages of fruit ripeness. ²⁷ According to previous studies and the transcriptomic analysis, we identified ten reference genes. There were five actin-related genes (actin-related protein 2/3 complex subunit 2A (actin 1), actin-7 (actin 2), actin-depolymerizing factor 2 (actin 3), actin-97 (actin 4), and actin-related protein 3 (actin 5)), 50S ribosomal protein L27 (RP), ubiquitin-conjugating enzyme E25 (UBE), ubiquitin domain-containing protein DSK2b (UBQ), NADH dehydrogenase (NADH), and adenylate kinase 4 (ADK), all of which have been widely used for gene validation in plants and fruits (Table S1). Simultaneously, RNA from the green, red, and blue fruits was extracted, and pure total RNA was quantified. Furthermore, the primers for RT–qPCR were designed (Table S2), and later these primers were confirmed using melting curves. The results showed that a single peak for each primer pair was present, which indicates exceptionally high primer specificity (Fig. 1).

Expression levels of the candidate reference genes

To evaluate the expression levels of the candidate reference genes, the cycle threshold (Ct) values of the ten candidate reference genes in green, red, and blue fruits were acquired with three replicates. The expression levels of these ten reference genes were dramatically varied, with Ct values ranging from 17.61 to 29.73 cycles. According to the Ct values, actin 2 expression levels in the three fruit samples were the highest, and actin 5 was the lowest (Fig. 2). In the green fruits, all the gene expression levels ranged from actin2 > actin3 > NADH > actin4 > UBE > ADK > UBQ > actin1 > RP > actin5 (Table S3). In thered and blue fruits, the arrangement changed in the order actin2 > actin3 > NADH > actin4 > ADK > UBE > UBQ > actin1 > RP > actin5 (Table S3). According to these data, we also found that the expression levels of actin3, NADH, ADK, and UBE had narrow variances across the three fruit development stages, which convincingly suggested that their expression levels were stable in the three fruits (Fig. 2). These data demonstrated that these genes could be used as reference genes.

Ct values variance among three fruit development stage. RT-qPCR was performed using Bio-Rad CFX Manager, and Ct values of each reference genes in the samples (S1, S2 and S3) were acquired later. The picture were drawn by R, using these Ct values.

Expression stability of candidate reference genes

geNorm analysis

The average expression stability (M value) of the ten genes was calculated by geNorm (version 3.5). The M values of all ten genes across green, red, and blue fruits were lower than the threshold of 1.5, which means that all the selected genes could be used as reference genes (Table 1). However, the M value of these genes showed a great variability; the M value of the RP gene was the highest, while those of the NADH and ADK genes were the smallest, suggesting that they had the highest expression stability (Table 1). These results showed that NADH and ADK are the best reference genes, followed by actin2. To further assess the optimal number of reference genes, pairwise variation (V_n/V_n + 1) was analysed using geNorm, and 0.15 was used as a cut-off value according to a previous study. The results showed that all the V_n/V_n + 1 values were less than 0.15. Moreover, the pairwise variation (V2/V3) value was significantly less than 0.15 (Fig. 3), suggesting that only two reference genes were necessary for the gene profile; thus, we surmised that no additional reference genes were required.

Pairwise variation (Vn/Vn + 1) values for three fruit ripen stage. Vn/Vn + 1 > 0.15 indicates another (n + 1) reference was required, while Vn/Vn + 1 ≤ 0.15 demontrates only n reference was enough.

NormFinder analysis

In addition to the geNorm analysis, NormFinder analysis is also a necessary method for confirming the best reference gene. The gene with the smallest S-value was the most suitable reference gene. According to the NormFinder analysis, the S-value of all the candidate genes was less than 1.0, which indicates that all the reference genes we evaluated could be used for RT–qPCR analysis. Of them, the S-value of the actin1 gene (0.03) was the smallest and that of RP (0.785) was the largest (Table 1), which indicates that the actin1 gene is the most stable reference gene, and the results were not uniform with geNorm analysis.

BestKeeper analysis

Finally, BestKeeper was used to evaluate the stability of the reference genes through the standard deviation (SD) and the coefficient of variance (CV) of the Ct values. If the SD value of the Ct values were > 1, the genecould not be set as a reference gene. According to BestKeeper analysis, both Ct values of all the genes were smaller than 1.0, which indicates that all the genes we selected met the standard of BestKeeper. Among them, however, actin3, UBE, and NADH are the three genes with lowest SD values. Like the SD value, these three genes had smaller CV values (Table 1). However, UBQ, actin5, and RP both had higher SD values and CV values (Table 1). Thus, while these results demonstrate that all the ten genes we selected could be used as reference genes, the actin3, UBE, and NADH genes are the best choices overall.

Recommended comprehensive ranking

Since the analysis by the three methods is not uniform, we compared and ranked the reference genes by RefFinder, which is a comprehensive, web-based tool combined with geNorm, Normfiner, BestKeeper, and Ct values. The overall final ranking based on the ranking from each program is shown in Table 2. These results showed that actin2, ADK, actin1, and NADH were the four genes with smaller values, and the overall ranking values (actin2, 2.59; ADK, 2.83; actin1, 3.15; NADH, 3.20) were very close (Table 2). However, partial Ct values of actin1 in red fruit were over 26, suggesting that its expression levels are slightly low and may not be suitable for RT–qPCR analysis of genes in all stages of fruit development. Considering all these results, we conclude that the three reference genes, actin2, ADK and NADH, are the best reference genes for RT–qPCR analysis across all stages of fruit development.

Identification of reference genes using marker genes in fruit development

Since ACO and PG genes are generally upregulated in fruit ripening programs, assessing their expression levels should be the best way to confirm the reliability of reference genes. Additionally, the ranking scores of the three best reference genes (actin2, ADK and NADH) is close, which means any one is suitable for gene expression analysis. To assess the reference gene, the expression levels of ACO1, ACO2, and PG were calculated using NADH, ADK, NADH + ADK, RP and actin5. We found that the transcript levels of ACO genes normalized by RP were not in accordance with their expression features (Fig. 4). Moreover, the transcript levels of PG normalized by actin5 in red fruit did not agree with the PG expression feature and the transcriptomic analysis (Fig. 4). However, the results calculated by NADH and ADK are like the expression features of these genes and the transcriptomic analysis (Fig. 4). Furthermore, NADH + ADK as the coreference gene is not consistent with NADH or ADK alone (Fig. 4). Taken together, either NADH or ADK are the most suitable reference genes for RT–qPCR analysis in VBT fruit development.

Evaluation of reference genes using the mark genes involved in fruit ripening. S1 indicates the green fruits, S2 indicates the red fruits, S3 indicates the blue fruits. Relative expression levels ot the three genes were normalized to reference genes NADH, ADK, NADH + ADK, RP, and actin5, respectively. Each RT-qPCR experiment was repeated at least three times. Labels for each reference genes as normalizer were shown in the picture. Pictures were drawn using GraghPad Prism 7.0.

Differential expressed genes from transcriptomic analysis validation using NADH

To further evaluate the reliability of the reference genes, the expression levels of 15 flavonoid synthesis-related genes were normalized to the expression level of NADH. These genes encode shikimate O-hydroxycinnamoyltransferase (SHT), flavonol synthase (FLS), flavonoid 3′-monooxygenase (F3MO), flavonoid 3′,5′-hydroxylase 2 (F35H2), dihydroflavonol 4-reductase (DFR), leucoanthocyanidin dioxygenase (LDOX), leucoanthocyanidin reductase (LAR), anthocyanidin reductase (ANR), UDP-glycosyltransferase 74B1 (UGT), caffeoyl-CoA O-methyltransferase (CcoA), probable caffeoyl-CoA O-methyltransferase (At4g), hydroxypalmitate O-feruloyl transferase (HOFT), transcription factor MYB113 (MYB113), acetylajmalan esterase (Ace), and 2-oxoglutarate-dependent dioxygenase (ODD) (Table S4). The results showed that the Pearson correlation coefficients of the 12 genes are above 0.9 (Fig. 5, Table S5), which indicates that their RT–qPCR results are highly consistent with the FPKM values in the transcriptome. In addition, the Pearson correlation coefficients of FLS and LAR2, and ANR2 arranged from 0.81 to 0.89 (Table S5), indicating a lower correlation. However, their expression patterns were also highly consistent with the transcriptomic data. Taken together, the reference gene NADH is highly reliable for RT–qPCR analysis in the fruit ripening stage.

Validation of the DEGs in three ripening stage using NADH. FLS, flavonol synthase; F3MO, flavonoid 3′-monooxygenase; F35H2, flavonoid 3′,5′-hydroxylase 2; DFR, dihydroflavonol 4-reductase; LDOX, leucoanthocyanidin dioxygenase; LAR, leucoanthocyanidin reductase; ANR, anthocyanidin reductase; UGT, UDP-glycosyltransferase; CcoA, caffeoyl-CoA O-methyltransferase; At4g, probable caffeoyl-CoA O-methyltransferase; HOFT, hydroxypalmitate O-feruloyl transferase; MYB113, transcription factor MYB113; Ace, acetylajmalan esterase; ODD, 2-oxoglutarate-dependent dioxygenase; SHT, shikimate O-hydroxycinnamoyl-transferase. The blue line indicates RT-qPCR data normalized by NADH, and the red lines represent FPKM values. Both RT-qPCR experiments and RNA-seq data were repeated at least three times. Pictures were drawn in the Microsoft excel based the RT-qPCR data and FPKM value.