Genome-Wide Differentiation of Various Melon Horticultural Groups for Use in GWAS for Fruit Firmness and Construction of a High Resolution Genetic Map

Melon (Cucumis melo L.) is a phenotypically diverse eudicot diploid (2n = 2x = 24) has climacteric and non-climacteric morphotypes and show wide variation for fruit firmness, an important trait for transportation and shelf life. We generated 13,789 SNP markers using genotyping-by-sequencing (GBS) and anchored them to chromosomes to understand genome-wide fixation indices (Fst) between various melon morphotypes and genomewide linkage disequilibrium (LD) decay. The F_ST between accessions of cantalupensis and inodorus was 0.23. The F_ST between cantalupensis and various agrestis accessions was in a range of 0.19–0.53 and between inodorus and agrestis accessions was in a range of 0.21–0.59 indicating sporadic to wide ranging introgression. The EM (Expectation Maximization) algorithm was used for estimation of 1436 haplotypes. Average genome-wide LD decay for the melon genome was noted to be 9.27 Kb. In the current research, we focused on the genome-wide divergence underlying diverse melon horticultural groups. A high-resolution genetic map with 7153 loci was constructed. Genome-wide segregation distortion and recombination rate across various chromosomes were characterized. Melon has climacteric and non-climacteric morphotypes and wide variation for fruit firmness, a very important trait for transportation and shelf life. Various levels of QTLs were identified with high to moderate stringency and linked to fruit firmness using both genome-wide association study (GWAS) and biparental mapping. Gene annotation revealed some of the SNPs are located in β-D-xylosidase, glyoxysomal malate synthase, chloroplastic anthranilate phosphoribosyltransferase, and histidine kinase, the genes that were previously characterized for fruit ripening and softening in other crops.

Introduction

Melon (Cucumis melo L.) is a phenotypically diverse eudicot diploid (2n = 2x = 12) which originated in Asia (Silberstein et al., 2003). According to the morphological observations of Jeffrey (1980) and Stepansky et al. (1999), varieties (vars.) cantalupensis (cantaloupe) and inodorus (honeydew) should be placed in subspecies melo and vars. momordica, conomon, dudaim, and chito in subspecies agrestis (Decker-Walters et al., 2002). Pitrat (2008) grouped melons into 15 widely accepted horticultural groups (cantalupensis, reticulatus, adana, chandalak, ameri, inodorus, chate, flexuosus, dudaim and tibish (in subsp. melo), and momordica, conomon, chinensis, makuwa, and acidulous in subsp. agrestis). Consumer demand for sweet melons has stimulated the selection and breeding of hundreds of cultivars belonging to numerous market types, with local, regional, and international distribution (Paris et al., 2012). Domestication of melons has not been intensively studied, the genetic control of domestication traits, and subsequent diversification and selection processes that led to various melon morphotypes is still poorly understood. A genome-wide sequence of melon of size 375 Mb (83.3% of estimated size) has been made available. This genome has enabled an exhaustive phylogenetic comparison of the melon genome with cucumber (Garcia-Mas et al., 2012). SNP discovery in diverse melon botanical groups will allow marker-anchoring to the whole genome sequence (WGS), thus giving researchers a better understanding of the genetic control of domestication and diversification as shown in the study of Argyris et al. (2015).

Advances in next-generation sequencing technologies have driven the costs of DNA sequencing down to the point that genotyping-by-sequencing (GBS) is now feasible for highly diverse species such as melons. This approach involves reduced-representation sequencing of multiplexed samples and is simple, quick, extremely specific, highly reproducible, and may reach important regions of the genome that are inaccessible to sequence capture approaches (Elshire et al., 2011; Poland and Rife, 2012). The flexibility and low cost of GBS makes this an excellent tool for building high density genetic maps and for use in genome-wide association studies (GWAS) (Poland et al., 2012; Nimmakayala et al., 2014; Reddy et al., 2014). A detailed understanding of population structure and linkage disequilibrium (LD) is paramount for association mapping of the QTLs that underlie various complex traits (Flint-Garcia et al., 2003; Wang et al., 2013; Rincent et al., 2014). The distribution pattern of LD across the genome directly depends on evolutionary forces such as genetic drift, population structure, levels of inbreeding across the genome, and map regions contributing genetic differentiation among the subpopulations. Emergence and maintenance of LD is based on these evolutionary forces and the associated pattern of selection (Ersoz et al., 2007). Esteras et al. (2013) developed a genotyping array for 768 SNPs from a collection of 74 melon accessions with the Illumina GoldenGate technology (Illumina Inc., San Diego, CA), identifying relatively low LD in melons. It is very important to precisely characterize LD blocks across various chromosomes for GWAS studies in diverse horticultural groups such as melon.

Melon is an important desert fruit with tremendous diversity that is a product of consumer preferences from different countries, ecologies, and cultures (Tomason et al., 2013). Understanding divergence and adaptation that underlie the formation of various morphotypes is very important to the development of disease resistant and high quality melons. Fruit quality is related to both internal variables such as fruit firmness, sugar content, acid content, and external variables including fruit shape, size, and texture (García-Ramos et al., 2005). Fruit firmness affects the quality of melon fruit, shelf life, and the ability to transport the fruit over long distances (Moreno et al., 2008; Dahmani-Mardas et al., 2010). Moreno et al. (2008) mapped important QTLs in the locations of candidate genes involved in ethylene regulation, biosynthesis and cell wall degradation using near-isogenic lines (NILs) derived from the non-climacteric melon parental lines PI 161375 and “Piel de Sapo.” Périn et al. (2002) performed genetic analysis for non-climacteric phenotype in fruit tissues on a population of recombinant cantaloupe Charentais × PI 161375 inbred lines to identify several QTLs for ethylene regulation. Dahmani-Mardas et al. (2010) constructed a mutant collection of 4023 melon M2 TILLING families to screen for 11 genes, of which four genes were involved in ethylene/fruit firmness/fruit ripening and identified a mutant for fruit firmness in the TILLING platform.

The current study is to resolve the genetic diversity and relatedness of melon germplasm with the melons of Asia and the western hemisphere using high density SNPs mapped to various chromosomes. Our objective is to compare LD across the chromosomes and to perform GWAS for fruit firmness. Other objectives of the current study were to construct a high-density genetic map for validation of QTLs and to understand genomic features such as colinearity with the publicly available melon genome sequence.

Materials and Methods

One hundred and twenty accessions of various melon horticultural groups representing a world-wide distribution were used for field evaluation. For generation of SNPs, 97 of the most diverse Plant Introductions from our collection were selected based on SSR data representing the botanical groups cantalupensis (51), inodorus (13), reticulatus (5), ameri (3), dudaim (6), flexuosus (5), conomon (9), makuwa (3), acidulus (1), and momordica (1) (Table S1). These PIs were self-pollinated, and the progeny were tested for two seasons (2013–2014) as three replications in two locations (Agriculture Experimental Station, West Virginia State University, Institute, WV and Alcorn State University, Lorman, MS) using a row to plant spacing of 180 × 70 cm. Ten plants per accession were grown per replication. Data were collected pertaining to fruit firmness by testing pressure in kg/cm2 were measured at full maturity using a FT 011 penetrometer (Model# FT 011, Effigy, Alfonsino, Italy) (Table S2). A cross was made between MR-1 (momordica) and “Hale's Best Jumbo” (P2), a western shipper cantaloupe. For building a genetic map, 103 F₂ progeny were generated from a single F₁ plant of P1 × P2 cross. Fruits of F₂ progeny were grown in greenhouse (Table S3). Fruit firmness was measured as pressure to compress (kg/cm2) using a FT 011 penetrometer employing a 0.15 cm tip.

SNP Discovery by GBS

Genomic DNA isolation from the seedlings involved the DNeasy plant mini kit (QIAGEN, Germany), and GBS was as described (Elshire et al., 2011; Reddy et al., 2014). DNA was treated with the restriction enzyme ApeKI, a type II restriction endonuclease, barcoded by accession, and sequenced on an Illumina HiSeq 2500 as described (Elshire et al., 2011). SNPs were extracted using TASSEL-GBS Discovery/Production pipeline (https://bitbucket.org/tasseladmin/tassel-5-source/wiki/Tassel5GBSv2Pipeline). Chromosomal assignment and position of SNPs on the physical map was deduced from the WGS draft of melon (version V3.5) at http://www.melonomics.net. SNPs are designated based on chromosome number and position (e.g., S10_172735351 meaning SNP located at 172735351th position on chromosome 10).

Population Structure Analysis and Divergence

Genetic diversity values were calculated by a neighbor-joining algorithm using TASSEL 5 (www.maizegenetics.net). To further validate the results of NJ-tree, we used principle component analysis (PCA) with the SNP and Variation Suite (SVS v8.1.5) (Golden Helix, Inc., Bozeman, MT, USA; www.goldenhelix.com). Estimation of F_ST was based on Wright's F statistic (Weir and Cockerham, 1984) with use of SVS v8.1.5.

Characterization of LD

For generating GBS data, we considered only SNPs that mapped to the melon whole-genome sequence draft, as the chromosome location of SNPs helps prevent spurious LD, thus reducing errors in GWAS. For haplotype estimation, we used “Minimize historical recombination,” a block-defining algorithm developed by Gabriel et al. (2002). This method is an iterative technique for obtaining maximum likelihood estimates of sample haplotype frequencies. The EM (Expectation Maximization) algorithm was used to estimate adjacent and pairwise measurements of linkage disequilibrium (LD) blocks using haplotype frequencies as formalized by Dempster et al. (1977).

Association Mapping

A set of markers, derived after removing minor allele frequencies, was used to estimate kinship (K) matrix using the software TASSEL 5.0 that uses the proportion of alleles shared between each pair of accessions in the study. PCA correction and method of stratification was followed as in Price et al. (2006). The mixed linear model (MLM) was used to reduce spurious marker trait associations (Type I error showing false positives) resulting from population structure as PCA vectors and K were used as covariance to adjust polygenic background in the analysis.

Construction of a High-Density Genetic Map Using SNP Markers

Linkage analysis and map construction of SNPs generated by GBS were performed using the MultiPoint package (http://www.multiqtl.com) (Mester et al., 2004; Korol et al., 2009; Reddy et al., 2014). GBS resulted in a disproportion between the high number of scored markers for the mapping population and population size. Multilocus ordering aims to pick the most informative markers for building a reliable skeletal map with additional markers being anchored to these framework markers using an algorithm based on evolutionary optimization strategy (Mester et al., 2003). MulitiPoint mapping is based on the maximum likelihood estimation to calculate pairwise recombination fractions (rf) for all marker pairs. In this study, preliminary clustering and assignment of markers to a linkage group (LG) was evaluated at an rf = 0.05 threshold. For example, marker mi may be assigned to an LGj if recombination between mi and at least one marker from LG j is lower than the threshold rf and is lowest compared to the distance to any other LG (Peleg et al., 2008). Selection of “delegates” (bin markers) with the highest information content and stability of their neighborhoods were tested by jackknife resampling, with repeated verification of marker order and removal of unreliable markers to increase the stability of multilocus ordering. SNP loci that mapped to the same location were binned and represented by a single delegate. Stable LGs were joined terminally by incrementally increasing the recombination threshold, with a final rf of 0.30. To avoid erroneous linkage groups based on incorrect marker phase, genotypes of unlinked loci or loci in fragment groups were converted to the alternate phase, reclustered, and assigned to linkage groups as published by Oliver et al. (2011).

Validation of GWAS Results Using Biparental QTL Mapping

Interval mapping and multiple QTL mapping (MQM) were performed using MapQTL5.0 (van Ooijen, 2004). For various QTLs, the genome-wide LOD significance threshold was calculated by the 1000 × permutation test, which restricted the occurrence of Type I statistical errors (false positives) to <5%.

Results

Genotyping-by-Sequencing

In total, 13,756 SNPs were called for uniquely aligned sequence tags and hence are singletons: 1391, 1188, 1151, 1467, 1170, 1314, 1121, 1290, 937, 706, 1072, and 949 SNPs were mapped to Chrs. 1–12, respectively. Following filtering for a minor allele frequency of 0.05%, total SNPs were reduced to 7609: 833, 637, 669, 754, 667, 746, 668, 679, 540, 284, 608, and 524 mapping to Chrs. 1–12, respectively. We found SNPs at average intervals of 30 kb across the melon genome.

Molecular Genetic Diversity among the Melon Collections

To investigate genetic differentiation due to population structure among melon horticultural groups as reflected by these 7609 genome-wide SNPs, we used Neighbor Joining (NJ) analysis. In a phylogenetic tree of all analyzed accessions, cantalupensis and inodorus groups clustered into two separate clusters (subspecies melo) along with some mixtures. The cluster of cantalupensis is grouped with one accession each of conomon, flexuosus, reticulatus, dudaim, and momordica. Cluster of inodorus had 5 cantalupensis, 3 reticulatus, and 3 ameri. The remaining melon horticultural groups vars. momordica, conomon, dudaim, flexousus, makuwa, and acidulous (subspecies agrestis) were grouped into a widely distributed third cluster (Figure 1). Principle component analysis (PCA) was carried out using 7609 SNPs. A PCA scatterplot of individuals on the first two dimensions corroborated the NJ clustering of samples cantalupensis, inodorus, and other groups, with some exceptions as described earlier (Figure 2 and Table S4).

FIGURE 1

Figure 1. Phylogenetic tree constructed with neighbor-joining.

FIGURE 2

Figure 2. First and second components of principal component analysis (PCA) of 7609 single nucleotide polymorphisms (SNPs) within a set of 97 melon accessions. See Table S4 for a list of accessions and respective eigen values for respective positions of individual accessions in the figure.

To further resolve differentiation between cantalupensis and inodorus, we estimated pairwise fixation index (F_ST) across all polymorphisms with MAF ≥ 0.05 (Table 1). All F_ST were highly significant (P < 0.001). The F_ST between accessions of cantalupensis and inodorus was 0.23. The F_ST between cantalupensis and various agrestis accessions was in a range of 0.19–0.53 and between inodorus and agrestis accessions was in a range of 0.21–0.59. Variety reticulatus is ancestral to the genomes of cantalupensis and inodorus, and very close and equidistant to both the groups with a F_ST value of 0.09.

TABLE 1

Table 1. Pairwise F_ST for combinations of various melon horticultural groups.

Haplotypes, LD Decay, and Chromosome-Wise Analysis of LD Blocks

Haplotype distribution is important in comparing common and unique patterns of genetic variation of melon gene pools and has a wide range of applications. The two major processes that shape haplotype structure are the divergence process and breeding history. We used “Minimize historical recombination,” a block-defining algorithm developed by Gabriel et al. (2002). The upper confidence bound was set to 0.98 and the lower bound was set to 0.70. SNPs below MAF of 0.05 were skipped. Maximum block length was set to 160 Kb. The EM (Expectation Maximization) algorithm was used for haplotype estimation with convergence tolerance 0.0001 and frequency threshold of 0.01. Maximum EM iterations were set to 50. We identified 4028 SNPs in 1436 haplotypes from the entire set of melon morphotypes studied (Table S5). We estimated LD by using an entire marker set with MAF ≥ 0.05 and identified 1937 associations in the entire melon collection used in the study (Table S6). The average genome-wide LD decay for the melon genome was noted to be 9.27 Kb, with the means across chromosomes 1–12 being 9.52, 10.11, 4.64, 8.97, 11.07, 7.15, 9.12, 8.87, 8.61, 8.92, 16.46, 7.85 Kb, respectively. Heat maps depicting individual LD blocks (r² for set of markers in LD) across the length of chromosomes are presented in Figures 3, 4.

FIGURE 3

Figure 3. Genome-wide distribution of LD (r²) across melon chromosomes 1–6. X-axis contains the SNP markers positioned on the physical map. Y-axis represents the intensity of LD as shown in the scale.

FIGURE 4

Figure 4. Genome-wide distribution of LD (r²) across melon chromosomes 7–12. X-axis contains the SNP markers positioned on the physical map. Y-axis represents the intensity of LD as shown in the scale.

High Density Genetic Map

A total of 7896 polymorphic SNPs were used to assemble the genetic linkage map using a mapping population that contained 91 F2 progeny, generated from a cross of MR-1 (momordica) and “Hales Best Jumbo” (P2), a Western Shipped Cantaloupe. We generated 807, 964, 628, 606, 724, 743, 520, 687, 631, 551, 450, and 585 SNPs for Chr-1 to Chr-12, respectively. Chromosome distribution of 431 skeletal markers is presented in Figure 5. The remaining “add on” or anchor markers are presented in Figure S1. Skeletal markers are the framework markers that have high confidence. In order to select skeletal markers, SNPs that are violating map stability upon mapping were removed and linkage groups were re-analyzed several times until the map attained complete stability. Chromosomes 1 through 12 consisted of 38, 39, 35, 36, 35, 42, 36, 28, 34, 32, 39, and 37 skeletal markers, respectively with the genetic lengths (cM) of 176, 185.5, 193, 236.1, 179.2, 249.6, 159, 192.6, 216, 132, 208.7, and 163.3, respectively. In addition, the current genetic map defines 1837 recombination events within the skeletal map. Each recombination bin or skeletal marker segregated with multiple add-on markers that aided in the development of the proposed high-density genetic map. This map consists of 79, 54, 79, 95, 91, 100, 57, 85, 65, 43, 63, and 54 markers on Chr-1 through 12, respectively.

FIGURE 5

Figure 5. Genetic map showing positions of skeletal markers on which a high-density genetic map is constructed. Blue is segregation distortion skewed toward female parent and green toward male parent. Red bars show the location of QTLs for fruit firmness.

We examined colinearity of genetic and physical maps for various chromosomes (Figure 6). Markers on Chr-2, -4, -7, -10, -11, and -12 were highly co-linear with respect to their physical locations. Chr-1 and -5 were moderately in agreement with the melon reference sequence. Chr-6, -8, and -9 showed the highest disagreement between the genetic and physical map having a large segment which was not in colinearity with the physical map.

FIGURE 6

Figure 6. Collinearity between genetic (vertical axis) and physical maps (horizontal axis) (markers that are distant from the “line of best fit” are not collinear).

Six hundred and thirty-two loci exhibited significant segregation distortion based on χ² test (P < 0.05), and those locations were mapped onto the final map (Figure 5). Sixteen segregation distortion regions were skewed toward the female parent MR-1 (momordica), and 15 regions showed segregation distortion toward male parent Hale's Best Jumbo. Chromosome-wide distribution of loci skewed toward either of the parent is shown in Table S7. Genome-wide Recombination Rate (GWRR) reflects recombination landscape and is one of the critical factors in shaping the cultivar divergence. GWRR was estimated using the formula cM/Mb. We observed wide variation of GWRRs within and among the chromosomes. A range of 0.01–60.2 noted across the genome (Figure 7). The highest GWRR was on chromosome 12 followed by chromosome 5.

FIGURE 7

Figure 7. Distribution of genome-wide recombination rate (GWRR) along chromosomes in the melon genome. In each plot, the horizontal axis (in Mb) represents the physical distance (PD) along the reference chromosomes and the vertical axis (cM/Mb) the genetic-to-physical distance ratio (green).

Implementation of a Medium-Resolution Genome-Wide Association Study for Fruit Firmness

A MLM was used in the current study to locate QTLs for fruit firmness two consecutive years in two locations. Lists of markers that were linked to various traits are listed along with their corresponding R², P-values, and allelic contributions across all the years are listed in Table 2. A majority of the QTLs detected are repeated across the years and locations indicating the robustness of the identified QTLs. In this study, we found 11 SNPs across five (6, 8, 9, 11, and 12) chromosomes tightly linked to fruit firmness, which is an important trait for ripening as well as shelf life and transportation of melon fruits. All SNPs except one (S11_14670800) were located on the coding regions of annotated genes (Table 3).

TABLE 2

Table 2. Common QTL for fruit firmness identified using GWAS across various years and locations.

TABLE 3

Table 3. Fruit firmness QTL identified using MQM mapping in a biparental cross.

QTL Validation in Biparental Mapping

F₂ progeny of MR-1 and Hale's Best Jumbo were evaluated in controlled conditions to validate various QTLs identified using GWAS. Entire list of QTLs identified using GWAS in the current study could be validated with the results of MQM mapping. Additional QTLs that are identified in MQM mapping but not in GWAS were not shown as these QTLs need further validation. Details pertaining to the QTLs that are validated in GWAS are presented in Table 4.

TABLE 4

Table 4. Location of SNPs on various chromosomes and annotation of corresponding genes.

Discussion

As genotypic data become easier to obtain, it is possible to analyse a more complete and accurate landscape of genetic diversity and linkage disequilibrium, reduce false positives arising from population structure, and target true biological associations (Lipka et al., 2015). In this research, we took advantage of whole-genome sequence data of melon to map SNPs generated by GBS to various chromosomes, thus estimating chromosome-wise divergence and characterized genome-wide LD. To date, the genetic basis of this diversity and the consequences of selection on genetic variation involving various horticultural groups have not yet been studied on a genome-wide basis (Blanca et al., 2012).

Several important studies have been conducted to generate SNPs for the melon genome for use in characterization of genetic variation, population structure, and LD (Blanca et al., 2012; Esteras et al., 2013; Leida et al., 2015). A true sense of genomic characterization is only possible after anchoring to the WGS (Argyris et al., 2015; Sanseverino et al., 2015). Sanseverino et al. (2015) resequenced seven melon varieties using a paired-end approach to generate 4,556,377 SNPs between melo and agrestis to study genome-wide nucleotide diversity. However, a large number of diverse accessions needs to be included, such as in the current study, to accurately understand genome-wide distribution of LD and to identify markers linked to various traits.

Our genetic diversity study clearly differentiated melo and agrestis into various clusters, this is in agreement with several previous studies (Stepansky et al., 1999; Garcia-Mas et al., 2000; Mliki et al., 2001; Staub et al., 2004; Nakata et al., 2005; Nimmakayala et al., 2009; Tomason et al., 2013). An admixture of melo and agrestis genomes by intentional and unintentional crossing is evident as shown in several previous genetic diversity studies in melon (Blanca et al., 2012; Esteras et al., 2013; Hu et al., 2015; Sanseverino et al., 2015). As cantalupensis and inodorus cultivars are fully cross-compatible with the groups of agrestis, the variability found in agrestis groups (such as conomon and momordica) has been used as a source of disease resistance for breeding purposes (Blanca et al., 2012; Sanseverino et al., 2015). In spite of conomon, flexuosus, and momordica mixtures in clusters belonging to cantalupensis and inodorus, our study clearly differentiates genetic boundaries as shown in NJ tree analysis and principle component analysis between cantalupensis and inodorus.

Pairwise F_ST indices for various SNP markers across the length of chromosomes could be used to identify important genomic regions contributing to genetic differentiation among horticultural groups in the study. Similar to our results, Sanseverino et al. (2015) compared melo and agrestis and identified genomic regions exhibiting extreme population differentiation regions on chromosome 1, 3, 7, 8, and 11, indicating that the process of genetic differentiation of melon subspecies is a genome-wide process.

In spite of rapid progress in sequencing technologies for creating affordable physical maps, high density genetic linkage maps are still indispensable for identification of genomic regions carrying quantitative trait loci (QTL) controlling agronomical traits. In addition, genetic linkage maps with high density SNPs and structural variants are a prerequisite for further map-based cloning and comparative genome analysis (Delourme et al., 2013; Reddy et al., 2014; Diaz et al., 2015; Ren et al., 2015). The high-resolution genetic map presented in the current study should prove useful in associating colinearity with the physical map. The disagreements between the genetic map and the physical map should be extremely useful for future melon genome sequencing endevours. Since the parents of the genetic mapping population are momordica and cantalupensis, map regions that are skewed to female or male parent will shed light on cryptic recombination sites between the morphotype genomes. A MLM was used in the current study to locate QTLs for two consecutive years in two locations for fruit firmness. The most intriguing part of the current study is that many of the QTLs identified based on GWAS could be validated with the QTL analysis using biparental mapping.

From our association study, we mapped S8_4742008, S12_22089329, S12_22581778, and S11_10348094 to the coding region of β-D-xylosidase, glyoxysomal malate synthase, chloroplastic anthranilate phosphoribosyltransferase, and histidine kinase respectively that are strongly associated with cell wall metabolism (Fischer and Bennett, 1991; Huber and O'Donoghue, 1993; Zablackis et al., 1995; Brummell and Harpster, 2001; Minorsky, 2002; Goujon et al., 2003; Itai et al., 2003; Tateishi et al., 2005; Di Santo et al., 2009) thereby affecting fruit firmness. β-xylosidase (AtBXL1) in Arabidopsis is thought to be involved in the organization and loosening cellulose deposition in the secondary cell wall (Goujon et al., 2003). During fruit ripening, several cell wall metabolizing enzymes contribute to changes in cell wall architecture (Fischer and Bennett, 1991). In addition, pectic and hemicellulosic polysaccharides become soluble causing depolymerization of the cells with the help of neutral sugars (Huber and O'Donoghue, 1993; Brummell and Harpster, 2001). The hemicellulosic components of the primary cell wall in dicots consist of well-characterized xyloglucans (Zablackis et al., 1995). The metabolism of xyloglucans in the cellulose microfibril network is believed to be important for cell-wall expansion. This loosening of the cell wall happens when xyloglucans tether the adjacent microfibrils for further modification (Minorsky, 2002). Various modifications of the cell-wall in developing and ripening fruits are thought to be mediated by cell-wall degrading enzymes such as β-D-xylosidase, which is involved in the breakdown of xylans. Earlier studies demonstrated that one of the functions of β-D-xylosidase is to control fruit development and ripening in tomato (Itai et al., 2003) and Japanese pear (Tateishi et al., 2005). This gene may play a role in melon fruit firmness by regulating rind thickness and pressure. Interestingly, ethylene production and accumulation of transcripts of β-D-xylosidase coincidentally occurred with increased fruit firmness in peach (Di Santo et al., 2009). Ethylene-stimulated malate synthase, a key enzyme responsible for malic acid synthesis in the glyoxylate cycle, converts lipids to carbohydrates and was found to be highly expressed during fruit ripening in banana (Pua et al., 2003). Moreover, malate has been shown to be involved in starch metabolism, ripening and postharvest softening in tomato (Centeno et al., 2011). Our GWAS identified another SNP (S12_22581778) located in the coding domain of a putative anthranilate phosphoribosyl-transferase, a gene that plays a role in cell wall metabolism in the presence of ethylene (Li et al., 2012). Histidine kinases (HKs) showed strong association with fruit firmness in the current research. HKs function in two-component regulation systems to transduce signals from hormones and external cues to multiple downstream functions in fruit development (Hwang et al., 2002; Singh and Kumar, 2012).

We characterized a GWAS panel offering the best operational representation of melon diversity so far. This study showed that the melon is one of the most diverse cultivars with considerably rapid decay of LD when assessed at the genome-wide scale. Rapid decay of LD in melon indicates a need for higher-density SNP panels for performing GWAS effectively. These GWAS results demonstrate that high-density SNP markers developed in the study provide an effective tool to dissect the genetic architecture of fruit firmness, although additional evidence is needed to support the identified loci and candidate genes.

Author Contributions

PN, UR, AL, WW, JM, and JG designed the study and drafted the manuscript. PN, YT, VA, AA, VV, GS, and GP generated field and fruit firmness phenotyping. PN, VA, AA, and VV extracted DNA and assisted to generate genomewide SNPs. AK, YR, and UR generated high resolution genetic map. JG provided whole genome sequence draft and mapped SNPs to the genome. UR, YT, PN, VA, and TS performed GWAS, biparental QTL analysis, and gene annotation.

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.

Acknowledgments

Funding support is provided by Gus R. Douglass Institute (Evans Allen Project to Nimmakayala) and USDA-NIFA (2010-02247 and 2012-02511).

Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fpls.2016.01437

Figure S1. High resolution genetic map of melon with 7153 loci.

Table S1. List of melon accessions used in the current study.

Table S2. Fruit pressure across various years and locations.

Table S3. Fruit pressure distribution in F2 progeny of a cross of MR-1 (momordica) and “Hale's Best Jumbo” (P2), a Western Shipped Cantaloupe.

Table S4. Eigen values for the first 10 principle components estimated for various accessions.

Table S5. Haplotype blocks across entire melon collection in the study.

Table S6. Linkage disequilibrium (LD) analysis for adjacent SNP pairs across melon genome.

Table S7. Linkage disequilibrium (LD) analysis for adjacent SNP pairs across melon genome.

References

Argyris, J. M., Ruiz-Herrera, A., Madriz-Masis, P., Sanseverino, W., Morata, J., Pujol, M., et al. (2015). Use of targeted SNP selection for an improved anchoring of the melon (Cucumis melo L.) scaffold genome assembly. BMC Genomics 16:4. doi: 10.1186/s12864-014-1196-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Blanca, J., Esteras, C., Ziarsolo, P., Pérez, D., Fernández-Pedrosa, V., Collado, C., et al. (2012). Transcriptome sequencing for SNP discovery across Cucumis melo. BMC Genomics 13:280. doi: 10.1186/1471-2164-13-280

PubMed Abstract | CrossRef Full Text | Google Scholar

Brummell, D. A., and Harpster, M. H. (2001). Cell wall metabolism in fruit softening and quality and its manipulation in transgenic plants. Plant Mol. Biol. 47, 311–340. doi: 10.1023/A:1010656104304

PubMed Abstract | CrossRef Full Text | Google Scholar

Centeno, D. C., Osorio, S., Nunes-Nesi, A., Bertolo, A. L., Carneiro, R. T., Araújo, W. L., et al. (2011). Malate plays a crucial role in starch metabolism, ripening, and soluble solid content of tomato fruit and affects postharvest softening. Plant Cell 23, 162–184. doi: 10.1105/tpc.109.072231

PubMed Abstract | CrossRef Full Text | Google Scholar

Dahmani-Mardas, F., Troadec, C., Boualem, A., Lévêque, S., Alsadon, A. A., Aldoss, A. A., et al. (2010). Engineering melon plants with improved fruit shelf life using the TILLING approach. PLoS ONE 5:e15776. doi: 10.1371/journal.pone.0015776

PubMed Abstract | CrossRef Full Text | Google Scholar

Decker-Walters, D., Chung, S.-M., Staub, J., Quemada, H., and López-Sesé, A. (2002). The origin and genetic affinities of wild populations of melon (Cucumis melo, Cucurbitaceae) in North America. Plant Syst. Evol. 233, 183–197. doi: 10.1007/s00606-002-0191-3

CrossRef Full Text | Google Scholar

Delourme, R., Falentin, C., Fomeju, B. F., Boillot, M., Lassalle, G., André, I., et al. (2013). High-density SNP-based genetic map development and linkage disequilibrium assessment in Brassica napus L. BMC Genomics 14:120. doi: 10.1186/1471-2164-14-120

PubMed Abstract | CrossRef Full Text

Dempster, A. P., Laird, N. M., and Rubin, D. B. (1977). Maximum likelihood from incomplete data via the EM algorithm. J. R. Stat. Soc. 39, 1–38.

Google Scholar

Diaz, A., Forment, J., Argyris, J. M., Fukino, N., Tzuri, G., Harel-Beja, R., et al. (2015). Anchoring the consensus ICuGI genetic map to the melon (Cucumis melo L.) genome. Mol. Breed. 35, 1–7. doi: 10.1007/s11032-015-0381-7

CrossRef Full Text | Google Scholar

Di Santo, M. C., Pagano, E. A., and Sozzi, G. O. (2009). Differential expression of α-L-arabinofuranosidase and α-L-arabinofuranosidase/β-D-xylosidase genes during peach growth and ripening. Plant Physiol. Biochem. 47, 562–569. doi: 10.1016/j.plaphy.2009.02.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Elshire, R. J., Glaubitz, J. C., Sun, Q., Poland, J. A., Kawamoto, K., Buckler, E. S., et al. (2011). A robust, simple genotyping-by-sequencing (GBS) approach for high diversity species. PLoS ONE 6:e19379. doi: 10.1371/journal.pone.0019379

PubMed Abstract | CrossRef Full Text | Google Scholar

Ersoz, E. S., Yu, J., and Buckler, E. S. (2007). “Applications of linkage disequilibrium and association mapping in crop plants,” in Genomics-Assisted Crop Improvement, eds R. K. Varshney and R. Tuberosa (Dordrecht: Springer), 97–119.

Esteras, C., Formisano, G., Roig, C., Díaz, A., Blanca, J., Garcia-Mas, J., et al. (2013). SNP genotyping in melons: genetic variation, population structure, and linkage disequilibrium. Theor. Appl. Genet. 126, 1285–1303. doi: 10.1007/s00122-013-2053-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Fischer, R. L., and Bennett, A. B. (1991). Role of cell wall hydrolases in fruit ripening. Annu. Rev. Plant Biol. 42, 675–703. doi: 10.1146/annurev.pp.42.060191.003331

CrossRef Full Text | Google Scholar

Flint-Garcia, S. A., Thornsberry, J. M., and Iv, B. (2003). Structure of linkage disequilibrium in plants*. Annu. Rev. Plant Biol. 54, 357–374. doi: 10.1146/annurev.arplant.54.031902.134907

PubMed Abstract | CrossRef Full Text | Google Scholar

Gabriel, S. B., Schaffner, S. F., Nguyen, H., Moore, J. M., Roy, J., Blumenstiel, B., et al. (2002). The structure of haplotype blocks in the human genome. Science 296, 2225–2229. doi: 10.1126/science.1069424

PubMed Abstract | CrossRef Full Text | Google Scholar

Garcia-Mas, J., Benjak, A., Sanseverino, W., Bourgeois, M., Mir, G., González, V. M., et al. (2012). The genome of melon (Cucumis melo L.). Proc. Natl. Acad. Sci. U.S.A. 109, 11872–11877. doi: 10.1073/pnas.1205415109

PubMed Abstract | CrossRef Full Text | Google Scholar

Garcia-Mas, J., Oliver, M., Gomez-Paniagua, H., and De Vicente, M. (2000). Comparing AFLP, RAPD and RFLP markers for measuring genetic diversity in melon. Theor. Appl. Genet. 101, 860–864. doi: 10.1007/s001220051553

CrossRef Full Text | Google Scholar

García-Ramos, F. J., Valero, C., Homer, I., Ortiz-Cañavate, J., and Ruiz-Altisent, M. (2005). Non-destructive fruit firmness sensors: a review. Span. J. Agric. Res. 3, 61–73. doi: 10.5424/sjar/2005031-125

CrossRef Full Text | Google Scholar

Goujon, T., Minic, Z., El Amrani, A., Lerouxel, O., Aletti, E., Lapierre, C., et al. (2003). AtBXL1, a novel higher plant (Arabidopsis thaliana) putative beta-xylosidase gene, is involved in secondary cell wall metabolism and plant development. Plant J. 33, 677–690. doi: 10.1046/j.1365-313X.2003.01654.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Huber, D. J., and O'Donoghue, E. M. (1993). Polyuronides in avocado (Persea americana) and tomato (Lycopersicon esculentum) fruits exhibit markedly different patterns of molecular weight downshifts during ripening. Plant Physiol. 102, 473–480.

PubMed Abstract | Google Scholar

Hu, J., Wang, P., Su, Y., Wang, R., Li, Q., and Sun, K. (2015). Microsatellite diversity, population structure, and core collection formation in melon germplasm. Plant Mol. Biol. Report. 33, 439–447. doi: 10.1007/s11105-014-0757-6

CrossRef Full Text | Google Scholar

Hwang, I., Chen, H.-C., and Sheen, J. (2002). Two-component signal transduction pathways in Arabidopsis. Plant Physiol. 129, 500–515. doi: 10.1104/pp.005504

PubMed Abstract | CrossRef Full Text | Google Scholar

Itai, A., Ishihara, K., and Bewley, J. D. (2003). Characterization of expression, and cloning, of β-d-xylosidase and α-l-arabinofuranosidase in developing and ripening tomato (Lycopersicon esculentum Mill.) fruit. J. Exp. Bot. 54, 2615–2622. doi: 10.1093/jxb/erg291

PubMed Abstract | CrossRef Full Text | Google Scholar

Jeffrey, C. (1980). A review of the Cucurbitaceae. Bot. J. Linn. Soc. 81, 233–247. doi: 10.1111/j.1095-8339.1980.tb01676.x

CrossRef Full Text | Google Scholar

Korol, A. M. D, Frenkel, Z., and Ronin, Y. (2009). Methods for Genetic Analysis in the Triticeae. New York, NY: Springer.

Google Scholar

Leida, C., Moser, C., Esteras, C., Sulpice, R., Lunn, J. E., De Langen, F., et al. (2015). Variability of candidate genes, genetic structure and association with sugar accumulation and climacteric behavior in a broad germplasm collection of melon (Cucumis melo L.). BMC Genet. 16:28. doi: 10.1186/s12863-015-0183-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Z., Xu, C., Li, K., Yan, S., Qu, X., and Zhang, J. (2012). Phosphate starvation of maize inhibits lateral root formation and alters gene expression in the lateral root primordium zone. BMC Plant Biol. 12:89. doi: 10.1186/1471-2229-12-89

PubMed Abstract | CrossRef Full Text | Google Scholar

Lipka, A. E., Kandianis, C. B., Hudson, M. E., Yu, J., Drnevich, J., Bradbury, P. J., et al. (2015). From association to prediction: statistical methods for the dissection and selection of complex traits in plants. Curr. Opin. Plant Biol. 24, 110–118. doi: 10.1016/j.pbi.2015.02.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Mester, D. R. Y., Minkov, D., Nevo, E., and Korol, A. (2003). Constructing large-scale genetic maps using an evolutionary strategy algorithm. Genetics 165, 2269–2282.

PubMed Abstract | Google Scholar

Mester, D. R. Y., Nevo, E., and Korol, A.b (2004). Fast and high precision algorithms for optimization in large-scale genomic problems. Comput. Biol. Chem. 28, 281–290. doi: 10.1016/j.compbiolchem.2004.08.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Minorsky, P. V. (2002). The wall becomes surmountable. Plant Physiol. 128, 345–353. doi: 10.1104/pp.900022

PubMed Abstract | CrossRef Full Text | Google Scholar

Mliki, A., Staub, J. E., Zhangyong, S., and Ghorbel, A. (2001). Genetic diversity in melon (Cucumis melo L.): anevaluation of African germplasm. Genet. Resour. Crop Evol. 48, 587–597. doi: 10.1023/A:1013840517032

CrossRef Full Text | Google Scholar

Moreno, E., Obando, J. M., Dos-Santos, N., Fernández-Trujillo, J. P., Monforte, A. J., and Garcia-Mas, J. (2008). Candidate genes and QTLs for fruit ripening and softening in melon. Theor. Appl. Genet. 116, 589–602. doi: 10.1007/s00122-007-0694-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Nakata, E., Staub, J. E., López-Sesé, A. I., and Katzir, N. (2005). Genetic diversity of Japanese melon cultivars (Cucumis melo L.) as assessed by random amplified polymorphic DNA and simple sequence repeat markers. Genet. Resour. Crop Evol. 52, 405–419. doi: 10.1007/s10722-005-2258-9

CrossRef Full Text | Google Scholar

Nimmakayala, P., Levi, A., Abburi, L., Abburi, V. L., Tomason, Y. R., Saminathan, T., et al. (2014). Single nucleotide polymorphisms generated by genotyping by sequencing to characterize genome-wide diversity, linkage disequilibrium, and selective sweeps in cultivated watermelon. BMC Genomics 15:767. doi: 10.1186/1471-2164-15-767

PubMed Abstract | CrossRef Full Text | Google Scholar

Nimmakayala, P., Tomason, Y. R., Jeong, J., Vajja, G., Levi, A., Gibson, P., et al. (2009). Molecular diversity in the Ukrainian melon collection as revealed by AFLPs and microsatellites. Plant Genet. Resour. 7, 127–134. doi: 10.1017/S1479262108098481

CrossRef Full Text | Google Scholar

Oliver, R. E., Jellen, E. N., Ladizinsky, G., Korol, A. B., Kilian, A., Beard, J. L., et al. (2011). New Diversity Arrays Technology (DArT) markers for tetraploid oat (Avena magna Murphy et Terrell) provide the first complete oat linkage map and markers linked to domestication genes from hexaploid A. sativa L. Theor. Appl. Genet. 123, 1159–1171. doi: 10.1007/s00122-011-1656-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Paris, H. S., Amar, Z., and Lev, E. (2012). Medieval emergence of sweet melons, Cucumis melo (Cucurbitaceae). Ann. Bot. 110, 23–33. doi: 10.1093/aob/mcs098

PubMed Abstract | CrossRef Full Text | Google Scholar

Peleg, Z., Saranga, Y., Suprunova, T., Ronin, Y., Röder, M., Kilian, A., et al. (2008). High-density genetic map of durum wheat × wild emmer wheat based on SSR and DArT markers. Theor. Appl. Genet. 117, 103–115. doi: 10.1007/s00122-008-0756-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Périn, C., Gomez-Jimenez, M., Hagen, L., Dogimont, C., Pech, J.-C., Latché, A., et al. (2002). Molecular and genetic characterization of a non-climacteric phenotype in melon reveals two loci conferring altered ethylene response in fruit. Plant Physiol. 129, 300–309. doi: 10.1104/pp.010613

PubMed Abstract | CrossRef Full Text | Google Scholar

Pitrat, M. (2008). “Melon (Cucumis melo L.),” in Handbook of Plant Breeding—Vegetables, Vol. 1, eds J. Prohens and F. Nuez (New York, NY: Springer), 283–315. doi: 10.1007/978-0-387-30443-4_9

CrossRef Full Text

Poland, J. A., Brown, P. J., Sorrells, M. E., and Jannink, J.-L. (2012). Development of high-density genetic maps for barley and wheat using a novel two-enzyme genotyping-by-sequencing approach. PLoS ONE 7:e32253. doi: 10.1371/journal.pone.0032253

PubMed Abstract | CrossRef Full Text | Google Scholar

Poland, J. A., and Rife, T. W. (2012). Genotyping-by-sequencing for plant breeding and genetics. Plant Genome 5, 92–102. doi: 10.3835/plantgenome2012.05.0005

CrossRef Full Text | Google Scholar

Price, A. L., Patterson, N. J., Plenge, R. M., Weinblatt, M. E., Shadick, N. A., and Reich, D. (2006). Principal components analysis corrects for stratification in genome-wide association studies. Nat. Genet. 38, 904–909. doi: 10.1038/ng1847

PubMed Abstract | CrossRef Full Text | Google Scholar

Pua, E. C., Chandramouli, S., Han, P., and Liu, P. (2003). Malate synthase gene expression during fruit ripening of Cavendish banana (Musa acuminata cv. Williams). J. Exp. Bot. 54, 309–316. doi: 10.1093/jxb/erg030

PubMed Abstract | CrossRef Full Text | Google Scholar

Reddy, U. K., Nimmakayala, P., Levi, A., Abburi, V. L., Saminathan, T., Tomason, Y. R., et al. (2014). High-resolution genetic map for understanding the effect of genome-wide recombination rate on nucleotide diversity in watermelon. G3 (Bethesda). 4, 2219–2230. doi: 10.1534/g3.114.012815

PubMed Abstract | CrossRef Full Text | Google Scholar

Ren, R., Ray, R., Li, P., Xu, J., Zhang, M., Liu, G., et al. (2015). Construction of a high-density DArTseq SNP-based genetic map and identification of genomic regions with segregation distortion in a genetic population derived from a cross between feral and cultivated-type watermelon. Mol. Genet. Genomics 290, 1457–1470. doi: 10.1007/s00438-015-0997-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Rincent, R., Moreau, L., Monod, H., Kuhn, E., Melchinger, A. E., Malvar, R. A., et al. (2014). Recovering power in association mapping panels with variable levels of linkage disequilibrium. Genetics 197, 375–387. doi: 10.1534/genetics.113.159731

PubMed Abstract | CrossRef Full Text | Google Scholar

Sanseverino, W., Hénaff, E., Vives, C., Pinosio, S., Burgos-Paz, W., Morgante, M., et al. (2015). Transposon insertions, structural variations, and SNPs contribute to the evolution of the melon genome. Mol. Biol. Evol. 32, 2760–2774. doi: 10.1093/molbev/msv152

PubMed Abstract | CrossRef Full Text | Google Scholar

Silberstein, L., Kovalski, I., Brotman, Y., Perin, C., Dogimont, C., Pitrat, M., et al. (2003). Linkage map of Cucumis melo including phenotypic traits and sequence-characterized genes. Genome 46, 761–773. doi: 10.1139/g03-060

PubMed Abstract | CrossRef Full Text | Google Scholar

Singh, G., and Kumar, R. (2012). Genome-wide in silico analysis of plant two component signaling system in woody model plant Populus trichocarpa. Res. Plant Biol. 2, 13–23.

Staub, J. E., López-Sesé, A. I., and Fanourakis, N. (2004). Diversity among melon landraces (Cucumis melo L.) from Greece and their genetic relationships with other melon germplasm of diverse origins. Euphytica 136, 151–166. doi: 10.1023/B:EUPH.0000030667.63614.bd

CrossRef Full Text | Google Scholar

Stepansky, A., Kovalski, I., and Perl-Treves, R. (1999). Intraspecific classification of melons (Cucumis melo L.) in view of their phenotypic and molecular variation. Plant Syst. Evol. 217, 313–332. doi: 10.1007/BF00984373

CrossRef Full Text | Google Scholar

Tateishi, A., Mori, H., Watari, J., Nagashima, K., Yamaki, S., and Inoue, H. (2005). Isolation, characterization, and cloning of α-L-arabinofuranosidase expressed during fruit ripening of Japanese pear. Plant Physiol. 138, 1653–1664. doi: 10.1104/pp.104.056655

PubMed Abstract | CrossRef Full Text | Google Scholar

Tomason, Y., Nimmakayala, P., Levi, A., and Reddy, U. K. (2013). Map-based molecular diversity, linkage disequilibrium and association mapping of fruit traits in melon. Mol. Breed. 31, 829–841. doi: 10.1007/s11032-013-9837-9

CrossRef Full Text | Google Scholar

van Ooijen, J. W. (2004). MapQTL® 5, Software for the Mapping of Quantitative Trait Loci in Experimental Populations. Wageningen: Kyazma B.V.

Wang, Y. H., Upadhyaya, H. D., Burrell, A. M., Sahraeian, S. M., Klein, R. R., and Klein, P. E. (2013). Genetic structure and linkage disequilibrium in a diverse, representative collection of the C4 model plant, Sorghum bicolor. G3 (Bethesda). 3, 783–793. doi: 10.1534/g3.112.004861

PubMed Abstract | CrossRef Full Text | Google Scholar

Weir, B. S., and Cockerham, C. C. (1984). Estimating F-statistics for the analysis of population structure. Evolution 38, 1358–1370. doi: 10.2307/2408641

CrossRef Full Text | Google Scholar

Zablackis, E., Huang, J., Muller, B., Darvill, A. G., and Albersheim, P. (1995). Characterization of the cell-wall polysaccharides of Arabidopsis thaliana leaves. Plant Physiol. 107, 1129–1138. doi: 10.1104/pp.107.4.1129

PubMed Abstract | CrossRef Full Text | Google Scholar