Combining QTL-seq and linkage mapping to uncover the genetic basis of single vs. paired spikelets in the advanced populations of two-ranked maize×teosinte

Background Teosinte ear bears single spikelet, whereas maize ear bears paired spikelets, doubling the number of grains in each cupulate during maize domestication. In the past 20 years, genetic analysis of single vs. paired spikelets (PEDS) has been stagnant. A better understanding of genetic basis of PEDS could help fine mapping of quantitative trait loci (QTL) and cloning of genes. Results In this study, the advanced mapping populations (BC3F2 and BC4F2) of maize × teosinte were developed by phenotypic recurrent selection. Four genomic regions associated with PEDS were detected using QTL-seq, located on 194.64–299.52 Mb, 0–162.80 Mb, 12.82–97.17 Mb, and 125.06–157.01 Mb of chromosomes 1, 3, 6, and 8, respectively. Five QTL for PEDS were identified in the regions of QTL-seq using traditional QTL mapping. Each QTL explained 1.12–38.05% of the phenotypic variance (PVE); notably, QTL qPEDS3.1 with the average PVE of 35.29% was identified in all tests. Moreover, 14 epistatic QTL were detected, with the total PVE of 47.57–66.81% in each test. The QTL qPEDS3.1 overlapped with, or was close to, one locus of 7 epistatic QTL. Near-isogenic lines (NILs) of QTL qPEDS1.1, qPEDS3.1, qPEDS6.1, and qPEDS8.1 were constructed. All individuals of NIL-qPEDS6.1(MT1) and NIL-qPEDS8.1(MT1) showed paired spikelets (PEDS = 0), but the flowering time was 7 days shorter in the NIL-qPEDS8.1(MT1). The ratio of plants with PEDS > 0 was low (1/18 to 3/18) in the NIL-qPEDS1.1(MT1) and NIL-qPEDS3.1(MT1), maybe due to the epistatic effect. Conclusion Our results suggested that major QTL, minor QTL, epistasis and photoperiod were associated with the variation of PEDS, which help us better understand the genetic basis of PEDS and provide a genetic resource for fine mapping of QTL. Supplementary Information The online version contains supplementary material available at 10.1186/s12870-021-03353-3.


Background
Improvement in grain yield was one of the main goals during maize domestication [1]. Among the many factors that affect the grain yield, grain number per ear is a major determinant. The number of kernels per maize ear had been raised to 200 or more from its closest wild relative teosinte that has 6 to 12 kernels per ear. The single spikelets in teosinte ear were completely transformed Open Access *Correspondence: liujian@sicau.edu.cn † Zhengjie Chen, Dengguo Tang and Kun Hu contributed equally to this work. 1 Maize Research Institute, Sichuan Agricultural University, No.211 Huiming Road, Wenjiang District, Chengdu 611130, Sichuan, China Full list of author information is available at the end of the article into paired spikelets in maize ears, with the two ranks of teosinte ear changed into multiple ranks of maize ears [2]; hence, the number of kernels per row or length of ear in maize was significantly increased compared to that in teosinte [3]. Notably, the transformation of single spikelet into paired spikelets could double the kernel number, which was one of the key steps in maize domestication. Dissecting the genetic architecture of single vs. paired spikelets could improve our understanding about procedures of maize domestication and the genetic mechanism of yield improvement. However, studying of single vs. paired spikelets has not progressed much in the past 20 years, and the genetic basis of single vs. paired spikelets remains poorly known.
Scanning electron microscopy (SEM) analysis of inflorescence development revealed that the inflorescence meristem (IM) transforms into spikelet pair meristems (SPM). Subsequently, each SPM produces two distinct spikelet meristems (SM), one sessile SM and one pedicellate SM [4]. Both sessile and pedicellate SM develop normally in maize ear, resulting in paired spikelets. In contrast, the pedicellate SM is aborted in teosinte ear, leading to single spikelets [5].
Several studies were performed to explore the inheritance of single vs. paired spikelets. As early as 1920, Collins and Kempton showed for the first time that the distributions of single vs. paired spikelets were continuous, deviating from the Mendelian pattern of 3:1 in the maize-teosinte F 2 population [6,7]. Later, Langham suggested that single vs. paired spikelets trait was controlled by a single gene [8], but Szabó and Burr inferred that two independent genes on chromosomes 4 and 8 were involved [9]. In contrast, Mangelsdorf [10], Rogers [11] and Doebley et al. [7] indicated that multiple genes plus epistasis governed single vs. paired spikelets. The differences in these results might be due to that different parents were used to develop the segregating populations and different analysis methods were used in these studies. Currently, the number of studies of QTL mapping for single vs. paired spikelets (PEDS) is less than 10. Doebley et al. mapped four QTL related to PEDS in the F 2 population of maize race chapalote (Sin 2) × Z. mays Ssp. Mexicana (Doebley 643), located on chromosomes 1 to 4, with the proportion of phenotypic variance (R 2 ) explained ranging from 0.05 to 0.24 [7]. One year later, three additional QTL for PEDS located on chromosomes 5, 6 and 7 were identified in the same dataset by a modified statistical analysis [12]. Doebley et al. constructed another F 2 population of maize race reventador (Nay 15) × Z. mays ssp. pamiglumis (Iltis Cochrane 81), and detected five QTL for PEDS on chromosomes 1, 2, 3, and 10, with the R 2 of each QTL was 0.07-0.16 [13]. Moreover, Doebley et al. identified four QTL associated with PEDS in the testcross (TC1) population, which were located on chromosomes 1, 3, 9, and 10, explaining 3.0-7.4% of the phenotypic variance [14]. The QTL (QTL-1 L and QTL-3 L) located on chromosomes 1 L and 3 L were repeatedly identified in different studies, explaining an average of 26.85 and 32.70% of the phenotypic variance, respectively [7,12,13]. In addition, a significant epistatic interaction between QTL-lL (umc107) and QTL-3 L (umc92) that increased the expression of PEDS was identified in the F 2 population, and the combined effect of these two QTL was 60% [12,15]. Nevertheless, the combined value was 7.3% when the maize alleles of these two QTL were introduced into teosinte background, which was much less than that in the F 2 population, suggesting that a number of epistatic interactions may be involved in the PEDS variation [15].
Traditional QTL mapping typically involves genotyping of dozens to hundreds of individuals in segregating populations, which is time and labor intensive and thus expensive [16,17]. In order to overcome these disadvantages, the selective genotyping was proposed to enhance the efficiency of QTL mapping through selectively genotyping individuals with extreme phenotypes [18]. Further, the procedures were simplified by bulked sample analysis (BSA), only genotyping four DNA pools, including a paired bulked DNA pools with opposite or extreme phenotypes and a pair of parents [19]. The rapid development of next-generation sequencing (NGS) technologies accelerated the application of BSA in QTL analysis and fine mapping because it overcame the limitations of marker density [20]. Combining with the NGS and BSA, several methods have been developed for identification of QTL/genes for quantitative traits, including SHOREmap [21], X-QTL [22], NGM [23], MutMap [24], QTL-seq [25,26], and QTG-Seq [27]. In maize, these new methods have been utilized for detecting candidate genomic regions, QTL and genes for both quality and quantitative traits such as plant height [28], maize lethal necrosis [29], resistance to gibberella stalk rot [30], and fertility restoration of maize CMS-C [31].
In this study, a two-ranked maize (SICAU1212) and teosinte (Z. mays ssp. mexicana) were used as parents to develop the advanced mapping populations (BC 3 F 2 and BC 4 F 2 ) using phenotypic recurrent selection. Firstly, QTL-seq was performed to rapidly identify candidate genomic regions associated with PEDS in five pairs of high PEDS bulks and the low PEDS bulk in the BC 3 F 2 populations. Subsequently, traditional QTL mapping was carried out to validate the genomic regions detected by QTL-seq and identify the epistatic interactions for PEDS in the BC 3 F 2 and BC 4 F 2 populations. Results from the study will enhance our understanding of the genetic basis of PEDS.

Phenotypic diversity in the BC 3 F 2 populations
The ear of maize parent SICAU1212 has four rows (two ranks, two spikelets per cupule), whereas teosinte MT1 has two rows (two ranks, one spikelet per cupule). Therefore, the only difference is paired spikelets in SICAU1212 ear versus single spikelet in MT1 ear. The phenotypic variability of PEDS was high, ranging from PEDS = 0 (same as ear of SICAU1212) to PEDS = 100% (same as ear of MT1) in the BC 3 F 2 populations (Table 1; Fig. 1 (Table 1). Notably, the frequency distributions were continuous variation (except in the 15JH environment) and did not follow normal distribution. Instead, they were dramatically skewed towards the maize phenotype (PEDS = 0%) ( Table 1), suggesting that PEDS may be a quantitatively inherited trait.

Whole-genome resequencing, mapping of reads and identification of SNPs
The six bulks and SICAU1212 were sequenced, and the whole-genome resequencing (WGRS) data were generated.

Candidate genomic regions associated with PEDS
To identify candidate genomic regions associated with PEDS, the SNP-index was estimated individually in each bulk and Δ (SNP-index) was calculated, and the corresponding graphs were plotted against the genome positions (Fig. S2). Generally speaking, the SNP-index graphs between high and low bulks would be identical for the genomic regions that are not relevant to the phenotypic difference, while they should exhibit unequal contributions from the two parental genomes if the genomic region harboring QTL contributed to the difference in the phenotype [25,32]. Based on the graphs of SNP-index of bulks in Fig. 2 Table 2). In all three regions, the alleles from teosinte parent MT1 increased PEDS value, whereas those from maize parent SICAU1212 decreased PEDS value. However, no significant candidate genome region was detected in other bulks tests.

Validation of the genomic regions
To validate the genomic regions detected by QTL-seq, the traditional QTL mapping was performed in the BC 3 F 2 population consisting of 280 and 333 plants from the 14JH and 18WJ environments, and the BC 4 F 2 population containing 651 plants from the 15WJ environment. In these populations, the frequency distributions did not follow a normal distribution, but were dramatically skewed to maize parent SICAU1212 (Fig. S3). Firstly, 74 markers (1 SSR and 73 InDel) were selected to check polymorphism between parents. Then, the polymorphic markers were utilized to genotype 36 F 2 plants. Finally, 31 polymorphic and no distorted segregation markers [including 8 markers on chromosome 1, 14 markers on chromosome 3, 4 markers on chromosome 6, and 5 markers on chromosome 8 (Table S2)] were used to genotype the BC 3 F 2 and BC 4 F 2 populations. The local genetic linkage maps were constructed for the candidate genomic regions on chromosomes 1, 3, 6, and 8 ( Fig. 3).
Using the local genetic linkage maps and the phenotyping data, two and two QTL related to PEDS were identified in the BC 3 F 2 population in 14WJ and 18WJ environments, and five QTL in the BC 4 F 2 population in 15WJ environment, located on chromosomes 1, 3, 6 and 8 (Table 3; Fig. 3). The phenotypic variance explained (PVE) of individual QTL ranged from 1.12 to 38.05%. Notably, the QTL qPEDS3.1 (a major and stably expressed QTL) explained 33.91-38.05% of the phenotypic variance identified in all three environments. The additive effect of qPEDS3.1 was negative, suggesting that the alleles from teosinte parent MT1 increased the value of PEDS. The consistent QTL qPEDS1.1 and qPEDS6.1 accounted for 1.12-7.15% of the phenotypic variance in two environments. The QTL qPEDS1.2 and qPEDS8.1 with PVE = 1.12-1.82% were only mapped in the BC 4 F 2 population in 15WJ environment. In summary, the QTL mapping results supported that there were QTL for PEDS in the candidate genomic regions on chromosomes 1, 3, 6 and 8 identified by QTL-seq.
Five, seven and five pairs of epistatic QTL were detected in the BC 3 F 2 population in the 14WJ and 18WJ environments, and the BC 4 F 2 population in the 15WJ environment, respectively (Table 4; Fig. 4). The PVE of individual epistatic QTL ranged from 3.85 to 21.31%; the total PVE of all the epistatic QTL varied from 45.81 to 66.81%, being larger than those of all the QTL mentioned above in each environment. However, there were two kinds of epistatic QTL [33]. One was that additive by additive effect was positive, increasing the expression of PEDS and explaining 56.85, 34.57 and 18.71% of phenotypic variation in each test, whereas the other was negative, decreasing the expression of PEDS and explaining 9.96, 11.23 and 28.85% of phenotypic variation.

Investigation of PEDS in the near-isogenic lines (NILs)
For each QTL related to PEDS identified by the traditional QTL mapping, NILs were developed through

Discussion
The two-ranked maize parent SICAU1212 was helpful for accurate identification of PEDS Single vs. paired spikelets (PEDS) is one of the key domesticated traits. The single spikelets of teosinte ears were totally transformed into paired spikelets of maize ears during maize domestication. However, little is known about the genetic basis of this transformation of single into paired spikelets. Understanding the genetic basis of PEDS is helpful for the map-based gene cloning. In the present study, a primitive four rows (two ranks) waxy maize inbred line (SICAU1212) was used as a female parent to develop the advanced mapping populations (BC 3 F 2 and BC 4 F 2 ). SICAU1212 is similar to the primitive maize of 6000 years ago in phenotypic traits, including small plant, multiple ears, rachis from the ears, 4-row ears and tens of grains [34,35]. Notably, both the SICAU1212 and teosinte MT1 ears were two ranks, and the ears of all plants in the segregating population of SICAU1212 × MT1 showed only two ranks, thus avoiding a possible confusion about the importance of the single spikelet vs. a decrease in the rank number in multi-rank ears. Therefore, the phenotype of PEDS was investigated more accurately, especially at the silking stage, which was essential for the QTL mapping.

Comparison of QTL for PEDS in this study with those in previous studies
Four candidate genomic regions associated with PEDS located on chromosomes 1, 3, 6, and 8 were identified using QTL-seq, followed by detecting the QTL for PEDS in each candidate genomic region by traditional QTL mapping. The QTL qPEDS1.1 identified in this study overlapped with the QTL flanked by umc11-umc83 and umc37b mapped by Doebley et al. [12,13], because they were all located at 184-258 Mb of chromosome 1. The QTL (QTL-1 L) accounting for ~ 19.5% of phenotypic variance was identified repeatedly in the previous studies [7,13,15], near to the gene teosinte branched1 (tb1) [15,36], which was included in the QTL qPEDS1.2. According to the results of QTL-seq, the region of 0-173 Mb on chromosome 3 was associated with PEDS variation; however, only one QTL qPEDS3.1 explaining 33.91-38.05% of phenotypic variance was identified in the BC 3 F 2 and BC 4 F 2 populations by traditional QTL mapping, possibly near to QTL (umc121-umc92 and umc32) mapped by Doebley et al. [7,14]. The QTL qPEDS3.1 may be a novel major and stably expressed QTL, and it is worthy of further study. In addition, QTL (QTL-3 L) was identified repeatedly in the previous studies [12,13,15], located at ~ 181 Mb on chromosome 3. However, no QTL for PEDS was identified near this region, probably the InDel markers used in this study did not contain this region, with the farthest InDel marker (PM16) being located at ~ 157 Mb on chromosome 3. It is interesting that the region of 0-173 Mb on chromosome 3 was completely retained in the BC 3 F 2 population developed by phenotypic recurrent selection for PEDS, inferring that several QTL related to PEDS located in this region may not be identified. Although one QTL (umc85-umc65) for PEDS located on chromosome 6 was mapped [12], we do not know how far away it is from the minor QTL qPEDS6.1 that explained 1.12-3.26% of phenotypic variance because the physical positions of markers umc85 and umc65 could not be found. The QTL qPEDS8.1 on chromosome  [9].

Epistasis played a vital role in PEDS variation
Five to seven epistatic QTL with the total PVE of 47.57-66.81% were identified in each test, which was more than the effect of main QTL (35.02-47.12%). Deobley at el. (1993,1995) identified significant epistasis between QTL-1 L (umc107) and QTL-3 L (umc92); the effect of each QTL was much lower in the advanced background population than the F 2 population, with the combined effect of two QTL decreasing about 10-fold [13,15]. No significant epistatic QTL between QTL-1 L and QTL-3 L was detected in this study; comparatively, one locus of EPqpeds-1 was included the region of QTL-1 L, one locus of EPqpeds-6 was close to QTL-1 L, and one locus of EPqpeds-11 was near QTL-3 L. Notably, one locus of epistatic QTL EPqpeds-2, EPqpeds-3, EPqpeds-4, EPqpeds-5, EPqpeds-7, and EPqpeds-13 overlapped with, or was close to, the major QTL qPEDS3.1 (Table 5); another locus was located on chromosomes 1, 3, 6, and 8, suggesting that the epistatic interactions between qPEDS3.1 and other loci were very strong, partly explaining the very low ratio of PEDS > 0 in the NIL-qPEDS3.1(MT1) families as well as in the other NILs of QTL. Moreover, 50 pairs of epistatic QTL were identified on 10 chromosomes, explaining 94.40% of the total phenotypic variance in the F 2 population of SICAU1212 × MT1 (unpublished data). The number and effect of epistatic QTL identified in the BC 3 F 2 and BC 4 F 2 populations decreased significantly, inferring that the advanced mapping population developed by phenotypic recurrent selection for PEDS eliminated partly the effect of epistasis. Taken together, the epistasis played a key role in PEDS variation. Unfortunately, the impact of epistasis was neglected in fine mapping or positional cloning of QTL/gene for complex traits such as PEDS [37]. For fine mapping of the major QTL qPEDS3.1, the two, three and four QTL will be combined to evaluate the effect; then, the residual heterozygous lines of qPEDS3.1 will be developed in future research by selecting the heterozygous alleles in the region of qPEDS3.1 and the homozygous alleles in other regions of QTL [38,39].

Photoperiod significantly affect the expression of single spikelet
The photoperiod-related trait is a primary domestication trait in maize. The geographic range of maize has rapidly expanded from tropical to temperate regions due to the loss of photoperiod sensitivity [40,41]. The effects of photoperiod on tassel branches and the number of anthers were reported by Bechoux et al. [42]. Boden et al. suggested that Photoperiod-1 (Ppd-1) responsible for flowering time has a significant inhibitory effect on paired spikelet formation [43]. In the present study, the expression of PEDS decreased significantly, and no genomic region for PEDS was identified in the 15JH environment using QTL-seq (Fig. S2E). Compared to the other three environments, the day length was 11.20 h in the seedling stage and 11.73 h 1 month later in the 15JH environment, which was shorter than 12.53 h in the 14WJ environment, 14.97 h in the 14EEDS environment, and 12.85 h in the 14JH environment. The QTL qPEDS8.1 was mapped on chromosome 8, which contained the gene ZCN8 involved in photoperiod sensitivity of tropical maize [44]. Investigation of NIL-qPEDS8.1(MT1) planted at Jinghong city, Yunnan province, in September 2015, revealed that the teosinte allele had a shorter flowering time (by about 7 days) than that of SICAU1212. However, the presence of single spikelet was not confirmed. In

Conclusion
In the present study, to accurately investigate single and paired spikelets, a two-ranked maize were chosen to cross to teosinte to develop the advanced mapping populations (BC3F2 and BC4F2). Four candidate genomic regions associated with PEDS were identified on chromosomes 1, 3, 6 and 8 using QTL-seq. One major QTL, four minor QTL and 14 epistatic QTL for PEDS were detected in these four genomic regions by linkage map analysis. The expression of single spikelet was very low (PEDS < 10%) in the 15JH environment, maybe due to the shorter day length compared to that in the other environments. The flowering time of plants of NIL-qPEDS8.1(MT1) was ~ 7 days shorter than that of "parent SICUAN1212". Taken together, these findings suggested that major QTL, minor QTL, epistasis, and photoperiod were associated with the transformation of single spikelet in teosinte ears into paired spikelets in maize ears during maize domestication, providing insights into the genetic basis of PEDS and an opportunity to fine map major QTL for PEDS in the future work.

Plant materials and development of the advanced mapping populations
The maize inbred line SICAU1212 and teosinte MT1 (Z. mays ssp. mexicana) were used as parental lines to develop the advanced segregating populations for PEDS. SICAU1212 was derived from a waxy maize landrace Silunuo [45,46] and is owned by our lab; it has two ranks with two spikelets per cupule. MT1 was obtained from the International Maize and Wheat Improvement Center (CIMMYT) (accession number 1 l394); it has two ranks with one spikelet per cupule. A cross was made between SICAU1212 and MT1 to create F 1 in April 2010; F 1 was self-pollinated to develop the F 2 population in October 2010. The plants with PEDS = 100% (exclusively singe spikelet) were selected from the F 2 population to backcross to SICAU1212 to create BC 1 F 1 in April 2011; then, BC 1 F 1 was self-pollinated to develop the BC 1 F 2 population in October 2011. Again, the plants with PEDS = 100% were chosen from the BC 1 F 2 to backcross to SICAU1212 to create BC 2 F 1 in April 2012; BC 2 F 1 was self-pollinated to develop the BC 2 F 2 population in October 2012. The procedure mentioned above had been performed repeatedly, until the BC 4 F 2 population was developed in October 2014 (Fig. S1). The presence (paired spikelets)/absence (single spikelet) of the pedicellate spikelet in each cupule can vary among cupules within a single ear; therefore, PEDS was recorded as percentage of cupules without a pedicellate spikelet [7,12]. For ease and accuracy, single and paired spikelets were counted on the basal-most secondary lateral inflorescence when the plant had silks visible. The cubic root transformations of phenotype data were carried out to reduce skewness and kurtosis because the PEDS distribution of BC 3 F 2 populations did not follow a normal distribution [12]. The descriptive statistics of the populations were analyzed using SPSS19.0 software (http:// www. spss. com).

Construction of bulks
Extreme bulks were prepared for PEDS based on precise phenotyping. Six extreme bulks were constructed as

Construction of sequencing libraries and Illumina sequencing
Genomic DNAs were extracted from the fresh leaves of individuals from the six bulks mentioned above and parent SICAU1212 by using the CTAB method [47]. DNA quality was evaluated using 1% agarose gel, and DNA concentration was quantified using a NanoDrop ® . The equimolar concentrations of DNA from individuals of each bulk were pooled together. To construct a library, two micrograms of DNA from each sample were first sheared using diagenode Bioruptor ® NGS (Diogenode, Liege, Belgium) and then were subjected to end repair and adapter ligation. The size selection of libraries was performed using 2% agarose gel to get a target insert size of 500-600 bp. The libraries were purified first and then enriched using the adaptor-compatible PCR primers. The size distribution of amplified DNA libraries was checked on an Agilent Technologies 2100 Bioanalyzer using a High Sensitivity chip. The DNA libraries were sequenced on Illumina HiSeqTM 2000 to generate 150-base pairend reads [48][49][50].

Calculation of SNP-index and Δ (SNP-index)
After generating the sequences of all seven samples, the statistical analysis of the sequencing reads was conducted using the raspberry tool of NGS-QCbox [51]. The maize reference genome (B73-RefGen_v4) was downloaded from MaizeGDB (http:// www. maize gdb. org/). The cleaned reads of each sample were aligned to B73-RefGen_v4 genome using the inbuilt BWA aligner [52], and the Coval software was used for post-processing and filtering of the alignment files [53]. The SNP-index as the proportion of short reads harboring SNPs that are different from the reference sequence [25,54], was calculated for each bulk using the formula [24]: SNP-index (at a position) = (count of alternate base) / (count of reads aligned). To improve accuracy of SNPs, the SNP positions with read depth < 4× or > 32× were filtered out because these SNPs may be false positives caused by sequencing, alignment errors, or genomic repeat sequences [25,32]. An average of SNP-indices of SNPs located in a given genomic interval was calculated using the sliding window analysis with a 1 Mb window size and 10 kb increment [25,54]. In this study, SNP-index = 1 if all the short reads came from the parent MT1, and SNP-index =0 if all the short reads came from the parent SICAU1212. Δ (SNPindex) was then calculated by subtracting SNP-index of low bulk from SNP index of high bulk. The SNP-index graphs for high bulks and low bulk, and corresponding Δ (SNP-index) Manhattan graphs were plotted against the position of each sliding window in the reference genome (B73-RefGen_v4). The Δ (SNP-index) value of candidate genomic regions associated with a target trait should be significantly different from 0; therefore, candidate genomic regions for PEDS were identified using the hypothesis proposed by Takagi et al. [25].

Validation of QTL-seq-derived genomic regions by traditional QTL mapping
The candidate genomic regions for PEDS identified by QTL-seq were validated using traditional QTL analysis in three environments. Briefly, InDel markers were identified by aligning high bulks and low bulk, and the primers were designed based on the methods described in the previous study [55]. The InDel markers located on chromosomes 1, 3, 6, and 8 were selected for polymorphism screening between the parents SICAU1212 and MT1. Then, the polymorphic markers were used in genotyping the F 2 population with 36 plants derived from a cross between SICAU1212 and MT1 to calculate the allele ratios, excluding the markers with serious segregation distortion as determined by chi-square analysis. Subsequently, the markers that followed the expected Mendelian segregation ratio were utilized for genotyping the BC 3 F 2 populations in the 14JH and 18WJ environments, and the BC 4 F 2 population in the 15WJ environment. The genotyping data were used to construct the local genetic linkage map using MAPMAKER/EXP 3.0 [56], and the genetic map distances (cM) were calculated following Kosambi mapping function [57]. Combining the genetic linkage information and phenotyping data, QTL analysis for PEDS was conducted with QTL IciMapping 4.1 (http:// www. isbre eding. net/) using the inclusive composite interval mapping of additive (ICIM-ADD) and the two locus epistatic QTL (ICIM-EPI) modules [33].