Original Research ARTICLE
Transcriptomic Analysis Identifies Differentially Expressed Genes (DEGs) Associated with Bolting and Flowering in Radish (Raphanus sativus L.)
- National Key Laboratory of Crop Genetics and Germplasm Enhancement, Key Laboratory of Biology and Genetic Improvement of Horticultural Crops (East China), Ministry of Agriculture of China, College of Horticulture, Nanjing Agricultural University, Nanjing, China
The transition of vegetative growth to bolting and flowering is an important process in the life cycle of plants, which is determined by numerous genes forming an intricate network of bolting and flowering. However, no comprehensive identification and profiling of bolting and flowering-related genes have been carried out in radish. In this study, RNA-Seq technology was applied to analyze the differential gene expressions during the transition from vegetative stage to reproductive stage in radish. A total of 5922 differentially expressed genes (DEGs) including 779 up-regulated and 5143 down-regulated genes were isolated. Functional enrichment analysis suggested that some DEGs were involved in hormone signaling pathways and the transcriptional regulation of bolting and flowering. KEGG-based analysis identified 37 DEGs being involved in phytohormone signaling pathways. Moreover, 95 DEGs related to bolting and flowering were identified and integrated into various flowering pathways. Several critical genes including FT, CO, SOC1, FLC, and LFY were characterized and profiled by RT-qPCR analysis. Correlation analysis indicated that 24 miRNA-DEG pairs were involved in radish bolting and flowering. Finally, a miRNA-DEG-based schematic model of bolting and flowering regulatory network was proposed in radish. These outcomes provided significant insights into genetic control of radish bolting and flowering, and would facilitate unraveling molecular regulatory mechanism underlying bolting and flowering in root vegetable crops.
The developmental transition from vegetative growth to bolting and flowering is one of the most important traits in plant life cycle. Bolting and flowering time must be appropriately determined to ensure reproductive success under most favorable conditions (Amasino and Michaels, 2010; Srikanth and Schmid, 2011). Plants have evolved an intricate bolting and flowering genetic circuitry in response to various endogenous and environmental signals including development, age, plant hormones, photoperiod, and temperature (Fornara et al., 2010; Capovilla et al., 2015; Kazan and Lyons, 2015). Molecular and genetic regulation of flowering has been extensively studied in the model plant Arabidopsis thaliana. Five major flowering pathways including vernalization, photoperiod, autonomous, aging and gibberellin (GA) pathways have been identified to govern bolting and flowering time (Amasino and Michaels, 2010; Fornara et al., 2010; Srikanth and Schmid, 2011), and a number of flowering-related genes involved in these pathways have been isolated and characterized in Arabidopsis (Fornara et al., 2010; Srikanth and Schmid, 2011).
The signals from flowering pathways converge on several floral pathway integrators such as FLOWERING LOCUS T (FT), SUPPRESSOR OF OVEREXPRESSION OF CO1 (SOC1) and LEAFY (LFY), which are integrated into the genetic networks of flowering (Moon et al., 2005; Parcy, 2005). Among these integrators, the florigen gene FT is a central node of floral transition, whose transcriptional expression is positively regulated by CONSTANS (CO) encoding a putative zinc finger transcription factor (Suárez-López et al., 2001), while it is negatively regulated by FLOWERING LOCUS C (FLC), a flowering repressor encoding a MADS-box transcription factor (Michaels and Amasino, 1999). Different environmental factors affect plant flowering by modulating the expression of floral integrators and stimulating changes in plant hormone levels (Yaish et al., 2011; Riboni et al., 2014; Kazan and Lyons, 2015). Increasing evidences have revealed the connections between flowering time and plant hormones including salicylic acid (SA), jasmonic acid (JA), GA, abscisic acid (ABA) and auxin (Davis, 2009; Kazan and Lyons, 2015). The effects of phytohormone signaling on flowering, particularly GA pathway, have been extensively described in Arabidopsis (Mutasa-Göttgens and Hedden, 2009). Therefore, understanding the roles of flowering-related genes and crosstalk between diverse genetic pathways is fundamental for elucidating the regulatory mechanisms underlying bolting and flowering in plants.
RNA sequencing (RNA-Seq), a powerful strategy for global discovery of functional genes, has provided a better qualitative and quantitative description of gene expressions under certain conditions in many plant species (Lister et al., 2009; Wang et al., 2009a). Digital gene expression (DGE) tag profiling is a revolutionary approach for identifying differentially expressed genes (DEGs) in diverse plant tissues, organs and developmental stages (Bai et al., 2013; Zhang et al., 2014a; Zhu et al., 2015). Moreover, RNA-Seq combined with DGE profiling has been employed for flowering-related gene discovery and expression analysis in some species such as bamboo (Gao et al., 2014), Lagerstroemia indica (Zhang et al., 2014b), sweetpotato (Tao et al., 2013) and litchi (Zhang et al., 2014c). However, to our knowledge, there are no studies on global expression profile analysis of bolting and flowering-related genes in radish (Raphanus sativus L.).
Radish (2n = 2x = 18), belonging to Brassicaceae family, is an important annual or biennial root vegetable crop worldwide. Premature bolting is a seriously destructive problem and results in poor root growth and the reduced harvest during radish production, especially in spring. Appropriate timing of bolting and flowering is significant for reproductive success at suitable conditions, as well as preventing the premature bolting in radish. Progress on bolting and flowering time control (Fornara et al., 2010; Srikanth and Schmid, 2011), especially in Arabidopsis, has provided a solid foundation and reference for identifying numerous functional genes during radish bolting and flowering. Recently, the transcriptomes from radish roots and leaves have been assembled and analyzed (Wang et al., 2013; Zhang et al., 2013; Xu et al., 2015). Moreover, a list of microRNAs (miRNAs) and functional genes related to bolting and flowering were successfully isolated from late-bolting radish based on transcriptomic datasets (Nie et al., 2015). Therefore, to further identify the DEGs involved in bolting and flowering regulation is of importance for understanding the genetic regulatory network of bolting and flowering in radish.
In this study, to investigate the gene expression patterns during the transition of vegetative growth to bolting stage in radish, using the late bolting radish advanced inbred line as material, two DGE libraries were constructed and sequenced with RNA-Seq technology. The aims were to comprehensively identify DEGs involved in bolting and flowering regulatory network and to explore their roles in determining radish bolting and flowering time. Expression patterns of several critical DEGs related to bolting and flowering were validated by quantitative real-time PCR (RT-qPCR) analysis. Finally, to characterize the bolting and flowering-related genes and miRNAs in flowering pathways, a putative miRNA-DEG-based model of bolting and flowering regulatory network was put forward in radish. These results could provide significant insights into the molecular mechanism underlying bolting and flowering regulation in radish and other root vegetable crops.
Materials and Methods
The late bolting radish advanced inbred line ‘NAU-LU127’, which was self-pollinated for more than 20 generations, was used in this study. The genetic background and structure of this line are stable and highly homozygous. After surface-sterilization, the seeds were sowed and grew in a growth chamber with 16 h light at 25°C and 8 h darkness at 16°C. The radish leaves used for DGE sequencing and RT-qPCR analysis were separately collected at two different developmental stages: vegetative stage (VS) and reproductive stage (RS), with three biological replicates. Each sample was collected at two developmental stages from three randomly selected individual plants, respectively. All the samples were immediately frozen in liquid nitrogen and stored at −80°C until use.
DGE Library Construction and Illumina Sequencing
Total RNA from radish leaves at vegetative stage and reproductive stage was individually extracted using Trizol® Reagent (Invitrogen) following the manufacturer's instructions. The equivalent quantity of total RNA from three replicates was pooled and used for library preparation and sequencing. Two cDNA libraries named NAU-VS and NAU-RS were constructed and sequenced according to the previously reported method (Xu et al., 2015). The library construction and Illumina sequencing were performed using HiSeq™ 2500 platform at Beijing Genomics Institute (BGI, Shenzhen, China). The RNA-Seq data were deposited in NCBI Sequence Read Archive (SRA, http://www.ncbi.nlm.nih.gov/Traces/sra/) with accession numbers of SRX1671036 (NAU-VS) and SRX1671054 (NAU-RS).
Data Processing and Expression Analysis of DEGs
The raw reads were primarily produced for data processing. After filtering low quality reads, adaptor sequences and reads containing ploy-N, the clean reads were obtained. These clean reads were then matched to the radish reference sequences which contained the public radish genomic survey sequences (GSS) and expressed sequence tag (EST) sequences and leaf transcriptome sequences from ‘NAU-LU127’ (Nie et al., in press) with no more than two mismatches. These sequences from radish leaf transcriptome have been deposited in NCBI Transcriptome Shotgun Assembly (TSA, http://www.ncbi.nlm.nih.gov/genbank/tsa/) database under the accession number GEMG00000000.
To screen the DEGs between two DGE libraries, the expression level of each transcript is calculated using RPKM (Reads Per kb per Million reads) method (Mortazavi et al., 2008). Prior to differential gene expression analysis, the read counts of each transcript were adjusted by edgeR program package (Robinson et al., 2010) through one scaling normalized factor. Trimmed Mean of M values (TMM), an appropriate normalization method implemented in the edgeR package (Robinson and Oshlack, 2010; Robinson et al., 2010), was employed to obtain the normalized read counts. The differential expression analysis of two libraries was performed using the DEGSeq R package 1.20.0 (Wang et al., 2010). Subsequently, the false discovery rate (FDR) was used to determine P-value threshold in multiple testing (Benjamini et al., 2001). A strict algorithm was used to further perform DEG identification according to the previous reports (Audic and Claverie, 1997). The absolute value of log2Ratio (NAU-RS/NAU-VS) ≥ 1, P < 0.05 and FDR ≤ 0.001 were used as threshold for judging the significance of gene expression difference. The cluster analysis of gene expression patterns was performed with cluster software and Java Treeview software (Saldanha, 2004).
Functional Annotation and Enrichment Analysis of DEGs
To investigate the biological function and involvement in functional pathways, all the identified transcripts were mapped to Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) database. For GO annotation, the unique transcripts were subjected to BLASTX searching against the NCBI Nr database using the E < 10−5. Then the Blast2GO (Conesa et al., 2005) and WEGO software (Ye et al., 2006) were used to obtain GO annotations and functional classifications. GO enrichment analysis of DEGs was implemented by the GOseq R package (Young et al., 2010). KOBAS software (Xie et al., 2011) was used to test the statistical enrichment of DEGs in KEGG pathways. The significantly enriched functional terms and pathways were identified using the criterion of a Bonferroni-corrected P ≤ 0.05.
Total RNAs from radish leaves were isolated and obtained as described above. RT-qPCR was performed on an iCycler Real-Time PCR Detection System (Bio-Rad, USA) with three replications according to previous reports (Nie et al., 2015; Xu et al., 2015). All PCR reactions were carried out in a total volume of 20 μL with RsActin gene as the internal control (Xu et al., 2012). The relative gene expression levels were calculated using 2 method (Livak and Schmittgen, 2001). The specific PCR primers were designed using Beacon Designer 7.0 (Premier Biosoft International, USA) and listed in Table S7.
DGE Library Sequencing and Data Analysis
To obtain global unique sequences from radish leaves, de novo assembly and analysis of transcriptome prepared from radish leaves of ‘NAU-LU127’ were carried out using Illumina RNA-Seq technology. Totally 111,167 contigs and 53,642 unigenes were generated from the radish leaf transcriptome (Nie et al., in press). The available dataset of radish leaf transcriptome integrating with the available radish GSS and EST sequences released in NCBI database enriched the radish reference sequences for DEG identification during radish bolting and flowering.
In this study, two DGE libraries from leaves of radish advanced inbred line ‘NAU-LU127’ at vegetative and reproductive stages were constructed and sequenced by Illumina HiSeq™ 2500 platform, respectively. As a result, 8,825,790 and 12,382,793 raw reads were obtained in NAU-VS and NAU-RS libraries, respectively (Table 1; Figure 1A). After removing adaptor sequences and low quality reads, 8,541,912 and 10,154,256 clean reads were generated in the two libraries (Table 1). These clean reads were then mapped to the radish reference sequences, resulting in the generation of 87.83% (7,502,455 reads) and 68.50% (6,955,386 reads) matched reads in NAU-VS and NAU-RS libraries, respectively. For the variation of clean reads mapping percentage between two libraries, it may arise from the sample differences and the specific pre-processing of obtained reads (Oshlack et al., 2010). The more mapped reads from NAU-VS library implied that some stage-specific genes may be expressed only at vegetative stage and differentially expressed between these two DGE libraries. Further analysis revealed that 4,124,632 reads (48.29 %) in NAU-VS library and 5,200,880 reads (51.22 %) in NAU-RS library were uniquely matched (Table 1).
Figure 1. Overview of the sequencing and mapping (A) and the identified DEGs (B) in NAU-VS and NAU-RS libraries.
Identification and Functional Enrichment Analysis of DEGs
The transcript abundance of each gene from two DGE libraries was calculated and analyzed by RPKM method. The threshold of |log2Ratio| ≥ 1 and FDR ≤ 0.001 were further used to determine the significantly DEGs. A total of 5922 significantly DEGs including 779 up-regulated and 5143 down-regulated genes were obtained from NAU-VS and NAU-RS libraries (Table S1; Figure 1B).
To better classify the functions of these identified DEGs, GO classification and enrichment analysis were carried out in this study. These DEGs were categorized into three main GO categories including 23 biological processes, 14 cellular components and 13 molecular functions (Figure 2). Functional enrichment analysis revealed that 140 GO terms were significantly enriched with a Bonferroni-corrected P ≤ 0.05 (Table S2). The terms of “metabolic process” (GO: 0008152) and “organic substance metabolic process” (GO: 0071704) were the dominant groups in biological processes; “cell” (GO: 0005623) and “cell part” (GO: 0044464) were the highly represented groups in the cellular components. For the molecular functions, a large proportion of genes were significantly enriched in “organic cyclic compound binding” (GO: 0097159) and “heterocyclic compound binding” (GO: 1901363) categories. Moreover, some enriched GO terms were related to plant flowering and meristem development, including “regulation of photoperiodism, flowering” (GO: 2000028), “regulation of timing of meristematic phase transition” (GO: 0048506), “meristem maintenance” (GO: 0010073), “meristem growth” (GO: 0035266), “meristem development” (GO: 0048507) and “flower development” (GO: 0009908) (Table S2).
To further understand the putative active biological pathways, all the identified DEGs were mapped to KEGG database by BLASTx with E ≤ 10−5 and Q ≤ 1. As a result, 5922 DEGs were successfully assigned to 128 KEGG pathways (Table S3). The dominant pathway was “Metabolic pathways,” followed by “Biosynthesis of secondary metabolites,” “Ribosome,” and “Plant hormone signal transduction.” Moreover, 17 pathways were significantly enriched (Q ≤ 0.05; Table 2), including “Circadian rhythm-plant” (ko04712), “Plant hormone signal transduction” (ko04075), “Photosynthesis” (ko00195), “Ribosome” (ko03010) and “Vitamin B6 metabolism” (ko00750).
DEGs Involved in Hormone Signal Transduction Pathway
In this study, pathway-based analysis showed that 37 DEGs representing 393 unique sequences were identified and involved in “Plant hormone signal transduction” (ko04075) pathway (Table 3; Table S4; Figure S1). These genes including AUX1, TRANSPORT INHIBITOR RESPONSE 1 (TIR1), AUXIN RESPONSE FACTOR (ARFs), GIBBERELLIN RECEPTOR 1 (GID1), and CORONITINE INSENSITIVE 1 (COI1), participated in the regulation of several hormone homeostasis and flowering time (Davis, 2009; Kazan and Lyons, 2015). Enrichment analysis revealed that most of genes were involved in auxin, GA, ABA, JA, and SA signaling pathways (Figure 3). In GA signaling pathway, one down- and three up-regulated transcripts were related to GID1 protein, while eight down-regulated transcripts encoded DELLA protein (Figure 3A). For the process of JA signaling, two down-regulated transcripts encoded JAR1 protein and one down-regulated transcript encoded COI1 protein (Figure 3D). In auxin signaling pathway, 10 down-regulated transcripts belonged to ARF genes, while one up- and seven down-regulated transcripts encoded auxin-responsive proteins (Figure 3E). In addition, some DEGs related to other phytohormone biosynthesis were also identified, including zeatin biosynthesis (ko00908, three DEGs), carotenoid biosynthesis (ko00906, four DEGs), cysteine and methionine metabolism (ko00270, seven DEGs), brassinosteroid biosynthesis (ko00905, eight DEGs), and phenylalanine metabolism (ko00360, three DEGs) (Table 3; Figure S1).
Figure 3. Heat map diagram of expression patterns for DEGs involved in some phytohormone signaling pathways, including GA (A), ABA (B), SA (C), JA (D), and auxin (E). Red and green colors indicate up- and down-regulated genes in NAU-RS library as compared with NAU-VS library, respectively.
DEGs Involved in the Transition of Vegetative Growth to Bolting in Radish
In this study, to identify DEGs during radish bolting and flowering, BLAST searching was performed and the putative functions of DEGs were assessed. A total of 95 DEGs representing 128 unique sequences related to bolting and flowering were identified (Table S5). The analysis of flowering pathways revealed that these genes were involved in five different flowering pathways including photoperiod, vernalization, autonomous, GA and aging pathways.
In the present study, some unigenes representing photoperiodic flowering genes were identified, including AGAMOUS-LIKE 24 (AGL24), APETALA2 (AP2), CO, CELL GROWTH DEFECT FACTOR 1 (CDF1), and CONSTITUTIVELY PHOTOMORPHOGENIC 1 (COP1; Table S5). Some genes related to circadian rhythm and light signaling pathway included CIRCADIAN CLOCK ASSOCIATED 1 (CCA1), CONSTANS-LIKE 1 (COL1), CRYPTOCHROME (CRY2), TIMING OF CAB EXPRESSION 1 (TOC1), and LATE ELONGATED HYPOCOTYL (LHY; Table S5).
For the vernalization pathway, one down-regulated transcript (CL1584.Contig3) belonging to FLC homolog was found in this study (Table S5), which is a major flowering repressor and integrates the autonomous and vernalization pathways (Michaels and Amasino, 1999). Many genes involved in vernalization pathway including FRIGIDA (FRI), FRIGIDA-like (FRL), FRIGIDA INTERACTING PROTEIN 2 (FIP2), EMBRYONIC FLOWER 2 (EMF2), VERNALIZATION 1 (VRN1), and VERNALIZATION 2 (VRN2) were also identified and implicated in regulating the expression of FLC. Furthermore, LUMINIDEPENDENS (LD), FPA, FVE, and FY involved in autonomous pathway were also identified (Table S5).
Moreover, some putative genes for GA and aging pathways were also found in the present study (Table S5). The candidate genes involved in GA pathway comprised GIGANTEA (GI), GNC, GA INSENSITIVE DWARF 1B (GID1B), DWARF AND DELAYED FLOWERING 1 (DDF1), and REPRESSOR OF GA 1-3 (RGA1-3). The candidate genes related to aging pathway included SQUAMOSA PROMOTER BINDING-LIKE PROTEIN 1 (SPL1), SPL2, SPL3, SPL9, SPL13, and SPL15. In addition, some floral integrators such as FT (FD571044), SOC1 (CL4258.Contig1) and LFY (Unigene29702), were also identified in this study (Table S5).
Expression Profile Analysis by RT-qPCR
To validate the differential expression patterns of DEGs during radish bolting and flowering, totally 21 functional genes were randomly selected and subjected to RT-qPCR analysis. These selected genes included six genes related to hormone signaling and 15 genes related to bolting and flowering regulation. The relative expression levels of these genes between vegetative growth and reproductive stage were analyzed and compared (Figure 4). Further, comparative analysis revealed that these gene expression trends except MYC2-CL4584.Contig1 were in agreement with the transcript abundance changes by RNA-Seq (Figure 4), indicating the highly accuracy and quality of DGE sequencing.
Figure 4. The validation of expression levels of selected DEGs related to radish bolting and flowering. (A) The relative expression levels of selected DEGs were compared with the transcript abundances from DGE sequencing. (B) Heat map diagram of expression patterns of DEGs in radish. Red and green colors in indicate up- and down-regulated genes in NAU-RS library as compared with NAU-VS library, respectively.
The Regulatory Network Underlying Bolting and Flowering in Radish
Considerable studies have revealed that some miRNAs regulating corresponding target genes played important roles in the transition from vegetative growth to bolting and flowering (Spanudakis and Jackson, 2014). In our recent study, several bolting and flowering-related miRNA-target gene pairs were identified and characterized in late-bolting radish (Nie et al., 2015). To better understand the genetic regulatory network of radish bolting and flowering, correlation analysis between the DEGs identified in the present study and bolting and flowering-related miRNAs previously reported (Nie et al., 2015) was performed. As expected, 24 miRNA-mRNA pairs including 16 miRNAs and 27 target DEGs were identified (Table S6). Among them, 19 miRNA-mRNA pairs showed negative correlations in expression patterns. Several DEGs including AP2 (targeted by miR172), VRN1 (targeted by miR5227), PRP39 (targeted by miR6273), and NF-YB3 (targeted by miR860), were found to be involved in bolting and flowering regulation (Wang et al., 2007; Kumimoto et al., 2008; Zhu and Helliwell, 2010).
To gain insights into the bolting and flowering regulatory network in radish, a putative model for summarizing the bolting and flowering-related DEGs and miRNAs was proposed (Figure 5). The critical genes involved in various flowering pathways and phytohormone signaling pathways were displayed in the schematic regulatory network of radish bolting and flowering. According to the known Arabidopsis flowering regulatory network (Fornara et al., 2010; Srikanth and Schmid, 2011), we speculated that the transcriptional regulations of several floral integrators including FT, CO, SOC1, FLC, and LFY, could integrate the signals from various pathways and modulate the radish bolting and flowering (Figure 5). Moreover, the models of miR172-AP2 and miR5227-VRN1 have been shown to be important participants in the regulatory network of bolting and flowering (Wang et al., 2007; Zhu and Helliwell, 2010; Nie et al., 2015).
Figure 5. The putative model of flowering regulatory network integrating various flowering pathways and transcriptional processes in radish. The FT, SOC1, LFY, FLC, and CO genes were flowering pathway integrators.
Radish bolting and flowering are integral stages in its complete life cycle. The timing of bolting and flowering is coordinately regulated by various endogenous and environmental signals integrating into a complexity of flowering regulation (Amasino and Michaels, 2010; Srikanth and Schmid, 2011). Recent advances in flowering genes and regulatory networks have greatly enhanced our knowledge of molecular basis underlying bolting and flowering-time control in Brassicaceae crops. However, no studies on comprehensive identification of DEGs related to radish bolting and flowering have been reported, and the regulatory mechanism of bolting and flowering-time control remains largely unexplored in radish. In this study, two cDNA libraries from leaves of radish advanced inbred line ‘NAU-LU127’ at vegetative and reproductive stages were constructed, respectively. A list of DEGs related to phytohormone signaling and transition from vegetative growth to bolting and flowering were identified and comprehensively profiled.
The Roles of Plant Hormone Signaling in Bolting and Flowering
Plant hormones are endogenously occurring compounds that regulate multiple aspects of plant growth and development including flowering time (Davis, 2009; Santner and Estelle, 2009). Various phytohormones have been implicated in the developmental transition of flowering (Davis, 2009; Domagalska et al., 2010). The pathways of several hormones including auxin, GA, ABA, SA, and JA signaling were significantly enriched by pathway-based analysis in our study (Table 3).
GA pathway is one of the genetic flowering pathways, which could interact with several pathways and is integrated into the flowering regulatory complexity (Srikanth and Schmid, 2011). The role of GA pathway in flowering time has been thoroughly investigated in Arabidopsis and several fruit trees (Wilkie et al., 2008; Mutasa-Göttgens and Hedden, 2009). Many genes related to GA metabolism and signaling were involved in GA-mediated regulatory process of flowering (Mutasa-Göttgens and Hedden, 2009; Domagalska et al., 2010). GA exerting its biological functions on floral transition and development is mainly dependent on the growth inhibitor DELLA proteins (Mutasa-Göttgens and Hedden, 2009). GA signaling promotes flowering through initiating the degradation of transcriptional regulator DELLA and activating the expression of SOC1, AGL24 and LFY (Davis, 2009; Mutasa-Göttgens and Hedden, 2009). As expected, both the decreased transcript abundance and expression level of DELLA (CL7599.Contig1) were detected in radish reproductive stage compared with vegetative phase (Figure 4). The ABA pathway, which is antagonistic to GA, has been demonstrated to delay flowering through modulating DELLA activity and affecting the transcriptional expression of floral repressor FLC (Achard et al., 2006; Domagalska et al., 2010). In the current study, unique transcripts annotated as PP2C and ABF, ABA signaling components, were identified and differentially regulated during radish bolting and flowering (Figure 3), which is consistent with the results in litchi (Zhang et al., 2014c) and soybean (Wong et al., 2013). These findings suggested that the differential expressions of ABA signaling-related genes may be associated with the timing of radish transition to bolting and flowering.
Function of SA in accelerating transition to flowering is pronounced by SA-deficient mutants of Arabidopsis (Martínez et al., 2004). SA could negatively regulate the floral repressor FLC and activate the flowering promoter FT which strongly highlights the positive role of SA in flowering transition (Martínez et al., 2004). SA promotes the activation of NON-EXPRESSOR OF PR-1 (NPR1) proteins, whose interaction with TGA transcription factors could induce the expression of PR genes (Wu et al., 2012). Moreover, JA is also implicated in flowering regulatory process and delays flowering in Arabidopsis (Krajnčič et al., 2006; Riboni et al., 2014). JA signaling pathway has been involved in three molecular elements including JA receptor gene COI1, transcriptional repressor JAZ protein and some transcription factors, e.g., the bHLH family (Krajnčič et al., 2006). Notably, recent studies have demonstrated the regulatory role of COI1 in delaying flowering mediating the repressed expression of FT (Zhai et al., 2015). In this study, some transcripts belonging to the main components of SA and JA signaling were found, including NPR1, TGA, PR, JAZ, COI1, and MYC2 (Table 3). In addition, previous studies reveal that auxin is necessary for flower initiation and floral organ identity (Cheng and Zhao, 2007). We also found the critical genes related to auxin signaling such as AUX1, SAUR, TIR1, and ARFs (Table 3; Table S4). Overall, these results reveal that phytohormone-mediated transcriptional reprogramming are crucial to the transition of bolting and flowering and participate in its regulatory network of radish. The characterization of critical genes in plant hormone signaling pathways would greatly help to illuminate the complex genetic network of bolting and flowering in radish.
The Complex Bolting and Flowering Regulatory Network in Radish
Multiple genetic flowering pathways integrating endogenous and environmental signals determine the transition from vegetative growth to reproductive development. Studies in Arabidopsis have revealed the participation of more than 200 flowering-related genes in the intricate regulatory network (Fornara et al., 2010; Srikanth and Schmid, 2011). In this study, 95 candidate genes related to bolting and flowering were isolated and involved in five major flowering pathways within genetic regulatory network (Table S5; Figure 5). It is inferred that known genetic pathways and critical flowering genes may conservatively present in radish, being consistent with the reports in maize (Dong et al., 2012), soybean (Jung et al., 2012), and citrus (Zhang et al., 2011). Gene expression profiling revealed that these genes were differentially expressed between NAU-VS and NAU-RS libraries, suggesting their putative important roles in radish bolting and flowering.
The complex regulatory network of Arabidopsis is composed of five major converging pathways (Fornara et al., 2010; Srikanth and Schmid, 2011). It is believed that endogenous developmental signals such as developmental stages of plants and phytohormones monitor flowering time through age, autonomous and GA pathways, while environmental cues regulate flowering time through the photoperiod and vernalization pathways in response to day length or temperature (Srikanth and Schmid, 2011; Capovilla et al., 2015). The signals from photoperiodic process are converted into the transcriptional regulation of key genes such as FT, CO, AP1, and AP2 to affect flowering time (Kikuchi and Handa, 2009; Amasino, 2010; Srikanth and Schmid, 2011). The florigen gene FT as a floral integrator is central for the photoperiodic flowering pathway of long-day plant Arabidopsis, which is perceived in leaves and transported to the shoot apex initiating floral transition (Huang et al., 2005; Parcy, 2005). The role of FT in promoting flowering has been proven by mutants and overexpressed transgenic analysis in Arabidopsis (Amasino, 2010; Srikanth and Schmid, 2011). As expected, the homolog of FT (FD571044) was up-regulated in reproductive stage of radish (Figure 4), indicating that the RsFT gene could positively regulate the development transition of bolting and flowering (Figure 5). Under long-day condition, the FT expression is activated by CO, which is a floral activator and modulated by the circadian clock and day length (Suárez-López et al., 2001; Amasino, 2010; Johansson and Staiger, 2015). The link between circadian clock and flowering control may be mainly mediated by the transcriptional expression of CO (Fujiwara et al., 2008; Johansson and Staiger, 2015). In Arabidopsis, two essential circadian clock components LATE ELONGATED HYPOCOTYL (LHY) and CIRCADIAN CLOCK ASSOCIATED1 (CCA1) function in photoperiodic flowering and regulate flowering pathway by controlling the rhythmic expression of CO and FT (Fujiwara et al., 2008). In this study, some transcripts belonging to CO, CCA1 and LHY homologs were found to be up-expressed in reproductive phase with DGE sequencing and RT-qPCR analysis (Figure 4), suggesting the critical roles of these genes in the transition of radish bolting and flowering.
The vernalization and autonomous pathways converge on the flowering repressor FLC, and many genes involved in these two pathways could control flowering time through affecting FLC expression (Amasino, 2010). The high level of FLC delays flowering and requires its activator FRI (Michaels and Amasino, 1999; Choi et al., 2011). Recently, several naturally occurring spliced transcripts of FLC were found and isolated from B. rapa (Yuan et al., 2009) and orange (Zhang et al., 2009), which were proven to be associated with variations in flowering time. The transcriptional co-expression analysis in B. rapa indicated that BrFLC2 may be the major regulator of flowering time in genetic flowering network (Xiao et al., 2013). In this study, we found putative homologs of FLC (CL1584.Contig3) from late-bolting radish, which was down-regulated in reproductive stage compared with vegetative stage, with similar patterns being detected in FRI and FRL (Figure 4). In addition, similar results were found in other homologous genes in vernalization pathway, including FIP2, EMF2, VRN1, and VRN2 (Table S5). These results indicate that the genetic elements of the vernalization pathway may be of importance for the manipulation of radish bolting and flowering time.
Furthermore, miRNAs as central regulators of gene expression have been shown to be implicated in multiple genetic pathways governing flowering time (Spanudakis and Jackson, 2014; Wang, 2014). The newly defined age pathway of flowering, which is controlled by miR156 and its target SPL transcription factors (Wang et al., 2009b), regulates flowering time and interacts with vernalization, photoperiodic and GA pathways (Zhou et al., 2013; Spanudakis and Jackson, 2014; Wang, 2014). Several members of SPL family were identified in this study, including SPL1, SPL2, SPL3, SPL9, SPL13, and SPL15 (Table S5). It was known that miR172 is down-regulated by the age-dependent expression of SPL9 (Wu et al., 2009; Spanudakis and Jackson, 2014). The target genes of miR172 are a class of AP2-like transcription factors including AP2, TARGET OF EAT 1-3 (TOE1-3), SCHLAFMÜTZE (SMZ), and SCHNARCHZAPFEN (SNZ), which act as floral repressors (Zhu and Helliwell, 2010). The levels of these AP2-like genes are relatively high during plant seedling stage and decline with plant development, ultimately relieving the repression of flowering to trigger flowering (Aukerman and Sakai, 2003; Zhu and Helliwell, 2010). Consistent with these evidences, the down-expressed patterns of AP2 (CL1275.Contig1), SMZ (Rsa#S42015352), and RAP2.7 (CL2600.Contig3) were detected at reproductive stage in this study (Figure 4). Moreover, correlation analysis revealed that some bolting and flowering-related DEGs were targeted by specific miRNAs forming the transcriptional model of miRNA-mRNA pairs (Table S6). These findings reveal that some miRNA-DEG models including miR5227-VRN1, miR6273-PRP39, and miR860-NF-YB3 are crucial participators and integrated into the intricate genetic networks of bolting and flowering in radish (Figure 5).
In summary, RNA-Seq technology was employed to systematically identify DEGs at transcriptome-wide level during radish transition from vegetative growth to bolting and flowering in this study. To our knowledge, this is the first investigation to illustrate the expression profiles of bolting-related genes and dissect the bolting and flowering regulatory network in radish. In this study, a total of 5922 DEGs were identified from late-bolting radish leaves. Several candidate genes related to plant hormone signal and bolting and flowering regulatory pathways were characterized and implicated in the complex networks of bolting and flowering regulation. Correlation analysis suggested that the miRNA-mRNA regulatory models played pivotal roles in determining bolting and flowering time. Moreover, a schematic regulatory network of radish bolting and flowering was put forward for characterization of DEGs and miRNAs. These results provided essential information for genetic control of radish bolting and flowering, and would facilitate unraveling the molecular regulatory mechanism underlying bolting and flowering in radish and other root vegetable crops.
SN, CL, and LL designed the research. SN, XS, and MT conducted experiments. SN, LX, and YW participated in the design of the study and performed the statistical analysis. SN analyzed data and wrote the manuscript. LL and EM helped with the revision of manuscript. All authors read and approved the manuscript.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
This work was partially supported by grants from the National Natural Science Foundation of China (31171956, 31372064), the National Key Technologies R & D Program of China (2012BAD02B01) and Key Technologies R & D Program of Jiangsu Province (BE2013429).
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/article/10.3389/fpls.2016.00682
Table S1. All the identified DEGs in NAU-VS and NAU-RS libraries.
Table S2. GO enrichment analysis for differentially expressed transcripts with corrected P ≤ 0.05.
Table S3. KEGG pathway analysis of differentially expressed transcripts in radish.
Table S4. The DEGs involved in plant hormone signal transduction pathway.
Table S5. The DEGs involved in the transition of vegetative growth to bolting in radish.
Table S6. The identified DEG and miRNA pairs during radish bolting and flowering.
Table S7. The primers of DEGs for RT-qPCR in radish.
Figure S1. The identified genes involved in plant hormone signal transduction by KEGG analysis.
Achard, P., Cheng, H., Grauwe, L., Decat, J., Schoutteten, H., Moritz, T., et al. (2006). Integration of plant responses to environmentally activated phytohormonal signals. Science 311, 91–94. doi: 10.1126/science.1118642
Bai, S., Saito, T., Sakamoto, D., Ito, A., Fujii, H., and Moriguchi, T. (2013). Transcriptome analysis of japanese pear (Pyrus pyrifolia Nakai) flower buds transitioning through endodormancy. Plant Cell Physiol. 54, 1132–1151. doi: 10.1093/pcp/pct067
Benjamini, Y., Drai, D., Elmer, G., Kafkafi, N., and Golani, I. (2001). Controlling the false discovery rate in behavior genetics research. Behav. Brain. Res. 125, 279–284. doi: 10.1016/S0166-4328(01)00297-2
Choi, K., Kim, J., Hwang, H., Kim, S., Park, C., Kim, S., et al. (2011). The FRIGIDA complex activates transcription of FLC, a strong flowering repressor in Arabidopsis, by recruiting chromatin modification factors. Plant Cell 23, 289–303. doi: 10.1105/tpc.110.075911
Conesa, A., Götz, S., García-Gómez, J., Terol, J., Talón, M., and Robles, M. (2005). Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics 21, 3674–3676. doi: 10.1093/bioinformatics/bti610
Domagalska, M. A., Sarnowska, E., Nagy, F., and Davis, S. J. (2010). Genetic analyses of interactions among gibberellin, abscisic acid, and brassinosteroids in the control of flowering time in Arabidopsis thaliana. PLoS ONE 5:e14012. doi: 10.1371/journal.pone.0014012
Dong, Z., Danilevskaya, O., Abadie, T., Messina, C., Coles, N., and Cooper, M. (2012). A gene regulatory network model for floral transition of the shoot apex in maize and its dynamic modeling. PLoS ONE 7:e43450. doi: 10.1371/journal.pone.0043450
Fujiwara, S., Oda, A., Yoshida, R., Niinuma, K., Miyata, K., Tomozoe, Y., et al. (2008). Circadian clock proteins LHY and CCA1 regulate SVP protein accumulation to control flowering in Arabidopsis. Plant Cell 20, 2960–2971. doi: 10.1105/tpc.108.061531
Gao, J., Zhang, Y., Zhang, C., Qi, F., Li, X., Mu, S., et al. (2014). Characterization of the floral transcriptome of moso bamboo (Phyllostachys edulis) at different flowering developmental stages by transcriptome sequencing and RNA-Seq analysis. PLoS ONE 9:e98910. doi: 10.1371/journal.pone.0098910
Huang, T., Böhlenius, H., Eriksson, S., Parcy, F., and Nilsson, O. (2005). The mRNA of the Arabidopsis gene FT moves from leaf to shoot apex and induces flowering. Science 309, 1694–1696. doi: 10.1126/science.1117768
Krajnčič, B., Kristl, J., and Janžekovič, I. (2006). Possible role of jasmonic acid in the regulation of floral induction, evocation and floral differentiation in Lemna minor L. Plant Physiol. Bioch. 44, 752–758. doi: 10.1016/j.plaphy.2006.10.029
Kumimoto, R. W., Adam, L., Hymus, G. J., Repetti, P. P., Reuber, T. L., Marion, C. M., et al. (2008). The Nuclear Factor Y subunits NF-YB2 and NF-YB3 play additive roles in the promotion of flowering by inductive long-day photoperiods in Arabidopsis. Planta 228, 709–723. doi: 10.1007/s00425-008-0773-6
Lister, R., Gregory, B. D., and Ecker, J. R. (2009). Next is now: new technologies for sequencing of genomes, transcriptomes, and beyond. Curr. Opin. Plant Biol. 12, 107–118. doi: 10.1016/j.pbi.2008.11.004
Martínez, C., Pons, E., Prats, G., and León, J. (2004). Salicylic acid regulates flowering time and links defence responses and reproductive development. Plant J. 37, 209–217. doi: 10.1046/j.1365-313X.2003.01954.x
Nie, S., Xu, L., Wang, Y., Huang, D., Muleke, E. M., Sun, X., et al. (2015). Identification of bolting-related microRNAs and their targets reveals complex miRNA-mediated flowering-time regulatory networks in radish (Raphanus sativus L.). Sci. Rep. 5:14034. doi: 10.1038/srep14034
Nie, S., Li, C., Xu, L., Wang, Y., Huang, D., Muleke, E. M., et al. (in press). De novo transcriptome analysis and identification of critical genes involved in bolting and flowering in radish (Raphanus sativus L.). BMC Genomics 17. doi: 10.1186/s12864-016-2633-2
Riboni, M., Robustelli, A., Galbiati, M., Tonelli, C., and Conti, L. (2014). Environmental stress and flowering time: the photoperiodic connection. Plant Signal. Behav. 9:e29036. doi: 10.4161/psb.29036
Robinson, M. D., McCarthy, D. J., and Smyth, G. K. (2010). edgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26, 139–140. doi: 10.1093/bioinformatics/btp616
Suárez-López, P., Wheatley, K., Robson, F., Onouchi, H., Valverde, F., and Coupland, G. (2001). CONSTANS mediates between the circadian clock and the control of flowering in Arabidopsis. Nature 410, 1116–1120. doi: 10.1038/35074138
Tao, X., Gu, Y., Jiang, Y., Zhang, Y., and Wang, H. (2013). Transcriptome analysis to identify putative floral-specific genes and flowering regulatory-related genes of sweet potato. Biosci. Biotech. Bioch. 77, 2169–2174. doi: 10.1271/bbb.130218
Wang, C., Tian, Q., Hou, Z., Mucha, M., Aukerman, M., and Olsen, O. (2007). The Arabidopsis thaliana ATPRP39-1 gene, encoding a tetratricopeptide repeat protein with similarity to the yeast pre-mRNA processing protein PRP39, affects flowering time. Plant Cell Rep. 26, 1357–1366. doi: 10.1007/s00299-007-0336-5
Wang, L., Feng, Z., Wang, X., and Zhang, X. (2010). DEGseq: an R package for identifying differentially expressed genes from RNA-Seq data. Bioinformatics 26, 136–138. doi: 10.1093/bioinformatics/btp612
Wang, J. W., Czech, B., and Weigel, D. (2009b). miR156-regulated SPL transcription factors define an endogenous flowering pathway in Arabidopsis thaliana. Cell 138, 738–749. doi: 10.1016/j.cell.2009.06.014
Wang, Y., Pan, Y., Liu, Z., Zhu, X., Zhai, L., Xu, L., et al. (2013). De novo transcriptome sequencing of radish (Raphanus sativus L.) and analysis of major genes involved in glucosinolate metabolism. BMC Genomics 14:836. doi: 10.1186/1471-2164-14-836
Wong, C. E., Singh, M. B., and Bhalla, P. L. (2013). The dynamics of soybean leaf and shoot apical meristem transcriptome undergoing floral initiation process. PLoS ONE 8:65319. doi: 10.1371/journal.pone.0065319
Wu, G., Park, M. Y., Conway, S. R., Wang, J., Weigel, D., and Poethig, R. S. (2009). The sequential action of miR156 and miR172 regulates developmental timing in Arabidopsis. Cell 138, 750–759. doi: 10.1016/j.cell.2009.06.031
Wu, Y., Zhang, D., Chu, J. Y., Boyle, P., Wang, Y., Brindle, I. D., et al. (2012). The Arabidopsis NPR1 protein is a receptor for the plant defense hormone salicylic acid. Cell Rep. 1, 639–647. doi: 10.1016/j.celrep.2012.05.008
Xiao, D., Zhao, J., Hou, X., Basnet, R. K., Carpio, D., Zhang, N., et al. (2013). The Brassica rapa FLC homologue FLC2 is a key regulator of flowering time, identified through transcriptional co-expression networks. J. Exp. Bot. 64, 4503–4516. doi: 10.1093/jxb/ert264
Xie, C., Mao, X., Huang, J., Ding, Y., Wu, J., Dong, S., et al. (2011). KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic Acids Res. 39, W316–W322. doi: 10.1093/nar/gkr483
Xu, L., Wang, Y., Liu, W., Wang, J., Zhu, X., Zhang, K., et al. (2015). De novo sequencing of root transcriptome reveals complex cadmium-responsive regulatory networks in radish (Raphanus sativus L.). Plant Sci. 236, 313–323. doi: 10.1016/j.plantsci.2015.04.015
Xu, Y., Zhu, X., Gong, Y., Xu, L., Wang, Y., and Liu, L. (2012). Evaluation of reference genes for gene expression studies in radish (Raphanus sativus L.) using quantitative real-time PCR. Biochem. Biophys. Res. Commun. 424, 398–403. doi: 10.1016/j.bbrc.2012.06.119
Yaish, M. W., Colasanti, J., and Rothstein, S. J. (2011). The role of epigenetic processes in controlling flowering time in plants exposed to stress. J. Exp. Bot. 62, 3727–3735. doi: 10.1093/jxb/err177
Yuan, Y., Wu, J., Sun, R., Zhang, X., Xu, D., Bonnema, G., et al. (2009). A naturally occurring splicing site mutation in the Brassica rapa FLC1 gene is associated with variation in flowering time. J. Exp. Bot. 60, 1299–1308. doi: 10.1093/jxb/erp010
Zhai, Q., Zhang, X., Wu, F., Feng, H., Deng, L., Xu, L., et al. (2015). Transcriptional mechanism of Jasmonate receptor COI1-mediated delay of flowering time in Arabidopsis. Plant Cell 27, 2814–2828. doi: 10.1105/tpc.15.00619
Zhang, J., Li, Z., Mei, L., Yao, J., and Hu, C. (2009). PtFLC homolog from trifoliate orange (Poncirus trifoliata) is regulated by alternative splicing and experiences seasonal fluctuation in expression level. Planta 229, 847–859. doi: 10.1007/s00425-008-0885-z
Zhang, J., Ai, X., Sun, L., Zhang, D., Guo, W., Deng, X., et al. (2011). Transcriptome profile analysis of flowering molecular processes of early flowering trifoliate orange mutant and the wild-type [Poncirus trifoliate (L.) Raf.] by massively parallel signature sequencing. BMC Genomics 12:63. doi: 10.1186/1471-2164-12-63
Zhang, Y., Peng, L., Wu, Y., Shen, Y., Wu, X., and Wang, J. (2014a). Analysis of global gene expression profiles to identify differentially expressed genes critical for embryo development in Brassica rapa. Plant Mol. Biol. 86, 425–442. doi: 10.1007/s11103-014-0238-1
Zhang, Z., Wang, P., Li, Y., Ma, L., Li, L., Yang, R., et al. (2014b). Global transcriptome analysis and identification of the flowering regulatory genes expressed in leaves of Lagerstroemia indica. DNA and Cell Biol. 33, 680–688. doi: 10.1089/dna.2014.2469
Zhang, H., Wei, Y., Shen, J., Lai, B., Huang, X., Ding, F., et al. (2014c). Transcriptomic analysis of floral initiation in litchi (Litchi chinensis Sonn.) based on de novo RNA sequencing. Plant Cell Rep. 33, 1723–1735. doi: 10.1007/s00299-014-1650-3
Zhang, L., Jia, H., Yin, Y., Wu, G., Xia, H., Wang, X., et al. (2013). Transcriptome analysis of leaf tissue of Raphanus sativus by RNA sequencing. PLoS ONE 8:e80350. doi: 10.1371/journal.pone.0080350
Keywords: Raphanus sativus L., bolting and flowering, RNA-Seq, hormone signaling, differentially expressed genes (DEGs)
Citation: Nie S, Li C, Wang Y, Xu L, Muleke EM, Tang M, Sun X and Liu L (2016) Transcriptomic Analysis Identifies Differentially Expressed Genes (DEGs) Associated with Bolting and Flowering in Radish (Raphanus sativus L.). Front. Plant Sci. 7:682. doi: 10.3389/fpls.2016.00682
Received: 08 January 2016; Accepted: 03 May 2016;
Published: 24 May 2016.
Edited by:Sarvajeet Singh Gill, Maharshi Dayanand University, India
Reviewed by:Aashish Ranjan, National Institute of Plant Genome Research, India
Krishna Kant Sharma, Maharshi Dayanand University, India
Narsingh Chauhan, Maharshi Dayanand University, India
Copyright © 2016 Nie, Li, Wang, Xu, Muleke, Tang, Sun and Liu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Liwang Liu, email@example.com