Transcriptome analysis of ovary culture-induced embryogenesis in cucumber (Cucumis sativus L.)

View article
Plant Biology

Introduction

Gametophyte cultures involve the stress-induced reprogramming of male or female gametophytes that develop into embryo-like structures and can be directly regenerated into completely homozygous, doubled haploid (DH) plants. The study of male gametophyte (microspore) embryogenesis began in the 1960s (Guha & Maheshwari, 1964). In vitro microspore embryogenesis can be easily monitored and provide a convenient experimental platform for large-scale physiological and biochemical analyses (Malik et al., 2007). However, the study of gynogenesis is more complicated. The embryo sac is small and embedded in the surrounding tissues, making early embryos difficult to observe and isolate. For plants that are recalcitrant to androgenesis (e.g., male-sterile and dioecious plants), gynogenesis is is a worthwhile investigation (Pazuki et al., 2018). In viro gynogenesis has mainly been used on the Cucurbitaceae, and previous studies have revealed that this technique is still imperfect with a low embryogenesis rate (Metwally et al., 1998; Gémes-Juhász et al., 2002; Shalaby, 2007; Diao et al., 2009; Li et al., 2013; Plapung, Khumsukdee & Smitamana, 2014; Tantasawat et al., 2015). The study of the embryogenesis mechanism has therefore remained unexplored, and the understanding of the mechanisms underlying gynogenesis is extremely limited.

In general, embryogenesis induction using in vitro culture is a stress-induced phenomenon. Microspores are cultured and develop into embryos in vitro, and stress treatments involving cold shock, heat shock or hormones are used as trigger factors to induce gametocytes to follow the sporophyte development pathway (Touraev, Vicente & Heberle-Bors, 1997). The genes associated with the reprogramming phase and the early stage of embryogenesis were successfully characterized and include AGAMOUS-LIKE15-related (AGL15-related) genes and the AGAMOUS-LIKE15 (AGL15) genes, the BABY BOOM (BBM) genes, and the HEAT SHOCK PROTEIN (HSP) genes (Perry & Fernandez, 1999; Boutilier et al., 2002; Maraschin, 2005). Recently, functional genomics has made it possible to identify more genes associated with different microspore embryogenesis stages. When analyzing in vitro microspore cultures of Brassica napus, it was found that genes related to embryogenesis, such as BBM1, LEC1 and LEC2, were expressed as early as 2–3 days after microspore culture. While embryogenesis was clearly established by 7 d of culture, the genes that were specifically expressed included ABSCISIC ACID INSENSITIVE3, ATS1, LEC1, LEC2 and FUSCA3 (Malik et al., 2007). Genes expression associated with metabolism, chromosome remodeling, transcription, and translation signaling was up-regulated during the stress treatment stage in tobacco (Hosp et al., 2007). Similarly, genes related to metabolism, cell wall, cell membrane, cell tissue control, cell communication, and signal transduction were detected in the early stage of rape embryogenesis (Rosa, Ana & María, 2013). High levels of BBM and LEC gene expression were confirmed in the early embryonic development of sweet pepper anther cultures (Irikova, Grozeva & Denev, 2012). The transition from microspores to developing embryos is mainly manifested in the induction of transcription factor genes that play an important role in early embryogenesis. Many genes involved in hormone biosynthesis and plant hormone signal transduction are additionally involved in secondary metabolism (Bélanger et al., 2018). Studies have shown that differentially expressed genes (DEGs) can be detected by gene chip technology in an early stage of the cucumber’s in vitro gynogenesis, suggesting that phenylalanine metabolism and synthesis may play an important role in the early in vitro development of cucumber in vitro (Zhang et al., 2013).

Cucumber (Cucumis sativus L.) is the fourth most important vegetable worldwide (Lv et al., 2012), and accounted for more than 2.2 million hectares and a total production of approximately 83 million tons in 2017 (http://www.fao.org/faostat). Cucumber is one of the oldest vegetable crops and was domesticated in China approximately 2,000 years ago (Golabadi, Golkar & Eghtedary, 2012). To meet the needs of production, breeders are constantly looking for valuable germplasm characteristics, particularly resistance to disease and environmental stresses such as cold, drought and salt (Wang et al., 2018). Although collecting extensive germplasm resources for variety improvements can greatly facilitate these breeding efforts, cucumber is a cross-pollinated plant with clear heterosis. Therefore, almost all of the cucumber varieties currently in production are hybrids with slow breeding processes. Breeding parents take 6–8 years of artificial self-breeding in order to develop a stable inbred line. To accelerate the purification of cucumber parents and improve breeding efficiency, haploid gametophyte cultures can be used to induce embryoids, since homozygous DH plants can be obtained in 1–2 years. DH technology is a powerful tool used to speed up plant breeding. However, additional research is needed to improve our understanding of genes and their roles in embryogenesis and the haploid induction mechanism in embryo sacs (Chen et al., 2011).

This study focused on the cucumber’s highly efficient in vitro ovary culture technology system. We divided the process from acquiring embryogenic potential to plant regeneration into three stages: early embryo development, embryo morphogenesis (from pro-embryos to cotyledon embryos), and shoot formation. Via ovary culture, we evaluated the embryo morphology and transcriptomes across different developmental stages during cucumber embryogenesis. The metabolic and biological process of embryogenesis in the cucumber gynogenesis were discussed, and several key genes that regulated embryogenesis were identified.

Materials and Methods

Plant materials

We used the F1 cultivar ‘SG033’ (http://kunming056278.11467.com/) for ovary culture because this material has a high embryo generation rate even after many screenings. Cucumber plants from southern China were cultivated and stored in our lab. The plant growth protocol was slightly adapted from Li et al. (2014). Specifically, the plants were grown in 12-h photoperiod greenhouse at the Horticulture Institute of Guizhou Academy of Agricultural Sciences, Guizhou, China, with an average daily air temperature of 25/15 °C (day/night), a relative humidity of 85%, and a photosynthetic luminous flux density of 500 μmol·m−2·s−1. We collected ovaries from ‘SG033’ after the first anthesis flower appeared.

Culture medium composition

For MS medium, we used a basal medium that was supplemented with 0.06 mg·L−1 thidiazuron (TDZ, Solarbio, Beijing, China) and 3% (w/v) sucrose, and was solidified with 7% (w/v) agar. The medium pH was adjusted to 5.9 before autoclaving at 116 °C and 1.1 kg/cm2 for 30 min.

Ovary culture

Ovaries were harvested 6 h before anthesis, and were treated for 24 h in a refrigerator at 4 °C. We removed their fruit tumor, sterilized with 75% ethanol for 30–60 s, and then rinsed in sterile distilled water three times, soaked in a 0.5% sodium hypochlorite for 20 min, and then rinsed them in sterile distilled water three times. The ovaries were cut into small round 2-mm slices under sterile conditions and placed on 30 mL of solid media in cylindrical flasks. All materials were treated with a high temperature of 33 °C for 2 d, after which recovery growth was allowed until shoot formation occurred at 25 °C under a 16/8 h (light/dark) photoperiod with a 4,000-lux light intensity.

Observation

We observed embryo development according to the methods previously published by our team (Deng et al., 2020). Specifically, SG033 explants (varian fragments) were collected before and after culture, and resin sections were obtained to observe early, embryo development. Fifty explants were collected separately at 0 d and 2 d fixed in formalin-acetic acid-alcohol (five ml of formalin with 38% formaldehyde, 5 ml of glacial acetic acid, and 90 ml of 70% ethanol) (FAA) at room temperature for 48 h, dehydrated with a graded ethanol series, and then embedded in the spur resin. Nine-micrometer-thick sections were stained with safranine and fast green, and were observed using a biological microscope (Leica DM2500, Wetzlar, Germany). After 10 d of culture, we evaluated the embryonic morphology using a stereo microscope (Leica M165 C, Wetzlar, Germany).

RNA isolation

The materials were harvested separately at the following time points: T0 (the fresh unpollinated ovaries were cultured for 0 d), T1 (the ovules were cultured for 2 d), T2 (the embryos were cultured for 10 d), T3 (the embryos were cultured for 20 d), T4 (the embryos were cultured for 30 d), and T5 (the shoots after 60 d of culture). We extracted ovules at T0 and T1 time points by hand under a stereo microscope, directly selected embryos at the time points from T2 to T4, and selected the shoots from T5. We used 200 mg for each sample as the starting material for the RNAseq and qPCR experiments.

These materials were stored at −80 °C for subsequent RNA extraction, and the remaining materials were maintained in culture in order to observe plantlet regeneration. All further experimental procedures, such as culture medium replacement and sample collection, were performed under similar conditions to minimize possible circadian effects.

As described above, we collected samples from six ovary culture time points for RNA-seq analysis. For each sample, the ovules or embryo-like structures were immersed in liquid nitrogen in a mortar and were ground into powder (three biological replicates per sample). Total RNA was isolated using TRIzol (Invitrogen, Carlsbad, CA, USA), and DNase I (Fermentas, Waltham, MA, USA); digestion was then performed for 30 min at 25 °C to remove DNA, according to the manufacturer’s instructions. We determined the integrity and quality of the total RNA using a NanoDrop 1,000 spectrophotometer and formaldehyde-agarose gel electrophoresis. RNA was used only when the OD260:OD280 ratio was above 1.8.

The total RNA concentration for all of the samples was ≥80 ng/μl, the OD: 260/280 > 2, 260/230 > 2, and the rin value was close to 8.0 (see Fig. 1 for a agarose gel electrophoresis), which is 28 s/18 s greater than 1.0. This indicated that the total concentration and quality of the prepared RNA was suitable.

Total RNA agarose gel electrophoresis of cucumber ovary embryos.

Figure 1: Total RNA agarose gel electrophoresis of cucumber ovary embryos.

1, 2, 3, 4, 5, 6……18 = T0-1, T0-2, T0-3, T1-1, T1-2, T1-3……T5-3

Library construction and sequencing

We performed library construction and sequencing as previously described by Liu et al. (2014). Specifically, the library was constructed with 3 micrograms of total RNA and enriched with magnetic beads coated with oligo (dT). After adding fragment buffer, the mRNA was converted to short fragments (200–300 bp). Next, we used mRNA fragments as templates to synthesize cDNA. The double-stranded cDNA was purified with a QiaQuick PCR kit, washed with EB buffer solution, and repaired at the end, and a single adenine (A) nucleotide was added. Finally, the sequencing adapters were connected to the fragments. The fragments were purified by using agarose gel electrophoresis and amplified by PCR. We then transferred the library products to an Illumina HiSeq 2000 for sequencing analysis. It was important to use this for base calling to convert the original image data into sequences in order to obtain clean reads before further analysis.

Transcriptome analysis

We analyzed the transcriptomes using the Illumina Hiseq 2000 platform. The raw data was uploaded to NCBI (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE147282). Using the Illumina GA Pipeline (v1.6) software, we removed more than 10% of unknown base and low-quality base readings from the raw sequencing data. Using TopHat2 (http://tophat.cbcb.umd.edu/) software, we aligned the filtered reads with the Cucumber (chinese long v2 Genome). Reference genome and gene database information were publicly available at the website http://cucurbitgenomics.org/ (Huang et al., 2009; Li et al., 2011). We performed differential expression analysis using DESseq (version 1.18.0). Benjamini and Hochberg’s methods were used to control the false discovery rate and adjust the obtained P values. Genes with an adjusted P value of <0.05 and an absolute value of |log2fold change ≥ 1| according to DEseq were considered differentially expressed. Gene ontology (GO) annotation was implemented by Blast2GO software, and GO association was performed via a BLASTX search of the NCBI NR database. We conducted the GO enrichment analysis of DEGs using the BiNGO plugin.

Over-represented GO terms were recognized using a hypergeometric test following Benjamini–Hochberg false discovery rate (FDR) correction with a significance threshold of 0.05. We used the KOBAS (2.0) software to carry out the KEGG enrichment analysis of theDEGs (Mao et al., 2005). A p value ≤ 0.05 was used as the threshold to judge significantly enriched DEG pathways.

Verifying RNA-seq data using quantitative real-time PCR (qRT-PCR)

To test the reliability of RNA-seq, we selected nine up-regulated genes for qRT-PCR. The total RNA of the RNA-seq samples were treated with DNase I enzyme, and were then transformed into single-strand cDNA using the GoScript Reverse

Transcription System (Promega, Madison, WI, USA) according to the manufacturer’s protocol. According to its cDNA sequence, Primer 5.0 was used to design specific gene primers (Table S1). Actin was used as an internal control to normalize minor differences in the number of templates. qRT-PCR was performed with GoTaq qPCR Master Mix (Promega, Madison, WI, USA) in a Bio-Rad iQ1 real-time PCR system (Bio-Rad, Hercules, CA, USA), the reaction conditions were as follows: pre-denaturation at 95 °C for 10 min, followed by denaturation at 95 °C for 15 s, and renaturation at 55 °C for 1 min for a total of 40 cycles; and at the end of the reaction, the system was kept at 95 °C for 15 s, followed by 60 °C for 1 min, and finally at last at 95 °C for 15 s. There were three biological replicates of each gene at each time point, qPCR experimental samples were the same as the RNA sequencing samples. The expression ratio was calculated using the 2−ΔΔCt method.

Results

Morphological and cytological characterization of in vitro-induced embryogenesis in cucumber

We conducted a negative and positive contrast of the cucumber ovary culture. After 2 days of culture, some ovules were clearly enlarged, the embryo sac was also significantly enlarged (Figs. 2A and 2C), but some ovules did not respond (Figs. 2A and 2B). After 10 days of culture, the non-embryogenic ovaries were presumed dead. In our study, we used the cells that underwent embryogenesis for sequencing analysis.(Fig. 2D).

Negative and positive contrast figure of cucumber ovary culture.

Figure 2: Negative and positive contrast figure of cucumber ovary culture.

(A), (B), (C) Ouvle of culture 2 d. (A) L and (B) non-responsive ovule. (A) R and (C) Enlarged ovule. (D) Ovaries of culture 10 days. Ov, ovule; Vb, vascular bundle; Es, embryo sac; Pl, placenta; In, intergument; Ie, inner epidermis. Enlarged ovule was marked with a white arrow.

We observed the early embryogenesis stage, embryo morphogenesis stage, and the shoot formation stage in order to determine the cytological changes that occur during embryogenesis in cucumber ovary culture. Ovaries were collected 6 h before anthesis and were pretreated at 4 °C. The ovaries were then sterilized, sliced, inoculated into the medium, incubated for 2 d at 33 °C in the dark, and then transferred to 25 °C conditions until plants formed. Notably, our morphological and cytological observations showed that critical developmental transitions in the embryos occurred at 2 d (T1), 10 d (T2), 20 d (T3), 30 d (T4) and 60 d (T5) of culture. Fresh unpollinated ovaries before culture were referred to as T0 (Fig. 3). At T0, we selected the fresh unpollinated ovaries with clear ovules (Fig. 3A) and mature embryo sacs (Fig. 3G). At T1, heat shock stress induced the embryogenesis of cucumber. The most obvious indicators of embryo development were enlarged cells in the embryo sac (Figs. 3B and 3H). Enlarged microspores have shown a correlation with embryogenic potential acquisition during androgenesis induction in many crop species (Hoekstra et al., 1996; Maraschin, 2005). After culturing for 2, 4 and 6 d, the expanded cells formed two (Fig. 3I), four (Fig. 3J), and multiple cell structures (Fig. 3K), respectively. At T2 (10 d), the ovules had successfully protruded from the surrounding tissues and a pro-embryo had formed (Fig. 3C). At T3 (20 d), the total number of embryos had increased in one ovary slice (Fig. 3D). From T3 to T4 (30 d culture), the whole embryo had formed with globular, heart-shaped, torpedo-shaped, and cotyledon-shaped (Fig. 3L). At the same time, the color of the embryonic cells appeared green, indicating the beginning of plant development from the embryo (Fig. 1E). At T5 (60 d after culture), shoots had formed (Fig. 3F).

Embryogenic process of ovary culture in cucumber.

Figure 3: Embryogenic process of ovary culture in cucumber.

(A)–(F) and (1) Observation under stereoscopic microscope. (A) Selected fresh unpollinated ovaries (0 d, T0). (B) Ovule were treated with the high temperature of 33 °C for 2 d; (C) Embryos cultured for 10 d, note the ovule enlargement and embryo initiation. (D) Embryos cultured for 10 d, with an average of six embryoids in each slice of ovary; (E) Embryos cultured for 30 d, embryo morphogenesis and cotyledon-embryo formation. (F) Cultured for 60 d, the embryo differentiated into shoot. (G–K) Histological observations: the histological sections were stained with Delafield’s hematoxylin. (G) Multicellular embryo sac before culturing, (H) One of the cells expands in the embryo sac, after 2 d of culture, (I) Cell mitosis occurred at 4 d of culture. (J) Cells continue to be divided into four cells at 6 d of culture. (K) A cell clumps structure was visible at 8 d of culture. (L) Globular-embryo switch into cotyledon shape-embryo. Ov, ovule; Ow, ovary wall; Es, embryo sac; Ec, egg cell; Sy, synergid; PN, polar nucleus; Ig, integument; Eb, embryo; Dc, developing cells; G-e, globular embryo; H-se, heart-shaped embryo; T-se, torpedo-shaped embryo; C-se, cotyledon-shaped embryo; S, shoot. The solid white lines represent 2 mm, and the black lines represent 200 μm.

At T1, ovule enlargement and cell enlargement were clearly observed, suggesting that the transformation of the developmental pathway from the gametophyte to the sporophyte was induced by stress, and that the stage between T0 and T1 was the key to inducing embryogenesis. After the switch in cell development, embryogenic potential was acquired. Mitosis then occurred to continuously form pro-embryos, globular embryos, heart-shaped embryos, torpedo-shaped embryos, and cotyledon-shaped embryos, similar to those observed during zygote development, and this process continued until T4. After 60 d of culture, the cotyledon-shaped embryos developed into shoots. We divided the six culture time points into three stages: T0 to T1 (early embryogenesis stage), T1 to T4 (embryo morphogenesis stage), and T4 to T5 (shoot formation stage).

Sequencing and aligning to the reference genome

Based on our observation of the morphological and cytological characteristics, we selected six time points and used three biological replicates for RNAseq, and ultimately represented a total of 18 libraries. The transcriptome data were generated using Illumina II HiSeq 2000 sequencing. Overall, more than 101.3 billion paired-end reads passing filter were generated, and we used these to align to the reference genes, which provided the general information for the project (Table S2). Reads counts were generated from BAM alignment files using HTseq 0.6.1p2 software (http://www-huber.embl.de/users/anders/HTseq). Data normalization was implemented with Reads Per Kilobase of transcript, per Million mapped reads (RPKM).

To analyze the ovary culture biochemical process on the cucumber transcriptome, we conducted transcriptomic analysis of T0 vs. T1 (the control group), and the experimental groups T1 vs. T2, T2 vs. T3, T3 vs. T4, and T4 vs. T5 to reveal significantly expressed genes. We determined which genes were differentially expressed using DESeq (version 1.18.0), and used P-adjust ≤ 0.05 and the absolute value of |log2foldchange| ≥ 1 as the threshold to judge the significance of the differential gene expression. Between T0 vs. T1; T1 vs. T2; T2 vs. T3; T3 vs. T4; and T4 vs. T5 there were 3,468, 2,341, 1,117, 1,478 and 1,383 up-regulated genes, and 3,065, 1,906, 574, 2,223 and 1,837 down-regulated genes, respectively (Fig. 4).

The total number of up-regulated and down-regulated genes.

Figure 4: The total number of up-regulated and down-regulated genes.

Black colored bars represent up-regulated genes, and gray ones represent down-regulated genes.

Validation of differentially expressed genes

Validation of the Illumina sequencing data and the DEG expression patterns revealed by RNA-seq was performed in order to examine the expression patterns of nine DEGs: including seven genes involved in plant hormone signal transduction and two plant-pathogen interaction genes. qRT-PCR results showed that seven genes were involved in plant hormone signal transduction, Csa7M452370.1, Csa1M006300.1, Csa5M223020.1, Csa5M434550.1, Csa5M623800.1, Csa6M383530.1 and Csa3M389850.1, with Pearson’s r = 0.77, 0.57, 0.74, 0.79, 0.69, 0.81 and 0.86, respectively. The two genes involved in plant-pathogen interaction, Csa7M420160.1 and Csa1M006320.1, had Pearson’s r values of 0.62 and 0.81, respectively. The results showed a correlation between the relative expression levels revealed by RNA-seq and qRT-PCR. The fold changes in the qRT-PCR analysis were different from those in the RNA-seq analysis, which may have been due to the different sensitivities between the qRT-PCR analysis and the RNA-seq technique. However, the qRT-PCR analysis showed that the up-regulation and down-regulation trends of the differential gene expression were consistent with those of the RNA-seq analysis (Fig. 5).

Validation of differentially expressed genes and correlation between RNA-seq and qRT-PCR.

Figure 5: Validation of differentially expressed genes and correlation between RNA-seq and qRT-PCR.

qRT-PCR of 9 up-regulated genes involved in plant hormone signal transduction and plant-pathogen interaction genes in whole ovary culture were analyzed. (A) Csa7M452370.1, (B) Csa1M006300.1, (C) Csa5M223020.1, (D) Csa5M434550.1, (E) Csa5M623800.1, (F) Csa6M383530.1, (G) Csa3M389850.1, (H) Csa7M420160.1, (I) Csa1M006320.1. Close correlations (Person’s r = 0.77, 0.57, 0.74, 0.79, 0.69, 0.81, 0.86, 0.62, 0.81) were observed between relative expression levels in RNA-seq and qRT-PCR, validating the RNA-Seq methodology described here for quantitative analysis of the cucumber transcriptome. Three biological replicates of each gene at each time point, and data were showed as mean ± SE (n = 3).

Embryogenesis-related genes expression in cucumber ovary culture

Ovule enlargement was noticeable at the early stage of embryogenesis, following heat shock stress (Fig. 3B). Cytological observations also showed that one of the cells in the embryo sac expanded (Fig. 3H), suggesting acquisition of embryogenic potential.

We determined that the differentially expressed genes at this stage (T0 vs. T1) were comprised of 3,468 up-regulated genes and 3,065 down-regulated genes (Fig. 4). We performed a GO enrichment analysis of these genes, and the results showed that the proteins encoded by these genes were assigned to four biological processes, three molecular functions, and nine cellular components (Table S3). The majority of DEGs were involved in the ‘single-organism process’ (GO:0044699) and the ‘oxidation-reduction process’ (GO:0055114). They were mainly distributed across the terms ‘membrane’ (GO:0016020), ‘membrane part’ (GO:0044425), ‘intrinsic component of membrane’ (GO:0031224), and ‘integral component of membrane’ (GO:0016021), with good molecular oxidoreductase activity (GO:0016491). Among the DEGs, the expression of related genes involved in functional classification might play an important role in early embryogenesis.

To compare and summarize the results of this stage, we performed a pathway analysis to identify potential target genes. Based on the KEGG database, the pathway enrichment analysis of these genes revealed significant involvement in eight distinct pathways (Fig. 6A): pathways involving plant hormone signal transduction, phenylpropanoid biosynthesis, plant-pathogen interaction, glutathione metabolism, cysteine and methionine metabolism, drug metabolism-cytochrome P450, phenylalanine metabolism, and metabolism of xenobiotics by cytochrome P450.

KEGG pathway enrichment analysis based on the differentially expressed genes.

Figure 6: KEGG pathway enrichment analysis based on the differentially expressed genes.

(A) In early embryo development, (B) in the stage of embryo morphogenesis, (C) in shoot stage. The longitudinal coordinates are the enriched KEGG Pathways, and the x-axis monitors the rich factor; the size of the dots in the graph represents the number of differentially expressed genes annotated to the mentioned pathway, and the color represents the significant P value of the pathway. The most significant pathways are shown in the picture.

Plant hormone signal transduction genes

The DEGs were significantly enriched in the plant hormone signal transduction pathway (Ko04075), We detected a large number of genes related to cytokinin and ethylene signalling pathways, histidine phosphotransferase protein encoders (Csa2M373410.1, Csa6M067360.1 and Csa7M452370.1), response regulators (RR) (Csa1M006300.1, Csa3M822100.1, Csa4M436980.1, Csa5M223020.1, Csa5M434550.1, Csa5M603910.1, Csa5M623800.1 and Csa6M383530.1), EIN3 binding (F-box) proteins (EBFs) (CsaUNG009930), and ethylene-responsive protein transcription factors (ERFs) (Csa3M389850.1) that were significantly up-regulated after heat treatment (Table S4). Moreover, the partially up-regulated DEGs encoding auxin influx carriers, auxin-responsive proteins, SAUR family proteins, abscisic acid receptor proteins, serine/threonine-protein kinases, ABA-responsive element binding factors, and brassinosteroid signal-responsive protein kinases were significantly enriched in eight pathways, some of which involved tryptophan metabolism as well as zeatin, carotene, and brassinosteroid biosynthesis, suggesting that hormone production and signal transduction were closely related to stress-induced embryogenesis. The up-regulated genes encoding major enzymes and the receptor proteins involved in the plant hormone signal transduction pathway were expressed during the whole ovary culture stage. Most of the genes were dynamically expressed after their up-regulated expression during embryogenesis. Six genes (Csa6M125240.1, Csa6M147590.1, Csa6M007440.1, Csa4M454150.1, CsaUNG009930 and Csa3M389850.1) were up-regulated only at the embryo initiation stage, and then appeared down-regulated in the later stages (Fig. 7). These genes mainly encoded auxin synthetase and induced proteins, ABA protein kinase, and ethylene-responsive transcription factors, indicating that they might be closely related to embryogenesis.

The matrix graph of up-regulated expression genes in the plant hormone signal transduction pathway in six time points.

Figure 7: The matrix graph of up-regulated expression genes in the plant hormone signal transduction pathway in six time points.

The absolute values of the relative expression of up-regulated genes were indicated by the depth of color. The x-axis of the matrix graph showed the six time points. The y-axis of the matrix graph showed the name of the up-regulated genes. The up-regulated genes were arrayed according to the hormonal signal transduction flows, which were classified into six groups: auxin, cytokinin, abscisic acid, ethylene, brassinosteroid, salicylic acid.

Stress response proteins related to plant defense mechanisms

DEGs were significantly enriched in the plant-pathogen interaction pathway, which is related to plant defense mechanisms and complex physiological responses to heat shock-induced embryogenesis. In this regulatory pathway, the expression of Pto-interacting protein 1 (Pti1, Csa7M420160.1), Heat shock protein 90 (HSP90, Csa3M183950.1), Enhanced disease susceptibility 1 (EDS1, Csa1M006320.1), and other related genes were up-regulated (Table S4). HSP90 was up-regulated under heat stress, had been highly conserved during evolution and has plant defense functions.

Embryogenesis-related transcription factors

Transcription regulation through transcription factors is a key step in a plant’s stress response. In the embryogenesis stage, a total of 32 transcription factor families were up-regulated after heat shock stress (Table S4). Among these, the AP2/ERF, WRKY, basic helix loop helix (bHLH), MYB, GATA, heat stress factor (HSF) and NAM/ATAF/CUC (NAC) families had the most transcription factors. There were 26 transcription factors up-regulated in the AP2/ERF family, which plays an important role in ethylene response, brassinolide response, and biological and abiotic stress responses. Two BBM genes (Csa3M827310.1 and Csa3M827320.1) were up-regulated. Nineteen transcription factors were up-regulated in the WRKY family, which are involved in physiological embryo development, metabolism, environmental stress, and defense response to pathogens. The bHLH family plays an active role in regulating environmental stress. The MYB family controls different stages of cell division and promoting embryogenesis. After heat shock stress, four HSF genes (Csa2M356690.1, Csa6M517310.1, Csa1M629180.1 and Csa3M822450.1) were specifically expressed that could activate the HSP90 gene in response to heat stress. The up-regulated expression of these transcription factors might be involved in embryogenesis and the active regulation of biological and abiotic stress processes.

Response of genes involved in microtubule organization in cucumber ovary culture

In the embryo morphogenesis stage, the gametophyte development pathway switched to the sporophyte development pathway, and the cells enlarged initiated mitosis, and formed two-cell, four-cell, and multicellular pro-embryos after heat stress (Figs. 3H3K). Afterward, globular embryos, heart-shaped embryos, torpedo-shaped embryos, and cotyledon embryos were formed (Fig. 3L). The DEG analyses of T1 vs. T2. T2 vs. T3, and T3 vs. T4 revealed 480 continuously expressed genes. Additionally, the numbers of uniquely expressed genes between T1 and T2, T2 and T3, and T3 and T4 were 2,236; 371 and 1,558, respectively (Fig. 8).

Venn diagram showed differentially expressed genes in the stage of embryo morphogenesis.

Figure 8: Venn diagram showed differentially expressed genes in the stage of embryo morphogenesis.

The DEG sets (T1 vs T2, T2 vs T3 and T3 vs T4) described in Fig. 4 were analyzed using the Venn method. The numbers marked in the diagram indicated the number of genes significantly expressed among the three DEGs sets.

The GO categories of the 480 continuously expressed genes were significantly assigned to two biological processes, cell movement or subcellular component (GO:0006928) and microtubule-based movement (GO:0007018), which are mainly associated with microtubules (GO:0005874), supramolecular fibers (GO:0099512), and polymeric cytoskeletal fibers (GO:0099513). Subsequently, the proteins encoded by these genes were classified into 11 functional categories: protein complex binding (GO:0032403), microtubule binding (GO:0008017), tetrapyrrole binding (GO:0046906), tubulin binding (GO:0015631), microtubule motor activity (GO:0003777), motor activity (GO:0003774), macromolecular complex binding (GO:0044877), heme binding (GO:0020037), oxidoreductase activity (GO:0016491), iron ion binding (GO:0005506), and oxidoreductase activity (GO:0016491), acting on paired donors or the molecular incorporation of oxygen (GO:0016705) (Table S5). These findings indicated that during certain stages, the cytoskeleton structure was continuously built and maintained by the action of microtubule binding proteins and enzyme modifications. This provided the basis for positioning the various organelles and implementation functions, in order to ensure the orderly activities of various cells in time and space.

To determine the involvement of these DEGs in embryo development, we performed a pathway analysis to identify the potential target genes. Seventeen significant pathways were obtained by mapping to the KEGG database (Fig. 6B): protein processing in the endoplasmic reticulum, phenylpropanoid biosynthesis, photosynthesis-antenna proteins, limonene and pinene degradation, meiosis-yeast, the estrogen signaling pathway, cell cycle, cell cycle-yeast, diterpenoid biosynthesis, chloroalkane and chloroalkene degradation, carotenoid biosynthesis, progesterone-mediated oocyte maturation, glycerolipid metabolism, histidine metabolism, fatty acid degradation, antigen processing and presentation, and ascorbate and aldarate metabolism. Among these, the protein processing pathways in the endoplasmic reticulum and phenylpropanoid biosynthesis were the most significant.

Main oxidation-reduction and metabolic process-related genes expression in cucumber ovary culture

In the shoot formation stage, the embryos further differentiated into shoots (Fig. 3F). A total of 3,320 genes were differentially expressed between T4, and T5, 1,383 of which genes were up-regulated and 1,837 were down-regulated in T5 (Fig. 4). Next, we performed an enrichment analysis of the genes using GO, which revealed eighteen biological processes, fifteen molecular functions, and five cellular components (Table S6). The terms annotated under the biological process category mainly included oxidation-reduction processes (GO:0055114), the regulation of primary metabolic processes (GO:0080090), the regulation of cellular metabolic processes (GO:0031323), and the regulation of macromolecule metabolic processes (GO:0060255). The molecular function category mainly included the terms oxidoreductase activity (GO:0016491) and DNA binding (GO:0003677). The cellular component category included mainly extracellular region (GO:0005576). Ovary culture development was mainly focused on the oxidation-reduction and metabolism processes.

The DEGs were subsequently annotated using the KEGG database in order to identify pathway enrichments. A variety of pathways were found to be significantly enriched (Fig. 6C), including those involved in phenylpropanoid biosynthesis; plant hormone signal transduction; stilbenoid, diarylheptanoid, and gingerol biosynthesis; metabolism of xenobiotics by cytochrome P450; drug metabolism-cytochrome P450; flavonoid biosynthesis; phenylalanine metabolism; starch and sucrose metabolism; and zeatin biosynthesis. Most of the DEGs were enriched in phenylpropanoid biosynthesis, plant hormone signal transduction, phenylalanine metabolism, and starch and sucrose metabolism. Pathways involving phenylpropanoid biosynthesis and phenylalanine are important for metabolizing secondary metabolites in plants and are closely related to cell differentiation and pigmentation in plant development. In the plant hormone signal transduction pathway, the expression of an auxin response factor, the gibberellin receptor GID1, ethylene-responsive transcription factor 1 (ERF1) and cyclin D3 were up-regulated, showing that they promoted plant development, morphogenesis, cell enlargement, division and stem developmen. Additionally, the DEGs involved in starch and sugar metabolism were also significantly enriched and mainly encoded 3-βD glucosidase, T6Ps, GlgB, endoglucanase, and alpha-trehalase, which satisfied the requirements of embryo growth and development.

Shoot formation-related transcription factors

In the shoot formation stage, the 62 transcription factors that were up-regulated belonged to 23 transcription factor families, with an up-regulation ratio between 2.03 and 214.13. The most up-regulated gene was the TCACG motif-binding factor (TGA) transcription factor (CSA4M625020.1) in the basic leucine zipper (bZIP) family. It which can specifically bind to the activation sequence with TGACG as the core in order to regulate the transcription level of target genes, and it plays an important role in organ development and the defense response against biological and abiotic stress. AP2/EFR, bHLH, MYB, WRKY, ZIP and WUSCHEL-related homeobox (WOX) are the families with a large number of transcription factors, which may be an important molecular indicator of shoot formation.

Discussion

In this present study, the whole process from acquisition of embryogenic potential to plant regeneration was investigated. We divided the ovary culture into three stages. In the early stage of embryogenesis, the acquisition of embryogenic potential by stress (e.g., low or high temperature or hormone induction) was observed together with the repression of gametophytic development, leading to cell dedifferentiation. In the embryo morphogenesis stage, cell divisions gave rise to the formation of pro-embryos (cell clumps), globular embryos, heart-shaped embryos, torpedo-shaped embryos, and cotyledon embryos. In the shoot formation stage, cotyledon embryos developed into shoots. There were molecular events that regulated embryo development at different stages of development.

Embryogenesis-related genes are expressed in cucumber ovary culture

In the early stage of embryogenesis, the female gametophyte was regulated by a variety of hormones, including auxin, abscisic acid, cytokinin and ethylene, similar to somatic cell embryogenesis and transformation (Ikeda-Iwai et al., 2003). A large number of DEGs were enriched in the plant hormone signal transduction pathway. Phytohormones such as cytokinin and ethylene were significantly up-regulated in multiple biosynthesis and metabolic pathways, and up-regulated and down-regulated genes related to auxin and abscisic acid were both observed. Generally, high-level endogenous cytokinin and low-level endogenous ethylene are beneficial for the acquisition of embryogenic ability. Embryogenesis is the comprehensive performance of interactions across different hormones, and all kinds of endogenous hormones show dynamic changes. The up-regulated expression of the main enzymes and receptor proteins in the plant hormone signal transduction pathway could promote embryogenesis, indicating that the high expression levels of related genes may play important roles in the process of switching from the gametophyte to the sporophyte development pathway. In this pathway, we found six genes related to hormone regulation (Csa6M125240.1, Csa6M147590.1, Csa6M007440.1, Csa4M454150.1, CsaUNG009930 and Csa3M389850.1) that were highly expressed only at the embryo initiation stage, then were down-regulated in the later stages. Additionally, we determined the functional categories and regulatory pathways of the specific functions that lay an important theoretical foundation for understanding this mechanism.

In in vitro culture, a certain stress condition is needed to block the development direction of the original gametophyte, which needs to be turned towards the direction of the sporophyte, and continue symmetry splitting in order to eventually lead to embryogenesis (Fan, Armstrong & Keller, 1988). In cucumber ovary culture, heat shock pre-treatment, silver nitrate, genotype and hormone combination factors can independently and simultaneously play key roles in embryo and callus production (Golabadi et al., 2017). In this study, heat shock stress was an important factor in the transformation from gametophyte to sporophyte. Previous studies have shown that HSP genes are specifically expressed in the spore stress-induction process of many crop species (Zarsky et al., 2010; Binarova et al., 1997; Sabehat & Weiss, 1998). HSP genes play important roles in cell cycle control, genomic silencing, protein transduction and signal transduction (Zhang et al., 2013). We found that 26 up-regulated genes encoded proteins and transcription factors related to heat shock, and that HSP90 (Csa3M183950.1) was involved in the regulation of the hypersensitive response in the plant and pathogen interaction pathway. The expression of these genes might play a regulatory role in the embryogenesis process. The BBM gene is preferentially expressed in developing embryos and seeds (Boutilier et al., 2002), and is a key regulator of plant cell totipotency because BBM ectopic overexpression induces the formation of somatic embryos in Arabidopsis seedlings without exogenous growth regulators or stress treatments (Horstman et al., 2017). We found that two BBM genes (Csa3M827310.1 and Csa3M827320.1) were up-regulated at the beginning of embryogenesis and were subsequently down-regulated. Additionally, AGL15, a member of the MADS family, is also considered a regulatory protein that acts at the start of cell division (Perry & Fernandez, 1999). We found that the AGL gene (Csa3M258140.1) was up-regulated in our study, indicating that all three of these genes (BBM, HSP90 and AGL) play important roles in the induction of embryogenesis by cucumber ovary culture.

Gametophyte embryogenesis and zygotic embryogenesis are similar processes, and the only difference is that the former occurs in sex cells while the latter occurs in fertilized egg cells during the development and differentiation of zygotic embryos. Microtubules play an important role in cells between the polarized diffuse growth and the acquisition stages (Kimata et al., 2016). Electron microscope-based observations of carrot culture showed that the appearance of microtubules was accompanied by the formation of somatic embryos (Halperin & Jensen, 1967). The different patterns of microtubule organization in the cells of the mature embryo sac cells reflect their structural adaptations for their future function (Huang & Sheridan, 1994). Similar to somatic embryos, a large number of genes participated in microtubule movement and bind to tubulin after embryonic potential is acquired, and the different embryo shapes form through active cell metabolism and rapid cell division. We concluded that the cellular reprogramming and morphological changes in embryos were controlled by microtubule organization genes.

Cotyledon embryos are further cultured to form shoots, and at this stage, a large number of transcription factors were expressed. The first constitutes transcription factors that contain a B3 domain (Stone et al., 2001). These transcription factors encode regulatory proteins involved in the embryonic development process, and maintaining embryo development during late embryonic development (Boutilier et al., 2002; Riechmann & Meyerowitz, 1998). In our study, a total of six transcription factors (Csa2M359980.1, Csa2M292240.1, Csa5M608380.1, Csa6M489980.1, Csa6M489940.1, and CsaUNG031640) containing a B3 domain were up-regulated. The second was the AIL gene from the AP2/ERF family, which is involved in key developmental processes throughout the whole plant life cycle. Some genes in the AP2/ERF gene family are expressed in many other tissues and participate in many plant developmental processes, such as embryogenesis and shoot development (Boutilier et al., 2002; Riechmann & Meyerowitz, 1998). Additionally, we found that the expression levels of nine genes (Csa1M423190.1, Csa3M114480.1, Csa3M652380.1, Csa3M114470.1, Csa4M644740.1, Csa5M175970.1, Csa5M608380.1, Csa6M496390.1 and CsaUNG031640) from the AP2/ERF family were up-regulated. Furthermore, WUS homologous domain proteins not only alter the cell fate of the shoot and flower meristem, but they also promote the development of somatic embryos into seedlings. WUS proteins have no direct connection to the embryo characteristics, but they will alter the development state of the tissue by maintaining cells in an un-differentiated state in response to different stimuli (Mayer et al., 1998; Gallois et al., 2004). WUS’s ability to transform the vegetative growth phase to the embryonic stage, as well as to eventually form somatic cells, indicates that this homologous domain protein also plays an important role in embryo morphogenesis and meristem (Palovaara & Hakman, 2008). In the shoot apical meristem, WUS defines the organizing center and is critical for the induction of these cells (Dai et al., 2017; Meng et al., 2017). During this stage, two genes (Csa6M301060.1 and Csa6M505860.1) encoding WUS homologous domain proteins were up-regulated. We concluded that those homologous domain proteins might directly or indirectly regulate shoot formation.

Because the ovary culture embryogenesis takes place in the surrounding tissue, it was difficult to observe the embryo development process and the transcriptome analysis sample was difficult to select. Additionally, this study sample contained other tissues, especially at the T0/T1 time points, and isolated ovules yielded sequencing materials that may have impacted the transcriptome data.

Conclusion

Our study revealed the molecular events of embryogenesis in cucumber ovary culture. The early stage of embryogenesis was the turning point for embryos formation, and a large number of hormone-related genes (Csa6M125240.1, Csa6M147590.1, Csa6M007440.1, Csa4M454150.1, CsaUNG009930, Csa3M389850.1. and BBM, HSP90 and AGL) were expressed at the stage of early embryo development. In the stage of embryo morphogenesis stage, GO analysis showed that continuously expressed genes participated in microtubule-based movement, movement of the cell, or subcellular component processes, indicating functional protein complex binding, microtubule binding, tetrapyrrole binding, tubulin binding and other microtubule activity, which are involved in protein processing in the endoplasmic reticulum and phenylpropanoid biosynthesis. Shoot formation was regulated by six transcription factors that contained a B3 domain, a total of nine genes in the AP2/ERF family, and two genes that encoded WUS homologous domain proteins. This evaluation of molecular events in gynogenesis may be useful for understanding the molecular mechanism of cucumber ovarian culture.

Supplemental Information

Primers for Real-time quantitative PCR.

DOI: 10.7717/peerj.12145/supp-1

Statistics of RNA-Seq alignment.

DOI: 10.7717/peerj.12145/supp-2

GO classification of common expressed genes in early stage of embryogenesis.

DOI: 10.7717/peerj.12145/supp-3

Main differentially expressed genes in T0 vs T1.

DOI: 10.7717/peerj.12145/supp-4

Supplemental Information 5

GO classification of common expressed genes in morphogenesis stage.

DOI: 10.7717/peerj.12145/supp-5

GO classification of common expressed genes in shoot formed stage.

DOI: 10.7717/peerj.12145/supp-6

qRT-PCR raw data.

DOI: 10.7717/peerj.12145/supp-7
2 Citations   Views   Downloads