A RAD-Based Genetic Map for Anchoring Scaffold Sequences and Identifying QTLs in Bitter Gourd (Momordica charantia)

Genetic mapping is a basic tool necessary for anchoring assembled scaffold sequences and for identifying QTLs controlling important traits. Though bitter gourd (Momordica charantia) is both consumed and used as a medicinal, research on its genomics and genetic mapping is severely limited. Here, we report the construction of a restriction site associated DNA (RAD)-based genetic map for bitter gourd using an F₂ mapping population comprising 423 individuals derived from two cultivated inbred lines, the gynoecious line ‘K44’ and the monoecious line ‘Dali-11.’ This map comprised 1,009 SNP markers and spanned a total genetic distance of 2,203.95 cM across the 11 linkage groups. It anchored a total of 113 assembled scaffolds that covered about 251.32 Mb (85.48%) of the 294.01 Mb assembled genome. In addition, three horticulturally important traits including sex expression, fruit epidermal structure, and immature fruit color were evaluated using a combination of qualitative and quantitative data. As a result, we identified three QTL/gene loci responsible for these traits in three environments. The QTL/gene gy/fffn/ffn, controlling sex expression involved in gynoecy, first female flower node, and female flower number was detected in the reported region. Particularly, two QTLs/genes, Fwa/Wr and w, were found to be responsible for fruit epidermal structure and white immature fruit color, respectively. This RAD-based genetic map promotes the assembly of the bitter gourd genome and the identified genetic loci will accelerate the cloning of relevant genes in the future.

Introduction

Bitter gourd (Momordica charantia; 2n = 2x = 22) is a tropical and subtropical vine in the family Cucurbitaceae, which is widely cultivated in Asia, Africa, and the Caribbean for its edible and medicinal fruit (Grover and Yadav, 2004; Marr et al., 2004; Van Wyk, 2015). Bitter gourd fruit provides a good source of phytonutrients like carbohydrates, minerals, and vitamins (Behera et al., 2010). The pharmacological composition and properties of bitter gourd have also been widely investigated (Tan et al., 2016). Although bitter gourd is consumed and cultivated, it remains an underutilized crop and scant efforts have been made for molecular improvement.

The construction of a genetic map is a common approach to detect QTLs and conduct gene mapping. The first genetic map of bitter gourd was constructed based on amplified fragment length polymorphism (AFLP) markers (Kole et al., 2012). A year later, the second genetic map was constructed, and it consisted of combined markers including simple sequence repeats (SSR), AFLP, and sequence-related amplified polymorphism (SRAP) markers (Wang and Xiang, 2013). Recently, two RAD-based genetic maps were constructed using two F₂ populations derived from the same parents (Matsumura et al., 2014; Urasaki et al., 2017), allowing for preliminary attempts at next-generation sequencing (NGS). Urasaki et al. (2017) reported a bitter gourd reference (OHB3-1 reference) and anchored 255 scaffolds to the genetic map. RAD markers are short fragments of DNA adjacent to each instance of a particular restriction enzyme recognition site (Baird et al., 2008). Combined with NGS, RAD-seq provides an inexpensive platform that allows high-density SNP discovery and genotyping in large populations (Davey et al., 2011).

In order to conduct marker-assisted selection (MAS) in bitter gourd breeding, it is necessary to target genetic loci underlying important horticultural traits. Similar to cucumber and melon, bitter gourd shows different sexual morphs. Gynoecy has potential applications in heterosis breeding and hybrid seed production. The gene that controls gynoecy (gy) was reported as a single recessive gene (Ram et al., 2006; Behera et al., 2009; Matsumura et al., 2014). An identified SNP marker linked to this gene has been mapped at a distance of 5.46 cM (Matsumura et al., 2014). In addition, bitter gourd fruits show extensive variation in epidermal color and structure, which directly effects exterior quality. Immature bitter gourd fruits have a broad color spectrum, ranging from white to dark green (Cheng et al., 2012). Several studies have shown that the green fruit color is monogenically dominant over the white fruit color (Srivastava and Nath, 1972; Hu et al., 2002; Dalamu et al., 2012). However, a more recent study indicated that this trait is controlled by quantitative genes (Huang and Hsieh, 2017). Bitter gourd fruit has either a smooth or a distinct warty exterior, and has broken or continuous ridges. The warty exterior is monogenically dominant over the smooth (Kole et al., 2012); however, research devoted to identifying genetic loci responsible for these important traits is limited.

As part of our ongoing efforts to sequence the bitter gourd genome, an efficient genetic map is instrumental for anchoring and orienting the assembled scaffolds. Here, we report a RAD-based linkage map of bitter gourd using an F₂ mapping population with 423 individuals derived from two cultivated inbred lines, the gynoecious line ‘K44’ and the monoecious line ‘Dali-11.’ Another major goal of the study was to map genetic loci underlying horticulturally important traits such as sex expression, fruit epidermal structure, and immature fruit color. To reach this objective, we carried out an elaborate investigation of phenotypic characters in three environments. The results of this study will be beneficial to research attempting to develop linked molecular markers in MAS breeding and to clone these genes in the future.

Materials and Methods

Mapping Population and Phenotype Evaluation

Two cultivated inbred lines of M. charantia were used as parents to produce the intra-specific F₂ mapping population. The female parent, ‘K44,’ is a gynoecious line, with white fruit color, dense round warty epidermis, and broken ridges; the male parent, ‘Dali-11,’ is a monoecious line, with green fruit color, non-warty epidermis, and continuous ridges (Figure 1). A total of 423 F₂ individuals were obtained from crosses of the two parental lines. The whole population was divided into three subsets: 205, 189, and 29 individuals, which were then grown in three different environments: Haikou (N 20.05°, E 110.20°) in spring 2014 (HKS), Guangzhou (N 23.13°, E 113.26°) in spring 2014 (GZS), and Guangzhou in autumn 2014 (GZA), respectively. In the present study, the horticultural traits of bitter gourd investigated included sex expression (gynoecy, gy; first female flower node, fffn; and female flower number, ffn), fruit epidermal structure (fruit wart, fwa; and width of ridge, wr), and immature fruit color (white fruit color, w; and two color measurements: lightness variable, L; and hue angle, H°, measured with a Minolta CR 300 Chroma Meter). In addition to gy, which was investigated in Haikou in spring 2014 (HKS) and in Guangzhou in autumn 2014 (GZA), and ffn, which was investigated in Haikou in spring 2014, the other traits were investigated in Haikou in spring 2014 and in Guangzhou in spring 2014 (GZS). The descriptions of each trait evaluated are shown in Table 1.

FIGURE 1

FIGURE 1. Phenotypic characterization of bitter gourd fruits. (A,C) Male parent inbred line ‘Dali-11;’ (B,D) female parent inbred line ‘K44;’ scale bar represents 1 cm.

TABLE 1

TABLE 1. Traits evaluated in the present study.

RAD Sequencing and SNP Identification

Following the protocol of a barcoding system for sample multiplexing, genomic DNA from ‘K44’ and 423 F₂ individuals were used to build RAD libraries (Baird et al., 2008). Briefly, each sample was digested with EcoRI (restriction enzyme cut site is 5′ G^∧AATTC 3′). This P1 adapter was ligated to the fragments. Then, the ligation products were pooled and randomly sheared with a Bioruptor (Diagenode, Belgium), and DNA fragments with a length of 350–550 bp were gel purified. The P2 adapter, a divergent “Y” adapter, was ligated to the chosen DNA fragments, and the ligation products were purified. The purified fragments were used for PCR amplification, products were gel purified again, and then DNA fragments of 350–550 bp were isolated. The prepared RAD libraries were sequenced on an Illumina Hiseq2000 platform at BGI-Shenzhen (Shenzhen, Guangdong, China) with both paired-end and single-end reads.

After performing standard quality control methods, including removal of the adapters, high missing reads (N > 5%), and low-quality reads, a clean data set was obtained and evaluated using the Q20 assessment. All of the clean reads were aligned against the Dali-11 reference genome (scaffold version, unpublished) using SOAP2 (Li et al., 2009). With the deep sequencing of the two parents (‘K44’ and ‘Dali-11’), more stringent quality control criteria, including sequencing depth of the site ≥10×, quality score of consensus genotype >20, rank sum test p-value < 1e-5, and average copy number of nearby region <1.1, were applied to identify the genotype of the two parents, finally generating their separated sites. Then, individuals from the F₂ population were genotyped at the separated sites of the two parents. To ensure the accuracy of genotyping, low coverage sequencing (<7× per site per individual) and low quality (quality score of consensus genotype Q < 20 per site per individual) sites were treated as missing information. After assigning accurate genotypes to F₂ individuals at the separated sites of the two parents, the genotypes of the individuals were recorded as ‘a’ and ‘b’ when it was the same as ‘K44’ and ‘Dali-11,’ respectively. Heterozygotes were recorded as ‘h’ and missing genotypes were recorded as ‘-.’ The final joined set of SNP genotyping was prepared to construct the genetic map.

Genetic Map Construction and Scaffold Anchoring

After filtering missing and heterozygous SNPs between two parental lines, and those SNPs with serious segregation distortion (p < 0.001, χ² > 13.82), the remaining SNPs were retained for genetic linkage analyses using JoinMap 4.0 (Van Ooijen, 2006). Then, an intercrossed filtering strategy was conducted between individuals and markers based on their missing rate. The logarithm of odds (LOD) threshold was set from three to 15. Finally, high-quality SNPs were split into different linkage groups that were estimated using the Kosambi mapping function. This procedure builds a map by adding loci one by one, starting from the most informative pair of loci, and continuing until all loci have been handled once. The first map node (Map 1), resulting from the first round, is perfect when measured by goodness of fit. Therefore, to further ensure the confidence of linkage groups, the resultant Map1 groups were adopted in this study. Scaffolds were assigned to linkage groups according to their linkage positions and the scaffold could be oriented when more than one marker existed.

QTL Detection

All traits investigated were used for QTL analysis; phenotypic data for each trait in different environments were analyzed separately. QTL analysis was performed using MapQTL 6 (Ooijen, 2009) software with the powerful MQM mapping (=composite interval mapping) method. Threshold values were calculated using 1,000 permutations and QTLs were designated when genome-wide LOD scores showed the presence of a significant peak at a level of p < 0.05. The QTL supported interval was identified by using a 1 LOD unit drop from the peak marker (if no marker exists, the nearest marker was selected). When the same QTL was detected in different environments, the overlapped region was considered the confidence interval for this QTL, and both the genetic interval and physical interval were delimited. Similarly, QTLs were merged if they did not overlap, but were close. The proportion of phenotypic variation explained by each QTL was estimated using the R² (%) value of the nearest markers.

Results

SNP Development From RAD Sequencing

We generated 193.72 Gb of clean data in total (∼22.28 billion reads) for 423 F₂ individuals and the female parent ‘K44’ (Supplementary Table S1) after using EcoRI to digest genomic DNA. The average genome coverage was approximately 18.00%, with an average GC content of 36.92%. The average Q20 of these samples was 98.60% (minimum 97.93%), indicating the high quality of the RAD data (Supplementary Table S1). Of these 22.28 billion reads, 93.00% were mapped to the reference genome and used to call SNPs. The number of called SNPs for each sample is presented in Supplementary Table S2. The resulting SNPs were filtered and formed a joined set containing 11,016 SNPs.

Genetic Map for Anchoring Scaffolds

After removing distorted loci and using an intercrossed filtering strategy between individuals and markers, 249 F₂ individuals (missing rate <35%) genotyped with 1,512 SNP markers (missing rate <50%) were included in downstream analyses. These SNP loci were assigned into 11 group nodes using JoinMap 4.0 at a LOD threshold of eight to 11, which correspond to the 11 chromosomes of bitter gourd. Then, we performed the map calculations, and the group node had a mapping node and two or three map nodes as child and grandchild nodes, respectively. Eleven Map 1 groups (first map node) containing 1,009 markers (Supplementary Table S3) were selected to represent 11 genetic linkage groups (MC01–MC11; Figure 2). This map spanned a total of 2,203.95 cM with a mean marker interval of 2.18 cM. Subsequently, a total of 113 scaffolds comprising 251.32 Mb (85.48% of the 294.01 Mb assembled genome) were anchored to the 11 linkage maps and 80 (∼70.80%) of the anchored scaffolds were oriented (Table 2). Based on this anchored assembly, the recombination rate of LGs ranged from 6.91 cM/Mb (MC04) to 11.04 cM/Mb (MC10), with an overall value of 8.77 cM/Mb.

FIGURE 2

FIGURE 2. The RAD-based genetic map of bitter gourd for ‘K44’ × ‘Dali-11.’

TABLE 2

TABLE 2. Anchoring of assembled scaffolds to the bitter gourd genetic map.

Inheritance of Traits Investigated

In the present study, three horticulturally important traits involved in sex expression (gynoecy, gy; first female flower node, fffn; and female flower number, ffn), fruit epidermal structure (fruit wart, fwa; and width of ridge, wr), and immature fruit color (white fruit color, w; lightness variable, L; and hue angle, H°) were investigated. The three characters of sex expression showed positive or negative correlations with each other in HKS, but gy and fffn investigated in GZS had no correlation with each other. In particular, the two fruit epidermal structure characters are highly correlated with each other, as are the three immature fruit color characters (Supplementary Tables S4, S5). The inheritance of these characters was also evaluated. Our data showed that gy, fwa, and w were qualitative characters. Interestingly, the F₂ plants exhibited monogenic segregation (3:1) for gy in GZA (χ² = 0.01), whereas two-genic segregation (15:1) was identified in HKS (χ² = 0.13). In addition, the F₂ plants showed monogenic segregation (3:1) for fwa in both HKS (χ² = 0.65) and GZS (χ² = 0.00) and monogenic segregation (3:1) for w in GZS (χ² = 0.02), but distorted segregation in HKS (χ² = 6.52; Figures 3A–C). Among the characters investigated, fffn, ffn, wr, H°, and L were determined to be quantitative characters, and fffn was more affected by the environment in which it was grown (Figures 3D–G and Supplementary Figure S1).

FIGURE 3

FIGURE 3. Distribution of phenotypic data for the seven characters investigated in this study. (A) gynoecy; (B) fruit wart; (C) fruit color; (D) first female flower node; (E) width of fruit ridge; (F) hue angle (H°); (G) lightness variable (L). HKS indicates the plants were grown in Haikou in spring 2014; GZS indicates the plants were grown in Guangzhou in spring 2014; and GZA indicates the plants were grown in Guangzhou in autumn 2014. N, the number of investigated plants; M, mean; S, standard deviation.

QTL Analysis

Sex Expression

Two major QTLs were identified on MC01 when bitter gourd was grown in repeat environments for gy, fffn, and ffn: gy1.1 and gy1.2; fffn1.1 and fffn1.2; and ffn1.1 and ffn1.2, respectively. We found that gy1.1, fffn1.1, and ffn1.1 were located in the same confidence interval, spanning a genetic interval and a physical interval of 17.74 cM and 976.53 Kb, respectively. Near this region, with an interval of ∼10 cM (∼712 Kb), fffn1.2, gy1.2, and ffn1.2 were detected, spanning 15.18 cM and 1,836.22 Kb for genetic and physical intervals, respectively. Therefore, we merged these two regions as the confidence interval of gy/fffn/ffn, which is responsible for sex expression, spanning 42.94 cM and 3,524.92 Kb for genetic and physical intervals, respectively (Figures 4A,B, Table 3 and Supplementary Figure S2 and Supplementary Table S6). After mapping the flanking RAD tags to the OHB3-1 reference, 602 predicted genes were found in this region (Supplementary Table S9).

FIGURE 4

FIGURE 4. QTLs detected for three horticultural traits. (A) The first female flower node (fffn); (B) gynoecy (gy); (C) fruit wart (Fwa); (D) width of fruit ridge (Wr); (E) white fruit color (w); (F) hue angle (H°); (G) lightness variable (L).

TABLE 3

TABLE 3. QTLs detected in ‘K44’ × ‘Dali-11’ F₂ population of bitter gourd.

Fruit Epidermal Structure

One QTL was detected on MC04 for each fruit epidermal structure character, Fwa and Wr for fwa and wr, respectively, explaining 17.70–56.70% of the phenotypic variation in the two environments (HKS and GZS). Both Fwa and Wr were located in the same confidence interval, indicating they are the same QTL, spanning 13.64 cM and 1,469.87 Kb for genetic and physical intervals, respectively (Figures 4C,D, Table 3 and Supplementary Table S7). Then, 218 predicted genes were found in this region (Supplementary Table S9).

Immature Fruit Color

One QTL (w/L/H°) was detected on MC10 for immature fruit color, which had a higher R² value (minimum 65.20%) in the two environments (HKS and GZS). We found that w, L, and H° occupied the same confidence interval, spanning 15.25 cM and 936.70 Kb for genetic and physical intervals, respectively (Figures 4E–G, Table 3 and Supplementary Table S8). Accordingly, 128 predicted genes were found in this region (Supplementary Table S9).

Discussion

One of the aims of this research was to construct a bitter gourd genetic map to anchor assembled scaffolds. RAD-seq has been widely used to construct genome wide genetic maps in ryegrass (Pfender et al., 2011), barley (Zhou et al., 2015), peanut (Gupta et al., 2015), apple (Sun et al., 2015), and soybean (Wang W. et al., 2016). In this study, we successfully created a RAD-based genetic map for bitter gourd, comprising 1,009 SNP markers that distributed on 11 linkage groups. The majority of assembled sequences (85.48%) can be anchored to this map, indicating a high genomic coverage of SNP markers. This map provided a basis for integrating assembled scaffolds into a draft genome of bitter gourd, which served as an important step in our ongoing work on sequencing the bitter gourd genome.

This RAD-based genetic map had a mean marker density of 0.46 (maker/cM), which is higher than the 0.30 mean marker density reported by Matsumura et al. (2014) and 0.42 by Urasaki et al. (2017). Nonetheless, it is far from being high resolution. Our study likely recovered a higher mean marker density due to the choice of restriction enzyme (Baird et al., 2008). We used EcoRI to digest genomic DNA for RAD-seq analysis, but Matsumura et al. (2014) and Urasaki et al. (2017) used PacI and AseI, respectively. In the future, we recommend the use of additional restriction enzymes when attempting to construct a high density or saturated RAD-based genetic map. Furthermore, the recombination rate is a critical consideration for fine mapping and breeding, and we estimated an overall value of 8.77 cM/Mb for bitter gourd. This recombination rate was significantly higher than 3.28 cM/Mb for cucumber (Huang et al., 2009) and 2.30 cM/Mb for watermelon (Ren et al., 2012), but lower than 13.19 cM/Mb for pumpkin (Zhang et al., 2015). Because a high recombination rate can indicate a large genetic distance over a narrow physical distance, our finding suggests that it is beneficial to use fine mapping methods to determine QTLs/genes for bitter gourd compared with the methods used for species with low recombination rates like cucumber and watermelon.

Another aim of this research was to map genetic loci of three horticulturally important traits involved in sex expression, fruit epidermal structure, and immature fruit color. These traits were evaluated using a combination of qualitative and quantitative data. The repeated QTLs/genes were identified in three environments, and both genetic and physical intervals for each locus were delimited. Currently, we are sequencing a bitter gourd genome dependent on this genetic map, and once the gene annotation is completed, fine mapping of these loci combined with candidate gene analysis will assist in cloning the relevant genes more efficiently.

The inheritance of gynoecy was not completely consistent with previous reports (Ram et al., 2006; Behera et al., 2009; Matsumura et al., 2014); it appeared to be controlled by two pairs of genes in the F₂ population in Haikou in spring 2014. Both correlation and inheritance analyses indicated that sex expression in bitter gourd is easily affected by the environmental factors in which it is grown, such as the temperature and photoperiod, which can regulate the sex of cucumber flowers (Yamasaki and Fujii, 2003). Matsumura et al. (2014) located the gy locus at the distal end of the linkage group, and a similar result was found in this study. The reported SNP marker (GTFL-1) that mapped at a distance of 5.46 cM to gy has been aligned at the site of scaffold44_7277273, consistent with our confidence interval (scaffold44_3793313–scaffold44_7318231). Because of its importance, this target region is necessary for further fine mapping. In melon and cucumber, the mechanism of sex determination has been studied in depth. Sexual morphs of melon plants are determined by a combination of alleles at gynoecious (G), andromonoecious (M), and androecious (A) sex loci, controlled by CmWIP1, CmACS7, and CmACS11, respectively (Boualem et al., 2015). Gynoecy in melons is governed by a combination of Cmwip1/Cmwip1, CmACS7/-, -/- (the minus symbol indicates any allele at the locus) at the G, M, and A loci, respectively. CmWIP1 encodes the zinc finger protein WIP1; CmACS7 and CmACS11 encode two 1-aminocyclopropane-1-carboxylic acid synthases. Interestingly, no homologs of these genes were annotated in this QTL (gy) region. Unlike for cucumber and melon, gibberellin (GA3) is the most effective growth regulator for increasing femaleness for bitter gourd (Ghosh and Basu, 1983). These findings imply that the gynoecy of bitter gourd may be controlled by other genetic mechanisms. Further work should be performed to unveil the genetic determination of sexual morphs in bitter gourd.

Fruit warts and ridges are distinct in shape and epidermal distribution, both of which act as identifying features of bitter gourd. The width of the fruit ridge is correlated with the presence of fruit warts, which provides a novel measurement to dissect the genetic loci of fruit epidermal structure. Our results indicated that the two characters of fruit epidermal structure were possibly controlled by a same dominant gene (Fwa/Wr), suggesting a pleiotropic effect. In cucumber, fruit warts consist of spines (non-glandular trichomes) and tubercules governed by csgl1/mict (Li et al., 2015; Zhao et al., 2015) and Tu (Yang et al., 2014), respectively, and csgl1/mict has an epistatic effect on Tu (Yang et al., 2014). Recently, tril, required for the initiation of trichomes, was reported to show epistasis on the csgl1/mict gene (Wang Y.L. et al., 2016). The fruit warts on bitter gourd have no spine, while trichomes can be produced on stems, leaves, fruit stalks, tendrils, and floral organs. Consequently, cloning the genes underlying both fruit warts and trichomes of bitter gourd is essential to understanding the molecular mechanisms that control the formation of fruit warts. Undoubtedly, our results will accelerate the cloning of the Fwa/Wr gene.

Immature fruit color is another important trait directly affecting exterior quality. We used data from our observations of fruit color (green and white) and two indirect measurements to delimit the determined gene, w, into a narrow physical interval (∼936.70 kb). In cucumber, the w gene was confirmed as an APRR2 (two-component response regulator-like) gene. The white color is caused by a single-base insertion in the w gene, which may disable its function in chlorophyll accumulation and chloroplast development (Liu et al., 2015, 2016). Similarly, the white color of bitter gourd fruit is a visible result of a lack of chlorophyll synthesis or disabled development of the chloroplast. The annotation of the bitter gourd genome anchored by the genetic map in the present study will help to directly analyze candidate genes in this physical interval.

Conclusion

In this paper, we constructed a bitter gourd genetic map using the RAD-seq method. This map successfully assisted in anchoring most of the assembled scaffold sequences. The genetic map and anchored assembly revealed very high recombination rates in bitter gourd. We further used this genetic map and component characters to detect the genetic loci controlling sex expression, fruit epidermal structure, and immature fruit color in bitter gourd. Three reliable QTL/gene loci responsible for these traits were identified through multiple phenotypic data and environments. This study provides the foundation for the assembly of the bitter gourd genome and for further dissecting these horticulture-associated loci.

Author Contributions

JuC, JiC, and KH conceived and designed the experiments. SL, YN, RH, QW, JS, and NM contributed to field management and conducted the DNA extractions. WH and ZD performed the bioinformatic analyses. JuC wrote the manuscript and KH revised the manuscript. All authors read and approved the final manuscript.

Funding

This research was supported by funding from the Modern Agricultural Industry Technology System of Guangdong Province, China (2016LM1108), Science and Technology Planning Project of Guangdong Province, China (2014B020202006), special fund for Agricultural Development and Rural Work of Guangdong Province, Guangdong Provincial Department of Science and Technology, China (2016B090918122), and the Science Technology and Innovation Committee of Shenzhen Municipality, China (JCYJ20160331190123578).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2018.00477/full#supplementary-material

FIGURE S1 | Distribution of bitter gourd female flower number investigated in an F₂ population in Haikou in spring 2014. N, the number of plants; M, mean; S, standard deviation.

FIGURE S2 | QTL location of traits involved in sex expression in bitter gourd.

TABLE S1 | The mapping rate and coverage of clean reads statistics for each sample.

TABLE S2 | The statistics of SNPs detected in each sample.

TABLE S3 | List of RAD tags as markers for genetic map construction in this study.

TABLE S4 | The correlation coefficient among bitter gourd characters investigated in Haikou in spring 2014 (HKS).

TABLE S5 | The correlation coefficient among bitter gourd characters investigated in Guangzhou in spring 2014 (GZS).

TABLE S6 | Genome-wide QTL mapping for six expression using multiple phenotypic data and environments.

TABLE S7 | Genome-wide QTL mapping for fruit epidermal structure using multiple phenotypic data and environments.

TABLE S8 | Genome-wide QTL mapping for immature fruit color using multiple phenotypic data and environments.

TABLE S9 | The candidate genes of three genetic loci controlling sex expression, fruit epidermal structure, and immature fruit color in bitter gourd.

References

Baird, N. A., Etter, P. D., Atwood, T. S., Currey, M. C., Shiver, A. L., Lewis, Z. A., et al. (2008). Rapid SNP discovery and genetic mapping using sequenced RAD markers. PLoS One 3:e3376. doi: 10.1371/journal.pone.0003376

PubMed Abstract | CrossRef Full Text | Google Scholar

Behera, T. K., Behera, S., Bharathi, L. K., John, K. J., Simon, P. W., and Staub, J. E. (2010). Bitter gourd: botany, horticulture, breeding. Hortic. Rev. 37, 100–141. doi: 10.1002/9780470543672.ch2

CrossRef Full Text | Google Scholar

Behera, T. K., Dey, S. S., Munshi, A. D., Gaikwad, A. B., Pal, A., and Singh, I. (2009). Sex inheritance and development of gynoecious hybrids in bitter gourd (Momordica charantia L.). Sci. Hortic. 120, 130–133. doi: 10.1016/j.scienta.2008.09.006

CrossRef Full Text | Google Scholar

Boualem, A., Troadec, C., Camps, C., Lemhemdi, A., Morin, H., Sari, M. A., et al. (2015). A cucurbit androecy gene reveals how unisexual flowers develop and dioecy emerges. Science 350, 688–691. doi: 10.1126/science.aac8370

PubMed Abstract | CrossRef Full Text | Google Scholar

Cheng, J. W., Cui, J. J., Tang, X., Wu, Z. M., and Hu, K. L. (2012). Correlation analysis of fruit color with SPAD reading and leaf color(Hue Angle) in bitter gourd. China Veg. 4, 54–57.

Google Scholar

Dalamu, T. K. B., Satyavati, C. T., and Pal, A. (2012). Generation mean analysis of yield related traits and inheritance of fruit colour and surface in bitter gourd. Indian J. Hortic. 69, 65–69.

Google Scholar

Davey, J. W., Hohenlohe, P. A., Etter, P. D., Boone, J. Q., Catchen, J. M., and Blaxter, M. L. (2011). Genome-wide genetic marker discovery and genotyping using next-generation sequencing. Nat. Rev. Genet. 12, 499–510. doi: 10.1038/nrg3012

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghosh, S., and Basu, P. S. (1983). Hormonal regulation of sex expression in Momordica charantia. Physiol. Plant. 57, 301–305. doi: 10.1111/j.1399-3054.1983.tb00915.x

CrossRef Full Text | Google Scholar

Grover, J. K., and Yadav, S. P. (2004). Pharmacological actions and potential uses of Momordica charantia: a review. J. Ethnopharmacol. 93, 123–132. doi: 10.1016/j.jep.2004.03.035

PubMed Abstract | CrossRef Full Text | Google Scholar

Gupta, S. K., Baek, J., Carrasquilla-Garcia, N., and Penmetsa, R. V. (2015). Genome-wide polymorphism detection in peanut using next-generation restriction-site-associated DNA (RAD) sequencing. Mol. Breed. 35:145. doi: 10.1186/1471-2229-12-148

PubMed Abstract | CrossRef Full Text | Google Scholar

Hu, K., Fu, Q., and Wang, G. (2002). Study on the heredity of fruit colour of Momordica charantia. China Veg. 1, 11–12. (In Chinese with English abstract),

Google Scholar

Huang, H. Y., and Hsieh, C. H. (2017). Genetic research on fruit color traits of the bitter gourd (Momordica charantia L.). Hortic. J. 86, 238–243. doi: 10.2503/hortj.MI-159

CrossRef Full Text | Google Scholar

Huang, S., Li, R., Zhang, Z., Li, L., Gu, X., Fan, W., et al. (2009). The genome of the cucumber, Cucumis sativus L. Nat. Genet. 41, 1275–1281. doi: 10.1038/ng.475

PubMed Abstract | CrossRef Full Text | Google Scholar

Kole, C., Olukolu, B. A., Kole, P., Rao, V. K., Bajpai, A., Backiyarani, S., et al. (2012). The first genetic map and positions of major fruit trait loci of bitter melon (Momordica charantia). J. Plant Sci. Mol. Breed. 1:1. doi: 10.7243/2050-2389-1-1

CrossRef Full Text | Google Scholar

Li, Q., Cao, C., Zhang, C., Zheng, S., Wang, Z., Wang, L., et al. (2015). The identification of Cucumis sativus Glabrous 1 (CsGL1) required for the formation of trichomes uncovers a novel function for the homeodomain-leucine zipper I gene. J. Exp. Bot. 66, 2515–2526. doi: 10.1093/jxb/erv046

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, R., Yu, C., Li, Y., Lam, T. W., Yiu, S. M., Kristiansen, K., et al. (2009). SOAP2: an improved ultrafast tool for short read alignment. Bioinformatics 25, 1966–1967. doi: 10.1093/bioinformatics/btp336

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, H., Jiao, J., Liang, X., Liu, J., Meng, H., Chen, S., et al. (2016). Map-based cloning, identification and characterization of the w gene controlling white immature fruit color in cucumber (Cucumis sativus L.). Theor. Appl. Genet. 129, 1–10. doi: 10.1007/s00122-016-2700-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, H., Meng, H., Pan, Y., Liang, X., Jiao, J., Li, Y., et al. (2015). Fine genetic mapping of the white immature fruit color gene w to a 33.0-kb region in cucumber (Cucumis sativus L.). Theor. Appl. Genet. 128, 2375–2385. doi: 10.1007/s00122-015-2592-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Marr, K. L., Mei, X. Y., and Bhattarai, N. K. (2004). Allozyme, morphological and nutritional analysis bearing on the domestication of Momordica charantia L.(Cucurbitaceae). Econ. Bot. 58, 435–455. doi: 10.1663/0013-0001(2004)058[0435:AMANAB]2.0.CO;2

CrossRef Full Text | Google Scholar

Matsumura, H., Miyagi, N., Taniai, N., Fukushima, M., Tarora, K., Shudo, A., et al. (2014). Mapping of the gynoecy in bitter gourd (Momordica charantia) using RAD-Seq analysis. PLoS One 9:e87138. doi: 10.1371/journal.pone.0087138

PubMed Abstract | CrossRef Full Text | Google Scholar

Ooijen, V. (2009). MapQTL 6, Software for the Mapping of Quantitative Trait Loci in Experimental Populations of Diploid Species. Wageningen: Kyazma BV.

Google Scholar

Pfender, W. F., Saha, M. C., Johnson, E. A., and Slabaugh, M. B. (2011). Mapping with RAD (restriction-site associated DNA) markers to rapidly identify QTL for stem rust resistance in Lolium perenne. Theor. Appl. Genet. 122, 1467–1480. doi: 10.1007/s00122-011-1546-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Ram, D., Kumar, S., Singh, M., Rai, M., and Kalloo, G. (2006). Inheritance of gynoecism in bitter gourd (Momordica charantia L.). J. Hered. 97, 294–295. doi: 10.1093/jhered/esj028

PubMed Abstract | CrossRef Full Text | Google Scholar

Ren, Y., Zhao, H., Kou, Q., Jiang, J., Guo, S., Zhang, H., et al. (2012). A high resolution genetic map anchoring scaffolds of the sequenced watermelon genome. PLoS One 7:e29453. doi: 10.1371/journal.pone.0029453

PubMed Abstract | CrossRef Full Text | Google Scholar

Srivastava, V. K., and Nath, P. (1972). Inheritance of some qualitative characters in Momordica charantia L. Indian J. Hortic. 29, 319–321.

Google Scholar

Sun, R., Chang, Y., Yang, F., Wang, Y., Li, H., Zhao, Y., et al. (2015). A dense SNP genetic map constructed using restriction site-associated DNA sequencing enables detection of QTLs controlling apple fruit quality. BMC Genomics 16:747. doi: 10.1186/s12864-015-1946-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Tan, S. P., Kha, T. C., Parks, S. E., and Roach, P. D. (2016). Bitter melon (Momordica charantia L.) bioactive composition and health benefits: a review. Food Rev. Int. 32, 181–202. doi: 10.1080/87559129.2015.1057843

CrossRef Full Text | Google Scholar

Urasaki, N., Takagi, H., Natsume, S., Uemura, A., Taniai, N., Miyagi, N., et al. (2017). Draft genome sequence of bitter gourd (Momordica charantia), a vegetable and medicinal plant in tropical and subtropical regions. DNA Res. 24, 51–58. doi: 10.1093/dnares/dsw047

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Ooijen, J. W. (2006). Joinmap4, Software for Calculation of Genetic Linkage Maps in Experimental Populations. Wageningen: Kyazma B.V.

Google Scholar

Van Wyk, B. E. (2015). A review of commercially important African medicinal plants. J. Ethnopharmacol. 176, 118–134. doi: 10.1016/j.jep.2015.10.031

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, W., Liu, M., Wang, Y., Li, X., Cheng, S., Shu, L., et al. (2016). Characterizing two inter-specific bin maps for the exploration of the QTLs/genes that confer three soybean evolutionary traits. Front. Plant Sci. 7:1248. doi: 10.3389/fpls.2016.01248

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Y. L., Nie, J. T., Chen, H. M., Guo, C. L., Pan, J., He, H. L., et al. (2016). Identification and mapping of Tril, a homeodomain-leucine zipper gene involved in multicellular trichome initiation in Cucumis sativus. Theor. Appl. Genet. 129, 305–316. doi: 10.1007/s00122-015-2628-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., and Xiang, C. (2013). Genetic mapping of QTLs for horticulture traits in a F2-3 population of bitter gourd (Momordica charantia L.). Euphytica 193, 235–250. doi: 10.1007/s10681-013-0932-0

CrossRef Full Text | Google Scholar

Yamasaki, S., and Fujii, N. H. (2003). Photoperiodic regulation of CS-ACS2, CS-ACS4 and CS-ERS gene expression contributes to the femaleness of cucumber flowers through diurnal ethylene production under short-day conditions. Plant Cell Environ. 26, 537–546. doi: 10.1046/j.1365-3040.2003.00984.x

CrossRef Full Text | Google Scholar

Yang, X., Zhang, W., He, H., Nie, J., Bie, B., Zhao, J., et al. (2014). Tuberculate fruit gene Tu encodes a C2H2 zinc finger protein that is required for the warty fruit phenotype in cucumber (Cucumis sativus L.). Plant J. 78, 1034–1046. doi: 10.1111/tpj.12531

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, G., Yi, R., Sun, H., Guo, S., Fan, Z., Jie, Z., et al. (2015). A high-density genetic map for anchoring genome sequences and identifying QTLs associated with dwarf vine in pumpkin (Cucurbita maxima Duch.). BMC Genomics 16:1101. doi: 10.1186/s12864-015-2312-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, J. L., Pan, J. S., Guan, Y., Zhang, W. W., Bie, B. B., Wang, Y. L., et al. (2015). Micro-trichome as a class I homeodomain-leucine zipper gene regulates multicellular trichome development in Cucumis sativus. J. Integr. Plant Biol. 57, 925–935. doi: 10.1111/jipb.12345

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, G., Zhang, Q., Zhang, X. Q., Tan, C., and Li, C. (2015). Construction of high-density genetic map in barley through restriction-site associated DNA sequencing. PLoS One 10:e0133161. doi: 10.1371/journal.pone.0133161

PubMed Abstract | CrossRef Full Text | Google Scholar