Liu et al. Horticulture Research (2019)6:127 https://doi.org/10.1038/s41438-019-0209-4 Horticulture Research www.nature.com/hortres ARTICLE Open Access Selection footprints reflect genomic changes associated with breeding efforts in 56 cucumber inbred lines Bin Liu 1, Dailu Guan 2 , Xuling Zhai1, Sen Yang1, Shudan Xue1, Shuying Chen1, Jing Huang3, Huazhong Ren1 and Xingwang Liu1 Abstract Cucumber selective breeding over recent decades has dramatically increased productivity and quality, but the genomic characterizations and changes associated with this breeding history remain unclear. Here, we analyzed the genome resequencing data of 56 artificially selected cucumber inbred lines that exhibit various phenotypes to detect trait-associated sequence variations that reflect breeding improvement. We found that the 56 cucumber lines could be assigned to group 1 and group 2, and the two groups formed a distinctive genetic structure due to the breeding history involving hybridization and selection. Differentially selected regions were identified between group 1 and group 2, with implications for genomic-selection breeding signatures. These regions included known quantitative trait loci or genes that were reported to be associated with agronomic traits. Our results advance knowledge of cucumber genomics, and the 56 selected inbred lines could be good germplasm resources for breeding. 1234567890():,; 1234567890():,; 1234567890():,; 1234567890():,; Introduction environmentally sustainable will probably still be chal- Cucumber (Cucumis sativus L.; 2n = 2x = 14) is gen- lenging in the future. erally considered to be an economically important vege- Traditionally, the breeding of cucumber germplasm was table crop worldwide1, as well as a well-characterized accomplished simply by selecting lines with desirable model for studying fleshy fruit development2. Cultivated characteristics for propagation, while the new methodol- cucumber varieties have evolved from their wild pro- ogy has accompanied the development of high-volume genitors under natural and artificial selection3,4. The parallel genotyping and sequencing technologies, for domestication of cucumber began ~3000 years ago, and instance, genomic selection6. However, selective cucum- since then, different breeding strategies have been used to ber breeding employing genomic markers was initiated produce desired characteristics in order to meet the only in the 1980s5,7–9, and the implementation of genomic demand for human nutrition and health5. However, it was selection in cucumber breeding is still difficult. Although reported that cucumber had a low species diversity due to artificial selection has greatly increased cucumber pro- narrow bottlenecks3; therefore, the breeding of new ductivity and quality with expected improvements in cucumber varieties that require low inputs and are traits, such as gynoecy, disease resistance, uniform ripening and bitterness10–13, genetic gains in genomic selection are pending. Next-generation DNA sequencing Correspondence: Huazhong Ren (renhuazhong@cau.edu.cn) or technologies now allow cost-effective genome sequencing Xingwang Liu (Liuxw01@cau.edu.cn) 1 Beijing Key Laboratory of Growth and Developmental Regulation for at a population scale, which has led to the construction of Protected Vegetable Crops, College of Horticulture, China Agricultural variation maps for crop plants such as maize14, rice15, University, Beijing 100193, P. R. China 2 soybean16, sorghum17, apple18, watermelon19, pepper20, Centre for Research in Agricultural Genomics (CRAG), CSIC-IRTA-UAB-UB, Campus Universitat Autònoma de Barcelona, Bellaterra 08193, Spain and cucumber3. Meanwhile, genome-wide association Full list of author information is available at the end of the article. © The Author(s) 2019 Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/. Liu et al. Horticulture Research (2019)6:127 Page 2 of 9 studies (GWASs) are used to detect genetic variations that We also included 115 lines sequenced by Qi et al. that underlie many important and complex traits in plants21. were publicly available from the National Center for To date, some loci or genes have been identified by Biotechnology Information (NCBI, https://www.ncbi.nlm. genomic studies in cucumber with large natural popula- nih.gov/)3. Resequencing of the 56 selected lines gener- tions3,22,23. However, very few studies have identified loci ated a total of 1.6 billion paired-end reads, with an average or genomic regions in populations exposed to long-term depth of ~18× and coverage of 98.4% (Supplementary artificial selection. Table 1). By aligning paired-end reads against the refer- In recent decades, we focused on breeding elite ence genome of the inbred cucumber line 99301, a total of cucumbers with outstanding phenotypes by crossing ~8.4 million single nucleotide polymorphism (SNPs) were excellent Asian and European germplasm and perform- identified (Fig. 1b). Among these SNPs, 67187 were ing artificial selection. On the one hand, we aimed to unique in our breeding lines, indicating that the artificial combine desirable agronomic traits from Asian and selection lines were enriched with cucumber genomic Eurasian lines in the progenies. On the other hand, we variations, which are likely important for cucumber always selected the seeds from healthy plants and breeding. removed the disease-infected plants; therefore, most of To examine the genetic structure of the 56 inbred lines, the lines performed well in the field and were free of principal component analysis (PCA)29 and ADMIXTURE disease infection. As a result of long-term selection, we analysis30 were employed. It was clear that our cucumber obtained 56 elite cucumber lines that retained phenotypic lines formed groups different from the previously pub- diversity. In this case, we would like to know the genomic lished EA, EU, Indian domesticated (IND), Indian wild characterizations and changes associated with this (INW) and Xishuangbanna (XSBN) groups (Fig. 1c, d, breeding history and the potential of our breeding Supplementary Figs. 2 and 3). Along the first component, materials. explaining 31% of the variation, our 56 selected lines Here, we resequenced the 56 cucumber lines and found clustered away from the others (Fig. 1c), which was sup- that their genetic background was significantly different ported by ADMIXTURE analysis when K = 2 (Fig. 1d). from that previously reported for 115 cucumber lines3. Within the 56 lines, separation was found between the The data from our 56 lines revealed specific selected two groups along the second component in the PCA and regions related to breeding efforts, and these selected at K = 6 in the ADMIXTURE analysis (Fig. 1c, d). When regions were associated with the agronomic performance the K-value with the lowest cross-validation error (K = 7) of cucumber varieties and harbored many reportedly was used, the G1, G2, EA, EU, IND, INW, and XSBN important genes. groups were assigned to independent clusters (Fig. 1d). However, the EA cluster was divided into two different Results and Discussion genetic backgrounds, which was somewhat different from Determination of cucumber genetic structure with the pattern observed with published data3 (Fig. 1d). In artificially selected inbred lines addition, this kind of clustering was also reinforced in the The genetic data of 56 lines subjected to artificial constructed neighbor-joining (NJ) tree31 (Supplementary selection represent both genomic variations and their Fig. 2), indicating that our 56 cucumber lines had a dif- usefulness in cucumber breeding improvement. Accord- ferent genetic structure than the 115 lines3 obtained by ing to our breeding aims, the 56 lines were assigned to selective breeding and suggesting that G1 and G2 repre- group 1 (G1), which has a background similar to that of sent ideal populations for cucumber breeding. East Asian (EA) cucumber lines, and group 2 (G2), which has the Eurasian (EU) background (Supplementary Table Evidence of introgression and genetic diversity in two 1; Supplementary Fig. 1). These lines were selected to groups of cucumber obtain fruit traits with commercial value, such as normal Using the TreeMix32 model and the seven populations, vs abnormal fruit shoulder (Fig. 1 a1), dull vs glossy fruit we detected that G1 and G2 were phylogenetically similar to skin24 (Fig. 1 a2), uniform vs nonuniform fruit color12 the EA population, and the migration edge signaled intro- (Fig. 1 a3), green vs yellow-green fruit color (Fig. 1 a4), gression of EU into G1 and G2, especially to G2 (Fig. 2a), trichome and tubercule presence vs trichome and tuber- when admixture was modeled to include two migrations, cule absence25,26 (Fig. 1 a5–6), large vs small tubercules27 suggesting that EU cucumbers were the donors used to (Fig. 1 a7, a7-B a7-S), and hard vs soft spines28 (Fig. 1 a8, improve G1 and G2 by breeding. This result perfectly a7-H a8-S). To generate the 56 lines, at least 13 founder matched our selection schemes (Supplementary Table 1). lines from the EA or EU market were selected for cross- When admixture was modeled to include five migrations, fertilization. Next, backcrossing or selfing was performed we detected that migration edges to EA, G1, and G2 lines to obtain recombinant inbred lines or near-isogenic lines were from XSBN. Surprisingly, we did not find any with the target traits. introgression from INW lines (Supplementary Fig. 4). Liu et al. Horticulture Research (2019)6:127 Page 3 of 9 a a7-B a7-S a1 a2 a3 a4 a8-H a8-S a5 a6 a7 a8 b c Principal Component 2 (11.62%) -0.10 -0.05 0.00 0.05 0.10 EA XSBN EU G1 IND G2 INW Resequencing -0.10 -0.05 0.00 0.05 0.10 Principal Component 1 (31.42%) d K=2 K=6 K=7 G1 G2 EA EU IND INW XSBN Fig. 1 (See legend on next page.) Liu et al. Horticulture Research (2019)6:127 Page 4 of 9 (see figure on previous page) Fig. 1 Determination of cucumber genetic structure and genomic variations in 56 artificial selection inbred lines. a Representatives of 56 cucumber lines with different phenotypes: (a1) Normal vs abnormal fruit shoulder (the part between the two white dotted lines). (a2) Dull vs glossy skin. (a3) Uniform vs nonuniform fruit color. (a4) Green vs yellow-green skin. (a5) Trichome and tubercule presence vs trichome and tubercule absence in yellow-green-skinned fruit. (a6) Trichome and tubercule presence vs trichome and tubercule absence in green-skinned fruit. (a7) Large tubercules (a7-B) vs small tubercules (a7-S). The tubercules are marked by red dotted lines. (a8) Hard spines (a8-H) vs soft spines (a8-S). b Venn diagram depicting unique and shared SNPs between the 56 resequenced lines and 115 previously reported lines. c Principal component analysis (PCA) of the 171 cucumber lines. The PCA considered principal components 1 (PC1) and 2 (PC2), which explained 31.42% and 11.62% of the variance, respectively. d ADMIXTURE analysis (K = 2, 6, and 7). Each bar represents one individual, and the length of the colored bar represents the proportion of the cucumber genome inherited from each ancestral population. East Asian (EA); Eurasian (EU); Indian domesticated (IND); Indian wild (INW); Xishuangbanna (XSBN); group 1 (G1); group 2 (G2) These signals of genetic introgression were supported by especially for G2. In addition, it is worth mentioning that D-statistic analyses (ABBA-BABA tests)33 (Supplemen- migration edges to EA, G1, and G2 lines were from XSBN tary Table 2). These results suggested that G1 and G2 (Supplementary Fig. 4). This can be explained by the fact have experienced gene flow from other cucumber popu- that XSBN is a semiwild landrace with some unique traits lations, which coincided with our breeding history in that are very useful for cucumber breeding35,36. However, which G1 was improved based on the EA background and wild species of Cucumis are important as potential genetic G2 was selected to change fruit traits mainly by crossing sources for breeding37; INW is the origin population of with the EU group. We also examined the residuals of the cucumber and thus harbors higher genomic diversity than model fit to identify the relationships that were not cap- XSBN (Fig. 2c) and can be used for breeding improve- tured by the maximum likelihood trees; the XSBN and EA ment in the future. groups stood out (Supplementary Fig. 5). Increased linkage disequilibrium (LD) and reduced Selection regions in the two groups of cucumber nucleotide diversity are expected in artificial selection To assess the extent of genetic differentiation in G1 and lines34. To understand the specific LD block patterns in G2, we used a method to calculate pairwise differences in our 56 bred cucumbers, the LD decay (r2) with increasing allele frequency (FST) between G1 and G2, as well as the physical distance between SNPs was calculated by com- composite likelihood ratio (CLR)38 to identify genomic bining the 115 publicly available lines. The results indi- regions differentially selected within each group. The FST cated that both G1 and G2 cucumbers exhibited higher value highlighted 89 signature sites that were distributed LD than previously reported for the 115 lines (Fig. 2b), on chromosomes 1–6 and were significantly selected suggesting that large introgressions were selected and between G1 and G2 (q < 0.05; Fig. 2d; Supplementary fixed in G1 and G2. We then evaluated nucleotide Table 3). However, it is intriguing that the CLR statistics diversity at the whole-genome level within different indicated selection signals over the threshold (top 1%) groups. The average values of genome-wide nucleotide mainly on chromosome 1 in G1 (Fig. 3a; Supplementary diversity (π) for the G1, G2, EA, EU, IND, INW, and Table 3) and on chromosomes 2, 3, and 5 in G2 (Fig. 3b; XSBN lines were 1.56 × 10−3, 1.71 × 10−3, 1.37 × 10−3, Supplementary Table 3). These results suggested that 2.26 × 10−3, 4.25 × 10−3, 9.27 × 10−3, and 1.37 × 10−3, different target regions were under selection between the respectively (Fig. 2c). The π values of the INW and IND two groups during the breeding process and were con- groups were higher than those of the other groups, con- sistent with different selection pressures leading to dis- sistent with previously published data22, while those of tinct groups. The selected regions of G1 and G2 detected both G1 and G2 were significantly lower (t-test, P < 2 × by merging FST and CLR statistics contained 4074 and 10−16). The high LD and low π values confirmed that G1 3508 genes, respectively (Supplementary Table 4). Gene and G2 are highly artificially selected populations. Ontology (GO) analysis of these genes showed enrich- Resequencing the 56 artificially selected cucumber lines ment in binding processes, metabolic processes and cat- also provided opportunities to investigate selection foot- alytic activity (Supplementary Table 5). prints reflecting the complex genomic and genetic varia- Previous studies identified several domestication-related tions resulting from cucumber breeding. LD and genes or QTLs in cucumber, including those for disease nucleotide diversity analysis demonstrated that G1 and resistance39,40, flowering time41, fruit length42, sex deter- G2 were exposed to higher selection pressure, while CLR mination42, fruit bitterness43, and trichome develop- analysis revealed clear differentiation via selection ment28. These reported genes or QTLs overlapped with between the two groups and multiple targets of selection the putative selective sweep regions detected in our study, in G1 and G2. In addition, we found that the EU group indicating that a subset of domestication loci might have was most likely the founder cultivar for G1 and G2, undergone selection for continuous improvement of Liu et al. Horticulture Research (2019)6:127 Page 5 of 9 a b XSBN r2 G2 Drift parameter Distance(kb) c d Fig. 2 Evidence of introgression and genetic diversity in two groups of cucumber. a TreeMix analysis of the inferred relationships among seven cucumber groups. Migration edges are colored according to percent ancestry received from the donor population. Scale bar shows 10 times the average standard error of the estimated entries in the sample covariance matrix. b LD decay of seven groups of cucumber measured by r2 implemented in PopLDdecay software55. c Nucleotide diversity of seven groups of cucumber measured by the π value, taking into account a 50 kb window and 5 kb sliding steps. d Population divergence (FST) in the 56 artificial inbred lines important agronomic traits. For example, in the G1 was strongly selected, and two transcription factors, Bl population, we selected a soft trichome phenotype (Fig. 1 and Bt, on chromosome 5 that regulate the pathway of a8-S). As a result, Csa1G056960, which is involved in soft biosynthesis of cucurbitacin C were detected in the spine formation, was detected in the putative selected selection region (Supplementary Table 4). It is known that region on chromosome 128 (Supplementary Table 4). selection imposed on the Bt gene during domestication Another selected trait in G1 was flowering time; led to the derivation of nonbitter cucurbits from their Csa1G651710, which is responsible for flowering41, was bitter ancestors43; therefore, our results proved that the also detected in the putative selected region on chromo- target region for selection can be explained by the selec- some 1 (Supplementary Table 4). In G2, a bitterness trait tion schemes applied. Liu et al. Horticulture Research (2019)6:127 Page 6 of 9 Fig. 3 Genes under artificial selection. a Selection signal in cucumber group 1. Red line indicates the top 1%. b Selection signal in cucumber group 2. Red line indicates the top 1%. c Venn diagram depicting unique and shared SNPs between the G1, G2, EA, EU, and INW populations in the powdery mildew resistance QTL interval Chr1:6815000–6858200. d Nucleotide diversity of 5 cucumber groups in the powdery mildew resistance QTL interval Chr1:6815000–6858200, as measured by the π value The different target selection regions in G1 and G2 powdery mildew resistance major-effect QTL Pm1.139 might reflect breeders’ preferences. In other words, most was located in a selection region. Then, we compared the artificial selection results from the intentional breeding of SNPs in the Pm1.1 interval among the G1, G2, EA, EU, target traits. However, selective breeding can also be and INW groups (Fig. 3c) and found that the INW, EA, unintentional44. For example, in recent decades, we did and EU groups had 934, 39, and 27 unique SNPs, not perform a disease resistance test to screen for resis- respectively, while there were no unique SNPs in G1 or tant cucumber plants, but we unintentionally kept the G2. In other words, the SNPs that existed in G1 and G2 healthy plants and removed disease-infected plants. To were perhaps derived from the INW, EA and EU groups date, most plants in the G1 and G2 populations have and fixed by artificial selection. Interestingly, G2 had a performed well in the field and have not shown any dis- smaller number of SNPs as well as a lower nucleotide ease infection. In this study, we surprisingly found that the diversity (π = 2.14 × 10−3, t-test P < 0.01) than the EA, EU, Liu et al. Horticulture Research (2019)6:127 Page 7 of 9 and INW groups (Fig. 3c, d), indicating that this interval indel regions. Generated files in Binary Alignment Map was homozygous as a result of intensive selection and (BAM) format were used to call single nucleotide poly- could be used to identify the causal region for powdery morphisms (SNPs) via the HaplotypeCaller function of mildew resistance. Genome Analysis Toolkit (GATK, version 3.8) with default parameters48. Conclusions To obtain high-quality SNPs, a hard filtering step was In this study, by using next-generation resequencing performed by following the GATK Best Practices technology, we identified highly informative SNPs, which recommendations48. The genotype calls were then were used to study the genetic relations between improved by imputing and phasing the polymorphism cucumber germplasm accessions. Here, we reported 56 dataset with Beagle 4.1 software based on genotype like- artificially selected cucumber inbred lines that exhibited lihoods49. Finally, the effects of SNPs were predicted with various breeding preference phenotypes and could be SnpEff 4.3 software50. divided into two groups. In addition, the two groups exhibited distinctive genetic structure due to their Investigating population genetic structure breeding history involving hybridization and selection, We further thinned the SNP dataset by using the fol- and they were significantly different from 115 previously lowing rules: (1) core chromosomes (“--chr 1–7”); (2) reported cucumber lines. Phylogenetic inference with the minor allele frequency > 0.05; (3) missing values <10%; (4) TreeMix model revealed a breeding history in which G1 a Hardy-Weinberg equilibrium test (“--hwe 0.001”); and was improved based on the EA background and G2 was (5) pruning SNPs by “--indep-pairwise 50 10 0.1”. Two selected to change fruit traits mainly by crosses with the parameters (“--pca” and “--distance-matrix”) in PLINK EU group. D-statistic values were consistent with positive software51 were used to calculate pairwise matrixes for the selection in the cucumber populations, possibly related to construction of a principal component analysis (PCA) plot domestication or selection on traits of interest. We also and neighbor-joining tree, respectively. In addition, the reported the possibility that selective breeding can be ADMIXTURE (version 1.3.0) tool52,53 was employed to unintentional. Our results revealed the molecular basis of estimate individual ancestry by setting the K-value from 2 selection during cucumber development, and the to 7. Nucleotide diversity measured by the π value for 56 selected inbred lines could be good germplasm each population was estimated by using VCFtools54 with a resources for cucumber breeding programs. For example, 50 kb window and 5 kb sliding steps. In addition, the based on specific traits, each of the inbred lines could be decay of linkage disequilibrium (LD) with physical dis- used to develop populations, and the highly informative tance between SNPs was calculated and visualized by SNPs are of great relevance to cucumber breeders. using PopLDdecay software55. Furthermore, TreeMix- 1.1232 was used to model the genetic drift of genome-wide Materials and methods allele frequency data to infer population splitting and Sample collection and sequencing mixing. The standard errors of migration proportions The cucumber population of 56 artificially selected were calculated by using the “-se” option. Migration edges inbred lines was collected, crossed and inbred at the Ren (with the “-m” option) were added gradually from 0 to 5. Laboratory of Vegetable Science, China Agricultural University. The population was grown in an experimental Identifying selective sweeps field in the Changping District, Beijing. Phenotyping was The difference in allele frequency between the G1 and conducted in Changping (N 40°13′, E 116°12′). G2 groups was estimated with BayeScan software (version Genomic DNA was extracted as previously described15. 2.0)56, which decomposed locus-population FST coeffi- At least 5 μg of genomic DNA was prepared from a single cients into a population-specific component (β value) and plant of each accession for sequencing, and a library was a locus-specific component (α value) using logistic constructed with an insert size of ~300–350 bp for all of regression. A positive α value suggests diversifying selec- the lines, following the manufacturer’s instructions (Illu- tion, whereas a negative value indicates balancing or mina HiSeq 2500, USA). purifying selection. A q-value <0.05, analogous to the p value but used under multiple testing, was used to Variant calling, filtering, and annotation indicate significance. BayeScan software56 was run under The paired-end reads were aligned against the 9930 the default parameters. To detect regions under complete cucumber reference genome using the MEM algorithm selection within the population, the SweeD program was implemented in BWA software45,46 after trimming low- employed to calculate the composite likelihood ratio quality reads by using Trimmomatic47. Picard tools (CLR)38,57, which identifies regions with significant (https://broadinstitute.github.io/picard/) was employed deviations from the neutral site frequency spectrum (SFS). for the removal of PCR duplicates and realignment of We set the top 1% as a significance threshold. Finally, the Liu et al. Horticulture Research (2019)6:127 Page 8 of 9 overlapping regions detected by FST and the CLR were 14. Hufford, M.-B. et al. Comparative population genomics of maize domestication identified as having signatures of selection. Gene Ontol- and improvement. Nat. Genet. 44, 808 (2012). 15. Huang, X. et al. Genome-wide association studies of 14 agronomic traits in rice ogy (GO) annotations were assigned using AgriGO for landraces. Nat. Genet. 42, 961 (2010). cucumber (http://systemsbiology.cau.edu.cn). 16. Lam, H. M. et al. Resequencing of 31 wild and cultivated soybean genomes identifies patterns of genetic diversity and selection. Nat. Genet. 42, 1053 Acknowledgements (2010). This study was supported by the National Key Research and Development 17. Mace, E. S. et al. Whole-genome sequencing reveals untapped genetic Program of China (2018YFD1000800), the National Natural Science Foundation potential in Africa’s indigenous cereal crop sorghum. Nat. Commun. 4, 2320 of China (31830080), the Beijing Municipal Natural Science Foundation (2013). (6184043), and the Project of Beijing Agricultural Innovation Consortium 18. Duan, N. et al. Genome re-sequencing reveals the history of apple and (BAIC01). We are grateful to members of the Ren Laboratory for technical supports a two-stage model for fruit enlargement. Nat. Commun. 8, 249 assistance and many discussions. (2017). 19. Guo, S. et al. The draft genome of watermelon (Citrullus lanatus) and rese- Author details quencing of 20 diverse accessions. Nat. Genet. 45, 51 (2013). 1 Beijing Key Laboratory of Growth and Developmental Regulation for 20. Pereira-Dias, L., Vilanova, S., Fita, A., Prohens, J. & Rodríguez-Burruezo, A. Genetic Protected Vegetable Crops, College of Horticulture, China Agricultural diversity, population structure, and relationships in a collection of pepper University, Beijing 100193, P. R. China. 2Centre for Research in Agricultural (Capsicum spp.) landraces from the Spanish centre of diversity revealed by Genomics (CRAG), CSIC-IRTA-UAB-UB, Campus Universitat Autònoma de genotyping-by-sequencing (GBS). Hort. Res. 6, 54 (2019). Barcelona, Bellaterra 08193, Spain. 3Department of Agronomy, College of 21. Huang, X. et al. Genome-wide association study of flowering time and grain Agriculture, Purdue University, West Lafayette, IN 47907, USA yield traits in a worldwide collection of rice germplasm. Nat. Genet. 44, 32 (2012). 22. Wang, X. et al. The USDA cucumber (Cucumis sativus L.) collection: genetic Author’s contributions diversity, population structure, genome-wide association studies, and core X.L. and H.R. proposed the project. B.L. and D.G. performed the data analysis collection development. Hort. Res. 5, 64 (2018). and designed the manuscript structure, with the participation of J.H., H.R., X.L., 23. Zheng, Y. et al. Cucurbit Genomics Database (CuGenDB): a central portal for B.L., X.Z., S.Y., S.X., and S.C. planted the population, scored the phenotypes and comparative and functional genomics of cucurbit crops. Nucleic Acids Res. 47, extracted DNA. B.L. wrote the paper, with help from D.G., X.L., and H.R. D1128–D1136 (2018). 24. Yang, X. et al. High-resolution mapping of the dull fruit skin gene D in Conflict of interest cucumber (Cucumis sativus L.). Mol. Breed. 33, 15–22 (2014). The authors declare that they have no conflict of interest. 25. Wang, Y.-L. et al. Identification and mapping of Tril, a homeodomain-leucine zipper gene involved in multicellular trichome initiation in Cucumis sativus. Supplementary Information accompanies this paper at (https://doi.org/ Theor. Appl. Genet. 129, 305–316 (2016). 10.1038/s41438-019-0209-4). 26. Yang, X. et al. Tuberculate fruit gene Tu encodes a C2H2 zinc finger protein that is required for the warty fruit phenotype in cucumber (Cucumis sativus L.). Received: 2 June 2019 Revised: 1 September 2019 Accepted: 13 October Plant J. 78, 1034–1046 (2014). 2019 27. Yang, S. et al. A CsTu‐TS 1 regulatory module promotes fruit tubercule for- mation in cucumber. Plant Biotechnol. J. 17, 289–301 (2018). 28. Guo, C. et al. Identification and mapping of ts (tender spines), a gene involved in soft spine development in Cucumis sativus. Theor. Appl. Genet. 131, 1–12 (2018). References 29. Patterson, N., Price, A. L. & Reich, D. Population structure and eigenanalysis. 1. Huang, S. et al. The genome of the cucumber, Cucumis sativus L. Nat. Genet. PLoS Genet. 2, e190 (2006). 41, 1275 (2009). 30. Zhou, H., Alexander, D. & Lange, K. A quasi-Newton acceleration for 2. Liu, B. et al. Silencing of the gibberellin receptor homolog, CsGID1a, affects high-dimensional optimization algorithms. Stat. Comput. 21, 261–273 locule formation in cucumber (Cucumis sativus) fruit. New Phytol. 210, 551–563 (2011). (2016). 31. Saitou, N. & Nei, M. The neighbor-joining method: a new method for 3. Qi, J. et al. A genomic variation map provides insights into the genetic basis of reconstructing phylogenetic trees. Mol. Biol. Evol. 4, 406–425 (1987). cucumber domestication and diversity. Nat. Genet. 45, 1510 (2013). 32. Pickrell, J. K. & Pritchard, J. K. Inference of population splits and mixtures from 4. Che, G. & Zhang, X. Molecular basis of cucumber fruit domestication. Curr. genome-wide allele frequency data. PLoS Genet. 8, e1002967 (2012). Opin. Plant Biol. 47, 38–46 (2019). 33. Durand, E. Y., Patterson, N., Reich, D. & Slatkin, M. Testing for ancient 5. Staub, J. E., Robbins, M. D. & Wehner, T. C. in Vegetables I (eds Prohens, J. & admixture between closely related populations. Mol. Biol. Evol. 28, Nuez, F.) 241–282 (Springer, 2008). 2239–2252 (2011). 6. Voss-Fels, K. P., Cooper, M. & Hayes, B. J. Accelerating crop genetic gains with 34. Hübner, S. et al. Sunflower pan-genome analysis shows that hybridization genomic selection. Theor. Appl. Genet. 132, 669–686 (2019). altered gene content and disease resistance. Nat. Plants 1, 54–62 (2018). 7. Kozik, E. U. in Cucurbitaceae 2016, XIth Eucarpia Meeting on Cucurbit Genetics & 35. Bo, K., Ma, Z., Chen, J. & Weng, Y. Molecular mapping reveals structural rear- Breeding, July 24–28, 2016, Warsaw, Poland. (Organizing Committee of rangements and quantitative trait loci underlying traits with local adaptation Cucurbitaceae, 2016). in semi-wild Xishuangbanna cucumber (Cucumis sativus L. var. xishuangban- 8. Xu, Y. & Crouch, J. H. Marker-assisted selection in plant breeding: from pub- nanesis Qi et Yuan). Theor. Appl. Genet. 128, 25–39 (2015). lications to practice. Crop Sci. 48, 391–407 (2008). 36. Pan, Y. et al. QTL mapping of domestication and diversifying selection related 9. Fan, Z., Robbins, M. D. & Staub, J. E. Population development by phenotypic traits in round-fruited semi-wild Xishuangbanna cucumber (Cucumis sativus L. selection with subsequent marker-assisted selection for line extraction in var. xishuangbannanesis). Theor. Appl. Genet. 130, 1531–1548 (2017). cucumber (Cucumis sativus L.). Theor. Appl. Genet. 112, 843–855 (2006). 37. Liu, B. et al. A new grafted rootstock against root-knot nematode for 10. Kooistra, E. Femaleness in breeding glasshouse cucumbers. Euphytica 16, 1–17 cucumber, melon, and watermelon. Agron. Sustain Dev. 35, 251–259 (2015). (1967). 38. Pavlidis, P., Živković, D., Stamatakis, A. & Alachiotis, N. SweeD: likelihood-based 11. Sitterly, W. R. Breeding for disease resistance in cucurbits. Annu. Rev. Phyto- detection of selective sweeps in thousands of genomes. Mol. Biol. Evol. 30, pathol. 10, 471–490 (1972). 2224–2234 (2013). 12. Yang, X. et al. Fine mapping of the uniform immature fruit color gene u in 39. Xu, X. et al. Fine mapping of a dominantly inherited powdery mildew resis- cucumber (Cucumis sativus L.). Euphytica 196, 341–348 (2014). tance major-effect QTL, Pm1. 1, in cucumber identifies a 41.1 kb region 13. Andeweg, J. & De, BruynJ. Breeding of non-bitter cucumbers. Euphytica 8, containing two tandemly arrayed cysteine-rich receptor-like protein kinase 13–20 (1959). genes. Theor. Appl. Genet. 129, 507–516 (2016). Liu et al. Horticulture Research (2019)6:127 Page 9 of 9 40. Win, K. T., Vegas, J., Zhang, C., Song, K. & Lee, S. QTL mapping for downy 49. Browning, B. L. & Browning, S. R. Genotype imputation with millions of mildew resistance in cucumber via bulked segregant analysis using next- reference samples. Am. J. Hum. Genet. 98, 116–126 (2016). generation sequencing and conventional methods. Theor. Appl. Genet. 130, 50. Cingolani, P. et al. A program for annotating and predicting the effects of 199–211 (2017). single nucleotide polymorphisms, SnpEff: SNPs in the genome of Drosophila 41. Lu, H. et al. QTL-seq identifies an early flowering QTL located near Flowering melanogaster strain w1118; iso-2; iso-3. Fly 6, 80–92 (2012). Locus T in cucumber. Theor. Appl. Genet. 127, 1491–1499 (2014). 51. Purcell, S. et al. PLINK: a tool set for whole-genome association and 42. Tan, J. et al. A novel allele of monoecious (m) locus is responsible for elon- population-based linkage analyses. Am. J. Hum. Genet. 81, 559–575 gated fruit shape and perfect flowers in cucumber (Cucumis sativus L.). Theor. (2007). Appl. Genet. 128, 2483–2493 (2015). 52. Alexander, D. H. & Lange, K. Enhancements to the ADMIXTURE algorithm for 43. Shang, Y. et al. Biosynthesis, regulation, and domestication of bitterness in individual ancestry estimation. BMC Bioinform. 12, 246 (2011). cucumber. Science 346, 1084–1088 (2014). 53. Alexander, D. H., Novembre, J. & Lange, K. Fast model-based estimation of 44. Purugganan, M. D. & Fuller, D. Q. The nature of selection during plant ancestry in unrelated individuals. Genome res. 9, 1655–1664 (2009). domestication. Nature 457, 843 (2009). 54. Danecek, P. et al. The variant call format and VCFtools. Bioinformatics 27, 45. Li, H. Aligning sequence reads, clone sequences and assembly contigs with 2156–2158 (2011). BWA-MEM. arXiv Prepr. arXiv 1303, 3997 (2013). 55. Zhang, C., Dong, S.-S., Xu, J.-Y., He, W.-M. & Yang, T.-L. PopLDdecay: a fast and 46. Li, H. & Durbin, R. Fast and accurate short read alignment with effective tool for linkage disequilibrium decay analysis based on variant call Burrows–Wheeler transform. Bioinformatics 25, 1754–1760 (2009). format files. Bioinformatics 10, 1786–1788 (2018). 47. Bolger, A. M., Lohse, M. & Usadel, B. Trimmomatic: a flexible trimmer for 56. Foll, M. & Gaggiotti, O. E. A genome scan method to identify selected loci Illumina sequence data. Bioinformatics 30, 2114–2120 (2014). appropriate for both dominant and codominant markers: a Bayesian per- 48. McKenna, A. et al. The Genome Analysis Toolkit: a MapReduce framework for spective. Genetics 2, 977–993 (2008). analyzing next-generation DNA sequencing data. Genome Res. 9, 1297–1303 57. Nielsen, R. et al. Genomic scans for selective sweeps using SNP data. Genome (2010). Res. 15, 1566–1575 (2005).