Phenotypic characterization of EP26 and EP28
In the present study, two dark purple eggplant advanced lines EP26 and EP28, with different fruit peel color intensity were compared (Fig. 1a). Based on the lightness (L*) of fruit skin, the fruit color intensity of EP26 was deeper than that of EP28 both at DDA15 and DDA25 stage (Fig. 1b). The total content of anthocyanin and chlorophyll from the two different stages (DDA15 and DDA25) of EP26 and EP28 were measured, respectively. The content of anthocyanin of EP26 at DDA15 and DDA25 were higher than those of EP28 (Fig. 1c). However, the level of chlorophyll content of EP26 was lower than that of EP28 at both stages (Fig. 1d). In addition, a total of eight anthocyanin compounds in fruit skin of EP26 and EP28 at DDA25 were detected, among which delphinidin 3-O-rutinoside is the most abundant both in EP26 and EP28 (Table 1).
Identification of differentially expressed genes (DEGs) between EP26 and EP28
In order to systematically study the phenotypic differences between EP26 and EP28, the transcriptomes of four eggplant fruit peel samples with three repeats (A: EP26 at DAA15, B: EP26 at DAA25, C: EP28 at DAA15, D: EP28 at DAA25) were obtained using RNA-Seq. In total, twelve libraries (two materials × two developmental stages × three biological replicates) were constructed. Differences in gene expression in the fruit peel of EP26 and EP28 at two developmental stages were assessed (A vs C and B vs D). DEGs were identified with the expression fold (log2Ratio ≥ 1) and false discovery rate ≤ 0.01 as the thresholds.
A total of 1456 DEGs (913 up-regulated and 543 down-regulated) and 1163 DEGs (821 up-regulated and 342 down-regulated), were detected in A vs C and B vs D, respectively (Fig. 2a). Overall, a total of 2218 DEGs were found between EP26 and EP28 at the two stages. Of all these DEGs, 401 genes with 231 up-regulated and 111 down-regulated were shared at both stages (Fig. 2b).
Gene ontology and KEGG pathway analyses of differentially expressed genes
The functions of the DEGs detected in the two developmental stages between EP26 and EP28 were classified according to the GO database respectively. Go terms were divided into three functional categories: biological process, cellular component, and molecular function. In the biological process category, the DEGs of the two stages were most enriched in the oxidation-reduction process (GO: 0055114). In the molecular functional category, the highly enriched terms among three stages were ATP binding (GO: 0005524) and structural constituent of ribosome (GO: 0003735). In the cellular component category, the GO terms cell periphery (GO: 0071944) were most enriched in A vs C comparison, while the GO terms including photosystem I (GO: 0009522), photosystem I reaction center (GO: 0009538), chloroplast thylakoid membrane (GO: 0009535) and photosystem II (GO: 0009523) were highly enriched terms in B vs D comparison (Table S1).
KEGG pathway enrichment analyses were performed according to the DEGs in the two stages respectively, to further explore the biological functions. As shown in Figure 3, the pathways of phenylpropanoid biosynthesis and flavonoid biosynthesis involved in anthocyanin biosynthesis were significantly enriched in the fruit skin of EP26 at DDA15. Particularly, the pathways of photosynthesis-antenna proteins, photosynthesis, and flavonoid biosynthesis were significantly enriched in the B vs D comparison.
Genes involved in anthocyanin biosynthesis
The anthocyanin biosynthetic pathway is an extension of the general flavonoid pathway. In the present study, structural genes participated in each step of the anthocyanin biosynthesis pathway were identified. Compared to EP28, the predicted encoding genes of phenylalanine ammonia-lyase (PAL), 4-coumarate-CoA ligase (4CL), chalcone isomerase (CHI), dihydroflavonol 4-reductase (DFR), flavanone 3-hydroxylase (F3H), and flavonoid 3-o-glucosyltransferase (UFGT) were up-regulated in EP26 at both two stages. At commercial maturity stage (DDA25), two unigenes encoding chalcone synthase (CHS) and one unigene encoding flavanone 3′, 5′-hydroxylase (F3′5′H) were also up-regulated in EP26 compared to EP28, but did not change expression in the comparison of EP26 vs EP28 at DDA15. Notably, all the identified structural genes including CHS, CHI, F3H, F3′5′H, DFR, ANS, and UFGT were up-regulated in EP26 vs EP28 comparison at DDA25 (Figure 4).
Genes involved in chlorophyll metabolism
In current study, DEGs related to chlorophyll biosynthesis were identified based on KEGG pathway annotations. Compared to EP28, only four genes involved in chlorophyll biosynthesis including HEMA (encoding glutamyl-tRNA reductase), CHLM (encoding Mg-protoporphyrin IX methyltransferase), CRD (encoding Mg-protoporphyrin IX monomethylester cyclase), and PORA (encoding protochlorophyllide oxidoreductase A) were up-regulated in EP26 at DDA25 stage. However, none of DEGs involved in chlorophyll biosynthesis were identified between EP26 and EP28 at DDA15 stage (Table 2). Moreover, compared to EP28, the expression levels of DEGs related to photosynthesis were mostly up-regulated in EP26 at DDA25 based on KEGG enrichment. Particularly, a total of 11 unigenes encoding chlorophyll a-b binding protein, which transfer excitation energy to the reaction centers in photosystem I and photosystem II, showed increased expression levels in EP26 at DDA25 (Table 2).
Identification of differentially expressed transcription factors
Transcription factors play important roles in regulating gene expression in anthocyanin synthesis. In the present study, a total of 131 TFs were found to be differentially expressed in EP26 at least one pairwise combination (A vs C and B vs D) (Table S2). These TFs could be classified into 29 families, among which bHLH constituted the largest group (20 genes), followed by MYB (16 genes), HD-ZIP (9 genes), ERF (8 genes), NAC (8 genes), and WRKY (7 genes) (Fig. 5a).
Among the 20 DEGs encoding bHLH TFs, only one DEG Sme2.5_03973.1_g0005 was shared by both stages, which was down-regulated in EP26 compared with EP28. The expression level of Sme2.5_00592.1_g00005, which is homologus to AtTT8 increased in the comparison of EP26 and EP28 at DDA25. Among the 16 unigenes encoding MYB TFs, two unigenes Sme2.5_05099.1_g00002 and Sme2.5_00492.1_g00008 were shared by both A vs C and B vs D comparison. Sme2.5_05099.1_g00002 was up-regulated and phylogenetic analysis showed that this unigene was homologous to Arabdopsis MYB113, which has been reported to be involved in the control of anthocyanin biosynthesis. Sme2.5_00492.1_g00008 was homologus to Arabdopsis MYB124, which was down-regulated in EP26 at both stages (Fig. 5b).
Validation of differential gene expressed genes by qRT-PCR
Based on the RNA-Seq results, gene expression of thirteen selected DEGs was further validated using qRT-PCR with gene-specific primers (Table S3). These genes included anthocyanin biosynthetic genes and regulatory genes, as well as genes involved in photosynthesis. A linear regression analysis showed an overall correlation coefficient of R2=0.7501, indicating the results of qRT-PCR were consistent with those of the transcriptome analysis (Figure 6).