In-silico Discovery of Genetic Diversity in Cucumis sativus var. hardwickii: A Wild Relative of Cultivated Cucumber

1. Introduction

Cucumis sativus var. hardwickii is one of the wild relatives of cucumbers. It belongs to the Cucurbitaceae family, primarily found in forests and mountain slopes, at altitudes of 700- 2000 m, in China, Myanmar, NE India, Nepal and Thailand. The origin of the cucumber is known to be India, where the wild ancestor C. sativus var. hardwickii round, and bitter fruit grows [1]. It is also found in South India, particularly in the state of Kerala. C. sativus var. hardwickii (wild cucumber) is a self-compatible, insect-pollinated taxon with a stable genome (2n = 14) and high reproductive output (85%) as compared to cultivated cucumber (73.3%). It is primarily composed of two commercially available vegetables: cucumber (C. sativus L. var. sativus; 2n=2x=14) and melon (C. melo L; 2n=2x=24) [2,3]. From India to China, the plant seems to have spread eastward, and from China westward, to Asia Minor, North Africa, and Southern Europe. As such, it represents an extreme of C. sativus variation and has the potential to expand the genetic diversity available for commercial cucumber breeding [4].

It possesses a number of advantageous characteristics, including a sequential multiple fruiting habit, multiple lateral branching, and a high overall fruit weight per plant [5], dark green-coloured fruits, early maturation, and fruit production during the rainy season. It includes several critical genes for the evolution of cultivated cucumber but is under represented in the world's central gene banks [6,7]. Genetic markers (morphological, biochemical, and molecular) have been used to characterize the genetic variation found in cucumber [8,9]. Thus, by analyzing a more robust selection of this species, we can better understand the variation present in this closely related population used in plant improvement programs [10].

Agricultural lands have been compromised by salination, air pollution, and global warming [11]. Studies of single nucleotide polymorphisms (SNPs) and significant advances in the use of genomic methods have created new, quick mapping procedures for identifying features and determining their causes. To expand multidisciplinary science and it is necessary to maintain traditional practices in the field of molecular genetics. Because of their low costs and high performance, microarrays and next-generation sequencing are commonly used [12]. Recent microarrays have been applied to biological samples.

Although Northern or quantitative PCR has previously been used in gene studies and the ability to assess a more significant number of genes in a single run improved sensitivity and resolution. The primary goal of next-generation sequencing is to read hundreds of thousands or millions of DNA molecules.

Comparative genomics serves as an approach to predict shared genes among different species and to uncover genes that are unique and novel to specific species. Homologous genes, including orthologues and paralogues, play a crucial role in tracing the evolutionary patterns of organisms, especially in relation to speciation or duplication events. Additionally, inferring orthology provides the most accurate method for describing similarities and differences in a genome [13]. Orthologous genes typically exhibit a syntenic relationship with closely related species and can be traced back to a common ancestor during speciation events. The presence of similar orthologous sequences in different species indicates their similar biological function and contrary if sequence show greater divergence then other species may perform distinct function in particular species [14]. In this study, we employed a graph-based approach using the OrthoVenn 3 tool [15] to conduct synteny analysis, comparing the genomes of C. melo and C. sativus var. hardwickii.

The most common type of DNA sequence variation between alleles is single nucleotide polymorphisms (SNPs), which can be revealed by using advanced technologies for genome analysis. SNPs are single base-pair differences (substitutions or deletions) among chromosome sequences that result from point mutations and account for much of the genetic variation in the genome. The binding affinity between a regulatory protein and a regulatory DNA binding site may change if there is a SNP in that site [16,17]. SNP markers have become increasingly popular in plant molecular genetics due to their abundance in the genome and their ease of detection at ultra-high-throughput. Despite these benefits, GBS continues to face challenges due to a lack of quick, reliable data imputation algorithms and powerful computers capable of processing and storing large amounts of sequencing data [18]. Even in the face of such disadvantages or threats, genomic and other next-generation advances increase the pace and benefits of plant breeding. SNPs are the first step toward unravelling evolutionary and genetic relationships and elucidating agronomic characteristics and plant diseases [19]. Therefore, we performed the genome-wide in-silico mining of SNPs and genomic variants from C. sativus var. hardwickii and predicted their effects in detail.

2. Materials and Methods

2.1. Raw reads and Mapping

In-Silico investigation has been done with the raw SRA data available in the public domain at NCBI (National Center for Biotechnology Information). The whole genomic data was used for the SNPs and variants calling using the whole genome sequencing data. Raw reads from the whole genome sequence of C. sativus var hardwickii were downloaded from the NCBI, SRA with an accession number SRR2096458 using the Illumina platform 92,921,066 raw reads were generated. These reads were mapped to the C. sativus genome using BWA, and the version of the Cucumber genome was ASM407v2 (var. sativus, cv. Chinese Long, inbred line 9930) [20].

2.2. Genome-Wide SNPs Calling

BWA-mem (v0.7.17) was used with its default parameters. The VCFstats feature of RTG-tools version 3.12 was used to calculate mapping statistics such as total mapped, unmapped, and paired mapped reads and perform additional analysis on the alignment files. We used Unified Genotyper from GATK (v.3.8-0) to call SNPs and indels for variant analysis. We separated the raw variant calling data into SNPs and indels using Select Variants from GATK (v. 4.0.1.2). Then we applied hard-filtering to these raw SNP calls using Variant Filtration from GATK (v. 4.0.1.2) based on these threshold criteria: Mapping Quality Rank Sum of < -12.5, QUAL < 30, QD < 3.0, FS > 30.0, MQ < 30.0, DP < 200, and genotype-filter-expression DP < 15. We then choose bi-allelic variants using VCFtools (version 0.1.15) [21,22].

2.3. Variant Effect Analysis

The mapped data was used for Variant effect prediction using Variant Effect Predictor (VEP) software (version (API), 102) to verify the effect of variants [23]. The script for Ensembl Variant Effect Predictor is available at (http://www.ensembl.org/info/docs/tools/vep/script/index.html) [24].

2.4. Synteny Analysis

In order to demonstrate the synteny relationship between the orthologous genes derived from C. melo and C. sativus var. hardwickii, syntenic analysis maps were generated using the OrthoVenn3 toolkit program, utilizing the default parameters.

3. Results

3.1. Genome-Wide Identification of SNPs

We identified a total 19, 74,213 polymorphic SNPs in C. sativus var. hardwickii in relation to the reference genome of cucumber with MNPs (0.00), Insertions (1, 33,468), deletions (1,43,237), Indels (75) from all the processed reads and 1,62,673 genotypes were missing (Table 1). Furthermore, the heterozygous/homozygous ratio was also predicted for SNPs (0.20%), insertion (0.99%), deletion (0.85%) with total Het/Hom ratio (0.26%) and not detectable for MNPs and indels (Table 1).

3.2. Genome-Wide Variant Analysis Using Ensembl Variant Effect Predictor

VEP analysis was performed to determine the effect of variants on genes, mRNA transcripts, protein-coding, and regulatory sequences that may help evaluate the genetic diversity within the species.

3.2.1. General Statistics Data Description of Analyzed Variants

A total of 22,28,224 variants were processed out of 6404727 input reads along with 23,777 number of overlapped genes and similar for overlapped transcripts with zero novel variants and overlapped regulatory features (Table 2).

3.2.2. Variants Classes

Single-nucleotide variants (SNV) are the most commonly present variants genome wide and detected with great value in present study. Total identified variants in present study with computational analysis are segregated mainly into four classes of genomic variants which includes 0.01% Sequence alteration (275), 5.94% Insertion (1,32,531), 6.37% Deletion (1,42,040) and 87.66% SNV (19,53,378) respectively (Table 3).

3.2.3. All Consequences and Most Severe Consequences of Variants

Most severe consequences were detected in 21 classes with significant difference in number, highest for upstream gene variants i.e., 1,029,734 and lowest for Start retained variant i. e., 02 only. Some major classes are frame shift variants (1,114), in frame insertion (588), in frame deletion (712), missense variant (38,641), protein altering variant (29), synonymous variant (42,823), coding sequence variant (21), 5’ prime UTR variant (20,079), 3’ prime UTR variant (27,691), intron variant (3,04,952), downstream gene variant (3,66,513) and intergenic variant (3,87,161) were also detected for most severe consequences. Similarly, all consequences were detected in 21 classes with great value change, highest for upstream gene variants i.e., 1,547,272 and lowest for start retained variant (13) (Table 4).

3.2.4. Coding Consequences and Variant Location on Chromosomes

A total of 86435 variants were recorded for Coding consequences (a variant that makes changes in coding sequence) and classified into 12 classes with variable count range and highly detected 50.55% (43,701) in Synonymous variant (A sequence variant changes one base of terminator codon but does not impact on encoded amino acids, the terminator remains), followed by 44.86% (38,779) Missense variant (A variant that make changes in one base of the transcript first codon) and least detected for 0.01% (13) Start retained variant class (Table 5). A total of 22,28,224 variants were detected by Ensembl VEP, although variants were identified in all chromosomes, the number of detected variants on chromosomes varied considerably. These variants distributed on all seven chromosomes of C. sativus var. hardwickii with high number 18.03% (401,920) of total detected variants on the chromosome 3 followed by 15.84% (352,961) chromosome 5, 15.46% (344,544) chromosome 1, 14.72% (328,039) chromosome 6, 12.92% (287,924) chromosome 2, 12.13% (270,449) chromosome 4 while less 10.87% (242,387) on chromosome 7 (Table 6).

3.3. Synteny Analysis

A synteny relationship can provide insight into the evolutionary history of a genome. To further investigate the synteny of genes in rep-presentative crop species, comparative syntenic maps were constructed for C. melo and C. sativus var. hardwickii genes (Figure 1). Interestingly, our findings indicate a significant correlation between the majority of genes on chromosome 10 and 11 of C. melo, and chromosome 5 of C. sativus var. hardwickii. Additionally, a majority of genes on chromosome 12 of C. melo exhibit synteny with chromosome 6 of C. sativus var. hardwickii. Similarly, it is observed that the majority of genes on chromosome 13 of C. melo display synteny with chromosome 1 of C. sativus var. hardwickii. In contrast, the genes located on chromosome 2 of C. melo did not exhibit clear synteny with C. sativus var. hardwickii. Instead, they were observed to be distributed across chromosomes 1, 2, 3, 4, 5, and 6 of C. sativus var. hardwickii.

According to the predicted summarized clustering result, the species are represented in the columns, and each cell in the row represents an orthologous cluster group. The program's prediction indicated the presence of a total of four cluster groups, with two in each Cucumis species, and the absence of a cluster group in two cells (Figure 2). The number of protein sequences predicted to be part of the cluster group varied across each Cucumis species, ranging from 183 to 15,677 in cumulative count (bar plot stacked at right in Figure 2). The C. melo and C. sativus var. hardwickii revealed to have 15677 number of protein sequences involved in cluster group.

4. Discussion

Cucumber is a vegetable crop, the economic importance of which varies according to the region of the world. Cucumber is widely considered an excellent model for studying sex determination in monoecious plants [25]. Cucumber plant breeding can benefit from research into the phenotypic diversity found in cucumber genotypes. Superior genotypes for fruit yield and component production could be used to increase fruit yield, and they are also recommended for hybridization in cucumber cultivation [13]. Next-generation sequencing (NGS) offers a new perspective on the genetic code, revealing new genes and regulatory sequences, as well as a wealth of molecular markers. When it comes to traits that are influenced by a large number of factors, expression studies provide breeders with a wealth of knowledge about the molecular basis of these traits [26]. Advances in next-generation sequencing technology have facilitated the widespread use of SNPs in crops. Computational approaches dominate SNP discovery methods because increasing sequence awareness in public databases makes it easier to identify insightful SNPs computationally; however, complex genomes complicate SNP discovery significantly more, necessitating the use of alternative strategies for some crops. Plant genomes, in their entirety, have been shown to be repetitive [27].

Genome-wide sequencing analysis for SNPs and Variants prediction, along with their genome-wide distribution, may help to improve our understanding of genetic diversity, evolutionary studies and superior economically important traits of crop plants [28]. In this study, we conducted genome-wide SNPs identification and variants analysis from whole-genome sequencing data of C. sativus var. hardwickii. SNPs are very useful for studying genetic diversity because they can reveal the connections between different varieties. SNPs have been used for several years to measure diversity in specific genes or genomic regions, and the results are used to infer phylogenetic relationships between species [29]. The successful mapping of 92.9M raw reads with the reference genome revealed 19, 74,213 SNPs and 22, 28,224 variants genome-wide. SNPs allowed us to better understand genetic diversity and agricultural improvement, which will result in the enhancement of critical crop traits.

Variants calling in cucumber mapped sequencing data resulted in four variant classes with high number of SNV (87.66%) followed by deletion (6.37%), insertion (5.94%) and negligible for sequence alteration (0.01%). Here we recorded 21 variants classes for the most severe consequences with a total count of 2228224, highest for upstream gene variants (1,029,734), which located in 2 Kb 5’ of a gene followed by Intergenic variant found in intergenic regions (387161) and lowest for Start retained variant (02). Similarly, 21 variant classes for all consequences type with total number 3741865, highly detected for upstream gene variants, i.e., 1,547,27, followed by Downstream gene variant, i.e., 1352690 and lowest for Start retained variant, i.e., 13. Based on the variant location in the protein-coding and non-coding regions, including non-synonymous, missense, nonsense or Frameshift variants [30].

Impact of variants for coding consequences category revealed high change is synonymous (50.55%) where one base of stop codon changed. Still, terminator remains and does not affect the encoding amino acid peptide, missense (44.86%), which change at least one base in start codon and affect the encoding protein, followed by frame shift (1.32%) which resulting in disruption of translating reading frame as insertion or deletion is not multiple of three and stop gained (1.00%) (Liu et al., 2019). The discovered Variants were distributed across all 7 chromosomes. Overall, 18.03% (401920) of the total discovered variants (2228224) were on chromosome 3, followed by 15.84% (352961) on chromosome 5 and at least 10.87% (242387) on chromosome 7. The above genome-wide analyzed SNPs and Variants may improve our understanding of genetic evolution and diversity in Cucumis sativus and significantly increase our knowledge about the variants' association with genomic diversity and plants' response to environmental factors.

Comparative genetic mapping proves valuable in uncovering the syntenic relationships present among closely related plant species [31,32]. Gaining insights into syntenic relationships among species enables the exploration of genome evolution and dynamics [33,34]. Moreover, it facilitates the utilization of genetic information from related species in gene isolation and molecular tagging experiments [35-37]. Our study's identification of synteny aligns with the earlier discoveries of Huang et al. [38], Fukino et al. [39], and Li et al. (2011) [40]. They revealed that following divergence from C. melo, five out of cucumber's seven chromosomes emerged through the fusion of ten ancestral chromosomes. Comparative mapping has proven successful in establishing syntenic relationships among closely related plant species across various families. For instance, in the Solanaceae family (pepper, tomato, and potato), studies by [41-43] Park et al. [41], Wu et al. [42], and Wu and Tanksley [43] have utilized comparative mapping to define these relationships. Similarly, in Gramineae grasses, Byrne et al. [44], Rohner et al. [45], Zhu et al. [46], Khahani et al. [47] and Shariatipour et al. [48] have applied this approach. Comparative mapping has also been employed in Fabaceae legumes by Mudge et al. [33], Kalo et al. [50], Yan et al. [51], Cannon et al. [52], Phan et al. [53], and McClean et al. [54]. Additionally, Wang et al. [55] explored comparative mapping in Brassicaceae, while Dirlewanger et al. [56] and Jung et al. [57] employed it in the Rosaceae family.

5. Conclusion

Genome-wide sequencing data have become a powerful tool for exploring the genetic diversity and functional significance of SNP) and variants in various organisms. These molecular markers can reveal the evolutionary history, population structure, and adaptive responses of different species and cultivars. In this study, we performed a comprehensive analysis of SNPs and variants in the re-sequenced genome of C. sativus var. hardwickii, a wild relative of cultivated cucumber. We detected a total of 1974213 SNPs and 2228224 variants across the seven chromosomes of C. sativus var. hardwickii, with chromosome 3 having the highest density and chromosome 7 having the lowest density. These genomic variations reflect the high level of genetic diversity and divergence between C. sativus var. hardwickii and other cucurbit species.

Furthermore, we annotated the functional effects of these SNPs and variants on genes and regulatory elements involved in various biological processes related to plant growth, development, stress tolerance, and disease resistance. Our results provide valuable resources for future studies on the molecular basis and phenotypic consequences of genome-wide variation in C. sativus var. hardwickii and its potential applications for cucumber breeding.