Source: https://hereditasjournal.biomedcentral.com/articles/10.1186/s41065-016-0021-1
Timestamp: 2019-04-20 02:58:08+00:00

Document:
The color of berry skin is an important economic trait for grape and is essentially determined by the components and content of anthocyanins. The fruit color of Chinese wild grapes is generally black, and the profile of anthocyanins in Chinese wild grapes is significantly different from that of Vitis vinifera. However, V. davidii is the only species that possesses white berry varieties among Chinese wild grape species. Thus, we performed a transcriptomic analysis to compare the difference of transcriptional level in black and white V. davidii, in order to find some key genes that are related to anthocyanins accumulation in V. davidii.
The results of anthocyanins detection revealed that 3,5-O-diglucoside anthocyanins is the predominant anthocyanins in V. davidii. It showed obvious differences from V. vinifera in the profile of the composition of anthocyanins. The transcriptome sequencing by Illumina mRNA-Seq technology generated an average of 57 million 100-base pair clean reads from each sample. Differential gene expression analysis revealed thousands of differential expression genes (DEGs) in the pairwise comparison of different fruit developmental stages between and within black and white V. davidii. After the analysis of functional category enrichment and differential expression patterns of DEGs, 46 genes were selected as the candidate genes. Some genes have been reported as being related to anthocyanins accumulation, and some genes were newly found in our study as probably being related to anthocyanins accumulation. We inferred that 3AT (VIT_03s0017g00870) played an important role in anthocyanin acylation, GST4 (VIT_04s0079g00690) and AM2 (VIT_16s0050g00910) played important roles in anthocyanins transport in V. davidii. The expression of some selected DEGs was further confirmed by quantitative real-time PCR (qRT-PCR).
The present study investigated the transcriptomic profiles of berry skin from black and white spine grapes at three fruit developmental stages by Illumina mRNA-Seq technology. It revealed the variety specificity of anthocyanins accumulation in V. davidi at the transcriptional level. The data reported here will provide a valuable resource for understanding anthocyanins accumulation in grapes, especially in V. davidii.
Grapevine (Vitis L.) is an important economic crop used as table fruit, dried raisins, and for wines or juice. The color of berry skin is an important economic trait for grapes. And it is closely related to the components and content of anthocyanins [1, 2]. In recent years, the profile of anthocyanins in grape berry skin has been widely studied. The anthocyanins are glycosides and acylglycosides of anthocyanidins. The main anthocyanidins found in grapes are pelargonidin, cyanidin, delphinidin, peonidin, petunidin as well as malvidin which is usually the predominant anthocyanidin in most red grapes [3, 4]. However, the kind of anthocyanins exhibits obvious differences in different grape species . In the V. vinifera cultivars, only 3-O-monoglucoside anthocyanins are detected, and just a few cultivars produce pelargonidin-based anthocyanins [6, 7]. In non-V. vinifera species, 3,5-O-disglucoside anthocyanins widely exist, even the pelargonidin-based anthocyanins are also detected [7, 8]. In addition, most V. vinifera cultivars possess acylated forms of the anthocyanins, but some other grape species do not produce acylated forms of anthocyanins, such as V. rotundifolia .
Anthocyanins in grape berry skin are synthesized via the flavonoid pathways, which has been extensively studied. In fact, the synthesis of anthocyanins shares the same upstream pathways with proanthocyanidins and flavonol derivatives . Phenylalanine is the precursor of flavonoid, which is used as substrate, phenylalanine ammonia lyase (PAL), cinnamate-4-hydroxylase (C4H) and 4-coumaroyl-CoA synthase (4CL) catalyze a series of reactions to produce 4-coumaroyl-CoA. The catalysis of chalcone synthase (CHS) is the first committed step in the flavonoid biosynthetic pathway, which can catalyze the synthesis of chalcones . Subsequently, after the action of chalcone isomerase (CHI), the basic three rings of the general C6-C3-C6 flavonoid skeleton is produced . The B ring of the naringenin flavanone can be further hydroxylated by flavonoid -3′-hydroxylase(F3′H) or flavonoid-3′5′-hydroxylase (F3′5′H) to form eriodictyol or pentahydroxyflavanone [12, 13]. The naringenin, eriodictyol and pentahydroxyflavanone can be modified by the catalysis of flavanone-3β-hydroxylase (F3H) to form the corresponding dihydrokaempferol, dihydroquercetin and dihydromyricetin, respectively. Besides, the dihydrokaempferol can also be catalyzed by F3′H or F3′5′H to produce other two dihydroflavonols, dihydroquercetin or dihydromyricetin . Then, dihydroflavonol-4-reductase (DFR) catalyzes these dihydroflavonols to form their corresponding leucoanthocyanidins . In the past, leucoanthocyanidin dioxygenase (LDOX) used to be considered as the first key enzyme that could catalyze the formation of anthocyanidins [14, 15] and lead the flavonoid flux into the anthocyanin branch. However, more studies showed that LDOX also play an important role in the biosynthesis of proanthocyanidins [16–19]. And the anthocyanidins which are the production of LDOX can also be catalyzed by anthocyanidin reductase (ANR) to produce the substrates for the proanthocyanidins synthesis [20, 21]. The glycosylation, methylation and acylation of anthocyanidins are very important for their stabilization. Glycosylation of the anthocyanidins is the key step to produce anthocyanins which is catalyzed by the UDP-glucose: anthocyanidin: flavonoid glucosyltransferase (UFGT) [1, 22–25]. In V. vinifera, the anthocyanidins can only be catalyzed by UFGT to glycosylate at C3 position . Thus, it is also called 3GT. The S-adenosyl-L-methionine (SAM) or O-methyltransferase (OMT) can catalyze the methylation of the hydroxyl groups at the C3 positions or both at the C3 and C5 positions on the B rings of the anthocyanins [26, 27]. Acylation can greatly enhance the structural diversity and stability of anthocyanins . It is catalyzed by the action of anthocyanin acyltransferases (AAT) [3, 29]. There are mainly two types of ACTs that are classified based on the acyl group donors: the BAHD family using acyl-CoA and the serine carboxypeptidase-like(SCPL) group using acyl-activated sugars [30, 31]. Anthocyanins are synthesized in the cytoplasm but accumulate in the vacuoles. Recent researches showed that glutathione S-transferase (GST), multidrug resistance-associated protein (MRP) and multidrug and toxic compound extrusion (MATE) were closely related to the transport of anthocyanins .
As part of the flavonoid pathway, the synthesis of anthocyanins are regulated by a complex regulation at the transcriptional level [4, 10, 15, 33–35]. Generally, the flavonoid pathway of anthocyanins biosynthesis is under the control of Myb transcriptional factors, basic helix–loop-helix proteins(bHLH) and WD40-like proteins, which also play crucial roles in the regulation of flavonols and proanthocyanidins [10, 17–19, 36–38]. The research showed that the bHLH gene, VvMYC1, was characterized as a component of the transcriptional complex regulating anthocyanins biosynthesis in grapevine . In grapes, a series of R2R3-Myb transcriptional factors are related to the synthesis of anthocyanins . The first Myb transcriptional factor Myba in grape was identified and isolated from V. labrusca hybrids , and the results suggested that VlmybA1-1, VlmybA1-2 and VlmybA2 transcriptional factors are involved in the regulation of anthocyanin biosynthesis in the grape by regulating the expression of the UFGT gene [41–43]. In red V. vinifera grapes, the functional VvmybA1 can regulate the expression of the UFGT gene to promote the synthesis of anthocyanin. However, in white V. vinifera grapes, a retrotransposon (Gret1) is inserted in the 5′-flanking region of VvmybA1 gene to form a non-functional VvmybA1a gene, resulting in the transcription factor losing its function . Subsequently, two other VvmybA regulator genes, VvmybA2 and VvmybA3 have been cloned and identified . The VvmybA1 and VvmybA2 are very similar, both of them can regulate the accumulation of anthocyanins in the grape berries, and the white-fruited grapes are caused by the mutation of these two similar and adjacent regulatory genes [46, 47]. More recent researches suggest that variation in anthocyanins content in grapes is involved with the VvmybA gene cluster [48–50].
China is the concentrated distribution area of East Asian Vitis species, which has more than 35 Vitis species [51–53]. In recent years, the research of anthocyanins in Chinese wild grapes has been carried out, and the profile of anthocyanins in Chinese wild grapes is significantly different from V. vinifera [5, 7, 54–56]. Generally, only 3-O-monoglucoside anthocyanins are detected in V. vinifera cultivars, and almost all of the cultivars of V. vinifera are devoid of pelargonidin-based anthocyanins. However, in Chinese wild grapes (such as V. amurensis, V. davidii and V. quinquangularis), the 3,5-O-disglucoside anthocyanins widely exist, and pelargonidin-based anthocyanins are detected in some wild species [5, 7, 55, 56]. The berry color of most Chinese wild grapes is black, except for the white-fruited varieties of V. davidii. The subsequent research found that VvmybA1a gene was detected in white-fruited varieties of V. davidii . The white fruit may be caused by this mutation of mybA1 gene. In order to find the causes that lead to the different profiles of anthocyanins between Chinese wild grapes and V. vinifera, we used the black and white varieties of V. davidii as the materials, and the method of transcriptome analysis to detect the differences at the transcriptional level that might be involved in anthocyanins accumulation of these two varieties.
A total of 24 kinds of anthocyanins were detected at three fruit developmental stages of black and white spine grapes (Table 1). Five categories of anthocyanins were detected: 11 kinds of anthocyanidin diglucosides and 13 kinds of anthocyanidin monoglucosides. And 14 kinds of coumaroylated anthocyanins were detected. The highest content of anthocyanins was detected in the third stage black spine grape (B3, DAF120), it was 1382.127 mg.kg−1. In black spine grape, the major anthocyanin was Malvidin-3,5-O-diglucoside, which accounted for 87% of the total anthocyanins. Besides, we also detected trace amount of anthocyanins in white V. davidii.
A total of 12 cDNA libraries were constructed from the total RNA of black and white spine grape berry skin in three fruit developmental stages (Fig. 1), two biological replicates were made at each stage. These cDNA libraries were subjected to pair-end reading with the Illumina HiSeq 2000 platform, generating from 56 million to 64 million pair-end raw reads of 100 bp in length, respectively. And the raw data have been submitted to SRA database of NCBI (Accession ID: SRP070860; Link: https://trace.ncbi.nlm.nih.gov/Traces/sra/?study=SRP070860). After removing the low-quality reads and trimming the adapter sequences, we obtained from 27 million to 31 million clean reads (Table 2). The subsequent analyses were based on the clean data. Though only two biological replicates were used for RNA-seq, the correlation analysis showed that the R2 of the two replicates of all the samples were greater than 0.95. It fully illustrated the consistency of two biological replicates and the reliability of the RNA-seq results.
We used the French-Italian Public Consortium for Grapevine Genome Characterization, publically accessible version of the complete V. vinifera genome at 12× coverage (ftp://ftp.ensemblgenomes.org/pub/release-23/plants/fasta/vitis_vinifera/dna/) , as the reference genome. After quality control, the reads from black spine grape (64–74%) and white spine grape (63–76%) successfully aligned to the reference genome. Most of the reads from each fruit developmental stage for black and white spine grape aligned to a single position. These uniquely mapped reads account on average for approximately 67 and 71% of the total number of sequenced reads for the black and white spine grape, respectively. The number of reads from each fruit developmental stage for both black and white spine grape mapped to ‘+’ and ‘–’ strand were mostly equal. The number of non-splice reads was approximately twice to splice reads (Table 3).
Approximately 80% of the reads from each fruit developmental stage for both black and white spine grape mapped to exons, 20% mapped to intergenic regions, and 1% mapped to introns (Fig. 2). The reads mapped to introns resulted from the residual of pre-mRNA and intron retention in the process of alternative splicing. However, the incomplete genome annotation led to the reads mapped to intergenic regions.
We used the RPKM value to show the gene expression level (Table 4). About 45–50% of the genes fall into the RPKM range of 0–1, 9–10% of the genes fall into the RPKM range of 1–3, 18–20% of the genes fall into the RPKM range of 3–15, 14–17% of the gens fall into the RPKM range of 15–60, and approximately 7% of the genes that RPKM value is greater than 60.
There were thousands of differential expression genes in three fruit developmental stages. But, we focused on the genes which were related to the change of berry skin color in fruit development. Thus, we compared the expression of genes in different stages, such as B2vsB1, B3vsB2, B3vsB1, B1vsW1, B2vsW2 and B3vsW3. Then, we used the differential expression genes which were gained from the comparison to make venn diagrams (Fig. 3). In the venn diagram of B2vsB1, B3vsB2 and B3vsB1, 1440 differential expression genes were found. And in the venn diagram of B1vsW1, B2vsW2 and B3vsW3, it was 1006. Accordingly, the scope of selecting the candidate differential expression genes is effectively narrowed.
Gene Ontology (GO, http://www.geneontology.org/) is an international standardized gene function classification system which can be used to classify the function of the predicted genes. GO include three parts: molecular function, biological process and cellular component. We used the differential expression genes between B2 and W2 to make GO enrichment analysis. Thirty significant terms were found, which belonged to biological process and molecular function (Fig. 4). The color of the berry skin is determined by the pigment, and the anthocyanins are the main pigment in grapes [1, 2]. The anthocyanins were synthesized in the endoplasmic reticulum, and transported into vacuole by a serious of anthocyanins transporters . Thus, in these 30 terms, transmembrane transport (GO:0055085) and pigment catabolic process(GO:0046149), which possibly contained the candidate genes, came into our attention. The transmembrane transport term contained 107 up-regulated expression genes and 110 down-regulated expression genes. The pigment catabolic process contained 6 up-regulated expression genes and 14 down-regulated expression genes (Additional file 1: Table S1).
To identify the biological pathways activated in grape fruit, pathway enrichment was quantified in the KEGG database . In grapes, anthocyanins synthesis is via the flavonoid pathways , which were significantly enriched in our study (Fig. 5). There were 33 background genes in the flavonoid biosynthesis pathway, 27 genes showed expressional differences (Table 5). In these 27 genes, VIT_00s0361g00040 was annotated as anthocyanidin reductase (ANR), VIT_01s0011g02960 and VIT_17s0000g04150 were annotated as leucoanthocyanidin reductase(LAR), which were related to the synthesis of proanthocyanidins [20, 21]. VIT_18s0001g03470 was annotated as flavonol synthase (FLS), which was related to the synthesis of flavonol aglycones [60, 61]. VIT_02s0033g00410 was annotated as transcription factor MYBA1, which can effectively regulate the anthocyanins biosynthesis [41–43]. VIT_02s0025g04720 and VIT_16s0039g02230 were respectively annotated as LDOX and UFGT, which were directly related to anthocyanins biosynthesis [14, 15, 22–25]. VIT_03s0063g00140, VIT_07s0031g00350, VIT_11s0016g02610 and VIT_12s0028g03110 were annotated as OMT, which participated in the methylation modification of anthocyanins [26, 27]. Other genes were related to the upstream pathways of the flavonoid biosynthesis.
We also used differential expression gene clustering methodology to find the genes that are related to the biosynthesis of anthocyanins in spine grape. The genes that showed same expression pattern were clustered together, and they may obtain the same function or belong to the same biological pathway. We used all of the differential expression genes to draw a heatmap (Fig. 6a). The genes that had the same or similar expression patterns were effectively clustered together. However, so many clusters were produced, and it was difficult to find the clusters that contained the genes closely related to the color of the berry skin. But, the previous studies have suggested that UFGT(VIT_16s0039g02230) [1, 24, 25] and MybA1(VIT_02s0033g00410)  were the key genes to regulate the synthesis of the anthocyanins. Thus, we selected the genes whose expression patterns were similar to these two genes. Then, we obtained 100 genes and used them to make a heatmap (Fig. 6b). The genes were divided into two groups that contained 68 and 32 genes respectively. And the expression pattern of the group that contained 68 genes was very similar to UFGT and MybA1. Then, we can select the genes that were related to the biosynthesis and transport of anthocyanins.
After these analyses, we finally obtained 41 differential expression genes, which were considered to be related to the biosynthesis and transport of anthocyanins in the berry skin of spine grape (Table 6). The expression patterns of these genes were very similar, mainly expressed in B2 and B3 (Fig. 7). In these candidate genes, just one gene (VIT_14s0128g00600) was not annotated. Yet its expression pattern was very similar to the other candidate genes, and we inferred that it was probably related to the anthocyanins accumulation.
Trace amount of anthocyanins were detected in white spine grape (Table 1). Thus, we inferred that some other UGTs are probably related to anthocyanins accumulation in white spine grape. We compared the third stage of white spine grape with the first stage, 62 significant differential expression genes were annotated as UGTs (Additional file 2: Table S2), containing 20 up-regulated and 42 down-regulated genes. From the 20 up-regulated genes, we selected the most 5 significant differential expression genes as the candidate genes (Table 7). All of these five candidate genes were above 4.5-fold higher in W3 than in W1. Especially, VIT_00s0324g00050 showed highly significant differential expression (P = 3.79E-129).
To confirm the results obtained from RNA-seq, relative expression profiles of 30 genes were analyzed by real-time RT-PCR at the three fruit developmental stages in black and white spine grapes. These 30 genes related to the biosynthesis and transport of anthocyanins were chosen for qRT-PCR analysis. The black spine grape showed much higher gene expression of the selected genes than that of the white spine grape. For most of the genes, the qRT-PCR results were consistent with those obtained from the expression profile determined from RNA-Seq data (Fig. 8).
In V. vinifera, only 3-O-monoglucoside anthocyanins were detected . However, in our research, 3,5-O-diglucoside anthocyanins was the predominant anthocyanins in V. davidii. The studies have shown that, the 3,5-O-diglucoside anthocyanins were usually detected in other Chinese wild grapes [5, 7, 55, 56]. Therefore, the V. davidii was closer to the other Chinese wild grapes from the perspective of the composition of anthocyanins. We inferred that the high 3-5-O-disglucoside anthocyanins concentration in black spine grape was closely related to the expression of VIT_09s0002g06590 (5GT).
Interestingly, we also detected trace amount of anthocyanins in white V. davidii. It was generally accepted that anthocyanins were present only in red grapes, and it was used as a standard to define the difference between red grape and white grape [46, 62]. But, a very recent research suggested that some white grapes contained measurable traces of anthocyanins . Previous researches showed that UFGT gene was not expressed in white grape berries, and the expression of UFGT was regulated by two very similar regulatory genes, VvMYBA1 and VvMYBA2, which were not transcribed in white grape berries [46, 62]. The present study revealed that UFGT (VIT_16s0039g02230) did not express in white V. davidii. However, why did we detect anthocyanins in white grape berries? A previous study showed that there were as many as 181 putative UDP-glycosyltranferases (UGTs) found in the genome of V. vinifera . Arapitsas et al.  detected measurable trace amounts of anthocyanins in some white grape varieties, and they inferred that the other UGTs expressed in berry skin accepts anthocyanidin as substrate and is therefore involved in the synthesis of trace amount anthocyanins detected in the berry skin of white grape. In our study, 5 candidate UGT genes were above 4.5-fold higher in W3 than in W1, especially VIT_00s0324g00050 showed highly significant differential expression (P = 3.79E-129). These genes are probably related to anthocyanins accumulation in white V. davidii.
The synthesis of anthocyanins via the phenylpropanoid and flavonoid pathway, and it is catalyzed by a serious of enzymes (Fig. 9) [1, 10, 37]. However, most of the genes coding for these enzymes are multi-copied in the grape genome . Phenylalanine ammonia lyase (PAL) catalyzes the first step of the phenylpropanoid pathway. In our study, two differential expression genes (VIT_06s0004g02620 and VIT_13s0019g04460) were annotated as PAL. VIT_06s0004g02620 was mainly expressed in B2 and B3, but just had a trace amount expression in W3. VIT_13s0019g04460 was also mainly expressed in B2 and B3, and had a high expression in W3. Thus, we inferred that VIT_06s0004g02620 played a major role in the pathway. Under the action of C4H and 4CL, 4-Coumaroyl-CoA was produced. We did not find any obvious differential expression genes annotated as C4H. However, VIT_16s0039g02040 and VIT_16s0050g00390 were annotated as 4CL. Chalcone synthase (CHS) is the first committed enzyme of the flavonoid pathway. Previous studies suggested that 3 CHS genes (CHS1, CHS2 and CHS3) played a role in this pathway. Of the CHSs, the mRNAs of CHS1and CHS2 were detected in both leaves and berry skins of white and red grape cultivars, whereas the mRNA of CHS3 was mainly accumulated in the berry skin of red cultivars during coloration [65–67]. In our study, only CHS2 (VIT_14s0068g00920) and CHS3 (VIT_05s0136g00260) showed significant expression differences between black and white spine grapes. Especially, CHS3 (VIT_05s0136g00260) was detected as having a very high expressional level in B2 and B3, and significantly higher than that in white spine grape. It was consistent with the previous studies [65–67]. Two CHI genes (VIT_13s0067g03820 and VIT_13s0067g02870) showed obvious differential expression levels between black and white spine grapes. These two CHI genes (CHI1 and CHI2) have been reported by previous researches [65, 67, 68]. The researches suggested that two F3H genes (F3H1 and F3H2) played a role in the flavonoid pathway, and the mRNA of F3H2 was detected at a high level in grape berry skin during coloration [65, 67]. However, in our study, only the expression of F3H1 (VIT_04s0023g03370) in black spine grape showed significantly higher than white spine grape. Maybe F3Hs play different roles in the flavonoid pathway between V. vinifera and V. davidii. F3′H and F3′5′H belong to the cytochrome P450 super family, and catalyze hydroxylation at the 3′ and 3′,5′ positions of the B-ring of the flavonoid to produce the precursors for cyanidin-based anthocyanins and delphinidin-based anthocyanins, respectively. Thus, components of anthocyanins in grape berry skins are closely related to the expression of F3′H and F3′5′H [62, 69]. The cDNAs of F3′H and F3′5′H were first isolated from petunia [70, 71]. The previous studies showed that grapevines contain two copies of F3′H and sixteen copies of F3′5′H, both of the F3′Hs are located on chr17, fifteen of F3′5′Hs located on chr6, and one of F3′5′Hs located on chr8 . The previous study showed that, the mRNA levels of F3′H and F3′5′H were high in grape berry skins at the harvest stage . In our study, two (VIT_11s0016g01020 and VIT_17s0000g07200) and four (VIT_06s0009g02830, VIT_06s0009g03010, VIT_06s0009g02810 and VIT_06s0009g02920) differential expression genes were annotated as F3′H and F3′5′H, respectively. Of the two F3′Hs, one F3′H gene (VIT_17s0000g07200) was located on chr17. However, the other F3′H gene (VIT_11s0016g01020) was located on chr11. In addition, VIT_11s0016g01020 showed more significantly differential expression level than VIT_17s0000g07200 between black and white grapes. We suggested that VIT_11s0016g01020 (F3′H) played a major role in anthocyanins accumulation in black spine grape. All of the four F3′5′H genes showed high expression level in B2 and B3, especially VIT_06s0009g02830 and VIT_06s0009g03010. But F3′5′H genes were almost not expressed at any of the three stages in white spine grape. Besides, we could see that the expression level of F3′5′H was obviously higher than F3′H in black spine grape. This led to more pentahydroxy-flavone and leucodelphinidin production, which were substrates of malvidin-based anthocyanins. Thus, we detected malvidin-based anthocyanins as the predominant kind of anthocyanins in black spine grape. The dihydroflavonols were catalyzed by DFR and LDOX to produce anthocyanidin. We found the differential expression genes (VIT_18s0001g12800 and VIT_02s0025g04720) annotated as DFR and LDOX, which was consistent with previous studies [67, 68, 74, 75].
It has been suggested that there was a mybA gene cluster on chromosome 2 of the grape, which contained MYBA1, MYBA2 and MYBA3 . The expression of MYBA1 and MYBA2 could promote the synthesis of anthocyanins through regulating the expression of UFGT [41, 46, 47]. In our study, MYBA1 (VIT_02s0033g00410) showed high expression levels in B2 and B3, but it did not express in white spine grape. MYBA2 (VIT_02s0033g00390) and MYBA3 (VIT_02s0033g00450) showed high expression levels in B2, B3 and W3. Jiao et al.  detected Gret1 in white-fruited varieties of V. davidii. Thus, it was reasonable that MYBA1 did not express in white spine grape. MYBA2 could also promote the expression of UFGT . Though MYBA2 showed high expression level in W3, UFGT (VIT_16s0039g02230) almost did not express and just trace amount anthocyanins were detected in white spine grape. We could infer from this observation that MYBA1 may play a major role in regulating the synthesis of anthocyanins in spine grapes. Besides, another R2R3-MYB gene VIT_01s0011g04760 was annotated as MYB4 expressed both in black and white spine grapes. But it showed significantly differential expression between them. The expression level of VIT_01s0011g04760 in black spine grape was obviously higher than in white spine grape. Thus, it can be regarded as a new candidate gene and needs further studies to confirm whether this gene regulates the biosynthesis of anthocyanins in grapes.
The flavonoid biosynthesis pathway is under the control of Myb transcriptional factors, basic helix–loop-helix proteins (bHLH) and WD40-like proteins [7, 10, 17–19, 36, 38]. The first bHLH was submitted as VvMYCA1 (accession number EF193002; gene ID VIT_15s0046g02560) [76, 77], which was considered to regulate the expression of UFGT. Subsequently, the second bHLH was called VvMYC1 (accession number EU447172; gene ID VIT_07s0104g00090) , which was characterized as a component of the transcriptional complex regulating anthocyanin biosynthesis in grapevine. However, both of these two genes did not show obvious differences at the transcriptional level between black and white spine grapes in our study. But, another gene (VIT_14s0060g01010) annotated as bHLH exhibited differential expression between black and white spine grapes. It was shown that WDR1 (Genbank accession number DQ517914) contributed positively to the accumulation of anthocyanins . However, we did not find any WD40 gene having obvious differential expression. In recent years, researches showed that NAC TFs were involved in the regulation of anthocyanins accumulation. A NAC TF has been proposed to be involved in the regulation of anthocyanins accumulation during the response of blood orange to cold exposure . In addition, PpNAC1 can activate the transcription of PpMYB10.1, resulting in anthocyanins pigmentation in blood-fleshed peach . In grape, there was no report that the NAC TFs were related to the anthocyanins accumulation. In our study, a candidate gene (VIT_14s0108g01070) was annotated as NAC, which may be related to anthocyanins accumulation in spine grapes.
Anthocyanidins need to be modified by glycosylation, methylation and acylation to form stabilized anthocyanins. In V. vinifera, the glycosylation was catalyzed by UDP-glucose: anthocyanidin: flavonoid glucosyltransferase (UFGT) at C3 position, and the gene encoding UFGT had been cloned [24, 80]. The action of UFGT was crucial for anthocyanin accumulation in grape berry skin [1, 22, 23, 25, 50, 67, 81, 82]. In our study, a differential expression gene (VIT_16s0039g02230) was annotated as UFGT, which had a high expression level in B2 and B3, and did not express in white spine grape. However, in non-V. vinifera species, the 3-5-O-disglucoside anthocyanins are widely present [7, 8]. There must be another UDP-glucose: anthocyanidin: flavonoid glucosyltransferase to catalyze glycosylation at C5 position. In recent years, a few studies reported the isolation of 5GT genes in grapes [83, 84]. Jánváry et al.  cloned functional Cha5GT and nonfunctional Dia5GT from the heterozygous hybrid cultivar ‘Regent’, a cross of V. vinifera cv. ‘Diana’ and the interspecific hybrid cv. ‘Chambourcin’. The functional analysis of Cha5GT and Dia5GT suggested that two mutations in the 5GT gene eliminated its enzymatic activity. Because of the absence of active 5GT, dis-glucosidic anthocyanins could not be produced in V. vinifera red grapes. He et al.  cloned the full-length cDNA of UDP-glucose: anthocyanin 5-O-glucosyltransferase (Va5GT) from V. amurensis Rupr. cv. ‘Zuoshanyi’. The results suggested that Va5GT was a key enzyme in the biosynthesis of dis-glucosidic anthocyanins in V. amurensis grape berries. In our study, a differential expression gene (VIT_09s0002g06590) showed a high expressional level just in B2 and B3, which was annotated as UGT. It has high homology with Va5GT. Thus, we inferred that the expression of VIT_09s0002g06590 (5GT) led to a high 3-5-O-disglucoside anthocyanins concentration in black spine grape.
In plants, methylated anthocyanidins accounted for a large proportion of the total reported anthocyanidins . Especially, three methylated anthocyanidins, peonidin, petunidin and malvidin, were commonly present in grape berry skin. The methylation was catalyzed by S-adenosyl-L-methionine (SAM) or O-methyltransferase (OMT) at the C3 positions or both at the C3 and C5 positions on the B rings of the anthocyanins in grape [26, 27]. And the cDNAs of several OMT had been cloned in grapes [85, 86]. Subsequently, a QTL for anthocyanins methylation variation was identified that was colocalized with a cluster of three putative OMT genes: VIT_01s0010g03470 (OMT3), VIT_01s0010g03490 (OMT2) and VIT_01s0010g03510 (OMT1) [87, 88]. Fournier-Level et al.  reported that OMT2 gene presented two SNPs associated with methylation level. It probably led to a structural change of the OMT2 protein, thus significantly affected the enzyme specific catalytic efficiency for the 3′-O-methylation of delphinidin-3-glucoside. In this study, three differential expression genes (VIT_07s0031g00350, VIT_01s0010g03510 and VIT_01s0010g03490) were annotated as OMTs. Of the three genes, VIT_01s0010g03510 (OMT1) and VIT_01s0010g03490 (OMT2) expressions were consistent with the previous researches [87, 88], which did not express in white spine grape. It is interesting that OMT1 showed higher activity than OMT2 in black spine grape. In addition, a new candidate gene (VIT_07s0031g00350) was also annotated as OMT, and showed high expressional level in B2 and B3. It probably played an important role in the methylation of anthocyanins in spine grape berry skin.
Acylation was a common modification of anthocyanins in grapes. It not only increased diversity of anthocyanins, but also improved the color stabilization and intensity for the anthocyanins . Two enzyme families (BAHD-ATs and SCPL-ATs) have been reported as related to the acylation of anthocyanins [30, 31]. A recent study identified a number of QTLs associated with variation in acylated anthocyanin levels in F1 progeny from a ‘Syrah’x‘Pinot Noir’ cross. The strongest candidate genes within these QTLs included those belonging to the BAHD and SCPL acyltransferase family . But, no QTL was found to cause the presence/absence of acylation in berries. In our study, two differential expression genes (VIT_14s0068g01440 and VIT_03s0017g00870) were annotated as acyltransferase (AT). Especially VIT_03s0017g00870 did not express in white spine grape and at the first stage of black spine grape. Yet it showed high expressional level in B2 and B3. Simultaneously, we detected acylated anthocyanins in black spine grape. Thus, this candidate gene can be very possibly related with the acylation of anthocyanins. It is known that V. vinifera cv. Pinot Noir does not synthesise acylated anthocyanins . After SNP analysis of the candidate genes, at the position 36 of the putative coding region of VIT_03s0017g00870, A is replaced by G, compared with the PN40024 grapevine genome , but it does not alter the predicted amino acid sequence. Rinaldo et al.  identified a gene, anthocyanin 3-O-glucoside-6-O-acyltransferase (Vv3AT), which encoded a BAHD acyltransferase protein. This protein can promote the synthesis of acylated anthocyanins, which was transcriptionally regulated by VvMYBA. This is consistent with our study.
Anthocyanins are synthesized in the endoplasmic reticulum and then transported into the vacuoles. Three kinds of anthocyanin transporters: glutathione S-transferase (GST), multidrug resistance-associated protein (MRP) and multidrug and toxic compound extrusion (MATE), were reported as related to the anthocyanins transport [32, 91]. In grapes, a few candidate anthocyanin transporters have been reported.
Conn et al.  cloned five GSTs (VvGST1,VvGST2,VvGST3,VvGST4 and VvGST5) from V. vinifera. cv. Gamay Fréaux and the study showed that VvGST1 and VvGST4 coded for the enzymes which had the function of anthocyanins transport. However, in our study, GST1 (VIT_19s0093g00320) was not detected, GST4 (VIT_04s0079g00690) and GST5 (VIT_19s0015g02730) showed significantly differential expression between black and white spine grape. Especially, GST4 just express in B2 and B3 (Fig. 10). GST2 (VIT_07s0005g00030) showed the highest expressional level in the first stage, and it possessed higher transcriptional level in white spine grape than black one. GST3 (VIT_12s0028g00920) also showed the highest expressional level at the first stage, but did not show obvious differential expression between black and white spine grapes. Thus, we inferred that GST4 (VIT_04s0079g00690) and GST5 (VIT_19s0015g02730) play an important role in anthocyanins transport in spine grapes. Besides, another two candidate genes (VIT_19s0015g02690 and VIT_19s0015g02880) were annotated as GSTs, which showed higher expressional level in black spine grape than in white one (Fig. 10). Hence, these two candidate genes are also probably related to anthocyanins transport in spine grapes.
Gomez et al.  identified three multidrug and toxic compound extrusion (MATE) genes as candidate anthoMATEs(AMs). In their study, just AM1 and AM3 were cloned from V. vinifera L. cv. Syrah, AM2 was not successfully cloned. They inferred that AM2 was probably not expressed at detectable levels in mature berry or was a pseudogene. Accordingly, their results revealed that AM1 and AM3 just could transport acylated anthocyanins in the presence of Mg ATP . Yet in our study, AM1 (VIT_16s0050g00900) was hardly expressed in spine grapes, the expressional level of AM3 (VIT_16s0050g00930) in white spine grape was higher than in black spine grape. However, AM2 (VIT_16s0050g00910) was mainly expressed in B2 and B3 (Fig. 11). Thus, we infer that AM2 plays an important role in spine grape anthocyanins transportation.
The plant ATP binding cassette (ABC) transporters, in particular from the ABCC subfamily (formerly named multidrug resistance proteins [MRPs]) had been reported as related to the accumulation of flavonoid in vacuoles . In grapes, an ABC transporter, ABCC1 was identified, which localizes to the tonoplast and is involved in the transport of glucosylated anthocyanidins . In our study, four candidate genes were annotated as ABC transporters: VIT_16s0050g02480 (ABCC1), VIT_10s0003g04390 (ABCC2), VIT_09s0002g02430 (ABCC8) and VIT_03s0017g01290 (ABCG11). Francisco et al.  were unsuccessful to clone the ABCC2 from the exocarp of grape, and inferred it might not be expressed in detectable amounts during the ripening stage. However, our study showed that both ABCC1 and ABCC2 were expressed in black and white spine grape berry skins. All of the four candidate ABCs were higher in expression in black spine grape than in white spine grape, but the differences were not very significant (Fig. 12). Thus, more studies are needed to confirm the key ABC transporters in grapes.
V. davidii is the only Chinese wild grape species which possesses white berry varieties. High levels of 3,5-O-diglucoside anthocyanins were detected in the black berry skin of V. davidii. The present study investigated the transcriptome profiles of the berry skin from black and white spine grapes at three fruit developmental stages by Illumina mRNA-Seq technology. The examination of absolute expression count for every gene has enabled us to carry out a global investigation of gene expression at these three key time-points in black and white spine grapes. The transcriptome analysis presented thousands of DEGs. We used gene clustering and the enrichment of GO and KEGG to describe the transcriptional patterns of genes involved in anthocyanins accumulation. We found 41 differentially expressed genes probably related to anthocyanins accumulation in V. davidii, including the genes that encode enzymes, transcription factors and transporters involved in anthocyanins biosynthesis, regulation and transport. Some genes were consistent with the previous studies in other grape species, some were newly found in our study. Trace amounts of anthocyanins were detected in berry skins of white V. davidii. Five candidate UGT genes were probably related to the biosynthesis of anthocyanins in white V. davidii. In conclusion, the present study provides new insights into the understanding of anthocyanins accumulation in grapes.
Berries of black and white spine grapes (V. davidii) were used in this study. Both of the black and white spine grapes were collected from the wild in Hunan, China. The plants were grown in China National Germplasm Resources Repository of Grape, Zhengzhou Fruit Research Institute, Chinese Academy of Agricultural Sciences. Skins of grape berries were collected at three different fruit developmental stages: before veraison (Veraison refers to the stage that the berries begin to color and be soft) (40 days after flowering [DAF]), at veraison (80 DAF) and fruit ripe period (120 DAF) (Fig. 1). Two biological replicates were made at each stage. At each stage, we harvested 6–8 clusters at each sampling date, and 5 berries were collected from each cluster. Thus, 30–40 berries were used to collect the skin at each sampling date. After collection, samples were flash-frozen in liquid nitrogen and stored at −80 °C until further processing.
Extraction of total anthocyanins was as described by He et al.  with some modifications. The peels were put into mortars and ground in liquid nitrogen. Aliquots of 0.5 g ground powder were added into 10 ml centrifuge tubes with 8 ml 2% formic acid-methanol solution. After ultrasonic oscillation for 10 min, the extracts were put on a shaker in the dark in the table concentrator at 25 °C, 200 rpm for 30 min, followed by centrifugation at 4 °C, 12000 rpm for 10 min. The supernatants were transferred into 50 ml centrifuge tubes. The residues were re-extracted 3 times. The organic fractions were pooled, evaporated by a vacuum rotary evaporator (BUCH, USA) at 40 °C. The residual parts were poured into activated solid phase extraction cartridges. The solid phase extraction cartridges were washed with 5 ml water for 2 times. After removing the leacheate, the solid phase extraction cartridges were eluted with 10 ml methanol for 2 times. The filtrates were collected and evaporated to dryness, then re-dissolved in 5 mL 0.5% hydrochloric acid- methanol solution. Finally, the solutions were filtered through a 0.22 μm Millipore filter for analysis.
The anthocyanins content was determined by an ACQUITY Ultra Performance Liquid Chromatography (UPLC) system (Waters, Milford, MA, USA) linked to both a PDAeλ photo diode array detector (Waters, Milford, MA, USA) and a Micromass Quattro microTM API benchtop triple quadrupole mass spectrometer (Waters MS Technologies, Manchester, UK), with a electrospray ionization (ESI) source operating in multiple reaction monitoring (MRM) mode. Sample solutions were injected into a ACQUITY UPLC®HSS T3 column (2.1 × 150 mm i.d,with 1.8 μm particle size, Waters, Milford, MA, USA), which was maintained at 30 °C. The mass spectrometric acquisition parameters were as follows: ESI source temperature 150 °C, desolvation gas temperature 400 °C, desolvation gas flow rate 800 L/h, cone gas flow rate 50 L/h, collision gas (high purity argon gas) flow rate 0.14 mL/min. The mobile phase had acetonitrile as solvent A, and 0.5% hydrochloric acid solution as solvent B. The gradient profile began with 5–10% A for 1 min, 10–25% A for 16 min, 25–40% A for 18 min, 40–100% A for 19 min, and then returned to initial conditions for 20 min and for 5 min. The flow rate was 1.0 ml min−1 and the column temperature was set at 40 °C. The injection volume was 2.0 μl. The detection wavelength was 520 nm.
Total RNA was extracted using TIANGEN RNAprep Pure Plant Kit (Tiangen Biotech Beijing, China). RNA degradation and contamination was monitored on 1% agarose gels. RNA purity was checked using the NanoPhotometer® spectrophotometer (IMPLEN, CA, USA). RNA concentration was measured using Qubit® RNA Assay Kit in Qubit® 2.0 Flurometer (Life Technologies, CA, USA). RNA integrity was assessed using the RNA Nano 6000 Assay Kit of the Bioanalyzer 2100 system (Agilent Technologies, CA, USA).
A total amount of 3 μg RNA per sample was used as the input material. Sequencing libraries were generated using NEBNext® Ultra™ Directional RNA Library Prep Kit for Illumina® (NEB, USA) following manufacturer’s recommendations and index codes were added to attribute sequences to each sample. Briefly, mRNAs were purified from the total RNAs using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in NEBNext First Strand Synthesis Reaction Buffer (5X). First strand cDNA was synthesized using random hexamer primer and M-MuLV Reverse Transcriptase (RNaseH−). Second strand cDNA synthesis was subsequently performed using DNA Polymerase I and RNase H. In the reaction buffer, dNTPs with dTTP were replaced by dUTP. The remaining overhangs were converted into blunt ends via exonuclease/polymerase activities. After adenylation of 3′ ends of DNA fragments, NEBNext Adaptor with hairpin loop structure were ligated for hybridization. In order to select cDNA fragments of preferentially 150–200 bp in length, the library fragments were purified with AMPure XP system (Beckman Coulter, Beverly, USA). Then 3 μl USER Enzyme (NEB, USA) was used with size-selected, adaptor-ligated cDNAs at 37 °C for 15 min followed by 95 °C for 5 min before PCR. Then PCR was performed with Phusion High-Fidelity DNA polymerase, universal PCR primers and Index (X) Primer. The PCR products were purified (AMPure XP system) and library quality was assessed on the Agilent Bioanalyzer 2100 system.
The clustering of the index-coded samples was performed on a cBot Cluster Generation System using TruSeq PE Cluster Kit v3-cBot-HS (Illumia) according to the manufacturer’s instructions. After cluster generation, the library preparations were sequenced by an Illumina Hiseq 2000 platform and 100-base paired-end reads were generated.
Raw data (raw reads) of fastq format were first processed through in-house PERL scripts. In this step, clean data (clean reads) were obtained by removing reads containing adapter, reads containing ploy-N and low quality reads from raw data. At the same time, Q20, Q30 and GC content of the clean data were calculated. All the downstream analyses were based on the clean data with high quality. Reference genome and gene model annotation files were downloaded from genome website directly (ftp://ftp.ensemblgenomes.org/pub/release-23/plants/fasta/vitis_vinifera/ dna/). Index of the reference genome was built using Bowtie v2.0.6 and base paired-end clean reads were aligned to the reference genome using TopHat v2.0.9. We selected TopHat as the mapping tool since TopHat can generate a database of splice junctions based on the gene model annotation file and thus a better mapping result than other non-splice mapping tools.
HTSeq v0.6.1 was used to count the reads numbers mapped to each gene. And then RPKM of each gene was calculated based on the length of the gene and reads count mapped to this gene. RPKM (Reads Per Kilobase of exon model per Million mapped reads) was used for the effect of sequencing depth and gene length for the reads count at the same time, and is currently the most commonly used method for estimating gene expression levels .
Differential expression analysis of two conditions/groups (two biological replicates per condition) was performed using the DESeq R package (1.10.1). DESeq provided statistical routines for determining differential expression in digital gene expression data using a model based on the negative binomial distribution. The resulting P-values were adjusted using the Benjamini and Hochberg’s approach for controlling the false discovery rate. Genes with an adjusted P-value < 0.05 selected by DESeq were assigned as differentially expressed.
We used Gene Ontology (GO) and KEGG enrichment analysis of differentially expressed genes to select candidate genes. GO enrichment analysis of differentially expressed genes was implemented by the GOseq R package, in which gene length bias was corrected. GO terms with corrected P-value less than 0.05 were considered significantly enriched by differential expressed genes. KEGG is a database resource for understanding high-level functions and utilities of the biological system, such as the cell, the organism and the ecosystem, from molecular-level information, especially large-scale molecular datasets generated by genome sequencing and other high-through put experimental technologies (http://www.genome.jp/kegg/). We used KOBAS software to test the statistical enrichment of differential expression genes in KEGG pathways. We also used the differentially expressed genes clustering methodology to find the candidate genes. The clustering analysis was made based on the expression pattern of the differentially expressed genes. We selected the genes that the expression patterns were very similar to the genes closely related to the anthocyanins biosynthesis as the candidate genes.
Thirty genes were chosen for validation using qRT-PCR. Specific primer pairs for selected genes used in qRT-PCR were designed as shown in Additional file 3: Table S3. The cDNA was transcribed from 1 μg of total RNA using the TOYOBO ReverTra Ace qPCR RT Master Mix (TOYOBO, Japan) in 10 μL of reaction mixture. The qRT-PCR was performed with the Roche LightCycler480 Real-Time Detection System (Roche) with Roche SYBR Green I (Roche, USA). The thermal profile for SYBR Green I RT-PCR was 95 °C for 15 min, followed by 40 cycles of 95 °C for 10s and 55 °C for 30s and 72 °C for 30s. Each plate was repeated three times in independent runs for all reference and selected genes. The reference gene (Actin) was used for normalization. The comparative CT method (2-ΔΔCT method) was used to analyze the expression levels of the different genes .
Picard-tools v1.96 and samtools v0.1.18 were used to sort, mark duplicated reads and reorder the bam alignment results of each sample. GATK2 software was used to perform SNP calling.
All of the experiments analyzed using data comparisons were repeated three times. Statistical analyses were performed using variance (ANOVA) followed by Duncan’s new multiple range tests with SPSS version 17.0 (SPSS, Chicago, IL, USA). A significance level of p < 0.01 was applied.
We thank Novogene Bioinformatics Technology Co., Ltd (Beijing, China) for its service in RNA-Seq and transcriptome analysis. We are grateful to Hongmin Zhang from Temple University for his assistance with English editing of the manuscript.
This work was supported by grants from the National Natural Science Foundation of China (31372024), China Agriculture Research System (CARS-30) and Agricultural Science and Technology Innovation Program (CAAS-ASTIP-2015-ZFRI).
LS, XF, YZ and CL conceived and designed the study. JJ and HS collected berry skin, prepared RNA samples. LS and XF extracted and detected anthocyanins from berry skin. LS and YZ performed the transcriptome analysis and qRT-PCR. LS and CL wrote and revised the manuscript. All authors read and approved the final manuscript.

References: V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 
 V. 

V. 
 V. 
 V. 
 V. 
 V.