Introduction
============

Flowering plants show extreme diversity in their reproductive systems. The determination of the environmental and genetic factors that control sexual diversification is a central problem in plant evolutionary biology ([@B8]). Plant sex determination that leads to the physical separation of male and female reproductive organs plays a crucial role in the evolutionary process from hermaphroditic to unisexual flowers ([@B21]; [@B64]). In order to facilitate outcrossing and more efficient allocation of resources to optimization of reproduction ([@B9]; [@B66]), approximately 10% of angiosperm species have produced unisexual flower reproductive systems in dioecious and monoecious plant taxa ([@B72]). Genetic, hormonal, and environmental cues are involved in the regulation of sex determination in flowering plants ([@B5]). Sex determination has evolved multiple times independently suggesting that various different genetic mechanisms underlie unisexual flower development ([@B2]; [@B18]; [@B24]).

In cucumber, sex determination is controlled by three genetic loci: *F*, *A*, and *M*. The *F* locus controls the degree of femaleness, the *A* locus increases maleness and the *M* locus is responsible for the selective abortion of stamens ([@B53]; [@B51]). The *F* locus is linked to the *1-aminocyclopropane-1 carboxylic acid synthase* (*ACS*) gene encoding the rate-limiting enzyme in ethylene (ETH) biosynthesis, which is present as a single copy in monoecious lines but as a duplicated copy in gynoecious lines ([@B39]; [@B48]). The *M* locus encodes the *CsACS2* gene, in which a conserved residue conversion (Gly33Cys) causes bisexual flower production ([@B44]). The ETH receptor *ETR1* is suppressed, causing female flower generation in *Arabidopsis*, suggesting that ETH perception is involved in stamen abortion ([@B69]). Sex determination in melon is controlled by two loci, *andromonoecious* (*a)* and *gynoecious* (*g)* ([@B55]), and their interaction results in a range of sexual types ([@B41]). The *a* locus is linked to the *GmACS-7* gene, in which a single missense mutation causes loss of enzymatic activity, resulting in the development of male organs to generate andromonoecious plants ([@B12]). CmACS-7 mainly represses stamen development and is not required for carpel development, although CmACS-7 is expressed in the carpel primordia during melon flower development ([@B12]). The *g* locus is links to a transcription factor, CmWIP1, which is a repressor of carpel development. The activation of *CmWIP1* by promoter methylation results in the transition from male to female flowers in gynoecious plants ([@B46]). *CmACS-11*, encoding a key enzyme of ETH biosynthesis, represses *CmWIP1* expression, and *CmACS-11* loss-of-function mutation lead to the transition from monoecious to androecious individuals ([@B13]). *CmACS-11, CmWIP1*, and *CmACS-7* coregulate sex determination in melon. In maize, seven mutants affecting male flower (tassel) development have been identified, including six *tasselseed* mutants and one *na1* mutant ([@B75]). *tasselseed1* (*ts1*) encodes a lipoxygensae involved in jasmonic acid (JA) biosynthesis. Loss of functional *ts1* causes a decrease in endogenous JA concentrations in developing inflorescences, and exogenous JA can rescue stamen development in *ts1* mutants ([@B1]). *tasselseed2* (*ts2*) encodes a short-chain alcohol dehydrogenase, required for the arrest of pistil primordia development in the maize tassel ([@B23]). *tasselseed4* (*ts4*) encodes a miR172 microRNA targeting *tasselseed6* (*ts6*)/*indeterminate spikelet1* (*ids1*), which is an *APETALA2* (*AP2*)-like gene. *ts4* and *ids1* are required for the carpel abortion in the tassel ([@B19]). *nana plant1* (*na1*) encodes a 5α-steroid reductase involved in brassinosteroid (BR) biosynthesis, and a mutation of *na1* causes a decrease in BR metabolites, leading to dwarf and feminizing phenotypes ([@B34]). In addition, *anther ear* (*an1*), *dwarf1* (*d1)*, *d2*, *d3*, *d5*, and *d8* mutants, which are deficient in gibberellic acid (GA) biosynthesis or perception, produce perfect florets in the tassel ([@B22]), and GA application to *an1* converts male to gynoecious florets in the tassel ([@B11]), indicating that GA is a negative regulator during stamen development in maize ([@B64]; [@B16]). However, GA is a positive regulator in the development of stamens and pollen in *Arabidopsis* and rice ([@B17]) and is considered to exert masculine effects in plant sex expression ([@B5]). In *Diospyros*, the autosomal *MeGI* (Japanese for "female tree") is a feminizing gene controlling anther fertility in a dosage-dependent fashion. A Y-chromosome *OGI* (Japanese for "male tree") encoding a small RNA suppresses *MeGI* expression, resulting in the generation of androecious individuals ([@B3]). These results suggest that at least three different regulation mechanisms of sex determination are present in flowering plants.

The exogenous application of auxin, BR, cytokinin (CK), ETH, GA, JA, and their inhibitors affects sex expression in flowering plants ([@B45]; [@B14]; [@B25]; [@B37]; [@B51]), and some of them have opposing effects on sex determination in different plant lineages ([@B71]). Environmental cues such as temperature, photoperiod, nutrition, drought, pH, and seasonality also effect sex determination in many species ([@B42]; [@B27]; [@B30]). However, the mechanism of environmental sex determination remains largely unknown, and epigenetic mechanisms are likely involved in the process ([@B54]).

*Jatropha curcas* L. is known as a potential biofuel plant, owing to its high seed oil content ([@B26]; [@B60]). Most *J. curcas* germplasms are monoecious ([@B57]), bearing male and female flowers on the same inflorescence. Low female:male flower ratio of approximately 1:13--29 in monoecious plants has been considered to be one of the factors causing low seed yield in *J. curcas* ([@B50]). Therefore, it is crucial to improve the female flower number per inflorescence for the increase in *J. curcas* seed yield. Previous research showed that the exogenous application of CK (6-benzyladenine, BA) on inflorescence buds significantly increases the female flower number and female:male flower ratio, resulting in an increase in seed yield ([@B50]; [@B15]; [@B28]). These results suggest that CKs are involved in *Jatropha* sex differentiation. In this study, a gynoecious genotype of *J. curcas*, whose male flowers are aborted at early stage of inflorescence development, was used for a comparative transcriptome analysis between gynoecious and monoecious inflorescences, to identify the regulation networks that participate in sex differentiation of *J. curcas.* The findings of this study will be helpful in understanding the sex determination mechanism in *J. curcas*.

Materials and Methods {#s1}
=====================

Application of Gibberellin (GA) and Collection of Samples for RNA Isolation
---------------------------------------------------------------------------

Gynoecious and monoecious plants were grown in the field in the Xishuangbanna Tropical Botanical Garden of the Chinese Academy of Sciences, located in Menglun town in Mengla County (21° N, 101° E), Yunnan province, China. Gynoecious inflorescence buds that occurred about 1--7 days were treated with 1.5 mM GA~3~ solution. Inflorescence buds that occurred about 3--4 days (stage I) and 8--9 days (stage II) were harvested from gynoecious and monoecious plants, respectively, for RNA isolation (**Figure [1](#F1){ref-type="fig"}**). Three biological replications were performed for each sample group. Morphologies of female and male flowers at developmental stages I and II were observed by a 3D super depth digital microscope (ZEISS Smartzoom 5).

![**Morphological comparison of inflorescences between gynoecious and monoecious plants at different developmental stages.** **(A--C)** Indicate gynoecious inflorescences at stages I, II, and III, respectively. **(D--F)** Indicate monoecious inflorescences at stages I, II, and III. Stage I spans approximately 1--7 days from the formation of visible inflorescence buds in which female and male flowers are normally develop in gynoecious plants. Stage II spans approximately 8--20 days from the formation of visible inflorescence buds in which female flowers are normally develop and male flowers are arrest in gynoecious plants. Stage III spans approximately 21--30 days from the formation of visible inflorescence buds in which gynoecious inflorescences have blooming female flowers. Red arrows indicate female flowers and blue arrows indicate male flowers; bar, 0.5 cm.](fpls-07-01953-g001){#F1}

RNA Isolation and Transcriptome Sequencing
------------------------------------------

Total RNA was extracted from inflorescences using the plant total RNA extraction reagent pBIOZOL (Bioflux, Germany). RNA purity was checked with a NanoPhotometer spectrophotometer (IMPLEN, Westlake Village, CA, USA). RNA concentration was determined with a Qubit 2.0 Fluorometer (Life Technologies, Carlsbad, CA, USA). RNA integrity was assessed with an Agilent Bioanalyzer 2100 system (Agilent Technologies, Santa Clara, CA, USA). Sequencing libraries were prepared with a NEBNext Ultra RNA Library Prep Kit for Illumina (NEB, USA) following the manufacturer's specifications. Library quality was assessed on the Agilent Bioanalyzer 2100 system. The clustering of the index-coded samples was performed on a cBot Cluster Generation System according to the manufacturer's instructions. After cluster generation, the libraries were sequenced on an Illumina HiSeq 2500 platform and 125-base paired-end reads were produced. Sequencing was performed by Novogene Bioinformatics Technology (Beijing, China).

*De novo* Transcriptome Assembly and Abundance Estimation
---------------------------------------------------------

Low-quality (with \< 20 Phred scores) reads in raw reads were removed with Fastq_clean ([@B76]) and assessed with FASTQC^[1](#fn01){ref-type="fn"}^. Clean reads were assembled using Trinity (version 2.0.6) with default parameters ([@B31]; [@B32]). The paired-end reads from each library were mapped to the *de novo* assemblies using Bowtie (version 1.1.1) ([@B43]), and abundance estimation of transcripts was performed using Corset (version 1.03) ([@B20]).

Identification of Differentially Expressed Transcripts
------------------------------------------------------

Cluster-level count data generated from Corset were processed with the edgeR package ([@B58]). Very low-count transcripts were removed and transcripts with at least one count per million reads (CPM) in at least three libraries were kept for analysis. A design matrix fitting generalized linear models (GLMs) was constructed and the Cox--Reid profile-adjusted likelihood (CR) method was used to estimate common, trended, and tagwise dispersions. Differential expressed transcripts were determined using GLM likelihood ratio test and transcripts with a false discovery rate (FDR) of \<0.001 were reported for annotation analysis. A Venn diagram of differentially expressed genes was generated using Venny diagrams^[2](#fn02){ref-type="fn"}^.

Transcript Annotation and Clustering Analysis
---------------------------------------------

Gene Ontology (GO) classification and enrichment analysis were performed with DAVID ([@B35]). Hierarchical clustering of genes was performed using the pheatmap R package (version 1.0.7)^[3](#fn03){ref-type="fn"}^.

Real Time Quantitative PCR (qPCR)
---------------------------------

Total RNA was extracted using the plant total RNA extraction reagent pBIOZOL (Bioflux, Germany). For each sample, total RNA (1 μg) was used for cDNA synthesis using PrimeScript kit (Takara Biotechnology, China). qPCR was performed on Roche Light Cycler 480II system (Roche) using LightCycler FastStart DNA MasterPLUS SYBR Green I Kit (Roche). Each sample has two or three biological replicates. Primers used for qPCR analysis were listed in **Supplementary Table [S7](#SM7){ref-type="supplementary-material"}**.

Availability of Supporting Data
-------------------------------

All RNA-Seq data for this project have been deposited in NCBI under the SRA accession SRP092157. For simplicity, only sequences of genes that appeared in the context of this paper are listed in **Supplementary Data Sheet [1](#SM8){ref-type="supplementary-material"}**.

Results
=======

Morphological Comparison between Gynoecious and Monoecious Plants
-----------------------------------------------------------------

To facilitate comparative research on inflorescence development between gynoecious and monoecious plants, the process of inflorescence development was roughly classified into three stages on the basis of inflorescence morphology. In stage I, spanning 1--7 days from the formation of visible inflorescence buds, both female and male flowers could develop in gynoecious plants that is similar to the monoecious plants (**Figure [1](#F1){ref-type="fig"}**). In this stage, 10 stamen primordia were found in an inner-outer whorl arrangement in male flowers, and three carpel primordia emerged in female flowers in both gynoecious and monoecious plants (**Figure [2](#F2){ref-type="fig"}**). Interestingly, stamen primordia also presented in monoecious female flowers, but were invisible in gynoecious ones (**Figure [2](#F2){ref-type="fig"}**). In stage II, spanning 8--20 days, the development of male flowers in gynoecious plants was clearly arrested and stamen primordia turned brown, whereas in monoecious male flowers stamen primordia developed into anthers (**Figures [1](#F1){ref-type="fig"}** and **[2](#F2){ref-type="fig"}**). In this stage, carpels closed up to develop into ovaries, and nectaries were formed in both gynoecious and monoecious female flowers; stamen primordia in monoecious female flowers were degenerated (**Figure [2](#F2){ref-type="fig"}**). In stage III, spanning 21--30 days, gynoecious inflorescences had blooming female flowers and vestigial male flowers, whereas monoecious inflorescences consisted of normal female and male flowers (**Figure [1](#F1){ref-type="fig"}**). Moreover, leaves of gynoecious plants were larger than ones of monoecious plants; gynoecious fruits were spindle-shaped, whereas monoecious fruits are round, but their seeds are similar in size and shape (**Supplementary Figures [S1](#SM9){ref-type="supplementary-material"}** and **[S2](#SM10){ref-type="supplementary-material"}**).

![**Morphological comparison of female and male flowers between gynoecious and monoecious plants at developmental stages I and II.** P, petal; St, stamen; Pi, pistil; N, nectary; red and blue arrows indicate stamen and carpel primordia, respectively; bar, 0.5 mm.](fpls-07-01953-g002){#F2}

Identification of Differentially Expressed Genes
------------------------------------------------

To investigate the mechanism of sex determination in *J. curcas*, we collected four groups of samples containing monoecious inflorescences at stages I and II (mI and mII) and gynoecious inflorescences at stages I and II (gI and gII) for RNA-Seq analysis to identify differentially expressed genes that might be involved in the transition between monoecy and gynoecy (**Supplementary Table [S1](#SM1){ref-type="supplementary-material"}**). A plot that shows the relationship between the four groups of samples based on multidimensional scaling with the biological coefficient of variation, displays that gynoecious and monoecious samples are well separated, indicating the presence of large numbers of differentially expressed genes (**Figure [3](#F3){ref-type="fig"}**). For simplicity, we will hereafter refer to the comparison of stage II vs. stage I in monoecious inflorescences as the pairwise mII vs. mI, the stage II vs. I in gynoecious inflorescences as pairwise gII vs. gI, the gynoecious vs. monoecious inflorescence at stage I as pairwise gI vs. mI, and the gynoecious vs. monoecious inflorescence at stage II as pairwise gII vs. mII. In all, 3,707 genes were significantly changed at an FDR of \<0.001 in four pairs (**Supplementary Table [S2](#SM2){ref-type="supplementary-material"}**). In pairwise mII vs. mI, only three differentially expressed genes were upregulated and 28 were downregulated, whereas 492 were upregulated and 759 were downregulated in pairwise gII vs. gI (**Figure [4](#F4){ref-type="fig"}**). In pairwise gI vs. mI, 1,037 genes were upregulated and 630 were downregulated, whereas 1,560 were upregulated and 1,072 were downregulated in pairwise gII vs. gI (**Figure [4](#F4){ref-type="fig"}**).

![**Relationships between 12 inflorescence samples based on multidimensional scaling.** In the plot, dimension 1 clearly separates monoecious and gynoecious samples. Dimension 2 unambiguously separates two different developmental phase samples in gynoecious plants, whereas two group monoecious samples cannot be clearly separated. mI containing m1, m2, and m3, indicates monoecious samples at stage I; mII containing m4, m5, and m6, indicates monoecious samples at stage II; gI containing g1, g2, and g3, indicates gynoecious samples at stage I; gII containing g4, g5, and g6, indicates gynoecious samples at stage II. BCV, biological coefficient of variation.](fpls-07-01953-g003){#F3}

![**Differentially expressed genes in the pairwise mII vs. mI, gII vs. gI, gI vs. mI and gII vs. mII in *J. curcas.* mII vs. mI indicates the pairwise stage II vs. I in monoecious inflorescences; gII vs. gI indicates the stage II vs. I in gynoecious inflorescences; gI vs. mI indicates the gynoecious vs. monoecious inflorescences at stage I and gII vs. mII indicates the gynoecious vs. monoecious inflorescences at stage II.** Blue lines indicate genes with a twofold change; red points indicate genes with a significantly different expression at a false discovery rate (FDR) of \<0.001. FC, fold change; CPM, counts per million mapped reads.](fpls-07-01953-g004){#F4}

The above results show that significantly more differentially expressed genes were found during inflorescence development from stage I to II in gynoecious plants than ones in monoecious plants, indicating that a very different physiological and biochemical process was occurring during this developmental transition between these two genotypes. Considering that the morphological differences between gynoecious and monoecious plants were mainly resulted from the arrest of stamen primordia that initiates at stage I (**Figure [2](#F2){ref-type="fig"}**) and the abortion of male flowers that occurs from stage I to stage II in gynoecious plants (**Figure [1](#F1){ref-type="fig"}**), we focused on the 171 coexpressed genes in pairwise gI vs. mI and gII vs. gI, excluding the 47 coexpressed genes in pairwise gII vs. gI, gI vs. mI and gII vs. mII (**Figure [5](#F5){ref-type="fig"}**) because these differentially expressed genes in pairwise gII vs. mII might be a consequence rather than a cause of stamen primordia arrest. GO analysis shows that the 171 coexpressed genes belong to 25 subset GO categories, and the "Membrane," "Regulation of biological process," "Regulation of cellular process," and "Ion binding" functional categories are predominant (**Figure [6](#F6){ref-type="fig"}**; **Supplementary Table [S3](#SM3){ref-type="supplementary-material"}**). The genes associated with reproductive (GO:0022414), reproductive developmental (GO:0003006), and developmental regulation (GO:0050793) processes are likely to play important roles during male flower abortion in gynoecious *Jatropha*.

![**Venn diagram shows the overlap of differentially expressed genes in the pairwise mII vs. mI, gII vs. gI, gI vs. mI, and gII vs. mII in *J. curcas.*** mII vs. mI, gII vs. gI, gI vs. mI, and gII vs. mII indicate the same pairs as shown in **Figure [4](#F4){ref-type="fig"}**.](fpls-07-01953-g005){#F5}

![**Gene Ontology (GO) categories of 171 coexpression genes in the pairwise gII vs. gI and gI vs. mI.** BP, biological process; CC, cellular component; MF, molecular function. gII vs. gI and gI vs. mI indicate the same pairs as shown in **Figure [4](#F4){ref-type="fig"}**.](fpls-07-01953-g006){#F6}

Differentially Expressed Genes Related to Floral Development
------------------------------------------------------------

Owing to that differentially expressed genes associated with the "Reproductive process" and "Reproductive developmental process" functional categories are enriched (*P* \< 0.05), we accordingly scanned our assembled data set and identified 32 genes that are homologous to floral development genes in model plants such as *Arabidopsis*, maize, and melon (**Supplementary Table [S4](#SM4){ref-type="supplementary-material"}**). In pairwise mII vs. mI, 11 genes were upregulated and 2 were downregulated, whereas 8 were upregulated and 11 were downregulated in pairwise gII vs. gI (Fold change ≥ 2.0 and FDR \< 0.05) (**Figure [7](#F7){ref-type="fig"}**; **Supplementary Table [S5](#SM5){ref-type="supplementary-material"}**). *AGAMOUS* (*AG*), *APETALA3* (*AP3*), *SHATTERPROOF 1* (*SHP1*), *SPOROCYTELESS* (*SPL*), *SEEDSTICK* (*STK*), *FLOWERING LOCUS D* (*FD*), and *UNUSUAL FLORAL ORGANS* (*UFO*) have similar expression patterns in pairwise mII vs. mI and gII vs. gI, suggesting that the changes in transcriptional levels may be caused by different development stages. In pairwise gI vs. mI, six genes were upregulated and three were downregulated, whereas six were upregulated and 16 were downregulated in pairwise gII vs. mII (**Figure [7](#F7){ref-type="fig"}**; **Supplementary Table [S5](#SM5){ref-type="supplementary-material"}**).

![**Hierarchical clustering of differentially expressed genes involved in floral development in the pairwise mII vs. mI, gII vs. gI, gI vs. mI, and gII vs. mII in *J. curcas.*** **(A,B)** Clustering of differentially expressed genes in the pairwise mII vs. mI and gII vs. gI, and in gI vs. mI and gII vs. mII, respectively. mII vs. mI, gII vs. gI, gI vs. mI, and gII vs. mII indicate the same pairs as shown in **Figure [4](#F4){ref-type="fig"}**.](fpls-07-01953-g007){#F7}

Differentially Expressed Genes Involved in Phytohormone Biosynthesis and Signaling Pathways
-------------------------------------------------------------------------------------------

Phytohormones play an important part in the process of sex expression ([@B5]), we monitored homologous genes involved in hormone metabolism and signaling pathways during the inflorescence development of *J. curcas*. In all, 70 differentially expressed genes were identified, in which eight were from abscisic acid (ABA), 16 from auxin, 9 from BR, 11 from CK, 4 from ETH, 7 from GA, 8 from JA, and 7 from salicylic acid (SA) pathways (**Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**). Among them, some genes encoding key rate-limiting enzymes in hormone biosynthesis pathways were significantly changed. In pairwise gII vs. gI and gII vs. mII, the expression of NINE-*CIS*-EPOXYCAROTENOID DIOXYGENASES 3 (NCED) and *NCED5* in ABA biosynthesis, TRYPTOPHAN AMINOTRANSFERASE OF ARABIDOPSIS 1 (TAA1) and YUCCA (YUC) in auxin, *1-AMINOCYCLOPROPANE-1-CARBOXYLIC ACID SYNTHASE* (*CmACS-7*) in ETH, *GA 20-oxidase 2* (*GA20ox2*) and *GA 3-oxidase 1* (*GA3ox1*) in GA, *DEFECTIVE ANTHER DEHISCENCE 1* (*DAD1*), *FATTY ACID DESATURASE 8* (*FAD8*), *LIPOXYGENASE 2* (*LOX2*), *ALLENE OXIDE SYNTHASE* (*AOS*), and *ALLENE OXIDE CYCLASE 3* (*AOC3*) in JA were downregulated, whereas *ATP/ADP isopentenyltransferases 5* (*IPT5*), *CYP735A1*, and *LONELY GUY 5* (*LOG5*) in CK were upregulated (**Figure [8](#F8){ref-type="fig"}**).

![**Expression patterns of differentially expressed genes that encode key enzymes in hormone biosynthesis pathways in the pairwise mII vs. mI, gII vs. gI, gI vs. mI, and gII vs.** mII in *J. curcas.* mII vs. mI, gII vs. gI, gI vs. mI, and gII vs. mII indicate the same pairs as shown in **Figure [4](#F4){ref-type="fig"}**.](fpls-07-01953-g008){#F8}

Validation of the Expression Profiles of Candidate Genes Involved in Sex Determination of *J. curcas* by Real Time qPCR
-----------------------------------------------------------------------------------------------------------------------

To validate the expression profiles of candidate genes obtained from the above comparative transcriptome analysis, we selected 10 candidate genes involved in sex determination of *J. curcas* to perform real time qPCR analysis. These 10 genes are *AOS*, *DAD1*, *KNOTTED1-LIKE HOMEOBOX GENE 6* (*KNAT6*), *MYC2*, *SHI-RELATED SEQUENCE 5* (*SRS5*), *SHORT VEGETATIVE PHASE* (*SVP*), *TAA1*, *TERMINAL FLOWER 1* (*TFL1*), *TASSELSEED2* (*TS2*) and *YUCCA 4* (*YUC4*). The results revealed by real time qPCR (**Figure [9](#F9){ref-type="fig"}**) show that the expression patterns of these genes were consistent with the ones obtained by transcriptome analysis (**Figure [7](#F7){ref-type="fig"}**; **Supplementary Tables [S4](#SM4){ref-type="supplementary-material"}** and **[S5](#SM5){ref-type="supplementary-material"}**), indicating transcriptome data in this study are reliable.

![**Validation of expression profiles of the candidate genes involved in sex determination of *J. curcas* by real time qPCR.** **(A)** *ALLENE OXIDE SYNTHASE* (*AOS*). **(B)** *DEFECTIVE ANTHER DEHISCENCE 1* (*DAD1*). **(C)** *KNOTTED1-LIKE HOMEOBOX GENE 6* (*KNAT6*). **(D)** *MYC2*. **(E)** *SHI-RELATED SEQUENCE 5* (*SRS5*). **(F)** *SHORT VEGETATIVE PHASE* (*SVP*). **(G)** *TRYPTOPHAN AMINOTRANSFERASE OF ARABIDOPSIS 1* (*TAA1*). **(H)** *TERMINAL FLOWER 1* (*TFL1*). **(I)** *TASSELSEED2* (*TS2*). **(J)** *YUCCA 4* (*YUC4*). **(K)** Expression patterns of the candidate genes in four pairs. *GAPDH* was used as the internal reference. The mean was calculated from three biological replicates. Error bars show the SD (*n* = 3). mI, mII, gI, and gII indicate the same samples as shown in **Figure [3](#F3){ref-type="fig"}**. mII vs. mI, gII vs. gI, gI vs. mI, and gII vs. mII indicate the same pairs as in **Figure [4](#F4){ref-type="fig"}**. FC, fold change.](fpls-07-01953-g009){#F9}

Effect of GA on Flower Organ Development in Gynoecious Plants
-------------------------------------------------------------

Since the expression of two GA biosynthesis genes *GA20ox2* and *GA3ox1* was found downregulated in gynoecious inflorescence (**Figure [8](#F8){ref-type="fig"}**; **Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**), we investigated the effect of GA on flower organ development in gynoecious plants. The result showed that exogenous application of GA~3~ promoted perianth development in male flowers and partly prevented pistil development in female flowers in gynoecious plants, which resulted in neutral flowers without pistils and stamens in gynoecious inflorescences (**Figure [10](#F10){ref-type="fig"}**).

![**Application of GA~3~ promotes perianth development in male flowers and prevents pistil development in female flowers in gynoecious inflorescences. (A)** Gynoecious inflorescence without treatment. Red arrows indicate female flowers and blue arrows indicate aborted male flowers. **(B)** Gynoecious inflorescence with GA~3~ treatment. Red arrows indicate neutral flowers originated from female floral meristem, and blue arrows indicate neutral flowers originated from male floral meristem. **(C)** Female flower as shown and indicated with red arrows in **(A)**. **(D)** Neutral flower originated from floral meristem generating female flowers after GA~3~ treatment as shown and indicated with red arrows in **(B)**. **(E)** Neutral flower originated from floral meristem generating male flowers after GA~3~ treatment as shown and indicated with blue arrows in **(B)**. Bar, 1.0 cm.](fpls-07-01953-g010){#F10}

Discussion
==========

Floral Development Genes May Be Key Regulators in the Transition from Monoecious to Gynoecious Plants
-----------------------------------------------------------------------------------------------------

At initiation of inflorescence buds, female and male flower primordia are emerged in both gynoecious and monoecious plants. However, in gynoecious plants, the development of male flowers is aborted, and only female flowers can normally grow; whereas the development of both female and male flowers is normal in monoecious plants. Transcriptome analysis showed that differentially expressed genes in pairwise stage II vs. I in monoecious inflorescence are far less than ones in gynoecious (**Figure [4](#F4){ref-type="fig"}**). This finding indicates that large physiological and biochemical changes occur in gynoecious inflorescences from stages I to II, in accord with the morphological change of male flower abortion. A large number of differentially expressed genes in pairwise gI vs. mI and gII vs. mII suggest that reproductive system differences between gynoecious and monoecious plants indeed result from inherent genetic variation.

In gynoecious plants, the abortion of male flowers was observed at stage II; however, the initiation of the process is expected to occur at stage I. Thus, the transition from stage I to stage II is a key period determining the sex expression of *J. curcas*. Genes differentially expressed at stage I contribute to the arrest of male flowers in gynoecious inflorescences. As shown in **Supplementary Table [S5](#SM5){ref-type="supplementary-material"}**, during the development of monoecious inflorescence from stage I to stage II, 11 of the 13 genes involved in the floral development were upregulated at transcriptional level. In contrast, the expression levels of 10 genes involved in the floral development were down-regulated during the development of gynoecious inflorescence from stage I to stage II, whereas the expression of these genes were not significantly changed in monoecious plants (**Supplementary Table [S5](#SM5){ref-type="supplementary-material"}**). These results suggest that different regulatory mechanisms of inflorescence development are present in monoecious and gynoecious *Jatropha*. The differentially expressed genes with similar expression patterns in pairwise gII vs. gI and gII vs. mII should be caused by male flower abortion, whereas ones with different expression patterns in pairwise gI vs. mI and gII vs. mII, are more likely to contribute to the arrest of stamen primordium in gynoecious plants (**Figure [7](#F7){ref-type="fig"}**). Among the differentially expressed genes shown in **Figure [7](#F7){ref-type="fig"}**, *KNAT6*, *MYC2*, *SRS5*, *SVP*, *TFL1*, and *TS2* were considered as strong candidate regulators of sex determination in *J. curcas*.

KNAT6 contributes to the maintenance of the shoot apical meristem and organ separation ([@B10]), and is an important regulator of *Arabidopsis* inflorescence architecture ([@B56]). MYC2 is a master regulator modulating various JA-dependent signaling in plant development and defense, and controlling the crosstalk between JA and the other hormone signaling pathways ([@B40]). SRS5 is a member of SHORT-INTERNODES/STYLISH (SHI/STY) family that may be involved in stamen and leaf development and flowering regulation ([@B61]). SVP and TFL1 are inhibitors of the transition from vegetative growth to flowering ([@B33]). SVP also represses the biosynthesis of GA at the shoot apex in *Arabidopsis* ([@B4]). TS2 is a short-chain alcohol dehydrogenase that represses pistil primordium in maize tassels and both *ts1* and *ts2* mutants result in conversion from staminate to pistillate tassels ([@B23]). Application of JA to developing tassels can rescue stamen development in *ts1* and *ts2* mutants ([@B1]). Our results suggest that some floral development genes and phytohormones (mainly GA and JA) likely coordinate to modulate sex expression in *J. curcas*. It is worth mentioning that the expression levels of homologs of several well-known sex determination genes, such as *CsACS* in cucumber ([@B44]), and *CmACS* and *CmWIP1* in melon ([@B13]), were not significantly different between monoecious and gynoecious plants, indicating that the mechanism of sex determination in cucumber and melon is different from that in *J. curcas*.

GA Plays Important Roles in the Development of *J. curcas* Flowers
------------------------------------------------------------------

GA20oxs and GA3oxs are key enzymes in GA biosynthesis, both being responsible for the conversion from GA~12~ to bioactive GAs ([@B70]). Expression of *GA20ox2* and *GA3ox1* are downregulated in gynoecious inflorescences (**Figure [8](#F8){ref-type="fig"}**; **Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**), indicating that the GA biosynthesis is lower in gynoecious than in monoecious inflorescences. Furthermore GA~3~ treatment promoted perianth development in male flowers and repressed pistil development in female flowers to produce neutral flowers in gynoecious inflorescences (**Figure [10](#F10){ref-type="fig"}**). Our previous studies showed that GA~3~ treatment induced the emergence of neutral and bisexual flowers in monoecious plants ([@B52]). These results suggest that GA plays essential roles in floral organ development in *J. curcas*. GA treatment did not resume the development of stamens in gynoecious male and female flowers, whereas GA promoted the development of stamens in monoecious female flowers to produce bisexual flowers, implying that the mechanisms controlling the arrest of stamen primordium between gynoecious and monoecious inflorescences are different. In maize, GA is a negative regulator in the development of stamen ([@B64]; [@B16]), but a positive regulator in *Arabidopsis* and rice ([@B17]). Our results showed that GA has multiple functions in floral organ development in *J. curcas*, and also suggested that regulatory mechanisms controlling sex expression are different in various plant species.

The Biosynthesis of Auxin Is Repressed in Gynoecious Inflorescences
-------------------------------------------------------------------

The TAA/YUC two-step pathway is the main auxin biosynthesis pathway ([@B77]) and is widely present in plants ([@B74]). TRYPTOPHAN AMINOTRANSFERASE OF ARABIDOPSIS 1 (TAA1) catalyzes the conversion from L-tryptophan (Trp) to indole-3-pyruvic acid (IPA) ([@B62]), and YUCCA (YUC) flavin-containing monooxygenases (FMOs) catalyze the oxidative decarboxylation of IPA to generate indole-3-acetic acid (IAA) ([@B47]). Both TAA1 and FMOs are rate-limiting enzymes in the regulation of auxin biosynthesis. The expression of *TAA1*, *YUC4*, and *YUC8* was downregulated in pairwise gI vs. mI and gII vs. mII, indicating that auxin biosynthesis is lower in gynoecious than in monoecious inflorescences. This finding is supported by the observation that expressions of auxin responsive genes *BTB AND TAZ DOMAIN PROTEIN 2* (*BT2*), *INDOLE-3-ACETIC ACID INDUCIBLE 14* (*IAA14*), *IAA29*, and *SMALL AUXIN UP RNA 20* (*SUR20*), auxin transport genes *PLEIOTROPIC DRUG RESISTANCE 9* (*PPDR9*), and *PIN6* were down-regulated (**Figure [8](#F8){ref-type="fig"}**; **Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**). In *Arabidopsis* and barley, high-temperature injury causes a specific decrease of endogenous auxin content in developing anthers, resulting in the arrest of pollen development, and exogenous application of IAA can completely reverse male sterility ([@B59]), indicating that auxin is an important regulator in the development of male flowers. This result is in agreement with the finding that auxin biosynthesis is suppressed in gynoecious inflorescences (**Figure [8](#F8){ref-type="fig"}**), which might result in male flower abortion in gynoecious *Jatropha*.

In addition, auxin plays an indispensable role in regulating inflorescence branching in *Arabidopsis* ([@B65]). A loss-of-function mutation in PIN-FORMED 1 (PIN1), an auxin carrier, abolishes the formation of an inflorescence meristem, resulting in naked inflorescence stems ([@B68]). In *J. curcas*, gynoecious inflorescences contain fewer branches than monoecious inflorescences (**Figures [1C,F](#F1){ref-type="fig"}**), which may be associated with the decrease of auxin biosynthesis in gynoecious inflorescence. This suggests that auxin may play a similar role in controlling inflorescence branching in *J. curcas*.

Two *Arabidopsis* cytochrome P450 monooxygenases, CYP79B2 and CYP79B3, catalyze the conversion of Trp to indole-3-acetaldoxime (IAOx), which is a common intermediate for the biosynthesis of IAA, 3-indolylmethyl-glucosinolate (IG) and camalexin ([@B36]; [@B49]; [@B79]; [@B29]). Overexpression of *CYP79B2* increases IAA level leading to auxin overproduction phenotypes of transgenic plants ([@B78], [@B79]). Correspondingly, *cyp79B2 cyp79B3* double mutants display IAA-deficient phenotypes in *Arabidopsis* ([@B79]). Thus, IAOx is an important intermediate in auxin biosynthesis distinct from the TAA/YUC pathway ([@B79]; [@B63]). However, IAOx has been detected only in *Arabidopsis* and related species, and the CYP79B family is also limited to Brassicaceae species, suggesting that the IAOx pathway is a Brassicaceae species-specific pathway in IAA biosynthesis ([@B47]). In *J. curcas*, the homologs *CYP79B2* and *CYP79B3* are upregulated in pairwise gI vs. mI, indicating that IOAx biosynthesis is higher in gynoecious than in monoecious inflorescences (**Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**); however, the expression of *NITRILASES 1* (*NIT1*), *NIT2*, and *NIT3*, which are involved in the conversion from IAN to IAA, was not significantly changed (**Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**). Moreover, the homolog *CYP83B1*, which converts IAOx to aci-nitro compounds in IG biosynthesis ([@B7]; [@B6]), is upregulated in pairwise gI vs. mI and gII vs. mII (**Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**). These results suggest that IAOx is more likely to be channeled into IG biosynthesis than into IAA, or even that IAOx is not involved in IAA biosynthesis, although an IAOx biosynthesis pathway is present in *J. curcas*.

JA Signaling Pathway May Participate in the Abortion of Male Flowers
--------------------------------------------------------------------

Previous study showed that JA plays an essential role in regulating the development of inflorescences, flowers, and male organs, and JA biosynthesis and signaling pathway during reproductive development are conserved in *Arabidopsis*, tomato, and rice ([@B73]). In *Arabidopsis*, *DEFECTIVE ANTHER DEHISCENCE 1* (*DAD1*) encodes a chloroplastic phospholipase A1 that catalyzes the initial step of JA biosynthesis, and the *dad1* mutant displays defects in anther dehiscence, pollen maturation, and flower opening ([@B38]). *FATTY ACID DESATURASE 8* (*FAD8*), *LIPOXYGENASE 2* (*LOX2*), *ALLENE OXIDE SYNTHASE* (*AOS*), and *ALLENE OXIDE CYCLASE 3* (*AOC3*) are involved in JA biosynthesis in plastids ([@B67]). In *J. curcas*, the homologs of these genes are all downregulated in pairwise gI vs. mI or/and gII vs. mII, implying that JA biosynthesis is lower in gynoecious than in monoecious inflorescences (**Figure [8](#F8){ref-type="fig"}**; **Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**). The expression of *FAD8*, *DAD1*, *AOS*, and *AOC3* is upregulated during inflorescence development from stage I to II in monoecious plants, whereas *DAD1* and *LOX2* are downregulated in gynoecious plants (**Figures [8](#F8){ref-type="fig"}** and **[9](#F9){ref-type="fig"}**; **Supplementary Table [S6](#SM6){ref-type="supplementary-material"}**). These results indicate that JA biosynthesis is also conserved in *Jatropha*, suggesting that the decrease of JA biosynthesis is likely to contribute to the abortion of male flowers in gynoecious plants. However, treatment with low concentration (6 μM) of MeJA did not resume the development of stamens of male and female flowers in gynoecious plants (data not shown). Because *Jatropha* is extremely sensitive to MeJA, and treatment with 10 μM of MeJA caused severe injury to inflorescences (data not shown), the concentration of MeJA used for treatment of *Jatropha* inflorescences will be optimized in future studies.

Regulatory Networks Participate in the Transition from Monoecy to Gynoecy in *J. curcas*
----------------------------------------------------------------------------------------

The phytohormones auxin, BR, CK, ETH, GA, and JA influence sex differentiation in flowering plants ([@B45]; [@B14]; [@B25]; [@B37]; [@B51]). Our results showed that the expression of some key genes involved in ABA, auxin, CK, ETH, GA, and JA biosynthesis pathways was significantly changed during inflorescence development and sex expression in *J. curcas*, suggesting that these hormones may participate in these processes. Thus, it is interesting to investigate how these hormones interact with one another to regulate the abortion of male flowers in gynoecious plants. However, only *MYC2*, a hub gene, was identified. MYC2 is a "master switch" in the crosstalk between JA and the other hormone signaling pathways, and has a wide range of functions in the regulation of developmental programs ([@B40]). Both GA and JA promote *MYC2* expression by a similar mechanism to modulate downstream gene expression ([@B40]). Moreover, many reverse genetics experiments support that GA and JA play crucial roles in the development of stamens and pollen ([@B17]; [@B73]). Based on the results, we speculate that the arrest of stamen primordium, which may result from one or more gene mutation in monoecious *Jatropha*, causes the abortion of male flowers to generate gynoecious individuals, in which some floral development genes and plant hormone signaling pathways participate in.

Author Contributions
====================

Z-FX, M-SC, and B-ZP designed experiments. M-SC and Z-FX wrote the manuscript. M-SC, LN, and YD analyzed and interpreted sequencing data. QF, Y-BT, JM-H, JN, and M-LZ carried out experiments and reviewed 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.

**Funding.** This work was supported by funding from the National Natural Science Foundation of China (31370595, 31300568, and 31500500), and the CAS 135 Program (XTBG-T02).

The authors gratefully acknowledge the Central Laboratory of the Xishuangbanna Tropical Botanical Garden for providing the high-performance computing and other research facilities.

<http://www.bioinformatics.babraham.ac.uk/projects/fastqc>

<http://bioinfogp.cnb.csic.es/tools/venny/index.html>

<https://github.com/cran/pheatmap>

Supplementary Material
======================

The Supplementary Material for this article can be found online at: <http://journal.frontiersin.org/article/10.3389/fpls.2016.01953/full#supplementary-material>

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

###### 

Click here for additional data file.

[^1]: Edited by: *Michael Deyholos, University of British Columbia, Canada*

[^2]: Reviewed by: *Fan Chen, Institute of Genetics and Developmental Biology (CAS), China; Ying Wang, South China Institute of Botany (CAS), China*

[^3]: This article was submitted to Plant Genetics and Genomics, a section of the journal Frontiers in Plant Science
