Carmenza Montoya1†
Fernan Santiago Mejia-Alvarado1†
David Botero-Rozo1
Ivan Mauricio Ayala-Diaz1
Hernan Mauricio Romero1,2*- 1Oil Palm Biology and Breeding Research Program, Colombian Oil Palm Research Center—Cenipalma, Bogotá, Colombia
- 2Department of Biology, Universidad Nacional de Colombia, Bogotá, Colombia
Parthenocarpy is the development without fertilization of seedless fruits. In the oil palm industry, the development of parthenocarpic fruits is considered an attractive option to increase palm oil production. Previous studies have shown the application of synthetic auxins in Elaeis guineensis, and interspecific O×G hybrids (Elaeis oleifera (Kunth) Cortés × E. guineensis Jacq.) induces parthenocarpy. The aim of this study was to identify the molecular mechanism through transcriptomics and biology system approach to responding to how the application of NAA induces parthenocarpic fruits in oil palm O×G hybrids. The transcriptome changes were studied in three phenological stages (PS) of the inflorescences: i) PS 603, pre-anthesis III, ii) PS 607, anthesis, and iii) PS 700, fertilized female flower. Each PS was treated with NAA, Pollen, and control (any application). The expression profile was studied at three separate times: five minutes (T0), 24 hours (T1), and 48 h post-treatment (T2). The RNA sequencing (RNA seq) approach was used with 27 oil palm O×G hybrids for a total of 81 raw samples. RNA-Seq showed around 445,920 genes. Numerous differentially expressed genes (DEGs) were involved in pollination, flowering, seed development, hormone biosynthesis, and signal transduction. The expression of the most relevant transcription factors (TF) families was variable and dependent on the stage and time post-treatment. In general, NAA treatment expressed differentially more genes than Pollen. Indeed, the gene co-expression network of Pollen was built with fewer nodes than the NAA treatment. The transcriptional profiles of Auxin-responsive protein and Gibberellin-regulated genes involved in parthenocarpy phenomena agreed with those previously reported in other species. The expression of 13 DEGs was validated by RT-qPCR analysis. This detailed knowledge about the molecular mechanisms involved in parthenocarpy could be used to facilitate the future development of genome editing techniques that enable the production of parthenocarpic O×G hybrid cultivars without growth regulator application.
Introduction
Parthenocarpy is the development of seedless fruits without fertilization. Parthenocarpy is commercially exploited because fruit seeds could be undesirable due to their hard texture, bitter taste, and toxic compounds and allergens that can affect palatability. Therefore, seedless fruits are desirable to improve the quality of fresh and processed fruits. Furthermore, replacing seeds and seed cavities with edible tissue attracts consumers (Dhatt and Kaur, 2016). Developing parthenocarpic fruits is considered an attractive option to increase palm oil production. It allows for higher oil content by eliminating the kernel and increasing the mesocarp, resulting in more oil stored in the fruit. Previous studies have shown the efficacy of auxins to induce parthenocarpy in oil palm interspecific O × G [E. oleifera (Kunth) Cortés × E. guineensis Jacq.] hybrids without interfering with fruit development, oil synthesis, fatty acid profile, oil quality, or bunch components (Daza et al., 2020). Using the synthetic auxin naphthalene acetic acid (NAA), it is possible to reach more than 10 tons per hectare per year of High-Oleic Palm Oil (HOPO) without increasing the planted oil palm area by the induction of parthenocarpic fruits in the O × G hybrids (Romero et al., 2021).
NAA is widely used in Colombia and Ecuador, where interspecific O × G hybrids are planted in large extensions to replace the African oil palm, E. guineensis. The use of O × G hybrids is increasing because the African oil palm is susceptible to bud rot disease (BR) (Navia et al., 2014; Avila-Diazgranados et al., 2016), the most limiting factor of oil palm production in America (Sundram and Intan-Nur, 2017). The O × G hybrids resist BR, have high productivity potential, and produce a high-quality oil rich in oleic acid (Romero et al., 2022). However, the oil extraction rate under natural pollination conditions is low due to limitations in bunch filling that may be associated with low pollen viability and germinability and indehiscent peduncular bracts in female inflorescences that make pollen entry to the flowers difficult (Hormaza et al., 2012). Under natural pollination, the O × G oil potential is much lower than the African oil palm. However, the oil potential is increased with assisted pollination, a mandatory practice in the O × G interspecific hybrids (Rincon et al., 2013). It is even higher when NAA is applied to induce parthenocarpic fruits (Romero et al., 2021).
The idea of applying growth regulators to oil palm began 50 years ago. However, a genuine interest was shown recently due to its potential applicability. In African oil palm, the exogenous application of the synthetic auxin 2,4,5-tri-chlorophenoxy propionic acid (2,4,5-TP) has been associated with a seedless phenotype. The expression of genes related to the synthesis of auxins and gibberellins like EgGH3.8 (indole -3-acetic acid (IAA)-amido synthetase GH3.8), EgGH3.1 (IAA-amido synthetase GH3.1), EgARG7 (IAA induced ARG7 like), EgTAA3 (tryptophan amino transferase-related protein 3-like) and EgFMO1 (flavin-containing monooxygenase) were associated with parthenocarpic fruit development (Somyong et al., 2018). On the other hand, Htwe et al. (2022) studied gene expression in E. guineensis under natural pollination conditions and boron applications. Gene expression results showed eight genes associated with programmed cell death (PCD) and three groups of genes related to PCD: 4-coumarate-CoA ligase (4CL), S-RNase, and MADS-box. These findings suggest that a synergic upregulation of 4CL and MADS-box TFs activated the expression of PCD genes in the pollen tube, leading to fertilization failure and the seedless phenotype. Although these genes exhibited differential expression between treated and untreated fruits, more research is needed to identify which genes are involved in parthenocarpy and their function in the metabolic pathways of embryo formation.
No transcriptome studies show gene expression patterns in the O × G hybrid when NAA induces parthenocarpic fruit formation. Thus, this study aims to identify genes that respond to the application of NAA, to understand the molecular and genetic basis of induced parthenocarpic fruits. In the long term, the goal is to use genetic editing methodologies that could facilitate the production of parthenocarpic O×G hybrid cultivars to produce fruits without the external application of growth regulators.
Materials and methods
This research was carried out at the Palmar de la Vizcaína Experimental Field (CEPV) of the Colombian Oil Palm Research Centre - Cenipalma, located in the rural area of Barrancabermeja, Santander, Colombia, at 6°58ʹ863ʹʹN and 73°42ʹ186ʹʹW. With an elevation of 102 m above sea level, an average annual rainfall of 3,472 mm, an average temperature of 29°C, and relative humidity (HR%) between 72% and 77%.
Inflorescences collection, RNA isolation, and library construction
The study included 27 O × G interspecific hybrid palms (Coari × La Mé). Female inflorescences were selected in the phenological stage (PS) 601 or pre-anthesis I (Hormaza et al., 2012). Each inflorescence was enclosed in a polyester bag (PBS International, UK) to prevent cross-contamination with pollen or pollinating insects. Treatments of NAA, Pollen, and controls were applied when the isolated female inflorescences reached three different PS: 1) pre-anthesis III or PS 603 (Figure 1A), 2) anthesis or PS 607 (Figure 1B), and 3) fertilized female flower or PS 700 (Figure 1C). The phenological stages were chosen because, in PS 607 and PS 700, farmers are applying NAA solution to the O × G hybrid to obtain parthenocarpic fruits. The PS 603 was used to detect the earliest pollen and NAA signaling pathways and elucidate early differences between treatments. Three oil palms of each PS were treated with 150 mL of NAA solution (1,200 ppm), another independent three oil palms of each PS were treated with a mixture of Pollen and inert talcum powder (1:9 ratio), and finally, another three oil palms of each PS were kept without any application as control. After each treatment, samples were collected at three different times for each PS: 5 min post-treatment (T0), 24 h post-treatment (T1), and 48 h post-treatment (T2). Eighty-one samples were collected from 27 treatments and three replicates per treatment. From each inflorescence, 20 spikes were collected and immediately immersed in liquid nitrogen. The samples were homogenized with liquid nitrogen before RNA isolation. Total RNA was extracted from 100 mg of inflorescence tissue using the NucleoSpin RNA Plant and Fungi Mini kit for RNA from plants and fungi (Macherey-Nagel SAS, Dueren, Germany) according to the manufacturer’s protocol. The concentration and quality of total RNA were examined using an Agilent 2100 Bioanalyzer Instrument (Santa Clara, CA, United States). RNA samples with RIN ≥7 were precipitated with 3M sodium acetate and ethanol and sent to RNA-Seq library construction by Genohub Inc. (Austin, TX, United States). The libraries were constructed using the KAPA mRNA HyperPrep kit (Roche Diagnostics Corporation, Indianapolis, IN, United States) and sequenced in the Illumina NovaSeq 6000 platform (Illumina, Inc. San Diego, CA, United States).
FIGURE 1. Inflorescences of oil palm interspecific O × G hybrids and single female flower by phenology states. (A) PS_603, (B) PS_607, and (C) PS_700.
Quality control and DEGs analyses
The quality of the 81 raw samples was explored by FastQC (Andrews, 2010) and MultiQC (Ewels et al., 2016). The ten first bases were trimmed, and reads shorter than 25 were discarded with Trimmomatic (Bolger et al., 2014). Transcriptome assembly and downstream analyses were performed using the Trinity protocol, according to Haas et al. (2013). The clean reads were assembled with Trinity using in silico normalization by sample and default parameters (K25, maxC200, minC1, maxCV10000) (Grabherr et al., 2011). Longest open reading frames (LongORF) were predicted with TransDecoder. The transcriptome was functionally annotated through Trinotate with BLAST 2.11.0+ against the Uniprot database (Consortium, 2020); protein domains were searched with HMMER 3.3.2 against the Pfam database (Eddy, 2010; Punta et al., 2012). To evaluate the quality of the transcriptome, clean reads were mapped using Burrows-Wheeler Alignment (BWA) tool mem (Li and Durbin, 2009). Based on the lengths of the assembled transcriptome contigs, N50 length statistics were calculated.
The abundance of genes and transcripts per sample was calculated using the free-alignment method, Kallisto (Bray et al., 2016). These estimations were used for the matrix counts construction and values normalization. DEGs were statistically assessed by DESeq2 (Love et al., 2014) using the specific contrasts of treatments, PS, and time (i.e., NAA_603_T0 vs. Control_603_T0). Gene ontology enrichment was performed with Trinotate and GOseq (Young et al., 2010). Finally, an annotation report was generated with Trinotate for downstream analysis.
Gene co-expression network reconstruction
Using the DEGs, co-expression networks were constructed for NAA and Pollen contrasts for each PS 603, 607, and 700. A p-value ≤0.1 and a Log2FoldChange ≥ |2| were set in the DESeq2 algorithm as thresholds for assessing significant genes. The igraph R package was used to construct the co-expression networks. Each network was built based on correlating all gene pair distances among normalized read counts (variance stabilization in DESeq2). Only sample columns of the corresponding phenomena were used (i.e., for the NAA_603 network, only columns of T0, T1, and T2 of NAA and control at PS 603 were used). The edge among genes was considered when the correlation was ≥0.8. The “edge betweenness” function identified communities or modules of genes correlated among them (Blondel et al., 2008). HUB’s genes were determined according to a HUB score using the algorithm developed by Kleinberg (1999). Additional metrics of the networks were estimated as density, diameter (Wasserman and Faust, 1994), weighted diameter, centralization degree, centralization closeness, centralization betweenness (Brandes et al., 2016), and average path length (West, 2001).
Transcriptome and gene co-expression network analysis
To study the DEGs involved in NAA-induced parthenocarpy compared to the pollination mechanism, Venn diagrams were constructed using the web-based tool InteractiVenn (Heberle et al., 2015) to identify common and unique genes among PS (603, 607, and 700), treatment (NAA and Pollen) and time post-treatment (T0, T1, and T2). Furthermore, topological centrality measurements, such as HUB score, degree, betweenness, and closeness, were analyzed to identify HUB genes and representative modules in the co-expression network. Finally, the GO terms and categories with p-value ≤0.05 were significantly enriched. The overrepresented and enriched DEGs were analyzed in the REVIGO web application tool (Supek et al., 2011) to visualize how DEGs were clustered in categories of GO terms and their significance.
Validation of DEGs by quantitative real-time PCR
To confirm RNA-Seq transcriptome, 13 genes were selected according to the differential expression among PS, highly HUB scored, and literature reports (target and reference genes sequences are listed in Supplementary Table S1). Primer sequences were designed with the program Primer3web V4.1.0 (https://primer3.ut.ee/). The RT-qPCR reaction was conducted in a 10 µL reaction of QuantiNova SYBR Green One-Step RT-qPCR kit (QIAGEN, Germany) according to the manufacturer’s instruction. Three µL of RNA (50 ng/μL) were added to each reaction as a template. The RT-qPCR was run in a ROCHE LightCycler® 480 Real Time PCR System (Roche Diagnostics International AG, Rotkreuz, Switzerland). The PCR conditions were set as follows: 10 min at 50°C for reverse transcription, followed by the initial denaturation at 95°C for 2 min, then 40 cycles of denaturation at 95°C for 10 s, annealing at 58°C for 30 s and extension at 72°C for 30 s. For the melting curve analysis, the ramp rate was set at 0.06°C/s from 60°C to 95°C. Relative expression for each gene was calculated using the delta-delta of Ct method (∆∆Ct), and CYP2 and GRAS genes were utilized as normalizers. The association between RNA-Seq and RT-qPCR results was established by a correlation coefficient in R and plotted with the ggplot package.
Results
Transcriptome sequencing (RNA-seq) analysis
On average, the 81 raw samples had 48 million reads per sample and 7.8 billion total reads. After cleaning with Trimmomatic, around 34 million reads per sample and a total of 5.5 billion reads were obtained. 445,920 genes and 797,051 transcripts (isoforms) were assembled with Trinity. Most of the paired reads mapped correctly with BWA tool mem. At least half of the assembled bases (N50) were found in contigs at least 1,765 bases in length, 1,206, when the longest isoforms per gene were used. Sequencing generated 150 bp paired-end reads.
Gene co-expression networks
Network metrics comparing NAA and Pollen at the three stages are summarized in Supplementary Table S2. The Pollen had a smaller number of nodes than NAA. Also, in Pollen, there was a decrease in the number of nodes from PS 603 to 607 and an increase in PS 700. In contrast, NAA networks predicted fewer nodes in PS 603 but increased around three times in PS 607 and 700; also, there was a considerable difference between these last stages (488 nodes). Accordingly, densities in the PS 607 and 700 NAA networks were smaller than in the other networks by one order of magnitude; thus, the proportion of present edges from all possible edges in the network was smaller than in different networks. The centralization closeness of NAA_603 was one order of magnitude greater than the other networks, so the inverse of the node’s average geodesic distance to other nodes in the network was the greatest in this network. This was expected since if the number of nodes is smaller, the geodesic distance among nodes is smaller; however, it is interesting to note that the Pollen_607 network had a similar number of nodes network, but it did not have a centralization closeness value in the same magnitude.
The general co-expression network analysis of NAA phenomena simultaneously activated several pathways and genes in PS 603, PS 607, and PS 700 through a modulated response (Figure 2). In the co-expression network for PS 603_NAA (Figure 2A), the under-expressed genes TRINITY_DN47848_c0_g3 and TRINITY_DN12574_c3_g1 were predicted as HUB genes. The co-expression network in PS 607_NAA (Figure 2B) was constructed from 3,145 DEGs, and the gene TRINITY_DN22399_c0_g1 (Doubtful hypothetical protein ∼ unknown_gene) was predicted as HUB gene. This gene was highly over-expressed in T0, T1, and T2. TRINITY_DN56636_c0_g2 (possible protein SCAR2-like) was found over-expressed and predicted as the HUB gene in the co-expression network for PS 700_ NAA (Figure 2C).
FIGURE 2. The co-expression network for each phenology stage under NAA (A–C) and the pollen (D–F) treatment. Each network was constructed with whole T0, T1, and T2 DEGs per treatment. Each node (sphere or bead-like shape) represents a gene, and groups of nodes with the same color indicate a module of genes. Black edges represent a direct correlation between genes, and red lines inverse correlation. The nodes’ size is proportional to the mean level expression of the gene represented by the node. The igraph R package was used to construct the phenology stage-specific co-expression networks under NAA and Pollen treatment.
The co-expression network to predict the pollination mechanism (Pollen treatment) was constructed with fewer DEGs than the NAA network. It displayed less modularity and diameter compared to the NAA treatment. The PS 603 Pollen co-expression network (Figure 2D) predicted the genes TRINITY_DN63896_c2_g1 and TRINITY_DN24805_c2_g1 as HUB genes. In the PS 607_Pollen gene co-expression network (Figure 2E), TRINITY_DN91842_c1_g1 was predicted as the most critical gene in this model. Concerning the PS 700_Pollen gene co-expression network (Figure 2F), TRINITY_DN99294_c0_g1 and TRINITY_DN97297_c0_g1 (Probable metal-nicotianamine transporter YSL6) were predicted as HUBs.
DEGs and gene ontology analysis
A comparative analysis was performed among Pollen and NAA samples in different combinations of time and phenology stages to identify DEGs. In general, NAA treatment expressed differentially more genes than Pollen (Figure 3A). PS 607 and PS 700 under NAA treatment showed a higher level of DEGs than the Pollen treatment; additionally, PS 700 presented a significant quantity of under-expressed genes.
FIGURE 3. (A) DEGs over-expressed (Up) and under-expressed (Down). (B) Venn diagrams for PS 603. 1. Commons genes shared in T0 between NAA and Pollen treatment. 2. Commons genes shared in T1 between NAA and Pollen treatment. 3. Commons genes shared in T2 between NAA and Pollen treatment. Conventions: 607_48_Pollen = Phenology stage_time post treatment_treatment.
The general comparison of PS 603 (603_T0, 603_T1, and 603_T2 of Pollen and NAA treatments) did not share common genes among times and treatments (Supplementary Figure S1). However, expression changes were observed when time and treatments were compared (i.e., NAA _603_T0 vs. Pollen _603_T0). For the PS 603_T0 comparison (Figure 3B1), 56 genes were shared, with six genes drastically changing their expression at 5 minutes (T0) post-NAA application compared to Pollen. The most relevant genes were TRINITY_DN57151_c0_g2, TRINITY_DN51871_c1_g1 (Probable receptor-like protein kinase At1g33260), and TRINITY_ DN35796_c1_g1 (Calcium uniporter protein 5, mitochondrial). Then, in the PS 603_T1 (Figure 3B2), 73 genes were shared between NAA and Pollen treatments. Here, ten genes presented changes in expression level between treatments: the most relevant genes were TRINITY_DN57151_c0_g2, TRINITY_DN57151_c0_g2 (Subtilisin-like protease SBT1.5), and TRINITY_DN46085_c0_g1. Finally, in the PS 603_T2 (Figure 3B3), 20 genes were shared between NAA and Pollen. Notably, three genes were under-expressed in the NAA treatment: TRINITY_DN23722_c2_g1, TRINITY_DN7013_c2_g1, and TRINITY_DN261446_c1_g2 (Putative uncharacterized protein ART2).
The Venn diagrams for PS 607 represented 11 shared genes among T0, T1, and T2 between NAA and Pollen treatment (Figure 4A). From those, two genes presented relevant changes at the expression level between NAA and Pollen: TRINITY_DN141433_c0_g1 and TRINITY_DN2888_c0_g2 (possible protein Lateral Root Primordium 1-like). The PS 700 shared 14 DEGs among T0, T1, and T2, and NAA and Pollen treatment (Figure 4B). From those, two genes presented a relevant change in expression level between NAA and Pollen: TRINITY_DN795_c0_g2 (possible Alpha-galactosidase-1) and TRINITY_DN1989_c2_g2.
FIGURE 4. (A) General Venn diagrams for PS 607. (B) General Venn diagrams for PS 700. Unique genes for the phenology stage are highlighted in bold, and shared genes among all phenology stages and times are highlighted in the purple fog. Conventions: 603_48_Pollen = Phenology stage_time post treatment_treatment.
Several TF families were identified between NAA and Pollen treatment. The expression of the most relevant TF families such as WRKY, Ethylene-responsive transcription factor (ERFs), Auxin response factor (ARFs), MYB, bHLH, and Heat stress transcription factors (HS) was variable and dependent on the stage and time post-treatment (Table 1). In PS 603, all TF families mentioned above were over-expressed. But, in PS 607 and PS 700, they changed their expression levels, even in the lapse of time of 24 h. The ERF family was under-expressed in PS 607 and PS 700 under NAA treatment. WRKY was found present in all PS and times under NAA treatment. The PS 700_T0 point showed the same TFs between Pollen and NAA. However, the expression was inverse among them.
TABLE 1. Dynamics of the main TF expressed under Pollen and NAA treatment in oil palm interspecific O × G hybrids.
Concerning GO terms enrichment (Supplementary Table S3), REVIGO categorized the PS 603 under NAA treatment into biological process terms such as organic and aromatic compound biosynthesis, response to stimulus, macromolecules biosynthesis process, and organic cyclic compound biosynthetic process. The general enrichment for PS 603 under Pollen treatment involved the cellular amide metabolic process, intracellular signal transduction, and transcription. In PS 607, the most frequent functional annotations for NAA treatment were response to stimulus, signaling, and the catabolic process cellular response to stress. The most frequent terms for Pollen treatment in PS 607 were ncRNA Processing, cellular protein metabolic process regulation, and aminoglycan metabolic process. The most frequent terms for the PS 700 under NAA were cellular biosynthetic process, cellular nitrogen compound biosynthetic process, and macromolecule biosynthetic process. This stage under Pollen treatment was categorized with terms such as catabolic process, DNA recombination, and cell division.
Validation of DEGs by RT-qPCR
The expression of 13 DEGs was validated by RT-qPCR analysis. RT-qPCR data highly correlated with the RNA-Seq values predicted for each gene (Supplementary Figure S2). Four genes were over-expressed under NAA treatment: Probable WRKY transcription factor 33 (TRINITY_DN143542_c1_g1), 1-aminocyclopropane-1-carboxylate oxidase (TRINITY_DN11220_c0_g2), NAC domain-containing protein 68 (TRINITY_DN4308_c0_g2), Acetyl-coenzyme A carboxylase carboxyl transferase subunit alpha, chloroplastic (TRINITY_DN38599_c0_g1). Meanwhile, the genes Flowering-promoting factor 1-like protein 1 (TRINITY_DN23510_c0_g1) and TRINITY_DN63896_c2_g1 were found over-expressed under Pollen treatment. TRINITY_DN63896_c2_g1 has an unknown function and presented the highest HUB score in the Pollen co-expression network but was not significant in NAA treatment. The genes Indole-3-acetic acid-induced ARG7-like (XM_010942329.2), Auxin-responsive protein SAUR71-like (XM_010943300.2), Flavin-containing monooxygenase 1 (XM_019853925.1), Auxin-induced protein 22D (TRINITY_DN12260_c0_g1), Gibberellin-regulated protein 14 (TRINITY_DN123905_c0_g1) and Auxin-responsive protein IAA9 (TRINITY_DN3523_c0_g1) were found under-expressed in NAA treatment. While in Pollen treatment, the transcription factor MYC2 (TRINITY_DN62356_c0_g1) was under-expressed.
Discussion
We studied the transcriptome of three flower phenology stages of the interspecific O×G hybrid Coari x La Me, involved in the parthenocarpic fruit formation induced by NAA. The PS 603, or pre-anthesis III is the phenology stage before anthesis. In this stage, we expected to find initial changes in the molecular mechanism of parthenocarpic fruit development. The PS 603 co-expression network for the NAA application showed a non-modular form (round shape) compared with the PS 603 co-expression network of the Pollen application (Figures 2A, D). A possible biological interpretation is the simultaneous activation of several pathways caused by the exogenous application of NAA in inflorescence. After pollen application, we did not observe critical gene expression changes across the times (T0, T1, and T2). However, we found over-expressed genes involved in plant growth, adaptation, flowering, and fruit ripening, such as class I Heat Shock proteins, Alcohol dehydrogenase 1, and Tubulin alpha-3 chain (Jin et al., 2016; ul Haq et al., 2019). In both Pollen and NAA treatments, the transcriptions factors Ethylene-responsive and Probable WRKY transcription factor 33 were over-expressed across times T0 to T2 in PS 603, suggesting a possible mediation of natural processes (Table 1).
An important fact to understand the parthenocarpic phenomena is to analyze the proteins involved in the reception of the stimuli. In this sense, we identified receptor-like protein kinases RLKs over-expressed in the NAA treatment compared to the Pollen treatment in PS 603_T0. We hypothesize that probable receptor-like protein kinase At1g11050 responds to the exogenous application of NAA and activates components of the calcium-conducting subunit (Calcium uniporter protein 5). A signal mediated by secondary messengers, such as calcium, ROS, or H2O2, induces response gene expression. Similar results were found in the eggplant transcriptome analysis, where two calcium-binding proteins, PBP1, were upregulated in the natural parthenocarpic line PP05 (Chen et al., 2017). Following this signal, we compared the NAA and Pollen treatments in the PS 603_T0 and found the gene 4-coumarate—CoA ligase 2 (4CL2) over-expressed. This gene has been reported in producing CoA thioesters and phenolic compounds in Arabidopsis (Ehlting et al., 1999). In E. guineensis, synergic over-expression of 4CL2 and MAD-box TF activate the expression of programmed cell death genes and induce the formation of seedless phenotype (Htwe et al., 2022). The exogenous application of NAA might be perceived as a chemical stress and immediately triggers a response against this stress. Thus, the identified stimuli receptors and early candidate response genes might play an essential role in parthenocarpic fruit development in the O×G hybrid.
At 24 h post-NAA treatment (PS 603_T1), genes related to stress (Stress-related protein) and plant-pathogen interactions (Subtilisin-like protease SBT1.5) were over-expressed. These findings reinforce the chemical stress hypothesis mentioned above. However, four non-annotated genes were dramatically under-expressed in this stage compared with the Pollen application: TRINITY_DN261449_c1_g1, TRINITY_DN47630_c1_g1, TRINITY_DN261946_c0_g1, and TRINITY_DN65993_c0_g1. The gene co-expression network grouped these genes by expression profile in common modules; additionally, NCBI BLAST® (https://blast.ncbi.nlm.nih.gov/Blast.cgi) and Interpro (https://www.ebi.ac.uk/interpro/) alignments did not show significative results. Considering the expression level, further studies will be necessary to identify the role of these genes in the NAA-induced parthenocarpic process in the interspecific O×G hybrid. After 48 h of Pollen application (PS 603_T2), we observed over-expressed genes involved in cell wall biogenesis, such as Beta-galactosidase 2 and Alcohol dehydrogenase 1. Their presence could be related to the natural formation of reproductive structures. In contrast, under NAA treatment, the mechanism of the cell wall formation might be affected by the under-expression of genes related to the xylem cell wall thickening, Cellulose synthase A catalytic subunit 7. In this stage, the response to stimuli and organic cyclic compound metabolic process were the most common GO terms (Supplementary Table S2).
During anthesis (PS 607), processes co-occur around fertilization. Pollen tube elongation, pollen deposition on the stigmas, nuclear fusion, gametophyte development, and anther dehiscence are some of the main processes (Drews and Yadegari, 2002). Under Pollen treatment in PS 607, we identified a subset of genes that control some plant developmental processes, such as the transition from flower to seed and the female gametophyte development. The APETALA2/Ethylene responsive factor (AP2/ERF), NUCLEAR FUSION DEFECTIVE1 (NFD1) gene, and SCARECROW gene were found over-expressed in T0. Those genes are required for gametophyte development and cell specification of surrounding stem cells (Wysocka-Diller et al., 2000; Portereiko et al., 2006). At 24 h (T1) of Pollen treatment in PS 607, the gene ABC transporter D family member 1 (ABCD1) was over-expressed. This gene has been found in the peroxisomal membrane importing substrates into the peroxisome for β-Oxidation and enhancing fertilization (Footitt et al., 2007). Finally, at 48 h (T2), the gene MADS-box transcription factor 6 (MADS6) involved in floral patterning and gametophyte development was found over-expressed, too (Ohmori et al., 2009).
In contrast, no genes related to fertilization and gametophyte development described above were significant under NAA treatment in PS 607. Instead, we found many stress-related genes and enzymes modulating oxidative stress and senescence mechanisms. Among them, the genes serine/threonine-protein kinase (At1g54610), trans-resveratrol di-O-methyltransferase (ROMT), subtilisin-like protease SBT1.1 (SBT1.1) and peroxidases 18 (PER18), formerly reported in biotic and abiotic stress response, were found over-expressed in T0. ROMT, SBT1.1, and PER18 were induced by UV light and AlCl3 treatments in grape and pathogen infections in Arabidopsis (Schmidlin et al., 2008; Figueiredo et al., 2014). In T1, the genes chalcone synthase 4 (CHS4), SBT1.1, superoxide dismutase [Fe] 2 (sodM), PER18, and Protein ACCELERATED CELL DEATH 6 (ACD6) were found over-expressed as well. The CHS4 genes have been reported in flavonoid/isoflavonoid biosynthesis in plant defense (Dao et al., 2011), while sodM and PXG4 catalyze radical anion scavenging and remove H2O2, respectively. The gene ACD6 plays a vital role in plant response to abiotic stresses and nitrogen remobilization (Jasinski et al., 2021). Remarkably, in T2, genes related to defense against pathogens such as G-protein coupled receptor (Protein CANDIDATE G-PROTEIN COUPLED RECEPTOR 2 CAND2), receptor-like protein EIX2 and nitrate transporters (Protein NRT1/PTR FAMILY 6.3 NPF6.3) were found over-expressed. Interestingly, in PS 607, some genes related to vascular development and functioning and abaxial cell fate in flower meristems were common and highly over-expressed between NAA and Pollen treatment. The gene Dof zinc finger protein 4 (DOF4) acts at the crosstalk of various developmental pathways in vascular plants, while the YABBY 3 (YAB3) gene acts redundantly to specify abaxial identity in Arabidopsis (Siegfried et al., 1999; Le Hir and Bellini, 2013). In this sense, the NAA treatment favors the transition from flower to seed, activating stress-related genes and nitrogen remobilization pathways. These processes could disturb the molecular mechanism involved in the development of reproductive organs but do not affect the mechanism involved in the functional and structural fruit process.
The post-anthesis, or PS 700, presented the most DEGs among all stages. The female flower is fertilized, the stigma turns black, and the fruit development begins. At this point, the control of cell division and expansion are critical factors for early fruit development. At 24- or 48 h post-pollination begins the seed filling and the exponential growth of the fruit until the next 45 days (Hormaza et al., 2012). Under Pollen treatment, we observed a subset of genes related to the biosynthesis of phospholipids, long-chain fatty acid catabolic pathway, senescence, and expression of several TF in T0. The genes Acyl-CoA-binding domain-containing protein 3 (ACBP3), Peroxisomal (S)-2-hydroxy-acid oxidase GLO4 (GLO4) Probable linoleate 9S-lipoxygenase 4 (LOX1.4) were found highly over-expressed. On the other hand, the TF probable WRKY transcription factor 72 (WRKY72), ethylene-responsive transcription factor 71 (ERF071), NAC domain-containing protein 58 (NAC58), and transcription factor bHLH162 were found over-expressed as well. Among them, the TF WRKY72 has been reported in the abscisic acid signal and auxin transport pathway in Arabidopsis inflorescence (Yu et al., 2010), the TF ERF071 has been related to sweet cherry fruit development (Alkio et al., 2014) and with the cell expansion during olive fruit development (Camarero et al., 2023). However, NAC58 and bHLH162 have not been reported in fruit development.
At 24 h post Pollen treatment, we observed genes involved in calcium signaling in the cytosol, such as Calcium-dependent protein kinase 5 (CDPKs). Curiously, all the TF expressed in T0 were not differentially expressed in T1. Instead, other TF such as NAC domain-containing protein 4 (NAC4), transcription factor MYB78, and heat stress transcription factor B-1 (HSFB1) were over-expressed. At 48 h, we found over-expressed the genes GATA transcription factor 16 (GATA16) and thaumatin-like protein (TLP1), both are associated with developmental processes, tissue differentiation, and induced systemic resistance (Reyes et al., 2004; Singh et al., 2013). The pollination process in the oil palm O × G in PS 700 triggers signaling processes, carbohydrate catabolic processes, cell differentiation, and cell division. These processes are involved in flower maturation and fruit development and could be modulated by TF WRKY72, ERF071, NAC58, and bHLH162 in the first 24 h post-pollination.
Regarding NAA treatment in PS 700_T0, a structural molecular activity associated with flower transition, senescence, and sulfate transport was observed. The gene family Squamosa promoter-binding-like protein (SPL), auxin response factor 6 (ARF6), and Cyclic dof factor 2 (CDF2) were identified in this stage. These genes have been associated with flower and fruit development, cereal grain yield improvement, and senescence (Nagpal et al., 2005; Chen et al., 2010). Additionally, the gene Sulfate transporter 2.1, (SULTR2.1), which plays a vital role in the uptake of sulfate into Arabidopsis seeds, was found over-expressed (Awazuhara et al., 2001). The expression of auxin response factor 7 (ARF7), auxin response factor 12 (ARF12), probable WRKY transcription factor 33 (WRKY33), transcription factor bHLH162, EPIDERMAL PATTERNING FACTOR-like protein 2 (EPFL), and protein ASPARTIC PROTEASE IN GUARD CELL 1 (ASPG1) genes In T1, suggest an essential role in the plant epidermal cell growth process, dormancy, germination of the seed and viability.
Finally, in PS NAA_700_T2, the genes GIGANTEA, FAR1-RELATED SEQUENCE 6 (FRF), Gibberellin 2-beta-dioxygenase 8 (GA2ox8), and MYB78 were the most relevant over-expressed. GIGANTEA and FRF have been involved in flowering time regulation, light signaling, sucrose signaling, and starch accumulation (Mishra and Panigrahi, 2015; Ma and Li, 2018). While GA2ox8 has been related to the degradation of bioactive gibberellin in seedless citrus fruit (Zhang et al., 2017). The exogenous application of NAA in this stage accelerates the transition from flower to fruit. The activation of genes that control flowering and light could trigger senescence, while the activation of the starch metabolism could induce cell growth and seed filling. The reduction of GAs could be co-related with the previously described process and promote the development of the oil palm hybrid seedless genotype.
Essential genes involved in parthenocarpic fruit development
The PS 607 is the phenological stage of flowering. Here, several hormones regulate fruit development and determine the path of senescence or fruit development depending on the signaling process mediated by auxins or ethylene (An et al., 2020). Consistent with the knowledge about signal pathways involved in parthenocarpy phenomena, we emphasize the role of ARF, SAUR, IAA, and GA genes. The expression of ARFs plays an essential role in the plant’s growth and control of different developmental processes (Pandolfini et al., 2009). Comparative transcriptome analyses for parthenocarpic fruit development in some vegetables, such as eggplant, tomato, and cucumber, show that genes related to auxin signaling are the most affected (Chen et al., 2017; Li et al., 2017; Zhang et al., 2021).
The gene ARF8 presented differences among the treatments. Under NAA treatment, ARF8 was under-expressed in the post-anthesis (PS 700) stage for T0 and T1, but in Pollen, the treatment was not significant for DESeq. ARF8 has been reported as a critical negative regulator involved in parthenocarpic fruits in eggplant. It inhibits the formation of the AUX/IAA9 complex, blocks the fertilization process, and can activate the fruit initiation genes (Du et al., 2016). In contrast, the ARF7 was strongly over-expressed at post-anthesis (PS 700) for T1 under NAA treatment, but in Pollen, the treatment was not significant for DESeq. This is a particular finding because ARF7 mediates cross-talk between auxin and gibberellin signaling during tomato fruit set and development (De Jong et al., 2011). Okabe et al. (2019) recently observed that the expression of SlARF7 in tomatoes, encoding a negative regulator of auxin signaling, sharply decreased after anthesis. In our case, the expression of ARF7 was different, and further analysis should be followed to identify the role of this gene in NAA-induced parthenocarpy.
The SAUR is a family of auxin-responsive genes critical to plant growth, development, stress responses, and cell elongation in Arabidopsis (van Mourik et al., 2017; Stortenbeker and Bemer, 2019). Our results showed that SAUR genes exhibit a negative regulation of gene expression under NAA treatment at anthesis (PS 603) across the three time periods evaluated and did not show any relevant performance for the other stages/treatments. Similar results were reported in oil palm (E. guineensis), where SAUR71-like genes were downregulated in the auxin treatment, compared with controls (Somyong et al., 2018). In zucchini (Cucurbita pepo L.) parthenocarpy transcriptome, SAUR gene expression decreased under pollination and auxin treatment (Pomares-Viciana et al., 2019). In tomatoes, the SAUR genes showed an upregulated expression in parthenocarpic tomato line “R35-P” and downregulated expression in non-parthenocarpic tomato lines “R35-N” “R35-N” at −2 days after anthesis (Zhang et al., 2021). Therefore, our results confirm that this gene family is involved in generating parthenocarpic fruits in the O×G hybrid, and unknown factors inhibit the positive regulation.
According to Luo et al. (2018), the Auxin/Indole-3-Acetic Acid (Aux/IAA) gene family has been identified as short-lived nuclear proteins that control the expression levels of genes activated by the auxin response factor (ARF) family and exhibit an essential role in parthenocarpic fruit development. This study showed a more significant response of Aux/IAA genes (at least 12) under-expressed in pre-anthesis (PS 603) under NAA treatment across the three sample times. In contrast, only two Aux/IAA genes were expressed at PS 700 for T0 and T1. In the Pollen treatment, Aux/IAA genes tended to be under-expressed at post-anthesis PS 700_T0. These results show that the downregulated expression for IAA is a typical response across different plant species (Chai et al., 2019; Molesini et al., 2020). Furthermore, the over-expression of the oil palm EgDREB1 gene in tomatoes increased production of LeAux/IAA and may inhibit the transcription of auxin-associated genes and prevent fruit set due to failure in the fertilization process (Azzeme et al., 2017). In the transcriptome analysis of oil palm inflorescences treated with the synthetic auxin, 2, 4, 5-tri-chlorophenoxy propionic acid (2, 4, 5-TP), the EgARG7 gene was under-expressed in auxin-treated samples (Somyong et al., 2018), which corresponded with our results. Also, for the zucchini (C. pepo L.) parthenocarpy transcriptome, 10 Aux/IAAs were expressed differently, whereas IAA26 and IAA18 were found to downregulated in parthenocarpic and non-parthenocarpic cultivars (Pomares-Viciana et al., 2019).
Gibberellins (GAs) are critical phytohormones for germination, shoot elongation, tuber formation, flowering, and fruit set process (Serrani et al., 2007). Exogenous applied GAs can induce parthenocarpic fruits in species such as tomatoes, eggplant, apples, and others (Pandolfini et al., 2009; Galimba et al., 2019). In addition, the expression of GA genes has been studied in artificially induced parthenocarpy to understand their role in signaling pathways. In this study, the transcriptional profile of five isoforms of GA genes was under-expressed in the NAA treatment at post-anthesis (PS 700) across T0, T1, and T2. In contrast, only one isoform was downregulated at anthesis (PS 607) for T2, and two isoforms were downregulated for PS 700 for T0 and T1.
GA genes responded late to NAA and Pollen treatments in the O×G hybrid. Galimba et al. (2019) observed similar results in' Honeycrisp' apples. RNA expression profiles (clustered expression patterns of DEGs) were highly similar in GA3-induced and hand-pollinated fruits. In zucchini, a subset of genes associated with GAs biosynthesis modulated during pollen treatment showed downregulation patterns (Pomares-Viciana et al., 2019). On the other hand, applications with auxins also displayed some shared GAs genes with similar under-regulation patterns. In tomatoes, Okabe et al. (2019) showed a crosstalk between GA and auxin and suggested that stamen development negatively regulates fruit set by repressing the GA biosynthesis. These results reinforce the relevant role of GA genes for fruit set and parthenocarpy, independently of the species or auxin treatment. Additionally, in our study, the high quantity of GA genes expressed at post-anthesis (PS 700) under NAA application compared with Pollen could explain the effectivity of auxin for parthenocarpy fruit development when the anthesis stage has been finished.
Conclusion
The gene expression comparison between Pollen and NAA-treated inflorescences at PS 603 to induce parthenocarpic fruits revealed candidate genes involved in the early response of parthenocarpic induced by exogenous auxin. In PS 607, the inhibition of genes involved in karyogamy, the regulation of TF, and the under-expression of auxin response genes in NAA treatment, might accentuate the parthenocarpic phenomena induced by NAA (Figure 5).
FIGURE 5. Possible model to understand de parthenocarpy fruit development induced by NAA in oil palm interspecific O × G hybrids at the anthesis (PS-607). Gene expression levels are indicated for treatments T0 (5 min post-treatment), T1 (24 h post-treatment) and T2 (48 h post-treatment). The color scale corresponds to the L2FC (upper right). Conventions: CAND2, Protein candidate G-protein coupled receptor 2; CPK4, Calcium-dependent protein kinase 4; PI3K/4K, 1-phosphatidylinositol-3-phosphate 5-kinase; ECA4, Calcium-transporting ATPase 4, endoplasmic reticulum-type, and LRR, Probable LRR receptor-like serine/threonine-protein kinase At5g10290. The role of genes such as WRYK, ERF, MYB, HS, bHLH, SAUR, IAA, YABBY, DOF, and ARF10 are detailed in the discussion section. The figure was partly generated using Servier Medical Art, provided by Servier, licensed under a Creative Commons Attribution 3.0 unported license (https://creativecommons.org/licenses/by/3.0/; https://bioicons.com).
Our findings agreed with the transcriptomic results published about parthenocarpy phenomena and the predicted signal pathways (Htwe et al., 2022; Sharif et al., 2022). Critical genes such as ARF, ERF, SAUR, IAA, and GA exhibited transcriptional profiles similar to those reported previously, which implies that the oil palm O × G hybrids have similar genetic mechanisms of response to exogenously applied auxins to induce parthenocarpy (Somyong et al., 2018; Htwe et al., 2022).
The validation of candidate genes and HUBs found in this study by RT-qPCR confirms that the transcriptome is reliable and supports the findings of Daza et al. (2020) and Romero et al. (2021) for induction of parthenocarpic fruits in interspecific O×G hybrids. It is especially worthwhile to pinpoint the role of the HUB genes in fruit set and development, which could be a target for genome editing techniques to facilitate the production of parthenocarpic fruits free of plant growth regulator applications.
Data availability statement
The original contributions presented in the study are publicly available. This data can be found here: https://www.ncbi.nlm.nih.gov/. Accession number: PRJNA901998.
Author contributions
Conceptualization, CM, IA-D, HMR; Methodology, CM, HMR, IA-D, FSM-A, DB-R; Formal analysis, CM, FSM-A, DB-R, HMR; Investigation, CM, FSM-A, DB-R, IA-D, HMR; Writing—original draft, CM, FSM-A, DB-R; Writing—review and editing, HMR, CM, IA-D, DB-R; All authors have read and agreed to the published version of the manuscript.
Funding
This research was funded by the Colombian Oil Palm Promotion Fund (FFP) administered by Fedepalma.
Acknowledgments
The authors would like to thank Arley Fernando Caicedo Zambrano, Wilson Diaz Castillo, and Pedro Parada Suarez for technical support in collecting samples and the Physiology Laboratory personnel for carrying out this research study. Also, to Molecular Biology Laboratory staff for processing the samples. To Stephany Guataquira for providing the single female flower photographs.
Conflict of interest
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.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2023.1099489/full#supplementary-material
SUPPLEMENTARY FIGURE S1 | General Venn diagrams for PS 603. Unique genes for the phenology stage are highlighted in bold, and shared genes among all phenology stages and times are highlighted in the purple fog. Conventions: 603_48_Pollen = Phenology stage_time post treatment_treatment.
SUPPLEMENTARY FIGURE S2 | Validation of 13 genes by RT-qPCR. Association between RNA-Seq and RT-qPCR.
SUPPLEMENTARY TABLE S1 | Primer sequences used for RT-qPCR validation.
SUPPLEMENTARY TABLE S2 | Summary metrics of co-expression networks.
SUPPLEMENTARY TABLE S3 | Top 20 GO terms by phenology stage under NAA and Pollen treatment.
References
Alkio, M., Jonas, U., Declercq, M., Van Nocker, S., and Knoche, M. (2014). Transcriptional dynamics of the developing sweet cherry (Prunus avium L.) fruit: Sequencing, annotation and expression profiling of exocarp-associated genes. Hortic. Res. 11, 11. doi:10.1038/hortres.2014.11
An, J., Almasaud, R. A., Bouzayen, M., Zouine, M., and Chervin, C. (2020). Auxin and ethylene regulation of fruit set. Plant Sci. 292, 110381. doi:10.1016/j.plantsci.2019.110381
Andrews, S. (2010). FastQC: A quality control tool for high throughput sequence data [Online]. Available at: http://www.bioinformatics.babraham.ac.uk/projects/fastqc/
Avila-Diazgranados, R. A., Daza, E. S., Navia, E., and Romero, H. M. (2016). Response of various oil palm materials (Elaeis guineensis and Elaeis oleifera x Elaeis guineensis interspecific hybrids) to bud rot disease in the southwestern oil palm-growing area of Colombia. Agron. Colomb. 34 (1), 74–81. doi:10.15446/agron.colomb.v34n1.53760
Awazuhara, M., Takahashi, H., Watanabe-Takahashi, A., Hayashi, H., Fujiwara, T., and Saito, K. (2001). “Function of the sulfate transporter Sultr2;1 in seeds of Arabidopsis thaliana,” in Plant nutrition. Developments in plant and soil sciences (Dordrecht: Springer).
Azzeme, A. M., Abdullah, S. N. A., Aziz, M. A., and Wahab, P. E. M. (2017). Oil palm drought inducible DREB1 induced expression of DRE/CRT- and non-DRE/CRT-containing genes in lowland transgenic tomato under cold and PEG treatments. Plant Physiology Biochem. 112, 129–151. doi:10.1016/j.plaphy.2016.12.025
Blondel, V. D., Guillaume, J.-L., Lambiotte, R., and Lefebvre, E. (2008). Fast unfolding of communities in large networks. J. Stat. Mech. theory Exp. 2008 (10), P10008. doi:10.1088/1742-5468/2008/10/p10008
Bolger, A. M., Lohse, M., and Usadel, B. (2014). Trimmomatic: A flexible trimmer for Illumina sequence data. Bioinformatics 30 (15), 2114–2120. doi:10.1093/bioinformatics/btu170
Brandes, U., Borgatti, S. P., and Freeman, L. C. (2016). Maintaining the duality of closeness and betweenness centrality. Soc. Netw. 44, 153–159. doi:10.1016/j.socnet.2015.08.003
Bray, N. L., Pimentel, H., Melsted, P., and Pachter, L. (2016). Near-optimal probabilistic RNA-seq quantification. Nat. Biotechnol. 34 (5), 525–527. doi:10.1038/nbt.3519
Camarero, M. C., Briegas, B., Corbacho, J., Labrador, J., Gallardo, M., and Gomez-Jimenez, M. C. (2023). Characterization of transcriptome dynamics during early fruit development in olive (olea europaea L.). Int. J. Mol. Sci. 24 (2), 961. doi:10.3390/ijms24020961
Chai, P., Dong, S., Chai, L., Chen, S., Flaishman, M., and Ma, H. (2019). Cytokinin-induced parthenocarpy of san Pedro type fig (Ficus carica L.) main crop: Explained by phytohormone assay and transcriptomic network comparison. Plant Mol. Biol. 99 (4), 329–346. doi:10.1007/s11103-019-00820-2
Chen, X., Zhang, M., Tan, J., Huang, S., Wang, C., Zhang, H., et al. (2017). Comparative transcriptome analysis provides insights into molecular mechanisms for parthenocarpic fruit development in eggplant (Solanum melongena L.). PLoS ONE 12 (6), 1–17. doi:10.1371/journal.pone.0179491
Chen, X., Zhang, Z., Liu, D., Zhang, K., Li, A., and Mao, L. (2010). SQUAMOSA promoter-binding protein-like transcription factors: Star players for plant growth and development. J. Integr. Plant Biol. 52, 946–951. doi:10.1111/j.1744-7909.2010.00987.x
Consortium, T. U. (2020). UniProt: The universal protein knowledgebase in 2021. Nucleic Acids Res. 49 (D1), D480–D489. doi:10.1093/nar/gkaa1100
Dao, T. T. H., Linthorst, H. J. M., and Verpoorte, R. (2011). Chalcone synthase and its functions in plant resistance. Phytochem. Rev. 10 (3), 397–412. doi:10.1007/s11101-011-9211-7
Daza, E., Ayala-Diaz, I., Ruiz-Romero, R., and Romero, H. M. (2020). Effect of the application of plant hormones on the formation of parthenocarpic fruits and oil production in oil palm interspecific hybrids (Elaeis oleifera Cortés × Elaeis guineensis Jacq.). Plant Prod. Sci. 24, 354–362. doi:10.1080/1343943X.2020.1862681
De Jong, M., Wolters-Arts, M., García-Martínez, J. L., Mariani, C., and Vriezen, W. H. (2011). The Solanum lycopersicum AUXIN RESPONSE FACTOR 7 (SlARF7) mediates cross-talk between auxin and gibberellin signalling during tomato fruit set and development. J. Exp. Bot. 62 (2), 617–626. doi:10.1093/jxb/erq293
Dhatt, A. S., and Kaur, G. (2016). Parthenocarpy: A potential trait to exploit in vegetable crops: A review. Agric. Rev. 37 (4), 300–308. doi:10.18805/ag.v37i4.6460
Drews, G. N., and Yadegari, R. (2002). Development and function of the angiosperm female gametophyte. Annu. Rev. Genet. 36, 99–124. doi:10.1146/annurev.genet.36.040102.131941
Du, L., Bao, C., Hu, T., Zhu, Q., Hu, H., He, Q., et al. (2016). SmARF8, a transcription factor involved in parthenocarpy in eggplant. Mol. Genet. genomics 291 (1), 93–105. doi:10.1007/s00438-015-1088-5
Ehlting, J., Büttner, D., Wang, Q., Douglas, C. J., Somssich, I. E., and Kombrink, E. (1999). Three 4-coumarate: coenzyme A ligases in Arabidopsis thaliana represent two evolutionarily divergent classes in angiosperms. plant J. 19 (1), 9–20. doi:10.1046/j.1365-313x.1999.00491.x
Ewels, P., Magnusson, M., Lundin, S., and Käller, M. (2016). MultiQC: Summarize analysis results for multiple tools and samples in a single report. Bioinformatics 32 (19), 3047–3048. doi:10.1093/bioinformatics/btw354
Figueiredo, A., Monteiro, F., and Sebastiana, M. (2014). Subtilisin-like proteases in plant–pathogen recognition and immune priming: A perspective. Front. Plant Sci. 5, 1–4. doi:10.3389/fpls.2014.00739
Footitt, S., Dietrich, D., Fait, A., Fernie, A. R., Holdsworth, M. J., Baker, A., et al. (2007). The COMATOSE ATP-binding cassette transporter is required for full fertility in Arabidopsis. Plant Physiol. 144 (3), 1467–1480. doi:10.1104/pp.107.099903
Galimba, K. D., Bullock, D. G., Dardick, C., Liu, Z., and Callahan, A. M. (2019). Gibberellic acid induced parthenocarpic ‘Honeycrisp’apples (Malus domestica) exhibit reduced ovary width and lower acidity. Hortic. Res. 6, 41. doi:10.1038/s41438-019-0124-8
Grabherr, M. G., Haas, B. J., Yassour, M., Levin, J. Z., Thompson, D., Amit, I., et al. (2011). Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat. Biotechnol. 29, 644–652. doi:10.1038/nbt.1883
Haas, B. J., Papanicolaou, A., Yassour, M., Grabherr, M., Blood, P. D., Bowden, J., et al. (2013). De novo transcript sequence reconstruction from RNA-seq: Reference generation and analysis with trinity. Nat. Protoc. 8 (8), 1494–1512. doi:10.1038/NPROT.2013.084
Heberle, H., Meirelles, G. V., da Silva, F. R., Telles, G. P., and Minghim, R. (2015). InteractiVenn: A web-based tool for the analysis of sets through Venn diagrams. BMC Bioinforma. 16 (1), 169–177. doi:10.1186/s12859-015-0611-3
Hormaza, P., Mesa Fuquen, E., and Romero, H. M. (2012). Phenology of the oil palm interspecific hybrid Elaeis oleifera × Elaeis guineensis. Sci. Agric. 69 (4), 275–280. doi:10.1590/S0103-90162012000400007
Htwe, Y. M., Shi, P., Zhang, D., Li, Z., Xiao, Y., Yang, Y., et al. (2022). Programmed cell death may be involved in the seedless phenotype formation of oil palm. Front. Plant Sci. 13, 832017. 1–17. doi:10.3389/fpls.2022.832017
Jasinski, S., Fabrissin, I., Masson, A., Marmagne, A., Lécureuil, A., Bill, L., et al. (2021). ACCELERATED CELL DEATH 6 acts on natural leaf senescence and nitrogen fluxes in Arabidopsis. Front. Plant Sci. 11 (611170), 1–15. doi:10.3389/fpls.2020.611170
Jin, Y., Zhang, C., Liu, W., Tang, Y., Qi, H., Chen, H., et al. (2016). The the alcohol dehydrogenase gene family in melon (Cucumis melo L.): Bioinformatic analysis and expression patterns. Front. Plant Sci. 7, 670. doi:10.3389/fpls.2016.00670
Kleinberg, J. M. (1999). Hubs, authorities, and communities. ACM Comput. Surv. (CSUR) 31, 5–es. doi:10.1145/345966.345982
Le Hir, R., and Bellini, C. (2013). The plant-specific Dof transcription factors family: New players involved in vascular system development and functioning in Arabidopsis. Front. Plant Sci. 4, 164. doi:10.3389/fpls.2013.00164
Li, H., and Durbin, R. (2009). Fast and accurate short read alignment with Burrows–Wheeler transform. bioinformatics 25 (14), 1754–1760. doi:10.1093/bioinformatics/btp324
Li, J., Xu, J., Guo, Q. W., Wu, Z., Zhang, T., Zhang, K. J., et al. (2017). Proteomic insight into fruit set of cucumber (Cucumis sativus L.) suggests the cues of hormone-independent parthenocarpy. BMC Genomics 18, 896. 1–18. doi:10.1186/s12864-017-4290-5
Love, M. I., Huber, W., and Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15 (12), 550–621. doi:10.1186/s13059-014-0550-8
Luo, J., Zhou, J.-J., and Zhang, J.-Z. (2018). Aux/IAA gene family in plants: Molecular structure, regulation, and function. Int. J. Mol. Sci. 19 (1), 259. doi:10.3390/ijms19010259
Ma, L., and Li, G. (2018). FAR1-RELATED SEQUENCE (FRS) and FRS-RELATED FACTOR (FRF) family proteins in Arabidopsis growth and development. Front. Plant Sci. 9, 692. doi:10.3389/fpls.2018.00692
Mishra, P., and Panigrahi, K. C. (2015). Gigantea - an emerging story. Front. Plant Sci. 6, 8–15. doi:10.3389/fpls.2015.00008
Molesini, B., Dusi, V., Pennisi, F., and Pandolfini, T. (2020). How hormones and mads-box transcription factors are involved in controlling fruit set and parthenocarpy in tomato. Genes 11 (12), 1441. doi:10.3390/genes11121441
Nagpal, P., Ellis, C. M., Weber, H., Ploense, S. E., Barkawi, L. S., Guilfoyle, T. J., et al. (2005). Auxin response factors ARF6 and ARF8 promote jasmonic acid production and flower maturation. Development 132 (18), 4107–4118. doi:10.1242/dev.01955
Navia, E. A., Restrepo, E. F., and Romero, H. M. (2014). Response of six sources of oil palm planting materials from Malaysia planted in the eastern plains of Colombia to bud rot. J. Oil Palm Res. 26 (1), 73–83.
Ohmori, S., Kimizu, M., Sugita, M., Miyao, A., Hirochika, H., Uchida, E., et al. (2009). MOSAIC FLORAL ORGANS1, an AGL6-like MADS box gene, regulates floral organ identity and meristem fate in rice. Plant Cell 21 (10), 3008–3025. doi:10.1105/tpc.109.068742
Okabe, Y., Yamaoka, T., Ariizumi, T., Ushijima, K., Kojima, M., Takebayashi, Y., et al. (2019). Aberrant stamen development is associated with parthenocarpic fruit set through up-regulation of gibberellin biosynthesis in tomato. Plant Cell Physiology 60 (1), 38–51. doi:10.1093/pcp/pcy184
Pandolfini, T., Molesini, B., and Spena, A. (2009). Parthenocarpy in crop plants. Annu. Plant Rev. 38, 326–345.
Pomares-Viciana, T., Río-Celestino, D., Román, B., Die, J., Pico, B., and Gómez, P. (2019). First RNA-seq approach to study fruit set and parthenocarpy in zucchini (Cucurbita pepo L.). BMC plant Biol. 19 (1), 61–20. doi:10.1186/s12870-019-1632-2
Portereiko, M. F., Sandaklie-Nikolova, L., Lloyd, A., Dever, C. A., Otsuga, D., and Drews, G. N. (2006). NUCLEAR FUSION DEFECTIVE1 encodes the Arabidopsis RPL21M protein and is required for karyogamy during female gametophyte development and fertilization. Plant physiol. 141 (3), 957–965. doi:10.1104/pp.106.079319
Punta, M., Coggill, P. C., Eberhardt, R. Y., Mistry, J., Tate, J., Boursnell, C., et al. (2012). The Pfam protein families database. Nucleic Acids Res. 40, D290–D301. doi:10.1093/nar/gkr1065
Reyes, J. C., Muro-Pastor, M. I., and Florencio, F. J. (2004). The GATA family of transcription factors in Arabidopsis and rice. Plant Physiol. 134 (4), 1718–1732. doi:10.1104/pp.103.037788
Rincón, S. M., Hormaza, P. A., Moreno, L. P., Prada, F., Portillo, D. J., García, J. A., et al. (2013). Use of phenological stages of the fruits and physicochemical characteristics of the oil to determine the optimal harvest time of oil palm interspecific OxG hybrid fruits. Industrial Crops Prod. 49, 204–210. doi:10.1016/j.indcrop.2013.04.035
Romero, H. M., Daza, E., Ayala-Díaz, I., and Ruiz-Romero, R. (2021). High-Oleic Palm Oil (HOPO) production from parthenocarpic fruits in oil palm interspecific hybrids using naphthalene acetic acid. Agronomy 11 (2), 290. doi:10.3390/agronomy11020290
Romero, H. M., Guataquira, S., and Forero, D. C. (2022). Light interception, photosynthetic performance, and yield of oil palm interspecific O×G hybrid (Elaeis oleifera (Kunth) Cortés × Elaeis guineensis Jacq.) under three planting densities. Plants 11 (9), 1166. doi:10.3390/plants11091166
Schmidlin, L., Poutaraud, A., Claudel, P., Mestre, P., Prado, E., Santos-Rosa, M., et al. (2008). A stress-inducible resveratrol O-methyltransferase involved in the biosynthesis of pterostilbene in grapevine. Plant Physiol. 148, 1630–1639. doi:10.1104/pp.108.126003
Serrani, J. C., Fos, M., Atarés, A., and García-Martínez, J. L. (2007). Effect of gibberellin and auxin on parthenocarpic fruit growth induction in the cv Micro-Tom of tomato. J. Plant Growth Regul. 26 (3), 211–221. doi:10.1007/s00344-007-9014-7
Sharif, R., Su, L., Chen, X., and Qi, X. (2022). Hormonal interactions underlying parthenocarpic fruit formation in horticultural crops. Hortic. Res. 9, uhab024. Horticulture research. doi:10.1093/hr/uhab024
Siegfried, K. R., Eshed, Y., Baum, S. F., Otsuga, D., Drews, G. N., and Bowman, J. L. (1999). Members of the YABBY gene family specify abaxial cell fate in Arabidopsis. Development 126 (18), 4117–4128. doi:10.1242/dev.126.18.4117
Singh, N. K., Kumar, K. R. R., Kumar, D., Shukla, P., and Kirti, P. B. (2013). Characterization of a pathogen induced thaumatin-like protein gene AdTLP from Arachis diogoi, a wild peanut. PLoS ONE 8 (12), 1–18. doi:10.1371/journal.pone.0083963
Somyong, S., Walayaporn, K., Jomchai, N., Naktang, C., Yodyingyong, T., Phumichai, C., et al. (2018). Transcriptome analysis of oil palm inflorescences revealed candidate genes for an auxin signaling pathway involved in parthenocarpy. PeerJ 6 (e5975), 1–22. doi:10.7717/peerj.5975
Stortenbeker, N., and Bemer, M. (2019). The SAUR gene family: The plant’s toolbox for adaptation of growth and development. J. Exp. Bot. 70 (1), 17–27. doi:10.1093/jxb/ery332
Sundram, S., and Intan-Nur, A. M. A. (2017). South American bud rot: A biosecurity threat to south east asian oil palm. Crop Prot. 101, 58–67. doi:10.1016/j.cropro.2017.07.010
Supek, F., Bošnjak, M., Škunca, N., and Šmuc, T. (2011). REVIGO summarizes and visualizes long lists of gene ontology terms. PloS one 6, e21800. doi:10.1371/journal.pone.0021800
ul Haq, S., Khan, A., Ali, M., Khattak, A. M., Gai, W. X., Zhang, H. X., et al. (2019). Heat shock proteins: Dynamic biomolecules to counter plant biotic and abiotic stresses. Int. J. Mol. Sci. 20, 1–31. doi:10.3390/ijms20215321
van Mourik, H., van Dijk, A. D. J., Stortenbeker, N., Angenent, G. C., and Bemer, M. (2017). Divergent regulation of Arabidopsis SAUR genes: A focus on the SAUR10-clade. BMC plant Biol. 17 (245), 1–14. doi:10.1186/s12870-017-1210-4
Wasserman, S., and Faust, K. (1994). Social network analysis: Methods and applications (structural analysis in the social sciences). Cambridge: Cambridge University Press.
Wysocka-Diller, J. W., Helariutta, Y., Fukaki, H., Malamy, J. E., and Benfey, P. N. (2000). Molecular analysis of SCARECROW function reveals a radial patterning mechanism common to root and shoot. Development 127 (3), 595–603. doi:10.1242/dev.127.3.595
Young, M. D., Wakefield, M. J., Smyth, G. K., and Oshlack, A. (2010). Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 11 (R14), 1–12. doi:10.1186/gb-2010-11-2-r14
Yu, S., Ligang, C., Liping, Z., and Diqiu, Y. (2010). Overexpression of OsWRKY72 gene interferes in the abscisic acid signal and auxin transport pathway of Arabidopsis. J. Biosci. 35 (3), 459–471. doi:10.1007/s12038-010-0051-1
Zhang, S., Gu, X., Shao, J., Hu, Z., Yang, W., Wang, L., et al. (2021). Auxin metabolism is involved in fruit set and early fruit development in the parthenocarpic tomato “R35-P”. Front. Plant Sci. 12, 671713. doi:10.3389/fpls.2021.671713
Keywords: oil palm, interspecific O × G hybrids, NAA, auxins, transcriptome, gene coexpression networks, parthenocarpy
Citation: Montoya C, Mejia-Alvarado FS, Botero-Rozo D, Ayala-Diaz IM and Romero HM (2023) Parthenocarpy-related genes induced by naphthalene acetic acid in oil palm interspecific O × G [Elaeis oleifera (Kunth) Cortés × Elaeis guineensis Jacq.] hybrids. Front. Genet. 14:1099489. doi: 10.3389/fgene.2023.1099489
Received: 15 November 2022; Accepted: 08 March 2023;
Published: 20 March 2023.
Edited by:
Karthikeyan Adhimoolam, Jeju National University, Republic of KoreaReviewed by:
Vikranth Kumar Chandrasekaran, University of Missouri, United StatesMeilina Ong-Abdullah, Malaysian Palm Oil Board, Malaysia
Renganathan V. G., Tamil Nadu Agricultural University, India
Copyright © 2023 Montoya, Mejia-Alvarado, Botero-Rozo, Ayala-Diaz and Romero. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Hernan Mauricio Romero, hmromeroa@unal.edu.co
†These authors have contributed equally to this work and share first authorship