- Research article
- Open Access
First RNA-seq approach to study fruit set and parthenocarpy in zucchini (Cucurbita pepo L.)
BMC Plant Biologyvolume 19, Article number: 61 (2019)
Zucchini fruit set can be limited due to unfavourable environmental conditions in off-seasons crops that caused ineffective pollination/fertilization. Parthenocarpy, the natural or artificial fruit development without fertilization, has been recognized as an important trait to avoid this problem, and is related to auxin signalling. Nevertheless, differences found in transcriptome analysis during early fruit development of zucchini suggest that other complementary pathways could regulate fruit formation in parthenocarpic cultivars of this species. The development of next-generation sequencing technologies (NGS) as RNA-sequencing (RNA-seq) opens a new horizon for mapping and quantifying transcriptome to understand the molecular basis of pathways that could regulate parthenocarpy in this species. The aim of the current study was to analyze fruit transcriptome of two cultivars of zucchini, a non-parthenocarpic cultivar and a parthenocarpic cultivar, in an attempt to identify key genes involved in parthenocarpy.
RNA-seq analysis of six libraries (unpollinated, pollinated and auxin treated fruit in a non-parthenocarpic and parthenocarpic cultivar) was performed mapping to a new version of C. pepo transcriptome, with a mean of 92% success rate of mapping. In the non-parthenocarpic cultivar, 6479 and 2186 genes were differentially expressed (DEGs) in pollinated fruit and auxin treated fruit, respectively. In the parthenocarpic cultivar, 10,497 in pollinated fruit and 5718 in auxin treated fruit. A comparison between transcriptome of the unpollinated fruit for each cultivar has been performed determining that 6120 genes were differentially expressed. Annotation analysis of these DEGs revealed that cell cycle, regulation of transcription, carbohydrate metabolism and coordination between auxin, ethylene and gibberellin were enriched biological processes during pollinated and parthenocarpic fruit set.
This analysis revealed the important role of hormones during fruit set, establishing the activating role of auxins and gibberellins against the inhibitory role of ethylene and different candidate genes that could be useful as markers for parthenocarpic selection in the current breeding programs of zucchini.
Fruit set is defined as the transition of an ovary to a growing young fruit, and depends on the successful completion of pollination and fertilization. Pollination is the transfer of a pollen grain from the anther to the stigma. Fertilization occurs in the ovule, and requires pollen tubes growth in the stylar tissue and fusion with the egg cell . Both processes are affected negatively by unfavourable environmental conditions such as low/ high temperature or inadequate humidity that prevents fruit set in the majority of flowering plants .
In the case of zucchini, one of the most important morphotypes of Cucurbita pepo, these harsh conditions especially occur in off-season crops, causing economic losses due to low fruit yield. Parthenocarpy, fruit development in the absence of pollination/fertilization, has been recognized as an important trait to avoid this problem in different fruit crops , considering than each species shows specific responses following pollination/fertilization, In Cucurbita, the most practical mean of increasing fruit set of zucchini when pollination/fertilization is inadequate would be the use of cultivars with vegetative parthenocarpy, innate ability to set parthenocarpic fruit , but this kind of parthenocarpy is limited to a few cultivars in C. pepo. For several cultivars of zucchini, we provided previously an exhaustive description of the specific changes in fruit set observed after successful pollination or parthenocarpy, when compared with non-pollination .
Parthenocarpy has been related to certain plant hormones as auxins, gibberellins, cytokinin and brassinosteroids . Exogenous application of these plant hormones induces parthenocarpy fruit set in cucumber and zucchini , and high levels of endogenous IAA were found in parthenocarpic fruit respect to pollinated fruit in cucumber . Furthermore, natural parthenocarpic tomato mutants pat and pat-2 accumulate high levels of gibberellin in the unpollinated ovaries .
Transcriptome analyses have led to the successful identification of some genes associated with parthenocarpy. Over-expression of SlGA20ox1 was found at high levels throughout fruit growth in the pat mutant in tomato . Over expression of SlTIR1 (an auxin receptor) give rise to parthenocarpic fruit set in cucumber treated with exogenous auxins . Down-regulated expression of SlARF7 (Auxin Response Factor 7) and SlIAA9 (Aux/IAA gene) also induced parthenocarpy in transgenic tomatoes [10, 11]. On the other hand, mutations in ARF8 and IAA9 induced parthenocarpy in Arabidopsis and tomato respectively [12, 13]. These genes, ARF8, IAA9 and TIR1, showed downregulation during fruit set in pollinated fruit of non-parthenocarpic cultivars of C. pepo. Nevertheless, differences found in transcriptome analysis suggest that other complementary pathways could regulate fruit formation in parthenocarpic cultivars of this species . Consequently, more knowledge about transcriptome responses of fruit set is required to exploit parthenocarpy in zucchini.
The reference transcriptome of Cucurbita pepo L. ssp. pepo is a valuable resource for identification of transcripts involved in specific biological processes, improving genome annotation, elucidating phylogenetic relationships, providing SSR and SNP markers and large-scale expression analysis . Large-scale expression analysis, as microarray gene chip experiments, could be a powerful tool for monitoring the expression level of thousands of genes involved in fruit set and parthenocarpy in zucchini, but detecting transcripts by microarray technology is limited to genomic sequencing. Moreover, it is a laborious and capital intensive task, and the statistical reproducibility of the data is relatively poor. Nowadays, the development of next-generation sequencing technologies as RNA-sequencing (RNA-seq) opens a new horizon for mapping and quantifying transcriptome, creating opportunities to better understand the molecular basis of pathways that could regulate these processes. These methods improved sequencing capabilities with respect to amount of data, time, and cost, and do not depend on prior gene identification or assembly onto microarrays [15, 16]. RNA-seq technology has been used to improve the genome annotations, to detect areas of alternative splicing to discover new genes and novel transcribed regions and to perform differential expression analysis. This method has many advantages such as lower background signal, higher suitability for both known transcripts and new genes, and ability to quantify a large dynamic range of expression levels, with absolute rather than relative values, even with species that lacked a reference genome . RNA-seq has been used for analyzing the transcriptome responses to fruit set and parthenocarpy in cucumber. . However this type of technology has not been yet applied for studying these processes in zucchini.
In the current study, differences in fruit transcriptome of two cultivars of zucchini, a non-parthenocarpic cultivar and a parthenocarpic cultivar has been analysed through RNA-seq in an attempt to identify genes that play important roles during pollination, fruit set and parthenocarpy.
Analysis of RNA-seq libraries
Two cultivars of Cucurbita pepo spp. pepo morphotype zucchini were used in this study, a non-parthenocarpic cultivar MUCU-16(Acronym MU16), and a parthenocarpic cultivar Whitaker (Acronym WHT). Samples of pollinated fruit (PF), non-pollinated fruit (UF) and auxin treated fruit (AF) of both cultivars were obtained to construct 6 libraries to be sequenced for RNA Seq analysis. Raw data was generated through Illumina HiSeq 2500 sequencing and was subjected to initial treatment, low quality regions and adapter sequences have been trimmed. Samples ranged from 42.4 to 84.8 million reads with mean read length of 101 pb, which is enough for the quantitative analysis of gene expression (Table 1).
The reference transcriptome of C. pepo v3 consists of 108,062 transcripts that represent 73,239 unigenes clusters with an average length of 1052 pb. RNA-seq reads obtained from six libraries were mapped to the reference transcriptome of C. pepo v3, with a mean of 87% success rate of mapping, but high percentage of multiple reads (36%) were obtained (Table 1).
Therefore, CD-hit program  was used to analyse the redundancy of the reference transcriptome and sequences with 85% of homology were clustered. The longest unigene was chosen to represent each cluster, resulting in a total of 73,239 unigenes with an average length of 818 pb. Sequencing reads were newly mapped to the processed transcriptome, obtaining higher rate of unique mapped reads and lower number of multiple mapped reads, with a mean of 92% success rate of mapping (Table 1). This improved RNA-seq using the processed transcriptome was used for subsequent analysis.
Differentially expressed genes during fruit set
RNA-seq mapping revealed the number of genes expressed in each sample (Fig. 1). Comparison between pollinated, auxin and unpollinated fruit treatments in both cultivars showed that 52,240 and 54,761 genes were expressed in common in MUCU-16 and Whitaker respectively independently of treatment applied. Genes with differential expression during fruit set were filtered using the transcriptome of unpollinated fruit for both cultivars. A FDR of 0.05 and an absolute value of fold change ≥2 were used as the threshold to judge the significance of the differential gene expression(Additional file 1) In the non-parthenocarpic cultivar, 6479 and 2186 genes were differentially expressed in pollinated fruit (PF MU16 vs UF MU16) and auxin treated fruit (AF MU16 vs UF MU16) respectively. In the parthenocarpic cultivar, 10,497 genes were differentially expressed in pollinated fruit (PF WHT vs UF WHT) and 5718 genes in auxin treated fruit (AF WHT vs UF WHT). A comparison between transcriptome of the unpollinated fruit for each cultivar has been performed (UF WHT vs UF MU16), determining that 6120 genes were differentially expressed. It was also found that pollination treatment generated the highest number of expressed genes and the highest number of differentially expressed genes (Fig. 1).
Previous C. pepo gene annotation was incomplete because about 50% of DEGs had previous annotation in C. pepo transcriptome v3. Structural and functional annotation of DEGs was performed by BLAST analysis (e-value of 1e-25) with public databases. DEGs were compared to Arabidopsis genome and Arabidopsis proteins from TAIR, Nucleotide collection nr from NCBI, Swiss-Prot database, Cucurbita maxima genome and Cucurbita pepo genome v4 (Additional file 2). About 30% of DEGs had one significant hit with Arabidopsis genome and more than 40% of DEGs blasted with Arabidopsis proteins, improving the previous annotation with this species (Table 2). Blast analysis revealed that 55% of unigenes had significant matches in Nucleotide collection nr database, where the most of hits found were with Cucumis melo and Cucumis sativus gene. Sequence homology was also performed with the accurate database Swiss-Prot, obtaining more than 30% had significant matches (Table 2). This new blast analysis allows the annotation of 20% of unigenes that did not have previously annotated in the reference C. pepo transcriptome. Additionally, blast analysis was carried out with the genomes of Cucurbita maxima and Cucurbita pepo (Additional file 2).
Gene Ontology (GO) terms were further assigned to screened DEGs (Additional file 3) based on their sequence similarities to domains they contain in Interpro database. Over 25% of DEGs were assigned at least one GO term in biological process, molecular function and cellular component categories (Table 3). These unigenes, were further classified into different functional categories using a set of plant-specific GO slims (Additional file 3). Functional classifications of DEGs into plant specific GO slims (level 2) within the biological process, molecular function and cellular component categories were carried out. Metabolic process and cellular process were the most highly represented groups in biological process (Table 3). Genes involved in other important biological processes such as response to stimulus were also identified in pollination in the non-parthenocarpic cultivar and single-organism process in parthenocarpic cultivar. Under molecular function category, assignments were to the binding and catalytic activities (Table 3). It is worth noting that GO annotations revealed high number of transferases, kinases and hydrolases in this category, suggesting that genes involved in the secondary metabolite synthesis pathways were induced during fruit set (Additional file 3). In the case of cellular component, GO terms related to cell, membrane and organelle were well-represented (Table 3) DEGs also were annotated using KEGG pathway database (Additional file 3).
DEGs between pollination and auxin treatment in both cultivars were compared. Cluster analysis was used to screen the common expressed genes between them. A total of 1570 and 4522 were induced in both treatments in MUCU-16 and Whitaker, respectively (Fig. 2a). Of these subsets of DEGs, 1007 genes were common in both cultivars, which represented 64% of common DEGs in MUCU-16 cultivar and 22% in Whitaker cultivar. To further understand the function of these DEGs, GO enrichment analysis in the category of biological process was performed (p < 0.05). DEGs involved in cell cycle, DNA replications, microtubule movement, regulation of transcription or biosynthetic process were highly enriched (Fig. 2b).
It was also compared pollination treatment between both cultivars. 3939 DEGs were induced during pollination fruit set in both cultivars (Fig. 3a). GO enrichment analysis revealed that DNA replication, microtubule movement, regulation of biosynthetic process, regulation of biological process, cell proliferation, auxin-activated signalling pathway and cellular response to auxin stimulus were enriched processes in the category of biological process (Fig. 3b). On the other hand, 1213 DEGs were induced during auxin treated fruit set in both cultivars (Fig. 3a). GO enrichment analysis showed differences in biological processes activated during auxin treatment in comparison to pollination (Fig. 3). Microtubule movement, regulation of biosynthetic process and regulation of biological process were also induced, but processes related to auxins were not enriched in this treatment (Fig. 3c).
Expression of genes associated with cell division during fruit set
Cell division related genes were analysed based on previous DEGs results, cyclins and expansins (EXPs) were filtered. During pollination, 18 and 32 cyclins were induced in MUCU-16 and Whitaker, respectively. Of those cyclins, 16 DEGs were expressed in common in both cultivars, showing the same regulation. It was found 6 cyclins of type A, 4 cyclins of type B, 2 cyclins of type D, 2 cyclins of type P, CSK1 and CDKB2 genes. All cyclins were up-regulated during pollination, except CYCP2–1.Pollination also induced differential expression of expansins, and the majority of these genes were up-regulated. 4 expansins were expressed in common during pollinated fruit set in both cultivars, EXPA1, EXPA8, EXPA4 were upregulated and EXPLB1 was strongly downregulated (Additional file 4). During auxin treatment, 2 cyclins were found in the non-parthenocarpic cultivar. These cyclins were also induced in pollination treatment in this cultivar. In the case of the parthenocarpic cultivar, 24 cyclins were differentially expressed, activating also cyclins of type A, B, D and P. All of these cyclins were also induced during pollination in this cultivar, but its fold change was lower in auxin treatment. In the case of expansins, EXPLB1 was also strongly downregulated in both cultivars and higher number of expansins were induced in Whitaker in comparison to MUCU-16 (Additional file 4). On the other hand, Cyclin D6–1 and EXPLA1showed the same regulation and pollination/auxin treatment fruit set and during parthenocarpic fruit set (UF WHT vs UF MU16) (Additional file 4).
Identification of transcription factors related to fruit set
RNA-seq data showed that 207 and 254 transcription factors (TFs) displayed differential expression during pollination in MUCU-16 and Whitaker, respectively. In the case of auxin treatment, 124 and 220 were induced, a lower number in comparison to pollination treatment (Additional file 5). TFs were classified into several families including NAC, ERFs (Ethylene-responsive transcription factor), WRKY, bHLH (Basic helix-loop-helix protein), MYB, HD-Zip (Homeobox-leucine zipper protein) and MADS-box, among which NAC, AP2/ERF and WRKY recruited the major members. The vast majority of them were downregulated after pollination and auxin treatment (Fig. 4, Additional file 5). In the case of UF WHT vs UF MU16, 158 TFs were differentially expressed. Expression of these genes was different in comparison to pollination and auxin treatment (Additional file 5).
MADS-box and HD-Zip families were reviewed in detail. MADS-box family proteins, including Agamous-like MADS box proteins (AGL), represented lower rate in comparison to NAC, AP2/ERF and WRKY families during fruit set (Fig. 4). DEGs analysis revealed differential expression of 6 MADS-box genes during pollinated fruit set, CMB1 (also called AGL3), SOC1 (also called AGL20), AGL8, AGL11, AGL16 and AGL18 in the parthenocarpic variety. CMB1, SOC1, SVP, FBP24, AGL12 and AGL16 were expressed in the non-parthenocarpic variety. CMB1, SOC1 and AGL16 were commonly downregulated during pollinated fruit set in both cultivars. During fruit set dependent of auxin treatment, AGL8, AGL16 and FBP24 were strongly downregulated in MUCU-16 cultivar. Downregulation was also observed in AGL19 and SOC1 in Whitaker cultivar (Additional file 5).
On the other hand, HD-zip family represented greater rate than MADS-box family (Fig. 4). The overwhelming majority of HD-zip genes were down-regulated in pollination and auxin treatment. HAT5, HOX11, HOX20, ATHB7, ATHB12, ATBH21 and ATHB40 were expressed in common in both cultivars during pollination and auxin treatment. Stronger downregulation was observed in ATHB7 and ATBH21, more than 20 fold of expression (Additional file 5).
Comparison of the carbohydrate-related responses during fruit set
Carbohydrate metabolism was strongly induced during fruit set in zucchini (Additional file 6). During pollinated fruit set, 45 genes were expressed in common in both cultivars, 17 genes were up-regulated and 28 genes were down-regulated. Over-expression was found in key enzymes involved in glycolisis as enolase1, glyceraldehyde-3-phosphate dehydrogenase and fructose-bisphosphate aldolase, and genes related to starch accumulation as phosphoglucan phosphatase amyloplastic, and starch biosynthesis as 1,4-alpha-glucan-branching enzyme. On the other hand, negative regulation was found in gluconokinase, phosphoenolpyruvate carboxykinase, sucrose synthase, amylases and alpha-trehalose-phosphate synthases (Additional file 6). High levels of glucose and fructose found during pollinated fruit set in both cultivars and low levels of sucrose are produced by activation of glycolisis via glyceraldehyde-3-phosphate dehydrogenase and fructose-bisphosphate aldolase and downregulation of sucrose synthase (Fig. 5). In addition, it was also found low levels of starch in spite of activation of enzymes related to starch biosinthesis (Fig. 5, Additional file 6). During auxin treatment, 26 DEGs were expressed in common in both cultivars. Of this subset of genes, 1,4-alpha-glucan-branching enzyme and polygalacturase were up-regulated in common in both cultivars in contrast to amylases and threalose phosphate synthases, which were down-regulated. A detailed analysis of auxin treatment in each cultivar was carried out. In MUCU-16, fructose 6-phosphate 1-phosphotransferase was up-regulated. On the other hand, pyruvate kinase and phosphoenolpyruvate carboxykinase 1 were up-regulated in Whitaker (Additional file 6). Over-expression of these enzymes corroborated high levels of glucose and fructose found during auxin treated fruit set in both cultivars (Fig. 5). Nevertheless, during parthenocarpic fruit set, glucose and fructose levels were similar to unpollinated fruit of MUCU-16 (Fig. 5). This is contrary to over-expression of fructose-1,6-bisphosphatase, pyruvate kinase, glyceraldehyde-3-phosphate dehydrogenase, sucrose synthase and 6-phosphofructokinase 5, that indicated the carbohydrate biosynthesis (Additional file 6).
Hormone metabolism and signalling during fruit set
Hormone related genes were filtered based on previous DEGs results, 225 and 255 DEGs related to hormone were screened in pollinated fruit in MUCU-16 and Whitaker, respectively. However, lower number of DEGs, 99 and 210, was expressed during auxin treatment. These screened genes were involved in biosynthesis and signalling of auxins, ethylene, gibberellins, abscisic acid, cytokinin, brassinosteroids, jasmonic and salicylic acid (Additional file 7). Percentages of DEGs related to these eight hormones, number and category, were similar between treatments except in the case of auxin treatment of MUCU-16 (Fig. 6). During fruit set and parthenocarpic fruit set, the most representative hormone-related DEGs were auxins, followed by ethylene, gibberellins and abscisic acid. Regulation of gene expression was also analysed showing that DEGs related to auxins displayed higher number of up-regulated genes in comparison with ethylene, gibberellins and abscisic acid (Fig. 6).
Cluster analysis was carried out to elucidate the regulation of common hormone-related genes during fruit set. Of the subsets of DEGs related to hormone analysed previously, 64 genes were common during pollination and auxin treatment in both cultivars. Regulation of gene expression was similar during pollination in both cultivars and during auxin treatment in both cultivars (Fig. 7). Most of these common genes were related to auxin (24 genes, 37.5%) and ethylene (25 genes, 39%). Auxin genes were involved in auxin signalling (Aux/IAA and WTR families), auxin polar transport (PIN Auxin efflux carriers), response to auxin (SAUR Auxin responsive proteins) and auxin homeostasis (indole-3-acetic acid-amido synthetase, GH3.1, and indole-3-acetate o-methyltransferase 1, IAMT1). Expression of IAA genes, IAA4, IAA14 and IAA16 clustered together with PIN6, WAT1 and SAUR50, showing over-expression. However, GH3.1, an auxin biosynthesis gene, and others members of SAUR family (SAUR32) and PIN family (PIN5) decreased their expression under pollination and auxin treatment. Remarkably, IAA19, AUX22 and WAT1 showed the strongest up-regulation during fruit set. On the other hand, DEG analysis showed a decrease in mRNA levels of ethylene related genes. Ethylene biosynthesis genes (ACC oxidases, ACO1 and ACO3) and ethylene signalling genes, including ERFs (12 ERFs), EIN3 (Ethylene insensitive 3) or ETR2 (Ethylene receptor 2) among others, displayed down-regulation during fruit set. Apart from these findings related to auxin and ethylene, a subset of genes associated with GAs biosynthesis and signalling has been also modulated during fruit set. Expression of GA20ox1 (Gibberellin 20 oxidase 1) was found to be strongly up-regulated in contrast to expression of GA2ox2 (Gibberellin 2-beta-dioxygenase 2) and GID1B (Gibberellin receptor). To a lesser extent, genes related to cytokinin (Cytokinin riboside 5′-monophosphate phosphoribohydrolases, LOG1 and LOG5, and Cytokinin synthase IPT5), abscisic acid (9-cis-epoxycarotenoid dioxygenases, NCED1 and NCDE3, and abscisic acid receptor PYL8), brassinosteroids (Cytochrome P450 90A1, CYP90A1) and salicylic acid (Salicylic acid-binding protein 2, SABP2) were differentially expressed in common. All genes were downregulated, except CYP90A1. CYP90A1, a BRs biosynthesis gene, was up-regulated and clustered together auxin genes (Fig. 7).
Auxin at the transcriptome level during zucchini fruit set
Genes related to auxin have been filtered during pollination and during auxin treatment in both cultivars, indicating that DEGs related to auxin biosynthesis, auxin metabolism and signalling were strongly induced during fruit set (Fig. 6, Fig. 7 and Additional file 7).
During pollinated fruit set, auxin biosynthesis genes were altered. YUCCA6 (indole-3-pyruvate monooxygenase) and IAMT1 were up-regulated and three genes of GH3 family, GH3.1, GH3.5 and GH3.6, were downregulated. Related to auxin signalling, ten Aux/IAAs were expressed differently, being only IAA26 and IAA18 down-regulated in both cultivars. Contrasting with Aux/IAA family, ARFs showed different regulation between cultivars, ARF5 and ARF16 were down-regulated in MUCU-16 cultivar, and ARF5, ARF7 and ARF18 were up-regulated in Whitaker cultivar. Over-expression was also found in auxin efflux carriers, PIN1A, PIN1B, PIN6, PIN8, and auxin influx carriers, AUX1, LAX2 and LAX4 (Additional file 7).
During auxin treatment, auxin biosynthesis genes were also differentially expressed, YUCCA8 was downregulated in MUCU-16, and YUCCA7 and IAMT1 was up-regulated in Whitaker. IAA and ARF genes were also modulated during auxin treatment and, most of them showed the same regulation than during pollinated fruit set. IAA4, IAA14, IAA16 and IAA19 were also up-regulated during auxin treated fruit set in both cultivars, but ARF5 and ARF18 were only up-regulated in Whitaker cultivar. Auxin treatment also produced changes in PIN proteins and AUX/LAX proteins. PIN1, PIN1D and PIN6 were over-expressed together two AUX/LAX transporters, indicating that efflux and influx of auxins also occurred in auxin treatment. In the case of parthenocarpic fruit set (UF WHT vs UF MU16 comparison), it was found DEGs related to auxin with similar expression to pollination or auxin treatment. ARF18, PIN1B, PIN8 and LAX1 were also upregulated as pollination/auxin treatment (Additional file 7).
Biosynthesis and signalling of ethylene and gibberellin during fruit set
Ethylene metabolism is one of the most represented hormones during pollinated fruit set. RNA-Seq data revealed a decrease in the mRNA levels of most of ethylene related gene (Fig. 6, Fig. 7, Additional file 7). ACO1 and ACO3, which catalyze the subsequent steps of ethylene biosynthesis via S-adenosyl-L-methionine, were strongly downregulated during pollinated fruit set in contrast to auxin biosynthesis genes. As occurred for ACO genes, ethylene signalling genes including ERFs, ETR2 and EIN3s appeared down-regulated during pollinated fruit set. Higher number of DEGs related to ethylene was also altered during auxin treatment, and their regulation was similar to found during pollinated fruit set (Fig. 7, Additional file 7).
During parthenocarpic fruit set, ACO1, ETR2 and five ERFs (ERF2, ERF11, ERF12, ERF17 and ERF60) were downregulated as pollination/ auxin treatment in both cultivars. Moreover, WIN1 were upregulated as pollination treatment in Whitaker cultivar (Additional file 7).
A subset of genes associated with GAs biosynthesis and signalling has also been modulated during pollinated fruit set (Table 4). Three GID1-like genes, that encode GA receptors, showed down-regulation in response to pollination as well as GID2, an F-box protein that is essential for GA-mediated DELLA protein degradation. GA2ox1 and GA2ox8, which play important role in gibberellins catabolism, displayed a decrease of gene expression. However, GA3OX4, involved in production of GA bioactive for reproductive development, and GASA4, gibberellin-regulated protein, was up-regulated during pollination. Remarkably, genes related to gibberellin biosynthesis, KAO1 (Ent-kaurenoic acid oxidase 1) and KAO2 (Ent-kaurenoic acid oxidase 2), which catalyze three successive oxidations of ent-kaurenoic acid giving gibberellin 12 (GA12), and GA20ox1 and GA20ox2, key enzymes of gibberellin biosynthesis that catalyzes the conversion of GA12 and GA53 to GA9 and GA20 respectively, were up-regulated during pollinated fruit set in zucchini. Five of these genes were also common between cultivars during auxin treatment (Table 4). In the case of parthenocarpic fruit set, 15 genes related to gibberellin were differentially expressed (Additional file 7). Of this subset of genes, four genes showed the same regulation as pollinated/ auxin treated fruit set. GA20ox1 and GASA4 were up-regulated in contrast to GA2ox8 and GID1B (Table 4).
A modification in C. pepo transcriptome allows more specific mapping during gene expression analysis in RNA-seq
RNA-seq is a high-throughput sequencing assay that combines transcript identification and quantification of gene expression, core activities to study many processes in the plant life cycle as fruit set. Two strategies are possible to read mapping and transcript identification when a reference sequence is available, genome mapping or transcriptome mapping . In the case of this assay, mapping was carried out against a reference transcriptome since the new version of C. pepo genome (https://bioinf.comav.upv.es/downloads/zucchini/) lacked of gene annotation .
In order to reduce the high percentage of multiple mapped reads and unmapped reads that could limit gene identification (Table 1), transcriptome sequenced from root, leaf and flower tissue was modified. This processed transcriptome improved mapping statistics of RNA-seq libraries. Total mapped reads were higher (Table 1) and similar to RNA-seq analysis using the last version of genome .
Previous gene annotation of the reference C. pepo transcriptome was also improved by a new blast analysis of DEGs with a more restrictive e-value (cut off e-value of 1e-25) carried out in this work, increasing the numbers of genes with functional annotation in the processed C. pepo transcriptome (Table 2, Additional file 2). In spite of this deeply blast analysis, there are still a number of non-annotated transcripts in C. pepo transcriptome (Additional file 2). These transcripts without annotation may be consequence to the short sequences generated by the 454 sequencing technology used to assembly the first transcriptome of C. pepo . Probably, these sequences lack the conserved functional domains or may be no-codings RNAs (Additional file 2).
Cluster analysis of DEGs elucidated differences between pollinated fruit set and auxin treated fruit set
Comparative analysis of gene expressed showed a high number of common expressed genes between pollinated, auxin treated and unpollinated fruit in zucchini, indicating that there are common pathways between treatments (Fig. 1). Thus, it was necessary to carry out DEG analysis to find key role genes induced by pollination and/or auxin treatment during fruit set . DEGs analysis revealed that pollination induced the highest number of DEGs and the highest values of fold change (Fig. 1, Additional file 1), although 1007 DEGs were expressed in common during pollination and auxin treatment (Fig. 2), which revealed that rapid activation of fruit metabolism is mediated by cell cycle, regulation of transcription, regulation of biosynthetic process or DNA replication (Fig. 2). However, DEGs cluster analysis and GO enrichment analysis revealed differences between pollinated fruit set and auxin treated fruit set. Pollination activated processes related to auxins as auxin-activated signalling pathway and cellular response to auxin stimulus, not found during auxin treated fruit set. On the other hand, auxin treatment activated many more processes related to cell division, not enriched during pollinated fruit set (Fig. 3). This corroborates the hypothesis that the addition of synthetic hormones does not initiate the same processes as natural pollination [11, 21].
Cell division is an enriched biological process during fruit set
During early fruit development in zucchini, fruit length and placental diameter grow rapidly suggesting that cell division is increased . Moreover, it was found that functional categories of DEGs by GO term enrichment analysis showed that cell cycle, cell cycle progress, mitotic cell cycle and mitotic cell cycle progress related genes were strongly induced (Fig. 2). Differential expression observed in cyclins, CDKs (cyclin dependent kinase enzymes) and cell wall genes during zucchini fruit set positively correlated with rapid cell division observed in zucchini, and also observed during early fruit development of tomato and cucumber [22, 23].
Cyclins are the most important cell cycle regulators and control the progression of cells through the cell cycle by activating CDKs . D-type cyclins (CYCDs) that regulate the G1–S transition and G2-M transition . A-type cyclins (CYCAs), that are mainly present from S phase to M phase, and B-type cyclins (CYCBs), that take part the G2–M transition and during M phase , were differentially expressed, evidencing a cell cycle progression during fruit set which induce the exponential growth phase observed previously . Of these subfamilies of cyclins, cyclin D proteins are often noted as sensors of environment, linking hormonal signal with external conditions, and reporting to cell on the favourable conditions for cell-cycle start . D-type cyclins as CYCD3;1, CYCD3;2 and CYCD5;1 were over-expressed during zucchini fruit set, (Additional file 4), and similar regulation had already been previously described in cucumber fruit  which would indicate their fundamental role in fruit set process. It is particularly important to point out that CYCD3;1 was reported to be a rate-limiting factor for G1 /S transition  and its expression is regulated by hormones signalling and the availability of nutrients as sugars Previous results reported that CYCD3;1 activated cell division at the G1-S cell cycle phase transition in response to cytokinin and auxin in Arabidopsis . Over-expression of this gene found in zucchini fruit could be related to over-expression of auxin and cytokinin biosynthesis and signalling genes also found (Fig. 7). Moreover, CYCD6;1, a cyclin involved G1-S cell cycle phase transition, were over-expressed during pollinated/auxin treatment and parthenocarpic fruit set (Additional file 4), indicating that its over expression is closely linked to parthenocarpy in zucchini. Together cyclins, CDKs regulate cell cycle through the association of these proteins to two key classes of regulatory molecules, CKS1 and CKS2 (Cyclin-dependent kinases regulatory subunits) . In the case of zucchini, CKS1 and CDKB2;2 were over-expressed with similar fold change during fruit set (Additional file 4),. These proteins have shown interaction between them in Arabidopsis, and this interaction may also occur in zucchini, regulating cell cycle .
Cell cycle progression gives rise to the formation of two daughter cells. Prior to this process, cell undergoes a rapid enlargement (cell expansion) that involves the wall loosening. Expansins have been recognized as the major wall-loosening agents and cause extension of plant cell walls by disrupting non-covalent binding between cellulose and hemicelluloses . DEG analysis revealed expasins of the four subfamilies identified in plants, A-expansins (EXPAs), B-expansins (EXPBs), expansin-like A (EXLA) and expansin-like B (EXLB) during zucchini fruit set (Additional file 4). Previous studies showed that expansins were also highly enriched in growing cucumber fruits with similar differential expression to found in zucchini fruits, EXPA1 and EXPA8 were up-regulated in contrast to the strongly downregulation showed by EXPLB1 . Interestingly, expansins are thought to linked cell enlargement and cell wall changes with auxin signalling . Over-expression of EXPA1 and EXPA8 found in zucchini fruit might be induced by auxin signalling in Arabidopsis . The relation between expression of expansins and auxin signalling in zucchini could be underline by downregulation of ARF7 during pollinated fruit set (Additional file 7). Arabidopsis mutants in ARF7 (loss of function mutant) showed less expression level of EXPA1 and EXPA8 .
Down-regulation of MADS-box and HD-ZIP transcription factors play key roles in zucchini fruit set
Gene expression at the level of transcription is mediated by transcription factors and is crucial for almost all biological processes, including fruit development. Different TFs were differentially expressed during fruit set, being NAC, ERF/AP2, WRKY and bhLH the most representative families (Fig. 4), similar to found during siliques development in Arabidopsis and fruit tomato [22, 32]. This indicates that rate between the different families of TFs could affect fruit development in zucchini.
A further focus was made on the MADS-box and HD-ZIP families, which play key roles in fruit, set, via down-regulating their transcription in tomato and Arabidopsis, and contain many members that are functionally well-characterized during fruit development [22, 32]. TF genes in this group are AGAMOUS like genes, MADS-box genes whose expression diminished in developing ovules when the integuments appear in Arabidopsis . In the case of zucchini, when the fruit is pollinated in both cultivars, three AGAMOUS like genes were downregulated, AGL16, SOC1 (AGL20), CMB1 (AGL3), indicating the development of fertilized ovules. Particularly, downregulation of AGL16, involved in repression of floral transition, suggested that the transition from floral stage to fruit stage occurred . However, there were no common genes between cultivars during auxin treatment suggesting that different pathways related to MADS-box proteins were induced (Additional file 5).
TFs of HD-Zip family reduced their expression during zucchini fruit set as previously reported in tomato fruit . Downregulation of HAT5 and strong downregulation of ATHB-40 were also occurred in tomato fruit (Additional file 5), showing that these genes could play key roles in fruit development regulation . In addition, this downregulation may be mediated by hormones accumulation in fruit tissues induced by pollination or auxin treatment . Increased levels of auxins and ABA during fruit set may produce downregulation of ATHB-40  and ATHB7 , respectively.
Carbohydrate metabolism is induced by pollination and auxin treatment
Carbohydrates are major category of compounds that include among others: reducing (glucose and fructose), non reducing (sucrose) sugars, starches and cellulose which play an important role in the structure, function of all cells and are crucial factors that influence fruit quality [35, 36]. High levels of glucose and fructose were found after pollination and auxin treatment in zucchini fruit, indicating that the imported photoassimilates becoming mainly glucose and fructose in these first phases of fruit development and to a lesser extend to sucrose and starch (Fig. 5). Increased of these carbohydrates could be a consequence of glycolisis activation by enolase 1, glyceraldehyde-3-phosphate dehydrogenase and fructose-bisphosphate aldolase during pollinated fruit set . Enolase 1 synthesizes pyruvate from D-glyceraldehyde 3-phosphate , glyceraldehyde-3-phosphate dehydrogenase (GAPC-1) catalyzes the first step of the pathway by converting D-glyceraldehyde 3-phosphate into 3-phospho-D-glyceroyl phosphate , and fructose-bisphosphate aldolase is involved in a subpathway that synthesizes D-glyceraldehyde 3-phosphate and glycerone phosphate from D-glucose (Additional file 6) . These enzymes are essential not only for carbohydrate metabolism but also for proper fruit development, because they can regulate seed and fruit development. Gapc-1 null mutant showed low seed number, altered embryo development and apical morphological alterations in siliques in Arabidopsis . Thus, its over-expression may be caused the proper fruit growth observed in pollinated fruit in zucchini .
In the case of auxin treatment, high levels of fructose and glucose may be due to the activation of fructose 6-phosphate 1-phosphotransferase, pyruvate kinase and phosphoenolpyruvate carboxykinase 1, enzymes also involved in glycolisis and gluconeogenesis . Fructose 6-phosphate 1-phosphotransferase catalyzes the phosphorylation of D-fructose 6-phosphate and can thus function both in glycolysis and gluconeogenesis , and pyruvate kinase catalyses the irreversible synthesis of pyruvate and ATP, which are both used in multiple biochemical pathways .
Complex carbohydrates as starch are formed by linkage of monosaccharides as glucose and fructose , however low levels of starch were found during pollinated/ auxin treated fruit set in spite of over expression of enzymes related to starch (Fig. 5, Additional file 6). These evidences suggest that fruit set induced the availability of glucose and fructose free, not linked in complex polymers. Nevertheless, glucose and fructose levels were lower during parthenocarpic fruit set (Fig. 5), despite the activation of pyruvate kinase or glyceraldehyde-3-phosphate (Additional file 6).
Auxin, ethylene and gibberellin interaction regulates fruit set in zucchini
The activation of auxin and gibberellin signalling pathways, associated with repression of ethylene signals, trigger fruit set in Arabidopsis or tomato, , but the role of these hormones is still poorly understood in zucchini fruit set, probably because only few components of these signalling pathways have been indentified and few researches on fruit set has been carried out [5, 43].
Auxins play a critical role in the development of zucchini fruit [42, 44, 45]. RNA-seq showed that genes related to this hormone were the most representative in pollinated fruit and auxin treated fruit. Processes related to auxin signalling were enriched during pollination (Fig. 7) and application of synthetic auxins induced fruit growth in zucchini, causing the development of the ovary into an induced parthenocarpic fruit, suggesting its key role in fruit set [21, 45]. The majority of DEGs related to auxin were up-regulated, and most of them belong to Aux/IAA and PIN families (Fig. 7). Current models suggest that Aux/IAA genes encode nuclear proteins which form heterodimers with ARFs, and these heterodimers restrain the transcription of the early auxin response genes [46,47,48]. However, it was found over-expression of IAA4, IAA14 and IAA16 during fruit set not only in zucchini but also in tomato . This particular regulation of these IAAs during fruit set, suggest that a minimum level of Aux/IAAs is required in order to create a negative feedback loop in the auxin signal transduction pathway, which enables the plant to fine-tune the strength of the auxin response. Moreover, over-expression of PIN1, PIN6 and WAT-1 found in zucchini fruit (Fig. 7), suggested that auxin response is also mediated by auxin polar transport in fruit tissues . This regulation is translated into cell division and into cross-talk with other hormones such as ethylene and gibberellins [31, 42], most representative hormones together auxins during zucchini fruit (Fig. 7). Probably, downregulation of ethylene related genes (Fig. 7, Additional file 7) in zucchini fruit is caused by auxin. Addition of synthetic auxin or pollination has significant effects on the expression of the auxin related genes, regulating negatively ethylene production and signalling during fruit set . ACO1, ACO3, ETR2 and EIN3-like genes were downregulated during fruit set (Fig. 7), and this downregulation has also reported previously in zucchini fruit , indicating that that ethylene production needs to be prevented for a proper fruit set. This crosstalk between auxin and ethylene has also been reported in tomato , displaying that the mRNA levels of several ethylene biosynthesis genes and genes involved in ethylene signalling decreased after pollination and recently, ETR2 inhibition resulted in earlier fruit set .
At the same time, downregulation of ABA biosynthesis related genes (Fig. 7), NCED1 and NCED3, occurred, suggested that auxin response also may imply the attenuation of ABA response . Conversely, cluster analysis of DEGs revealed the strongly over-expression of a key role in gibberellin signalling (Fig. 7), GA20ox1 and downregulation of genes related with gibberellin catabolism (GID family), indicating that fruit set depends not only on auxin signalling activation but also on induction of gibberellin signalling pathway (Fig. 7). Similar regulation of these genes has been described previously in Arabidopsis and tomato during fruit development [8, 42]. In the case of Arabidopsis, GA20ox1 expression was increased in the ovules, whereas GA2ox gene expression was downregulated. Moreover, in the case of tomato fruit, GID1B and GA2ox genes were downregulated in contrast to GA20ox1 during fruit set, which was strongly up-regulated. However, this regulation does not appear to be independent of auxin, since auxin stimulates gibberellin biosynthesis through the transcriptional activation of GA20ox1 [8, 42]. These findings suggested that zucchini fruit development depends on the successful activation of the auxin and gibberellin signalling pathways as occurred in tomato.
Key genes for parthenocarpic fruit set in zucchini
Coordination between auxin, ethylene and gibberellin signalling pathways has proven to be an essential process for fruit set also in zucchini. The role of these hormones is still unknown in parthenocarpic process in zucchini. Differential expression analysis carried out in this assay have revealed a group of genes that showed the same regulation during pollinated/auxin treated and parthenocarpic fruit set (Additional file 7). ARF18, PIN1B, PIN8 and LAX1 were up-regulated like GA20ox1 and GASA4 gibberellin related genes This expression profile was similar to found in tomato or Arabidopsis, indicating the important role of these genes in parthenocarpy [22, 42]. Remarkably, over-expression of GA20ox1 gene seems to be crucial for fruit set and parthenocarpy process in zucchini. Moreover, this gene regulated fruit set in coordination of auxin signalling through the transcriptional activation of ARFs/IAAs as ARF18 . The activation of auxin signaling repressed ethylene response in zucchini . This repression was also observed during parthenocarpic fruit set through downregulation of ACO1, ETR2, ERF11 or ERF17 (Additional file 7). Gene expression analysis of this subset of genes has represented an important advance in parthenocarpy study in zucchini, and makes possible the introduction of this knowledge in breeding programmes.
In the present work, transcriptomic changes that take place in zucchini during fruit set have been analysed and compared in two contrasting cultivars for this process, Whitaker, parthenocarpic fruit set, and MUCU-16, non-parthenocarpic fruit set. This study also highlights the crucial role of some pathways including cell cycle, regulation of transcription and carbohydrate metabolism during fruit set in zucchini. Moreover, it was elucidated the important role of hormones during fruit set, establishing the activating role of auxins and gibberellins against the inhibitory role of ethylene. Functional analysis of RNA-Seq data have revealed different candidate genes that could be useful as markers for parthenocarpic selection in the current breeding programs of zucchini.
Plants materials, growth conditions and treatments
Two cultivars of Cucurbita pepo spp. pepo morphotype zucchini were used in this study, the non-parthenocarpic cultivar MUCU-16 (held in COMAV-UPV Genbank, accession BGV004370, https://www.comav.upv.es/index.php/databasesgermplasm/databases/genebank-database) , and the parthenocarpic cultivar Whitaker (developed in CORNELL College of Agriculture and Life Science, Geneva, New York, USA) . Those varieties are results of independent traditional breeding process. Seed germination and plant cultivation were performed following standard local practices, with the appropriate permissions to growth plants in a greenhouse of the IFAPA research centre in Almeria (Spain). A mean of 14 h photoperiod, mean daily air temperature of 24 /15 C day/night and relative humidity of 75% were registered in the greenhouse during the experiment.
Female flowers of the above cultivars, at 100 days post transplanting, were protected with paper bags in order to prevent pollen contamination on the day before anthesis. This was followed by three kinds of treatments on each protected flower: non-pollination, pollination and auxin treatment. Paper bags were removed at 1 DPA (day post anthesis) to claim that ovaries were not pollinated in non-pollination treatment. In the case of pollination, paper bags were removed in anthesis and female flowers were self-pollinated by hand early in the morning. In the case of auxin treatment paper bags were removed in anthesis and female flowers were treated with 0.5 ml of 0.8% of synthetic auxin “fruitone” (0.45% 1-naphthalene acetic acid, 1.2% 1-naphthaleneacetamide) for pollinated treatment and auxin treatment, respectively. Ovaries were collected for each treatment from different plants at 2 days post anthesis.
Six fruit samples were collected for RNA-seq analysis, i.e. unpollinated fruit of MUCU-16 (UF MU16), pollinated fruit of MUCU-16 (PF MU16), auxin treated fruit of MUCU-16 (AF MU16), unpollinated fruit of Whitaker (UF WHI), pollinated fruit of Whitaker (PF WHI) and auxin treated fruit of Whitaker (AF WHI). For samples collection, three biological replicates were performed (20 ovaries in each replicate). All collected samples were immediately frozen in liquid nitrogen and stored a − 80 °C until further analysis.
Total RNA was extracted using the RNeasy Plant Mini Kit from (QIAGEN, Germany). 6 libraries for sequencing were constructed for 6 RNA samples, and each sample was the resultant mix of three RNA extraction, one of each sample replication. RNA samples were sent to STABvida (Caparica, Portugal) for RNA-Seq analysis. RNA quality control was performed to evaluate integrity and concentration using Agilen 2200 Tape Station system (Agilent technologies, CA, USA). All the samples were within suitable parameters, RNA amount > 10 μg and RIN > 6.9. Library construction of cDNA molecules was carried out using Illumina TruSeq Stranded mRNA Library Preparation Kit, following manufacter’s instructions. Generate DNA fragments (cDNA libraries) were sequenced in the Illumina HiSeq 2500 platform, using 100 bp paired-end sequencing reads.
Gene quantification and differential expression analysis
Generated sequence raw data was analyzed using CLC Genomics Worbench 10.0.1. High quality sequencing reads were mapped against the Curcubita pepo transcriptome v3 and a processed version of this transcriptome after removal of alternatives sequences, using the following parameters length fraction of 0.80 (at least 80% of the alignment must match the reference sequence before the read is included in the mapping) and similarity fraction of 0.80 (identity should be at least 80% for the read to be included in the mapping). Mapping served to determine the gene expression levels based on the TPM (Transcripts per Million) . Bioconductor package EdgeR was applied to identify differentially expressed genes with multi-factorial statistical analysis tool based on a negative binomial model . A p-value could denote its expression difference between two libraries, and false discovery rates (FDRs) were used to determine the threshold of p-value. Differentially expressed genes were filtered using standard conditions, fold change (≥ 2 or ≤ − 2) and FDR p-Value ≤0.05 [53, 54].
Gene annotation and gene ontology term enrichment analysis of DEGs
Gene annotation was performed by comparing sequences using algorithm blastn or blastx (cut off e-value of 1e-25) implemented by the program Blast2GO v.4.3.1  with public databases, Arabidopsis genome and Arabidopsis proteins from TAIR , Swissprot  and GenBank non redundant nucleotide database (nr) from NCBI . Also, blastn search comparison (cut off e value of 1e-25) was performed with Cucurbita pepo ssp. ovifera , and Cucurbita maxima , and with the last versions of Cucurbita pepo genome from Cucurbigene ( https://bioinf.comav.upv.es/downloads/zucchini/genome_v4.1 /). DEGs were scanned against InterPro database  to assign Gene Ontology terms to zucchini genes , and GO annotations were mapped to the plant-specific GO slim ontology using Blast2GO with default parameters. Enzyme code mapping was perfomed using GO annotations and pathways The GO term enrichment analysis was conducted using a cut-off p-value of 0.05 for significant represented GO Terms and corrections for p-value were performed using FDR (false discovery rate). KEGG  was used to identified metabolic pathways in DEGs for further understanding genes functions.
Total sugar measurement
Quantification of sugars in zucchini fruits was conducted using the same samples that have been applied to RNAseq analysis. 12 mg of fruit material were extracted twice with 80% ethanol solution (80% ethanol, 2.5 mM HEPES, pH 7.5) at 95 °C for 30 min, followed by one extraction with 50% ethanol solution (50% ethanol, 2.5 mM HEPES, pH 7.5). The supernatants were combined and used for assaying total soluble sugars, glucose, fructose and sucrose by measuring the difference of absorbance at 340 nm in buffer (75 mM HEPES/KOH pH 7, 2.3 mM ATP, 1 mM NADP and glucose-6-phosphate dehydrogenase) after sequencing adding of hexokinase, phosphoglucose, isomerase, and invertase. For starch determination, the pellets of the ethanol extraction were solubized by heating them to 95 °C in 0.1 M NaOH for 30 min. After acidification with an HCl/sodium-acetate mixture pH 4.9, part of the suspension was digested overnight with amyloglucosidase and α-amylase. The glucose content of the supernant was then used to assess the starch content of the sample by measuring the difference in absorbance at 340 nm after adding hexokinase in the same buffer mentioned above.
Auxin treated fruit
Auxin response factor
Basic helix-loop-helix protein
Basic local alignment search tool
Cluster database at high identity with tolerance
Cyclin dependet kinase
Differentially expressed gene
Ethylene response factor
Homeobox-leucine zipper protein
Auxin influx carriers
- PIN :
Auxin efflux carriers
- SAUR :
Auxin responsive proteins
Varga A, Bruinsma J. Tomato. In: Monselise SP, editor. CRC Handbook of Fruit Set and Development. Boca Raton: CRC Press; 1986. p. 461–80.
Nepi M, Cresti L, Guarnieri M, Pacini E. Effect of relative humidity on water content, viability and carbohydrate profile of Petunia hybrid and Cucurbita pepo pollen. Plant Syst Evol. 2010;284:57–64.
Gustafson FG. Parthenocarpy: natural and artificial. Bot Rev. 1942;8:599–654.
Robinson RW, Reiners S. Parthenocarpy in summer squash. Hortscience. 1999;34:715–7.
Pomares-Viciana T, Die J, Del Río-Celestino M, Román B, Gómez P. Auxin signalling regulation during induced and parthenocarpic fruit set in zucchini. Mol Breeding. 2017;37:56.
Ozga JA, Reinecke DM. Hormonal interactions in fruit development. J Plant Growth Regul. 2003;22:73–81.
Kim IS, Okubo H, Fujieda K. Endogenous levels of IAA in relation to parthenocarpy in cucumber (Cucumis sativus L). Sci Hortic. 1992;52:1–8.
Olimpieri I, Siligato F, Caccia R, Mariotti L, Ceccarelli N, Soressi GP, et al. Tomato fruit set driven by pollination or by the parthenocarpic fruit allele are mediated by transcriptionally regulated gibberellin biosynthesis. Planta. 2007;226:877–88.
Cui L, Zhang T, Li J, Lou Q, Chen J. Cloning and expression analysis of Cs-TIR1/AFB2: the fruit development-related genes of cucumber (Cucumis sativus L.). Acta Physiol Plant. 2014;36:139–49.
De Jong M, Wolters-Arts J, Feron R, Mariani C, Vriezen WH. The Solanum lycopersicum auxin response factor 7 (SlARF7) regulates auxin signalling during tomato fruit set and development. Plant J. 2009;57:160–70.
Wang H, Jones B, Li Z, Frasse P, Delalande C, Regad F, Chaabouni S, Latché A, Pech JC, Bouzayen M. The tomato aux/IAA transcription factor IAA9 is involved in fruit development and leaf morphogenesis. Plant Cell. 2005;17(10):2676–92.
Goetz M, Vivian-Smith A, Johnson SD, Koltunow AM. AUXIN RESPONSE FACTOR 8 is a negative regulator of fruit initiation in Arabidopsis. Plant Cell. 2006;18(8):1873–86.
Mazzucato A, Cellini F, Bouzayen M, Zouine M, Mila I, Minoia S et al. A TILLING allele of the tomato aux/IAA9 gene offers new insights into fruit set mechanisms and perspectives for breeding seedless tomatoes. Mol Breeding. 2015; 35(22):1-15.
Blanca J, Cañizares J, Roig C, Ziarsolo P, Nuez F, Picó B. Transcriptome characterization and high throughput SSRs and SNPs discovery in Cucurbita pepo (Cucurbitaceae). BMC Genomics. 2011;12:104.
Wang Z, Gerstein M, Snyder M. RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009;10(1):57–63.
Da Fonseca RR, Albrechtsen A, Themudo GE, Ramos-Madrigal J, Sibbesen JA, Maretty L, et al. Next-generation biology: sequencing and data analysis approaches for non-model organisms. Mar Genomics. 2016;30:3–13.
Conesa A, Madrigal P, Tarazona S, Gomez-Cabrero D, Cervera A, McPherson A, et al. A survey of best practices for RNA-seq data analysis. Genome Biol. 2016;17:13.
Li J, Cui ZWJ, Zhang T, Guo Q, Xu J, Li J, et al. Transcriptome comparison of global distinctive features between pollination and parthenocarpic fruit set reveals transcriptional phytohormone cross-talk in cucumber (Cucumis sativus L). Plant Cell Physiol. 2014;55(7):1325–42.
Fu L, Niu B, Zhu Z, Wu S, Li W. CD-HIT: accelerated for clustering the next-generation sequencing data. Bioinformatics. 2012;28(23):3150–2.
Montero-Pau J, Blanca J, Bombarely A, Ziarsolo P, Esteras C, Martí-Gómez C, et al. De novo assembly of the zucchini genome reveals a whole genome duplication associated with the origin of the Cucurbita genus. Plant Biotechnol J. 2017. https://doi.org/10.1111/pbi.12860.
Vriezen WH, Feron R, Maretto F, Keijman J, Mariani C. Changes in tomato ovary transcriptome demonstrate complex hormonal regulation of fruit set. New Phytol. 2008;177:60–76.
Tang N, Deng W, Hu G, Hu N, Li Z. Transcriptome profiling reveals the regulatory mechanism underlying pollination dependent and parthenocarpic fruit set mainly mediated by auxin and gibberellin. PLoS One. 2015;10(4):e0125355.
Li J, Yan S, Yang W, Li Y, Xia M, Chen Z, et al. Transcriptomic analysis reveals the roles of microtubule-related genes and transcription factors in fruit length regulation in cucumber (Cucumis sativus L.). Sci Rep. 2015;26(5):8031.
Mironov V, De Veylder L, Van Montagu M, Inze D. Cyclin-dependent kinases and cell division in plants- the nexus. Plant Cell. 1999;11(4):509–22.
Perrot-Rechenmann C. Cellular responses to auxin: division versus expansion. Cold Spring Harb Perspect Biol. 2010;2(5):a001446.
De Veylder L, Beeckman T, Beemster GT, Krols L, Terras F, Landrieu I, et al. Functional analysis of cyclin-dependent kinase inhibitors of Arabidopsis. Plant Cell. 2001;13:1653–68.
Nieuwland J, Menges M, Murray JAH. The plant cyclins. In: Inze D, editor. Cell cycle control and plant development, vol. 2007. Oxford: Wiley-Blackwell Publishing; 2007. p. 33–61.
Menges M, Samland AK, Planchais S, Murray JA. The D-type cyclin CYCD3;1 is limiting for the G1-to-S-phasetransition in Arabidopsis. Plant Cell. 2006;18:893–906.
Boruc J, Mylle E, Duda M, De Clercq R, Rombauts S, Geelen D, et al. Systematic localization of the Arabidopsis core cell cycle proteins reveals novel cell division complexes. Plant Physiol. 2010;152(2):553–65.
Sampedro J, Cosgrove DJ. The expansin superfamily. Genome Biol. 2005;6:242.
Esmon CA, Tinsley AG, Ljung K, Sandberg G, Hearne LB, Liscum E. A gradient of auxin and auxin-dependent transcription precedes tropic growth responses. Proc Natl Acad Sci. 2006;103:236–41.
De Folter S, Busscher J, Colombo L, Losa A, Angenent GC. Transcript profiling of transcription factors genes during siliques development in Arabidopsis. Plant Mol Bio. 2004;56:351–3662004.
Son O, Cho HY, Kim MR, Lee H, Lee MS, Song E, et al. Induction of a homeodomain-leucine zipper gene by auxin is inhibited by cytokinin in Arabidopsis roots. Biochem Biophys Res Commun. 2005;326(1):203–9.
Olsson ASB, Engstroem P, Seoderman E. The homeobox genes ATHB12 and ATHB7 encode potential regulators of growth in response to water deficit in Arabidopsis. Plant Mol Biol. 2004;55:663–77.
Merrow SB, Hopp RJ. Storage effects on winter squashes. Associations between the sugar and starch content of and the degree of preference for winter squashes. J Agric Food Chem. 1961;9:321–6.
Berg JM, Tymoczko JL, Stryer L. Carbohydrates. In: Freeman WH, editor. Biochemistry. 5th ed. New York: W H Freeman; 2002.
Prabhakar V, Löttgert T, Gigolashvili T, Bell K, Flügge UI, Häusler RE. Molecular and functional characterization of the plastid-localized phosphoenolpyruvate enolase (ENO1) from Arabidopsis thaliana. FEBS Lett. 2009;583(6):983–91.
Rius SP, Casati P, Iglesias AA, Gomez-Casati DF. Characterization of Arabidopsis lines deficient in GAPC-1, a cytosolic NAD-dependent glyceraldehyde-3-phosphate dehydrogenase. Plant Physiol. 2008;148(3):1655–67.
Van der Linde K, Gutsche N, Leffers HM, Lindermayr C, Müller B, Holtgrefe S, et al. Regulation of plant cytosolic aldolase functions by redox-modifications. Plant Physiol Biochem. 2011;49(9):946–57.
Lim H, Cho MH, Jeon JS, Bhoo SH, Kwon YK, Hahn TR. Altered expression of pyrophosphate: fructose-6-phosphate 1-phosphotransferase affects the growth of transgenic Arabidopsis plants. Mol Cells. 2009;27(6):641–9.
Baud S, Wuillème S, Dubreucq B, De Almeida A, Vuagnat C, Lepiniec L, et al. Function of plastidial pyruvate kinases in seeds of Arabidopsis thaliana. Plant J. 2007;52:405–19.
De Jong M, Mariani C, Vriezen WH. The role of auxin and gibberellin in tomato fruit set. J Exp Bot. 2009;60(5):1523–32.
Martínez C, Manzano S, Megías Z, Garrido D, Picó B, Jamilena M. Involvement of ethylene biosynthesis and signalling in fruit set and early fruit development in zucchini squash (Cucurbita pepo L.). BMC Plant Biol. 2013;13:139.
Serrani JC, Fos M, Atarés A, Garcia-martinez JL. Effect of gibberellin and auxin on parthenocarpic fruit growth induction in the cv. micro-tom of tomato. J Plant Growth Regul. 2007;26:211–21.
Mapelli S. Changes in cytokinin in the fruits of parthenocarpic and normal tomatoes. Plant Sci Lett. 1981;22:227–33.
Ulmasov T, Hagen G, Guilfoyle TJ. Activation and repression of transcription by auxin-response factors. Proc Natl Acad Sci U S A. 1999;96:5844–9.
Ulmasov T, Hagen G, Guilfoyle TJ. Dimerization and DNA binding of auxin response factors. Plant J. 1999;19:309–19.
Tiwari SB, Hagen G, Guilfoyle TJ. Aux/IAA proteins contain a potent transcriptional repression domain. Plant Cell. 2004;16:533–43.
Switzenberg JA, Beaudry RM, Grumet R. Effect of CRC:: etr1-1 transgene expression on ethylene production, sex expression, fruit set and fruit ripening in transgenic melon (Cucumis melo L.). Transgenic Res. 2015;24(3):497-507.
Nitsch LM, Oplaat C, Feron R, Ma Q, Wolters-Arts M, Hedden P, et al. Abscisic acid levels in tomato ovaries are regulated by LeNCED1 and SlCYP707A1. Planta. 2009;229(6):1335–46.
Mortazavi A, Williams BA, McCue K, Schaeffer L, Wold B. Mapping and quantifying mammalian transcriptomes by RNA-seq. Nat Methods. 2008;5(7):621–8.
Robinson MD, McCarthy DJ, Smyth GK. Edger: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2008;26(1):139–40.
Raza K, Mishra A. A novel anticlustering filtering algorithm for the prediction of genes as a drug target. Am J Bio Engi. 2012;2(5):206–11.
Van Iterson M, Boer JM, Menezes RX. Filtering, FDR and power. BMCBioinformatics. 2010;11:450.
Conesa A, Götz S, Garcia-Gomez JM, Terol J, Talon M, Robles M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005;21:3674–6.
Berardini TZ, Reiser L, Li D, Mezheritsky Y, Muller R, Strait E, Huala E. The Arabidopsis information resource: making and mining the “gold standard” annotated reference plant genome. Genesis. 2015. https://doi.org/10.1002/dvg.22877.
Bairoch A, Apweiler R. The SWISS-PROT protein sequence database and its supplement TrEMBL. Nucleic Acids Res. 2000;28(1):45–8.
Johnson M, Zaretskaya I, Raytselis Y, Merezhuk Y, McGinnis S, Madden TL. NCBI BLAST: a better web interface. Nucleic Acids Res. 2008;36:W5–9.
Wyatt LE, Strickler SR, Mueller LA, Mazourek M. An acorn squash (Cucurbita pepo ssp. ovifera) fruit and seed transcriptome as a resource for the study of fruit traits in Cucurbita. Hortic Res. 2015;2:14070. https://doi.org/10.1038/hortres.2014.70.
Zhang A, Ren GA, Sun YA, Guo H, Zhang SA, Zhang FA, et al. A high-density genetic map for anchoring genome sequences and identifying QTLs associated with dwarf vine in pumpkin (Cucurbita maxima Duch.). BMC Genomics. 2015;16:1101.
Finn RD, Attwood TK, Babbit AB, Bork P, Bridge AJ, Chang HY. InterPro in 2017-beyond protein family and domain annotations. Nucleic Acids Res. https://doi.org/10.1093/nar/gkw1107.
Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Sherlock G. Gene ontology: tool for the unification of biology. Nat Genet. 2000;25(1):25–9.
Kanehisa M, Araki M, Goto S, et al. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 2008;36:480–4.
The authors thank Nicholas Davies for his help in grammatical revision of the manuscript, and Jose Pérez for his help in bioinformatics analysis.
Research worked is supported by the project RTA2014–00078 from the Spanish Institute of Agronomy Research INIA (Instituto Nacional de Investigación y Tecnología Agraría y Alimentaria) and also PP.AVA.AVA201601.7, FEDER y FSE (Programa Operativo FSE de Andalucia 2007–2013 “Andalucía se mueve con Europa”). TPV is supported by a FPI scholarship from RTA2011–00044-C02–01/02 project of INIA. The funding agencies were not involved in the design of the study, collection, analysis, and interpretation of data and in writing the manuscript.
Availability of data and materials
All supporting data can be found within the manuscript and its additional files.
The sequencing data have been deposited in the NCBI Sequence Read Archive (http://0-www.ncbi.nlm.nih.gov.brum.beds.ac.uk/sra/) for data accessibility. Samples accessions are SAMN10345060, SAMN10345061, SAMN10345062, SAMN10345063, SAMN10345064, and SAMN10345065 (BioProject ID PRJNA499783).
Ethics approval and consent to participate
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Tables “PF MU16 vs UF MU16” (List of DEGs during pollinated fruit set in MUCU-16 with their respective accession, length, Log2 fold change, relative fold change, p-value, FDR p-value and Bonferroni), “AFMU16 vs UF MU16 (List of DEGs during auxin treated fruit set in MUCU-16 cultivar with their respective accession, length, Log2 fold change, relative fold change, p-value, FDR p-value and Bonferroni), “PF WHT vs UF WHT” (List of DEGs during pollinated fruit set in Whitaker cultivar with their respective accession, length, Log2 fold change, relative fold change, p-value, FDR p-value and Bonferroni), “AF WHT vs UF WHT (List of DEGs during auxin treated fruit set in Whitaker cultivar with their respective accession, length, Log2 fold change, relative fold change, p-value, FDR p-value and Bonferroni), “UF WHT vs UF MU16” (List of DEGs during parthenocarpic fruit set with their respective accession, length, Log2 fold change, relative fold change, p-value, FDR p-value and Bonferroni). (XLSX 3123 kb) (XLSX 3122 kb)
Tables “PF MU16 vs UF MU16” (Functional annotation of DEGs during pollinated fruit set in MUCU-16 cultivar based on blast search), “AF MU16 vs UF MU16 (Functional annotation of DEGs during auxin treated fruit set in MUCU-16 cultivar based on blast search), “PF WHT vs UF WHT” (Functional annotation of DEGs during pollinated fruit set in Whitaker cultivar based on blast search), “AF WHT vs UF WHT (Functional annotation of DEGs during auxin treated fruit set in Whitaker cultivar based on blast search), “UF WHT vs UF MU16” (Functional annotation of DEGs during parthenocarpic fruit set based on blast search). (XLSX 6796 kb) (XLSX 6880 kb)
Tables “PF MU16 vs UF MU16” (Ontology terms of protein products for DEGs during pollinated fruit set in MUCU-16 cultivar), “AF MU16 vs UF MU16 (Ontology terms of protein products for DEGs during auxin treated fruit set in MUCU-16 cultivar), “PF WHT vs UF WHT” (Ontology terms of protein products for DEGs during pollinated fruit set in Whitaker cultivar), “AF WHT vs UF WHT (Ontology terms of protein products for DEGs during auxin treated fruit set in Whitaker cultivar), “UF WHT vs UF MU16” (Ontology terms of protein products for DEGs during parthenocarpic fruit set). (XLSX 2238 kb) (XLSX 2237 kb)
List of DEGs related to cell division identified by accession, gene annotation, log2 fold change and fold change. (XLSX 20 kb) (XLSX 19 kb)
Tables “PF MU16 vs UF MU16” (List of DEGs related to TFs during pollinated fruit set in MUCU-16 cultivar identified by accession, gene annotation, log2 fold change and fold change), “AF MU16 vs UF MU16 (List of DEGs related to TFs during auxin treated fruit set in MUCU-16 cultivar by accession, gene annotation, log2 fold change and fold change), “PF WHT vs UF WHT” (List of DEGs related to TFs during pollinated fruit set in Whitaker cultivar by accession, gene annotation, log2 fold change and fold change), “AF WHT vs UF WHT (List of DEGs related to TFs during auxin treated fruit set in Whitaker cultivar by accession, gene annotation, log2 fold change and fold change), “UF WHT vs UF MU16” (List of DEGs related to TFs during parthenocarpic fruit set by accession, gene annotation, log2 fold change and fold change). (XLSX 73 kb) (XLSX 72 kb)
Tables “PF MU16 vs UF MU16” (List of DEGs related to carbohydrate metabolism during pollinated fruit set in MUCU-16 cultivar identified by accession, gene annotation, log2 fold change, fold change, ontology terms and KEGG pathway), “AF MU16 vs UF MU16 (List of DEGs related to carbohydrate metabolism during auxin treated fruit set in MUCU-16 cultivar by accession, gene annotation, log2 fold change, fold change, ontology terms and KEGG pathway), “PF WHT vs UF WHT” (List of DEGs related to carbohydrate metabolism during pollinated fruit set in Whitaker cultivar by accession, gene annotation, log2 fold change, fold change, ontology terms and KEGG pathway), “AF WHT vs UF WHT (List of DEGs related to carbohydrate metabolism during auxin treated fruit set in Whitaker cultivar by accession, gene annotation, log2 fold change, fold change, ontology terms and KEGG pathway), “UF WHT vs UF MU16” (List of DEGs related to carbohydrate metabolism during parthenocarpic fruit set by accession, gene annotation, log2 fold change, fold change, ontology terms and KEGG pathway). (XLSX 54 kb) (XLSX 53 kb)
Tables “PF MU16 vs UF MU16” (List of DEGs related to hormones during pollinated fruit set in MUCU-16 cultivar identified by accession, gene annotation, log2 fold change, fold change and hormone pathway), “AF MU16 vs UF MU16 (List of DEGs related to hormones during auxin treated fruit set in MUCU-16 cultivar by accession, gene annotation, log2 fold change, fold change and hormone pathway), “PF WHT vs UF WHT” (List of DEGs related to hormones during pollinated fruit set in Whitaker cultivar by accession, gene annotation, log2 fold change, fold change and hormone pathway), “AF WHT vs UF WHT (List of DEGs related to hormones during auxin treated fruit set in Whitaker cultivar by accession, gene annotation, log2 fold change, fold change and hormone pathway), “UF WHT vs UF MU16” (List of DEGs related to hormones during parthenocarpic fruit set by accession, gene annotation, log2 fold change, fold change and hormone pathway). (XLSX 72 kb) (XLSX 71 kb)