The inbreed lines (i.e., JEF and KWS) are obtained from the leading varieties “Joeun Baekdadagi” and “Gyeoulsal-i Cheongjang” from Fomer Heungnong Seeds Co. After selecting the individual that best characteristics represent of each group in the F₂ populations, two inbreeds were raised through self-fertilization. The resulted breed line i.e., JEF is gynoecious, which is semi-white fruit skin color with white spine and KWS is monoecious which is uniform dark green skin color with black spine (Figure 1A). The Cucumis sativus breeding line plants were directly harvested in June 2018 in a field in Wanju, Jeollabuk-do, South Korea (35°90′ N, 127°15′ E), near the National Institute of Agricultural Sciences. Sampled fruits are shown in Figure 1A, and the complete work flow followed in this study is given in Supplementary Figure S1.

Total DNA was isolated from the samples individually according to sequencing protocols. The isolated DNA was sequenced using two different sequencing systems, PacBio Sequel (Pacific Biosciences, Menlo Park, CA, United States) and Illumina HiSeq 2,500 (Illumina, San Diego, CA, United States), which are widely used in long- and short-read sequencing. For Illumina sequencing, DNA was prepared using the TruSeq Nano DNA Library Prep Kit (Illumina). For PacBio sequencing, DNA was prepared using the SMRTbell Express Template Prep Kit (Pacific Biosciences; catalog no. 101–357–000). The experimental procedures were fully conducted by DNA Link (Seoul, Korea), an authorized service provider in South Korea. The Illumina paired-end sequences were initially subjected to filtering of technical artifacts (i.e., base-calling errors [Phred quality score ≤ Q20]) and adapters using Trimmomatic v. 0.32 (Bolger et al., 2014). These Illumina reads were used for error correction of PacBio reads in CLC Assembly Cell v. 5.1.1 (Qiagen, Hilden, Germany). The corrected PacBio reads were used to prepare the initial draft version of the cucumber genomes in FALCON-Unzip v. 0.30, a haplotype assembler program (Chin et al., 2016). Finally, using the RaGOO method (Alonge et al., 2019), the genome contigs were clustered and reordered according to their alignment with chromosomal units in the reference genome (‘Chinese Long’ 9,930). The assembled genomes were assessed for completeness using BUSCO v. 4.1.4 with the Viridiplantae_odb10 reference dataset (Seppey et al., 2019).

To prepare a clean reference genome, it was necessary to remove bacterial contamination and organelle genomes from the database. The complete GenBank database, which contains draft and reference genomes of bacteria and organelles (mitochondria and plastids), was used as the reference to determine which reads should be removed from the raw sequences. All reference mapping of preprocessed reads was conducted using Bowtie 2 v. 2.2.8 (Langmead and Salzberg, 2012). Details regarding reference paths and sizes are given in Supplementary Table S1, and mapping statistics are given in Supplementary Table S2.

All the Illumina-preprocessed sequences from the paired-end library were subjected to genome size estimation based on k-mers. The k-mer frequencies (k-mer size = 17) were obtained using Jellyfish v. 2.0 (Marçais and Kingsford, 2011), and the genome size was calculated from the following formulas: genome coverage depth = (k-mer coverage depth × average read length)/(average read length–k-mer size +1); genome size = total base number/genome coverage depth. Here, the k-mer coverage depth is the major peak of the k-mer distribution.

Repeat regions in the cucumber genomes were predicted using RepeatModeler (www.repeatmasker.org/RepeatModeler/) and classified into subclasses using the repbase v. 20.08 reference database (www.girinst.org/repbase/) (Bao et al., 2015). Finally, the repeats were masked in the genome using RepeatMasker v. 4.0.5 (www.repeatmasker.org) with RMBlastn v. 2.2.27+. The results are shown in Supplementary Figure S3.

The mRNA library from the collected samples was prepared according to the TruSeq Stranded mRNA Prep Kit protocol (Illumina). The isolated mRNA was sequenced using the Illumina sequencer (Supplementary Tables S4 and Supplementary Tables S5).

The genes from the cucumber draft genomes were predicted using an in-house gene prediction tool that includes three modules: an evidence-based gene modeler (EVM), an ab-initio gene modeler, and a consensus gene modeler. The Illumina-sequenced transcriptomes were mapped to the respective repeat-masked draft genomes using TopHat, and Trinity v2.5.1 method was used to assemble the transcripts and mark gene structural boundaries (Trapnell et al., 2012). The ab-initio gene modeler and EVM, which included Exonerate (Slater and Birney, 2005), Geneid and AUGUSTUS (Stanke et al., 2006), were trained with several genomes. The final gene and transcript models were optimized using a consensus gene modeler and annotated using Trinotate v. 3.0.1 (Bryant et al., 2017).

Total proteins from the two cucumber genomes were subjected to ortholog analysis to provide insight into the differences between cucumber proteins and those of other plants. In total, 14 genomes from Cucurbitaceae (including the two assembled in this study) were used in the ortholog analysis, with Brassicaceae as outliers (Figure 1D and Supplementary Table S3). The complete proteins of the selected genomes were also subjected to ortholog analysis using OrthoMCL (Li et al., 2003). The single-copy genes from the given genomes were subjected to phylogenetic tree reconstruction using BEAST (Bayesian Evolutionary Analysis Sampling Trees) to assess the evolutionary time and the degree of similarity among the given genomes (Suchard et al., 2018). Furthermore, to assess the gain and loss of genes in the given genomes, the proteins were analyzed using CAFE v. 3.1 (Han et al., 2013).

Initially, the sizes of the cucumber genomes were estimated to be 267.7 (JEF) and 276.4 MB (KWS) (Figure 1B) based on ∼50 GB of short-read sequences (Table 1A and Supplementary Table S4), but 230.8 MB (JEF) and 231.1 MB (KWS) based on the representative scaffolds assembled from ∼30 GB of error-corrected long-read sequences (Table 1A,B). The N50s of the assembled genomes were 30.5 MB (JEF) and 31.3 MB (KWS), and 40% of the assembled contigs were covered by repeats, in which the long terminal repeat (LTR) elements dominated, accounting for 36% of contigs (Supplementary Figure S3). In total, 25,968 genes were predicted from the JEF genome and 26,011 from KWS, with average sizes of 4,111 and 4,114 bases respectively, and BUSCO scores of 97.88 and 98.35% completeness respectively. (Table 1C). Finally, 66.54% of JEF genes and 65.96% of KWS genes had homologous sequences in GenBank, while 60.25% of JEF genes and 59.82% of KWS genes had gene ontology descriptions (Table 1D). The two genomes were scaffolded onto the reference “Chinese Long” 9,930 genome using the RaGOO method. Overall, these genome assemblies have ∼5 MB of additional bases compared with the reference and similar BUSCO completeness scores, indicating that they are of good quality. Additionally, an average of 99% of both DNA and RNA sequences were mapped to the reference assembly as an additional measure to ensure the quality of the new assemblies (Supplementary Figure S2). The ortholog analysis revealed genome-specific genes, as well as gain and loss of genes, in the selected cucumber genomes (Figure 1C and Supplementary Figure S4). In addition, the RNA samples were collected from five different developmental stages, revealing that both genomes contain genes expressed differentially in different organs or at different stages (Supplementary Figures S5 and Supplementary Figures S6).