Analysis of Transcriptional Responses of the Inflorescence Meristems in Jatropha curcas Following Gibberellin Treatment

Jatropha curcas L. seeds an oilseed plant with great potential for biodiesel production. However, low seed yield, which was limited by its lower female flowers, was a major drawback for its utilization. Our previous study found that the flower number and female-to-male ratio were increased by gibberellin treatment. Here, we compared the transcriptomic profiles of inflorescence meristem at different time points after gibberellic acid A3 (GA3) treatment. The present study showed that 951 differentially expressed genes were obtained in response to gibberellin treatment, compared with control samples. The 6-h time point was an important phase in the response to exogenous gibberellin. Furthermore, the plant endogenous gibberellin, auxin, ethylene, abscisic acid, and brassinolide-signaling transduction pathways were repressed, whereas the genes associated with cytokinin and jasmonic acid signaling were upregulated for 24-h time point following GA3 treatment. In addition, the floral meristem determinacy genes (JcLFY, JcSOC1) and floral organ identity genes (JcAP3, JcPI, JcSEP1-3) were significantly upregulated, but their negative regulator (JcSVP) was downregulated after GA3 treatment. Moreover, the effects of phytohormone, which was induced by exogenous plant growth regulator, mainly acted on the female floral differentiation process. To the best of our knowledge, this data is the first comprehensive analysis of the underlying transcriptional response mechanism of floral differentiation following GA3 treatment in J. curcas, which helps in engineering high-yielding varieties of Jatropha.


Introduction
Jatropha curcas L. (Euphorbiaceae) is a native species of South America [1]. Depending on geographic and climatic conditions, up to 60% of oil in Jatropha seeds can be used directly or in transesterified form as a biodiesel [2,3]. However, its utilization as bio-energy plant was limited by its low seed yield, which was caused by low number female flowers [4,5]. Thus, seeking a useful way to improve the fruit yield and uncovering the underlying reasons associated with floral differentiation in J. curcas, especially female floral differentiation, were critical for the improvement of the cultivation of high-yielding J. curcas germplasms. The error bar means S.D., which was calculated by fifteen inflorescence meristems after the GA3 and control treatment, respectively. * is significant difference (p < 0.05), ** is very significant difference (p < 0.01).

Illumina Sequencing of Different cDNA Libraries
After the quality control, 6.27-10.46 G clean reads were produced from eighteen cDNA libraries in this study, respectively (Table S1). The error rate of RNA-seq was only 0.01%, all of the Q30 (Phred score) were more than 93%, and the GC content was more than 42% in each sample. In addition, when mapped to the genome data of J. curcas, the alignable reads ranged from 54.39 million (T24_2, 88.16%) to 36.15 million (T6_1, 86.50%) based on clean reads (Table S2). Among the mapped reads, 35.36-53.29 million were uniquely aligned reads, with T24_2 making up the highest percentage (86.39%). Moreover, 8532 novel transcripts were found in present RNA-seq analysis, which could contribute to obtain new genes and exons related to inflorescence meristem differentiation in J. curcas.

DEGs Annotation and Enrichment in Response to GA3 Treatment
The DEGs were detected in four comparisons, T6 vs. CK6, T12 vs. CK12, T24 vs. CK24, and T48 vs. CK48, to perform the transcriptional response analysis of the inflorescence meristems in J. curcas. In total, 951 DEGs were obtained from these four comparisons, and these DEGs were annotated to 43 orthologous terms (Figure 2), including 20 terms involved in biological processes that were mainly The error bar means S.D., which was calculated by fifteen inflorescence meristems after the GA 3 and control treatment, respectively. * is significant difference (p < 0.05), ** is very significant difference (p < 0.01).

Illumina Sequencing of Different cDNA Libraries
After the quality control, 6.27-10.46 G clean reads were produced from eighteen cDNA libraries in this study, respectively (Table S1). The error rate of RNA-seq was only 0.01%, all of the Q30 (Phred score) were more than 93%, and the GC content was more than 42% in each sample. In addition, when mapped to the genome data of J. curcas, the alignable reads ranged from 54.39 million (T24_2, 88.16%) to 36.15 million (T6_1, 86.50%) based on clean reads (Table S2). Among the mapped reads, 35.36-53.29 million were uniquely aligned reads, with T24_2 making up the highest percentage (86.39%). Moreover, 8532 novel transcripts were found in present RNA-seq analysis, which could contribute to obtain new genes and exons related to inflorescence meristem differentiation in J. curcas.

DEGs Annotation and Enrichment in Response to GA 3 Treatment
The DEGs were detected in four comparisons, T6 vs. CK6, T12 vs. CK12, T24 vs. CK24, and T48 vs. CK48, to perform the transcriptional response analysis of the inflorescence meristems in J. curcas. In total, 951 DEGs were obtained from these four comparisons, and these DEGs were annotated to 43 orthologous terms (Figure 2), including 20 terms involved in biological processes that were mainly associated with DNA replication and macromolecules biosynthesis process, 6 terms involved in cellular components that were significantly annotated to telomerase and thylakoid region, and 17 terms involved in molecular functions that were mainly related to the activity of peptidase and microtubule protein. Furthermore, the 951 DEGs were significantly enriched in protein processing in the endoplasmic reticulum, plant hormone signal transduction, photosynthesis, and flavonoid biosynthesis pathways (Figure 3). It indicated that the more actively biological synthesis and metabolism was showed in Jatropha inflorescence after GA 3 treatment. associated with DNA replication and macromolecules biosynthesis process, 6 terms involved in cellular components that were significantly annotated to telomerase and thylakoid region, and 17 terms involved in molecular functions that were mainly related to the activity of peptidase and microtubule protein. Furthermore, the 951 DEGs were significantly enriched in protein processing in the endoplasmic reticulum, plant hormone signal transduction, photosynthesis, and flavonoid biosynthesis pathways ( Figure 3). It indicated that the more actively biological synthesis and metabolism was showed in Jatropha inflorescence after GA3 treatment.    associated with DNA replication and macromolecules biosynthesis process, 6 terms involved in cellular components that were significantly annotated to telomerase and thylakoid region, and 17 terms involved in molecular functions that were mainly related to the activity of peptidase and microtubule protein. Furthermore, the 951 DEGs were significantly enriched in protein processing in the endoplasmic reticulum, plant hormone signal transduction, photosynthesis, and flavonoid biosynthesis pathways ( Figure 3). It indicated that the more actively biological synthesis and metabolism was showed in Jatropha inflorescence after GA3 treatment.   In addition, almost equal DEGs number was obtained between up-and down-regulation at each time point (Figure 4). However, the greatest overall number and larger fold change of DEGs were obtained at 6-h time point. Furthermore, GO (Gene Ontology) enrichment showed that the number of DEGs in each category, such as "metabolic process", "catalytic activity", "organic substance metabolic process", "cellular process", and "primary metabolic process", were highest at the 6-h time point ( Figure S1).All of these results indicated that 6-h time point may be an important phase in the response to exogenous gibberellins in J. curcas. In addition, almost equal DEGs number was obtained between up-and down-regulation at each time point ( Figure 4). However, the greatest overall number and larger fold change of DEGs were obtained at 6-h time point. Furthermore, GO (Gene Ontology) enrichment showed that the number of DEGs in each category, such as "metabolic process", "catalytic activity", "organic substance metabolic process", "cellular process", and "primary metabolic process", were highest at the 6-h time point ( Figure S1).All of these results indicated that 6-h time point may be an important phase in the response to exogenous gibberellins in J. curcas. Moreover, there were 23 DEGs co-detected in the four comparisons ( Figure 5a). Three DEGs of them, JC12650, Novel00988, and Novel00596 (Novel means a new transcript) were upregulated at 6 and 12-h time points, and the Novel00596 was still upregulated in the following time points ( Figure  5b), but no homologous sequence of these three genes was isolated in TAIR10 protein database (Table 1), which suggested that these sequences may be some specific genes involved in response to GA3 treatment in Jatropha inflorescence meristems. Furthermore, JC06413 (JcTT4), JC06540 (JcTT6), and JC22745 (JcTT7) were annotated to flavonoid biosynthetic process, and upregulated at 6-h and 12-h time points, respectively. An important transcription factor associated with the transitions from vegetative to reproductive growth, JcLFY (JC26359), was identified in these 23 DEGs library, and it was also upregulatedat 6-h and 12-h time points. Interestingly, eight DEGs were only significantly upregulated at 24-h time point after GA3 treatment, JC02534, JC08649, JC08881, JC08882, JC08884, JC08888, JC15089, and JC26020 ( Figure 5b). The Blastx with TAIR10 protein database showed that these genes were mainly related to phosphorylase superfamily proteins, which indicated that 24-h point may be the other important phase in the response to exogenous gibberellins. Moreover, there were 23 DEGs co-detected in the four comparisons ( Figure 5a). Three DEGs of them, JC12650, Novel00988, and Novel00596 (Novel means a new transcript) were upregulated at 6 and 12-h time points, and the Novel00596 was still upregulated in the following time points (Figure 5b), but no homologous sequence of these three genes was isolated in TAIR10 protein database (Table 1), which suggested that these sequences may be some specific genes involved in response to GA 3 treatment in Jatropha inflorescence meristems. Furthermore, JC06413 (JcTT4), JC06540 (JcTT6), and JC22745 (JcTT7) were annotated to flavonoid biosynthetic process, and upregulated at 6-h and 12-h time points, respectively. An important transcription factor associated with the transitions from vegetative to reproductive growth, JcLFY (JC26359), was identified in these 23 DEGs library, and it was also upregulatedat 6-h and 12-h time points. Interestingly, eight DEGs were only significantly upregulated at 24-h time point after GA 3 treatment, JC02534, JC08649, JC08881, JC08882, JC08884, JC08888, JC15089, and JC26020 ( Figure 5b). The Blastx with TAIR10 protein database showed that these genes were mainly related to phosphorylase superfamily proteins, which indicated that 24-h point may be the other important phase in the response to exogenous gibberellins.  - All of the genes were annotated with TAIR10 database. Gene ID is the gene number in RNA-seq database of Jatropha curcas. At. locus is the locus of homologous gene in Arabidopsis thaliana, At. name is the gene name of homologous gene in Arabidopsis thaliana. "-" indicated that no homologous genes was detected in TAIR10 protein database.

DEGs Involved in Gibberellin Biosynthesis, Metabolism and Signaling
Based on the above analysis, the DEGs library was significantly enriched in plant hormone signal transduction pathways ( Figure 3). In order to obtain more information about these details, 12 DEGs were identified to perform the further analysis involved in GA biosynthesis, metabolism and signaling pathways. Based on the Blastx with TAIR10 protein database, two DEGs related to bioactive GA formation ( Figure 6a, Table 2), JC04894 (JcGA20OX1) and JC04895 (JcGA20OX2), were not shown differential expression until 24-h time point after GA3 treatment ( Figure 6c). However, the DEGs encoding GA-inactivating enzymes, JcGA2OX2 and JcGA2OX8, were upregulated at  - All of the genes were annotated with TAIR10 database. Gene ID is the gene number in RNA-seq database of Jatropha curcas. At. locus is the locus of homologous gene in Arabidopsis thaliana, At. name is the gene name of homologous gene in Arabidopsis thaliana. "-" indicated that no homologous genes was detected in TAIR10 protein database.

DEGs Involved in Gibberellin Biosynthesis, Metabolism and Signaling
Based on the above analysis, the DEGs library was significantly enriched in plant hormone signal transduction pathways ( Figure 3). In order to obtain more information about these details, 12 DEGs were identified to perform the further analysis involved in GA biosynthesis, metabolism and signaling pathways. Based on the Blastx with TAIR10 protein database, two DEGs related to bioactive GA formation ( Figure 6a, Table 2), JC04894 (JcGA20OX1) and JC04895 (JcGA20OX2), were not shown differential expression until 24-h time point after GA 3 treatment (Figure 6c). However, the DEGs encoding GA-inactivating enzymes, JcGA2OX2 and JcGA2OX8, were upregulated at different time points after GA 3 treatment, especially JC21273 (JcGA2OX8). It indicated that the endogenous gibberellin production was reduced by exogenous GA 3 application in J. curcas.
The GA-GID1-DELLA complex was mediated via the binding of GA to GID1, which caused a rapid degradation of DELLAs through the ubiquitin-proteasome pathway (Figure 6b). The SLY1/GID2, a specific ubiquitin E3 ligase complex, was required to recruit DELLA protein for the subsequent degradation by the 26S proteasome. After that, the GA signaling could be transmitted into the following pathways. In the present study, two Jatropha ortholog genes encoding a GA receptor, JC22004 (JcGID1B) and JC23405 (JcGID1C), were screened in the DEGs dataset after GA 3 treatment (Table 2). Furthermore, JcGID1B was downregulated at each time point following GA 3 treatment, while the JcGID1C was not shown differential expression in T6 vs. CK6 and significantly downregulated at the 12 and 24-h time points (Figure 6c). Additionally, two DEGs involved in GA signal transduction pathways, JC01849 (JcZAR1) and JC23160 (JcPAP3), were also downregulated after GA 3 treatment. Meanwhile, JC20657 (JcRGL1), a gene encoding the DELLA protein, was upregulated at 12 and 48-h time points after GA 3 application. These results suggested that exogenous GA 3 reduced the endogenous gibberellin activity by downregulating GA receptor genes and upregulating suppressors of GA signaling. different time points after GA3 treatment, especially JC21273 (JcGA2OX8). It indicated that the endogenous gibberellin production was reduced by exogenous GA3 application in J. curcas. The GA-GID1-DELLA complex was mediated via the binding of GA to GID1, which caused a rapid degradation of DELLAs through the ubiquitin-proteasome pathway (Figure 6b). The SLY1/GID2, a specific ubiquitin E3 ligase complex, was required to recruit DELLA protein for the subsequent degradation by the 26S proteasome. After that, the GA signaling could be transmitted into the following pathways. In the present study, two Jatropha ortholog genes encoding a GA receptor, JC22004 (JcGID1B) and JC23405 (JcGID1C), were screened in the DEGs dataset after GA3 treatment (Table 2). Furthermore, JcGID1B was downregulated at each time point following GA3 treatment, while the JcGID1C was not shown differential expression in T6 vs. CK6 and significantly downregulated at the 12 and 24-h time points (Figure 6c). Additionally, two DEGs involved in GA signal transduction pathways, JC01849 (JcZAR1) and JC23160 (JcPAP3), were also downregulated after GA3 treatment. Meanwhile, JC20657 (JcRGL1), a gene encoding the DELLA protein, was upregulated at 12 and 48-h time points after GA3 application. These results suggested that exogenous GA3 reduced the endogenous gibberellin activity by downregulating GA receptor genes and upregulating suppressors of GA signaling.   Signaling pathways of phytohormones showed complex interaction networks to regulate floral differentiation in plants. In order to understand the roles of other phytohormones in inflorescence meristems in response to GA application, we identified homologous genes related to various hormonal transduction pathways in the 951 DEGs library by Blastx with TAIR10 database. In total, 20 DEGs were obtained to conduct the further analysis associated with auxin (IAA), cytokinin (CTK), jasmonic acid (JA), ethylene (ETH), abscisic acid (ABA), and brassinosteroid (BR)-signaling transduction pathways (Table 3). Nine DEGs were involved in auxin-activated signaling pathway (Table 3)      Two DEGs, JC23402 (JcAHP1) and JC17975 (JcHK3), were involved in cytokinin-activated signaling pathway (Table 3), and both of them were upregulated at 6 and 12-h time points (Figure 7a). JcAHP1, a gene encoding the histidine phosphotransfer protein 1, was the vital factor related to the CTK signaling transduction from cytoplasm to nucleus. JcHK3, a gene encoding the histidine kinases, acted as the cytokinin sensor in CTK signaling pathways. The upregulation of these genes could significantly enhance the cytokinin signaling after GA 3 application.
Meanwhile, three DEGs, JC14204 (JcMYC4), JC12669 (JcMYC2), and JC23114 (JcJAZ10), were associated with jasmonic acid-signaling pathway. Similarly with the CTK signaling after GA 3 application, JC14204 (JcMYC4) and JC12669 (JcMYC2) were upregulated at 6 and 12-h time points (Figure 7a), which were related to Jasmonic acid-activated signaling pathway (Table 3). However, JC23114 (JcJAZ10), a gene encoding jasmonate ZIM-domain (JAZ) proteins, which was a key negative regulator of JA signaling, was upregulated at 6-h time point after GA 3 treatment. These results indicated that the intricate agonistic interactions were shown between GA and JA in the differentiation of inflorescence meristems in J. curcas. The exogenous gibberellin could induce the endogenous jasmonicacid-signaling pathway within a certain range.
Based on the Blastx with TAIR10 database, three DEGs were related to ethylene-signaling transduction pathway (Figure 7a). JC13342 and JC07165 were annotated to Jatropha orthologs of EBF1, which was involved in the negative regulation of ethylene-activated signaling pathway ( Table 3). Both of them were upregulated at 6-h time point after GA 3 treatment (Figure 7a). Furthermore, JC17700 (JcERF1), a key ethylene-responsive transcription factor, was significantly downregulatedat 6-h time point, and then it was not detected differential expression at the following time points after GA 3 treatment. These results indicated that a negative crosstalk was shown between GA and ETH-signaling pathways in J. curcas.
Moreover, two DEGs, JC14274 (JcPYL9) and JC02934 (JcHAI2), were involved in abscisic acid-signaling transduction pathway (Table 3). JC14274 (JcPYL9), an ABA receptor belongs to PYR1-like, was associated with abscisic acid-activated signaling pathway. JC02934 (JcHAI2), a gene encoding the protein phosphatase 2C, was a negative regulator of ABA signaling. Both of them were downregulated at 6-h time point after GA 3 treatment (Figure 7a), whereas JC02934 (JcHAI2) was significantly upregulated at 12 and 48-h time points. In addition, one DEG, JC22124 (JcCYCD3), was detected and annotated to brassinosteroid-signaling transduction pathway. It was downregulated at 6 and 12-h time points after GA 3 application (Figure 7a).These results implied that ABA and BR signaling was repressed in the inflorescence differentiation process of J. curcas after exogenous GA 3 treatment.

MADS-Box Transcription Factor in Response to GA 3 Treatment
Gibberellin could regulate the floral differentiation by interacting with MADS-box transcription factors, which play vital regulatory roles in the floral organ differentiation process. 11 differentially expressed GA-responsive transcription factors were identified in the present study (Table 4). Based on the Blastx with TAIR10 database, two DEGs, JC15882 and Novel00351, were annotated to Jatropha orthologous SVP, a negative regulator of flowering. Both of them were downregulated at 6 and 48-h time points after GA 3 treatment (Figure 7b). However, JC26359 (JcLFY), a positive regulator of flowering, was significantly upregulated at 6, 12, and 48-h time points. Meanwhile, JcSOC1 (JC14482), an upstream gene of JcLFY, was also upregulated at 24-h time point after GA 3 treatment. It indicated that the flowering time could be accelerated by exogenous GA 3 treatment in J. curcas, which was consistent with the previous studies [28,29].
Additionally, floral organ identity genes played important roles in floral sex differentiation process. In our previous study, ABCDE model was formulated to regulate the floral organ identification in J. curcas [30]. In this study, 3 DEGs of E-function, JC11754, JC25593 and JC17987, were isolated and annotated to JcSEP1, JcSEP2, and JcSEP3, respectively ( Table 4). All of them were upregulated at 6 and 12-h time points, whereas they were downregulated at 24-h time point (Figure 7b). Moreover, 3 DEGs of B-function, JC04507 (JcAP3), JC13660 (JcAP3), and JC12153 (JcPI) were identified in the 951 DEGs library. Similarly with the DEGs in class E, JC13660 and JC12153 were significantly upregulated at 6 and 12-h time points and downregulated at 24-h time point (Figure 7b). It indicated that exogenous GA could promote the B-and E-function transcription factor, which might be a reason that the floral number was increased after GA 3 treatment.

DEGs Involved in Male and Female Floral Differentiation
In order to obtain more information about the DEGs involved in male and female floral differentiation after GA 3 treatment, the 951 DEGs detected in the present study was combined with our previous transcriptome dataset (Gene Expression Omnibus number: GSE102894) related to floral sex differentiation process [30], including male floral initiation stage (STD1 vs. IND, IND was the stage of inflorescence meristems, and STD1 was the stage of male floral initiation), male floral development stage (STD2 vs. STD1, STD2 was the stage of ten complete stamens formed), female floral initiation stage (PID1 vs. IND, PID1 was the stage of female floral initiation), and female floral development stage (PID2 vs. PID1, PID2 was the stage of complete carpel and ovary formed).
The results showed that 55 DEGs were isolated and associated with male floral initiation after GA 3 application ( Figure S2). Four DEGs of them (JC11710, JC18282, JC20688, JC21298) were significantly annotated to flavonoid biosynthetic process by KEGG enrichment analysis ( Figure S3a). Interestingly, these DEGs were mainly upregulated both in T6 vs. CK6 and STD1 vs. IND. Moreover, 32 DEGs were co-detected in the DEGs libraries of GA 3 treatment and STD1 vs. IND and STD2 vs. STD1, and 54 DEGs were identified both in the DEGs libraries of GA 3 treatment and STD2 vs. STD1 ( Figure S2). KEGG enrichment showed that these genes were mainly involved in DNA replication, protein biosynthesis, and sugar metabolism pathways ( Figure S3b), indicating that these genes were contributed to the male flower development after GA 3 treatment.
There were 16 DEGs co-detected in the DEGs libraries of GA 3 treatment and PID1 vs. IND ( Figure S2). Based on the KEGG enrichment analysis (Figure S4a), two DEGs of them, JC23402 (JcAHP1) and JC04692 (JcPIN3), were significantly associated with CTK and IAA signal transduction pathways, respectively (Table 3). In addition, eighteen DEGs were isolated both in the DEGs libraries of GA 3 treatment and PID1 vs. IND and PID2 vs. PID1, while one hundred and ten DEGs were obtained in the DEGs libraries of GA 3 treatment and PID2 vs. PID1. The number of these DEGs was about 1.5 times of the corresponding period in male floral differentiation process ( Figure S2). KEGG enrichment showed that these 128 DEGs were significantly related to plant hormone-signaling transduction pathways ( Figure S4b). It indicated that the female floral differentiation was significantly involved in plant hormone-signaling transduction pathways after exogenous GA 3 treatment.

qRT-PCR Validation
To validate the RNA-seq data, 48 DEGs detected by T6 vs. CK6, T12 vs. CK12, T24 vs. CK24, and T48 vs. CK48 were tested by qRT-PCR (Table S3). These genes were selected because of their important function identified in this study, including 20 downregulated genes and 28 upregulated genes. All of them were consistent with the same trend of upregulation or downregulation between the two different expression analysis platforms (Figure 8). The correlation of the two expression measurements was 0.83 between these 48 genes (R 2 = 0.83). In a word, the results of RNA-seq and qRT-PCR were consistent. significantly annotated to flavonoid biosynthetic process by KEGG enrichment analysis ( Figure  S3a). Interestingly, these DEGs were mainly upregulated both in T6 vs. CK6 and STD1 vs. IND. Moreover, 32 DEGs were co-detected in the DEGs libraries of GA3 treatment and STD1 vs. IND and STD2 vs. STD1, and 54 DEGs were identified both in the DEGs libraries of GA3 treatment and STD2 vs. STD1 ( Figure S2). KEGG enrichment showed that these genes were mainly involved in DNA replication, protein biosynthesis, and sugar metabolism pathways ( Figure S3b), indicating that these genes were contributed to the male flower development after GA3 treatment. There were 16 DEGs co-detected in the DEGs libraries of GA3 treatment and PID1 vs. IND ( Figure S2). Based on the KEGG enrichment analysis (Figure S4a), two DEGs of them, JC23402 (JcAHP1) and JC04692 (JcPIN3), were significantly associated with CTK and IAA signal transduction pathways, respectively (Table 3). In addition, eighteen DEGs were isolated both in the DEGs libraries of GA3 treatment and PID1 vs. IND and PID2 vs. PID1, while one hundred and ten DEGs were obtained in the DEGs libraries of GA3 treatment and PID2 vs. PID1. The number of these DEGs was about 1.5 times of the corresponding period in male floral differentiation process ( Figure S2). KEGG enrichment showed that these 128 DEGs were significantly related to plant hormone-signaling transduction pathways ( Figure S4b). It indicated that the female floral differentiation was significantly involved in plant hormone-signaling transduction pathways after exogenous GA3 treatment.

qRT-PCR Validation
To validate the RNA-seq data, 48 DEGs detected by T6 vs. CK6, T12 vs. CK12, T24 vs. CK24, and T48 vs. CK48 were tested by qRT-PCR (Table S3). These genes were selected because of their important function identified in this study, including 20 downregulated genes and 28 upregulated genes. All of them were consistent with the same trend of upregulation or downregulation between the two different expression analysis platforms (Figure 8). The correlation of the two expression measurements was 0.83 between these 48 genes (R 2 = 0.83). In a word, the results of RNA-seq and qRT-PCR were consistent.

Discussion
GA signaling is one of the important regulated factors in the network of floral induction pathways [9]. Recently, Chen et al. [31] found that the exogenous application of GA 3 could partially prevent pistil development to generate neutral flowers without stamens and pistils in gynoecious Jatropha. However, Makwana et al. [18] found that the female flower and fruit yield could be increased by exogenous GA 3 in wild J. curcas. Our previous study was consistent with the Makwana's research [8]. From the present study, the female flower number and female-to-male flower ratio were significantly increased after GA 3 treatment (Figure 1). Therefore, we supposed that whether the application of exogenous GA 3 to J. curcas inflorescences stimulates floral differentiation may be cultivar-dependent.
To advance our understanding of underlying GA-induced responses in J. curcas, the transcriptome analysis of Jatropha inflorescence meristems were carried out after GA 3 treatment. We also performed the control samples at the same conditions to remove some genes that were differential expression following the plant growth, but not the GA 3 treatment between different time points. Moreover, different plant species were given specific response times for various plant growth regulators [32]. In the previous study, GA content of the berry was substantially increased for 24 h following GA 3 application in grapevine flower [33]. Chen et al. [34] found that the DEGs were significantly increased at 12-h time point after exogenous cytokinin treatment in J. curcas. In order to determine the appropriate time point in the present study, we collected the five time point samples from 0 h to 48 h after 40 mg/L GA 3 treatment in J. curcas. Fortunately, 6.27-10.46 G clean data were obtained in GA 3 treated and untreated samples, respectively. 951 DEGs were isolated at the different time points, compared with control samples. In addition, 8532 novel transcripts were found in present RNA-seq analysis, which could contribute to obtain new genes and exons related to inflorescence meristem differentiation in J. curcas. Furthermore, we supported that 6-h time point was an important phase in the response to exogenous gibberellin in this study.
Gibberellins regulate various aspects of plant growth and development [35], it also play a major role in the network of floral induction pathways [36]. In the present study, the DEGs annotated to GAs biosynthesis and signaling were further analyzed. We found that GA20OX1 (JC04894) and GA20OX2 (JC04895), which catalyzed the formation of bioactive GAs [37], were not detected differential expression until 24-h time points after GA 3 treatment (Figure 6c). In contrast, the genes encoding GA2OX, a major GA inactive catabolic oxidase [11], were upregulated after GA 3 treatment. It indicated that exogenous GA could inhibit the endogenous GA biosynthesis process. This finding was consistent well with previous studies [33,38]. Therefore, we suggested that the content of bioactive GAs might show a feedback regulation to response exogenous GA 3 application in J. curcas. The other evidence supporting this suggestion was that the GA receptor genes, JcGID1B (JC22004) and JcGID1C (JC23405), were significantly downregulated at 6 and 24-h time points after GA 3 application in the present study, but JcRGL1 (JC20657), a gene encoding DELLA protein (Figure 6b,c), was upregulated at 12-h time point, which agreed with the previous researches [39,40]. In further study, a strategic approach is to investigate the concentration of endogenous gibberellins after exogenous GA 3 treatment in J. curcas.
In addition, gibberellins could jointly regulate the floral differentiation by crosstalk with various phytohormone-signaling pathways in J. curcas. In previous study, the exogenous CTK application could significantly downregulate the expression of the genes involved in ABA, and ETH signal transduction pathways in J. curcas [19]. The similar results in the present study were that the ABA, ETH, and BR signaling were also repressed after exogenous GA treatment. It indicated that these phytohormone-signalings were not important for floral differentiation process in J. curcas. Moreover, we detected that the antagonistic interactions were existed between GA and IAA in the development of inflorescence meristems in J. curcas. That might be the underlying reason why the female flowers were induced by spraying exogenous GA onto the inflorescence, but the highest concentration of GA resulted in the withering of the Jatropha inflorescence in previous studies [8,18].
However, the DEGs involved in the JA-signaling transduction pathway, such as JcMYC2 (JC12669) and JcMYC4 (JC14204), were upregulated after GA 3 treatment in the present study. The exogenous CTK could also promote the DEGs related to JA-signaling transduction pathway in J. curcas [34]. It has been reported that JA play an important role in pollen development process [41], which might be a reason that JA content was upregulated both in GA and CTK treatment. Interestingly, exogenous cytokinin could increase endogenous GA and CTK transcription in previous study [34], but exogenous gibberellins promoted the endogenous CTK signaling and repressed the GA production in this study (Figures 6c and 7a). It indicated that both CTK and GA were necessary for floral differentiation in J. curcas. The consistent results were that the floral transition was mediated by coordinate regulation of CTK and GAs activities in different species [42][43][44]. Furthermore, exogenous CTKs induced amounts of female flowers [5,7], but exogenous GAs just increased about 2-fold female flowers in the present study. Therefore, we suggested that the ratio of CTK and GA might significantly affect the female floral differentiation. The high percentage of CTK to GA may be appropriate for improving the female flower and increasing the fruit yield, and vice versa. A further study could focus on finding the best proportion of CTKs and GAs to obtain a high fruit yield in J. curcas.
Our previous study found that male floral initiation was associated with the flavonoid biosynthesis pathway, while female floral initiation was related to the phytohormone signal transduction pathway [30]. We combined the 951 DEGs obtained in this study with the previous DEGs dataset (GSE102894) detected in male or female floral initiation and development stages, respectively ( Figure S2). The consistent results were that the co-detected DEGs in the libraries of the male floral initiation stage and GA 3 treatment were also significantly enriched in flavonoid biosynthesis pathways, while the co-detected DEGs in the libraries of female floral differentiation and GA 3 treatment were also related to the phytohormone signal transduction pathways ( Figures S3 and S4). In previous studies, flavonoid biosynthesis inhibition resulted in the male sterility in Petunia [45], but not in Arabidopsis thaliana [46]. This suggested that flavonoids might act various roles in different species. Thus, a further study could focus on exploring the Jatropha flavonoid pathways to uncover the regulatory events associated with male floral differentiation. Additionally, the DEGs after exogenous CTK treatment were also significantly associated with phytohormone-signaling transduction pathways in J. curcas [19,34]. It indicated that the effects of phytohormone, which was induced by exogenous plant growth regulator, were mainly acted on female floral differentiation process in J. curcas.

Plant Material and GA Treatment
Three-year-old Jatropha clone was selected as the experimental material, which was planted in a field with latosolic red soil at a forestry trial base of South China Agricultural University (23.24 • N, 113.64 • E). The general climate has an average temperature of 20-22 • C and annual rainfall of 1720 mm in the rainy season from April to June. In our previous study, we found that spraying 40 mg/L GA 3 onto the inflorescence meristems of J. curcas could significantly increase the flower number and the female-to-male flower ratio [8]. In the present study, the inflorescence meristems (about 0.5 cm in diameter) were sprayed with 5 mL GA 3 working solution (40 mg/L) containing 0.05% (v/v) Tween-20, and the control inflorescence meristems were sprayed with distilled water containing 0.05% (v/v) Tween-20. For the RNA-seq analysis, the mixed samples of three inflorescence meristems were collected at 0, 6, 12, 24, and 48-h after the GA 3 and control treatment, respectively. Two biological replicates were performed for each time point as in other studies [47][48][49]. All of the samples were immediately frozen in liquid nitrogen and stored at −80 • C until further used for RNA extraction. Additionally, fifteen inflorescence meristems of the GA 3 and control treatment were used to investigate the effects of GA 3 on branching of inflorescence in J. curcas. The experiments were carried out in September 2016, and all of the GA 3 and control treatments were performed between 8:00 a.m. and 9:00 a.m.

Quantitative Real-Time PCR Validation
To validate the results of transcriptome, a total of 48 DEGs were selected for quantitative real-time PCR (qRT-PCR) analysis using the same plant materials used for RNA sequencing. The mixed samples of three inflorescence meristems were conducted for the DEG validation at different time points after the GA 3 and control treatment. In addition, two biological replicates were carried out for each time point with three technological replications for each gene. These genes were selected because of their important function after GA 3 treatment according to DEGs analysis in present study. The cDNA synthesis was conducted using PrimeScript ® II first Strand cDNA Synthesis Kit (TaKaRa, Kyoto, Japan). The specific primers were designed by Primer Premier 5.0 (PREMIER Biosoft, Palo Alto, CA, USA) (Table S3). qRT-PCR was performed on Roche LightCyler480 system (Roche, Basel, Switzerland) with SYBR Premix Ex Taq TM II (TaKaRa) [30]. The qRT-PCR procedure was as follows: 2 µL of cDNA dilution in H 2 O (about 50 ng/µL) was added to 10 µL of 2× SYBR ® buffer, with 1 µL (10 µM) of each primer and H 2 O to a final volume of 20 µL. The cycling reaction was 94 • C for 2 min, followed by 40 cycles of 94 • C for 10 s, 55 • C for 10 s and 72 • C for 20 s. JcGAPDH and β-actin (Jcactin) were used as internal controls [59]. The 2 −∆∆Ct method was used to calculate the relative expression level of the DEGs [60].

Conclusions
The present study performed the transcriptional response analysis of the inflorescence meristems in Jatropha curcas following gibberellin treatment. Our results showed that 6.27-10.46 G clean data were obtained in GA 3 treated and untreated samples, respectively. 951 DEGs were isolated at the different time points, compared with control samples. The 6-h time point was an important phase in the response to exogenous gibberellins in the present study. In addition, 8532 novel transcripts were found in present RNA-seq analysis, which could contribute to obtain new genes and exons related to inflorescence meristem differentiation in J. curcas. Furthermore, the plant endogenous IAA, ETH, ABA, and BR signaling were repressed after exogenous gibberellins treatment. The exogenous GA could also inhibit the endogenous GA biosynthesis and signaling pathways. However, the DEGs associated with a JA signal transduction pathway were upregulated following GA 3 treatment to contribute to pollen development process. Both CTK and GA might be necessary for floral differentiation in J. curcas; the ratio of CTK to GA significantly affected the female floral differentiation and fruit yield. Moreover, the floral meristem determinacy genes (JcLFY and JcSOC1) and floral organ identity genes (JcAP3, JcPI, and JcSEP1-3) were significantly upregulated, but their negative regulator (JcSVP) was downregulated after GA 3 treatment. Additionally, the effects of phytohormone, which was induced by exogenous plant growth regulator, mainly acted on female floral differentiation process in J. curcas.
This data will create a reference transcriptome for the genomics database of J. curcas for future studies.Our study will contribute to understanding the underlying transcriptional response mechanism of floral differentiation following GA 3 treatment. We expect that the dataset will serve as a foundation to study the genes' functions, which help in engineering high-yielding varieties in J. curcas.