Genome-wide analysis of transcription factors related to anthocyanin biosynthesis in carmine radish (Raphanus sativus L.) fleshy roots

Quantification of anthocyanin levels in carmine radish

Five local cultivars of carmine radish (‘Hongxin 1’, ‘Guanguan’, ‘Longquan 1’, ‘Yanzhi 1’ and ‘Yanzhi 2’) containing natural red pigment (red radish pigment) collected from Fuling were selected as experimental materials. To identify the pigment contents of the five local cultivars of carmine radish, the dynamics of anthocyanin in the development of fleshy roots in carmine radish were investigated by HPLC analysis. The study groups included fleshy roots from the seedling stage (SS, 15 days after planting), initially expanded fleshy roots (IE, 40 days after planting), fully expanded fleshy roots (FE, 70 days after planting), fleshy roots from the bolting stage (BS, 120 days after planting), fleshly roots from the initial flowering stage (IFS, 140 days after planting), fleshy roots from the full-bloom stage (FBS, 160 days after planting), and fleshy roots from the podding stage (PS, 200 days after planting). Briefly, fleshy roots were collected from three homozygous individuals from the five local cultivars of carmine radish and pooled together; three replicates were included. The fleshy root tissues were ground in liquid nitrogen, and the red pigment was then extracted with a solvent mixture containing methanol (40%, v/v), formic acid (0.1%, v/v) and acetone (40%, v/v). We used a 10 μL injection volume on a VDS C-18 column (4.6 × 250 mm, 5 μm, VDS Optilab, Berlin, Germany) with a 0.8 mL min⁻¹ flow rate. Based on the results from HPLC analysis, RNA-Seq was used to investigate young fleshy roots in the carmine radish ‘Hongxin 1’ cultivar obtained from the dynamic fleshy root developmental stages, including fleshy roots from the SS, the initial expansion (IE) stage, the full expansion (FE) stage, the BS, the IFS, the FBS and the PS (Fig. S1). In this study, based on our private RNA sequencing data (PRJNA565866), transcriptome analysis and functional validation of putative DETFs related to the development of fleshy roots in ‘Hongxin’ carmine radish were conducted.

Sample collection and RNA isolation

The carmine radish ‘Hongxin 1’ was cultivated in a greenhouse at the experimental farm of Yihe (Yangtze Normal University experiment base) in 2018. First, we sowed seeds of ‘Hongxin 1’ in sterilized soil for 2 weeks under normal growth conditions (23 °C, 16 h light/8 h dark). Next, 2-week-old plants were transferred and kept for 15 days in the cold room (5 ± 1 °C, 12 h light/12 h dark) for vernalization treatment. After the vernalization period, the plants were grown in a normal growth room under normal growth conditions (23 °C, 16 h light/8 h dark). At least three independent biological replicates of fleshy roots obtained from the development stages of carmine radish ‘Hongxin 1’ were collected for qRT-PCR analysis of miRNAs and their related targets. All harvested tissues were immediately frozen in liquid nitrogen and stored at −80 °C for qRT-PCR. Subsequently, using the mirVana™ miRNA Isolation Kit (Ambion) and Trizol Reagent (Invitrogen, Nottingham, UK) kit according to the manufacturers’ instructions, small and total RNAs were isolated from each sample.

Gene annotation, enrichment analysis and cluster analysis of DETFs

Transcription factors were obtained from the available ‘Hongxin’ carmine radish transcription profile data. Subsequently, the TFs were subjected to GO annotation analysis using the Gene Ontology database (http://www.geneontology.org/) and to KEGG pathway enrichment analysis using KOBAS software (Xie et al., 2011). DETFs were then screened by noiseqbio (Tarazona et al., 2012) and identified using a corrected P-value < 0.05 between each set of compared samples (fold changes in the expression levels of genes in six radish cultivars, including ‘IE_root’, ‘FE_root’, ‘BS_root’, ‘IFS_root’, ‘FBS_root’ and ‘PS_root’, were identified by comparison with ‘SS_root’—‘IE_root’ vs ‘SS_root’, ‘FE_root’ vs ‘SS_root’, ‘BS_root’ vs ‘SS_root’, ‘IFS_root’ vs ‘SS_root’, ‘FBS_root’ vs ‘SS_root’ and ‘PS_root’ vs ‘SS_root’). DETFs were then clustered according to their expression patterns in the dynamic growth stages of carmine radish and plotted using the neighbor-joining cluster method through a homemade R script.

Expression pattern analysis and regulatory pathway identification

To better understand the regulatory networks involving the DETFs, we searched putative DETFs targeted by miRNAs collected from the literature published by psRNATarget (Sun et al., 2015), and the parameters were set as follows: a maximum expectation of 3.5 and a target accessibility (UPE) of 50. Negative correlations between DETs and miRNAs were validated by quantitative real-time (Dai, Zhuang & Zhao, 2018) PCR.

Validation of DETFs related to anthocyanin biosynthesis using qRT-PCR

To monitor all DETFs related to anthocyanin biosynthesis, eight differentially expressed transcription factors (DETFs, cluster_24061_AP2/ERF-ERF, cluster_24166_AP2/ERF-ERF, cluster_59539_GRAS, cluster_44972_B3, cluster_22927_NAC, cluster_10013_C2C2_Dof, cluster_5503_bHLH, cluster_12215_bHLH, cluster_8116_WRKY and cluster_16741_C3HH) were selected for qRT-PCR (ABI 7500 real-time PCR System, United States) based on their expression patterns (Table S1). In addition, eight DETFs related to miRNAs (cluster_10013_C2C2_Dof, cluster_23883_C2C2_Dof, cluster_17168_bHLH, cluster_16851_AP2/ERF-ERF, cluster_11321_AP2/ERF-ERF, cluster_19966_AP2/ERF-ERF, Cluster_6812_AP2/ERF-ERF, cluster_16203_AP2/ERF-ERF, and cluster_16851_AP2/ERF-ERF) were also selected for validation using qRT-PCR. The primers for the qRT-PCR experiments were designed using Primer 5.0 software, and the radish actin gene was used as the reference gene (Table S2). First-strand cDNA synthesis was performed with 1 μg of total RNA from ‘Hongxin 1’ carmine radish using an M-MLV reverse transcriptase (Promega, Madison, WI, USA). The amplification programs were run according to the standard protocol of the ABI7500 system in triplicate as described by Gao et al. (2015). The relative quantitative method (2^−ΔΔCT) was used to calculate the fold changes in the expression levels of target genes (Schefe et al., 2006).

Validation of the candidate miRNAs by qRT-PCR

To validate the putative miRNAs related to anthocyanin biosynthesis, eight miRNAs, including miR408a, miR172d, miR1432-5p, miR1425-5p, miR166g, miR827a, miR172b and miR165a-5p, were selected based on their corresponding DETFs. The primers for qRT-PCR experiments were designed using Primer 5.0 software, and the radish gene actin was used as a standard control (Table S2). First, using the One Step miRNA 1st cDNA Synthesis Kit (Shenggong, Chengdu, China), microRNA reverse transcription reactions were performed in an Eppendorf Mastercycler (Eppendorf North America, Westbury, NY, USA) for 60 min at 37 °C and for 5 min at 95 °C; the samples were then stored at 4 °C until further use. The RT-PCRs were performed in a 10 μL volume containing 1 μL of the diluted reverse transcription product, 1 µL of PCR buffer, 0.2 mM dNTPs, 2.0 U of EasyTaq DNA polymerase (TransGen Biotech, Beijing, China), 0.5 μM specific miRNA primer and 0.5 μM universal primer (5-TTACCTAGCGTATCGTTGAC-3) on an Eppendorf Mastercycler. The amplification programs were performed according to the standard protocol of the ABI7500 system in triplicate as described by Gao et al. (2015). The relative quantitative method (2^−ΔΔCT) was used to calculate the fold changes in the expression levels of target genes (Schefe et al., 2006).

Transcription factors related to anthocyanin biosynthesis in the development of fleshy roots of carmine radish: a summary

Based on GO and KEGG pathway analyses, TFs related to anthocyanin biosynthesis were identified in ‘Hongxin 1’ carmine radish. A total of 1,747 TFs in 64 TF families related to anthocyanin biosynthesis in the dynamic developmental stages of fleshy roots were obtained from RNA-Seq analysis based on their DNA-binding domains. Of these TF families, the AP2/ERF-ERF, WRKY, NAC, bHLH, bZIP, MYB-related, C3H, C2H2, GARP-G2-like, MYB and GRAS families were found to be the most represented, with 125, 114, 84, 80, 78, 73, 63, 55, 48, 48 and 43 members, respectively, and were expressed in at least one of the seven libraries (Fig. 1A).

The results showed that 1,682, 1,677, 1,568, 1,632, 1,633, 1,656 and 1,586 TFs related to anthocyanin biosynthesis were expressed in fleshy roots from the SS, the IE stage, the FE stage, the BS, the IFS, the FBS and the PS, respectively (Table S3). Of those, 1,415 (80.99%) TFs were expressed in each of the seven libraries. Of those, five candidate TFs belonging to four TF families (AP2/ERF-ERF, GRAS, NAC and B3) were found to exhibit opposite dynamic anthocyanin profiles in the fleshy roots of carmine radish, but five candidate TFs belonging to four TF families (C2C2-Dof, bHLH, WRKY and C3H) were consistently found. Subsequently, those 10 putative TFs were selected and validated by qRT-PCR analysis. As shown in Fig. 1B, the results of qRT-PCR analysis were highly consistent with the expression profiles determined by RNASeq.

Differentially expressed TFs (DETFs) related to anthocyanin biosynthesis involved in the development of fleshy roots in carmine radish

A total of 161 DETFs were identified to be related to anthocyanin biosynthesis during the development of fleshy roots in carmine radish when compared with the SS_root group (P value <0.05, fold change <-1 or >1). In total, 116, 130, 136, 144, 141 and 136 DETFs were identified for ‘IE_root_vs_SS_root’, ‘FE_root’ vs ‘SS_root’, ‘BS_root’ vs ‘SS_root’, ‘IFS_root’ vs ‘SS_root’, ‘FBS_root’ vs ‘SS_root’ and ‘PS_root’ vs ‘SS_root’, respectively (Fig. 2A). In addition, these 161 DETFs were categorized into 43 families based on their DNA-binding domains (Table S4). Of those, the AP2/ERF-ERF (31 members), WRKY (nine members), bHLH (11 members), NAC (seven members), MYB (seven members) and HB-HD-zip (five members) families were found to be the most represented in the development of fleshy roots in carmine radish. Among the DETFs, approximately 71 were commonly differentially expressed in any one stage compared with the SS_root group using a venny graph (Fig. 2B). Subsequently, 161 DETFs were categorized into 4 distinct clusters using the K-means clustering algorithm (Fig. 2C). Among them, class I contained 24 substantially altered DETFs that were identified in the almost-whole stage and included the AP2/ERF-ERF (4 numbers), AUX/IAA (1 number), bHLH (1 number), HD-ZIP (3), NAC (2 numbers) and WRKY (four numbers) transcripts. A total of 45 genes belonging to class II showed slight dynamic changes during the development of fleshy roots. The expression patterns of 67 TFs in class IV were found to oppose those in class I, as they were substantially changed in the almost-whole stage and included AP2/ERF-ERF (22 numbers), bHLH (six numbers), MYB (eight numbers) and GRAS (two numbers) transcripts. Class III contained 25 genes that did not show significant trends throughout the development of fleshy roots. We further performed KOBAS analysis to annotate the DETFs in classes I and IV that are involved in the development of fleshy roots (Fig. 2C).

Validation of DETFs targeted by miRNAs related to anthocyanin biosynthesis in the development of fleshy roots in carmine radish

To investigate the functions of DETFs in carmine radish, the web-based program psRNATarget was selected to predict their related miRNAs. A total of 26 DETF transcripts, including MYB, WRKY, bHLH, ERF, GRAS, NF-YA, C2H2-Dof, and HD-ZIP, were predicted to be targets of 74 miRNAs belonging to 25 miRNA families (Table S5). In addition, we found that some miRNAs were associated with more than one DETF transcript. For example, miR395b could target ERF (Cluster_7379) and WRKY (Cluster_8116). In addition, some transcripts were regulated by more than one miRNA. For instance, ERF (Cluster_7379) was targeted by miR395a, miR395b and miR399b-5p, and WRKY (Cluster_8116) was targeted by miR319-5p, miR395b and miR535b (Table S5). Of those, eight DETF transcripts belonging to the C2C2-Dof, bHLH and ERF families and their eight corresponding miRNAs were selected for qRT-PCR to verify their functions in the development of fleshy roots in carmine radish (Fig. 3A). The expression levels of miR408a, miR1432-5p and miR166g were down-regulated in the development of fleshy roots in carmine radish ‘Hongxin 1’, whereas those of their corresponding target genes C2C2-Dof, bHLH and ERF, respectively, were upregulated in the ‘Hongxin 1’ cultivar (Figs. 3B, 3D, 3G). The expression levels of miR172b, miR165-5p, miR172d, miR827a and miR1425-5p were increased, while those of their corresponding target genes C2C2-Dof and ERF, respectively, were decreased in the ‘Hongxin 1’ cultivar (Figs. 3C, E, F, H and I). The results showed that inversely correlated expression patterns were found between the miRNAs and their corresponding targets.

Putative miRNA-target model associated with anthocyanin biosynthesis in carmine radish

Some TFs were shown to function as important regulators of anthocyanin biosynthesis, including MYB, bHLH, WRKY, ERF and HD-Zip. In addition, sucrose synthase, sugar/inositol transporter, and the WMBW complex were also found to potentially participate in anthocyanin biosynthesis. Here, based on the annotations of DETFs, we identified 26 transcripts of DETFs targeted by 74 miRNAs belonging to 25 miRNA families that might participate in anthocyanin biosynthesis in carmine radish (Table S5). Here, we propose a putative miRNA-target model associated with anthocyanin biosynthesis in carmine radish (Fig. 4). The miR395, miR319, miR5720, miR1432 and miR408 families potentially targeted WRKY, MYB and b-HLH through the WMBW complex; the miR171, miR408 and miR172, miR167 families possibly targeted GRAS, C2C2-Dof and HD-ZIP to play roles in auxin-mediated signaling; and miR827, miR165/166, miR172 and miR395/miR396 targeted ERF to activate sucrose synthase. All three biological processes (involving the WMBW complex, auxin-mediated signaling and sucrose synthase) influence the expression of structural genes (PAL, C4H, 4CL, CHS, F3H, DFR, ANS and UFGT) and subsequently modulate anthocyanin biosynthesis by forming a complex regulatory network.