Characterization of a major QTL and Genome-Wide Epistatic Interactions for the Transformation of Single Spikelet in Teosinte Ears into Paired Spikelets in Maize Ears During Maize Domestication

Maize ear carries paired spikelets, whereas the ear of its wild ancestor, teosinte, bears single spikelets. However, little is known about the genetic basis of the processes of transformation of single spikelets in teosinte ear to paired spikelets in maize ear. In this study, a two-ranked, paired-spikelets primitive maize and a two-ranked, single-spikelet teosinte were utilized to develop an F 2 population, and QTL mapping for single vs. paired spikelets (PEDS) was performed. Two QTL (qPEDS1.1 and qPEDS3.1) for PEDS located on chromosomes 1L and 3S were identied in the 162 F 2 plants using the inclusive composite interval mapping of additive (ICIM-ADD) module, explaining 1.93% and 23.79% of the phenotypic variance, respectively. Out of the 409 F 2 plants, 43 plants with PEDS = 0% and 43 plants with PEDS > 20% were selected for selective genotyping; the QTL (qPEDS3.1) accounting for 64.01% of the phenotypic variance for PEDS was also detected. Moreover, the QTL (qPEDS3.1) was validated in three environments, which explained 31.05%, 38.94% and 23.16% of the phenotypic variance, respectively. In addition, 50 epistatic QTLs were detected in 162 F 2 plants using the two-locus epistatic QTL (ICIM-EPI) module; they were distributed on all 10 chromosomes and explained 94.40% of the total phenotypic variance. The results contribute to a better understanding of the genetic basis of domestication of paired spikelets and provide a genetic resource for future map-based cloning; in addition, the systematic dissection of epistatic interactions underlies a theoretical framework for overcoming epistatic effects on QTL ne mapping.


Introduction
During maize domestication, the main goals were to improve the yield and ease of harvesting (Wang et al. 2005;Doebley 2004). The morphology of maize ear and that of its ancestral teosinte exhibit extreme differences, with taxonomists classifying them in different genera initially (Doebley 2004). There are four major traits distinguishing teosinte ear from maize ear (Doebley 1992(Doebley , 2004. First, compared to the teosinte kernels protected by a hard outer glume, the maize kernels are naked, enabling the kernels to be easily harvested and digested. Second, the teosinte ear has two ranks with one rank on each side, whereas the maize ear has two or more ranks, resulting in maize ear bearing hundreds of kernels. Third, teosinte has a shattering rachis and maize has a non-shattering one; loss of shattering facilitates easy collection of kernels, and has been selected for in domestication of all cereals. Fourth, each cupulate fruitcase of the teosinte ear has single mature spikelet (one kernel per spikelet), whereas paired mature spikelets are in each cupulate fruitcase of the maize ear, which could double the kernel number. In addition to the single vs. paired spikelets, the key genes affecting the other three traits have been cloned, including Teosinte glume architecture 1 (tga1) involved in the origin of naked kernels (Dorweiler et al. 1993;Wang et al. 2005); ZmSh1-1, ZmSh1-5.1 and ZmSh1-5.2 related to seed shattering (Lin et al. 2012); and FLORICAULA/LEAFY homologs z 1 and z 2, unbranched2 (ub2) and unbranched3 (ub3) and its regulatory locus KRN4 responsible for rank number of ear determining number of kernel rows (Bomblies et al. 2003;Bomblies and Doebley 2006).
Observing the development of the female in orescence by scanning electron microscopy, it was found that the transformation of the spikelet pair meristems (SPM) to the spikelet meristems (SM) was the crucial step involved in single vs. paired spikelets. In maize, the SPM produced two morphologically distinct spikelet meristems, one sessile and one pedicellate. Both of them developed into oral meristems (FM) that ultimately produced the oral organs (Upadyayula et al. 2006). In contrast, the pedicellate SM was aborted in teosinte (Doebley et al. 1995b). In the past 100 years, few studies were conducted to dissect the genetic basis of single vs. paired spikelets through analysis of Mendelian inheritance and linkages with a known gene or molecular marker loci. Colins and Kempton (1920) rst reported that a ratio of paired:single spikelets approaching 3:1 was obtained in teosinte-maize hybrid populations, but the inheritance of single vs. paired spikelets deviated from the Mendelian pattern due to continuous phenotypic variation. Conversely, Langham (1940) inferred that single spikelet (named pd) was controlled by a single Mendelian factor. Mangelsdorf (1947) and Szabó et al. (1996) suggested that the segregation of single vs. paired spikelets tted a model with two independent genes. Nevertheless, Rogers (1950) found that unifactorial inheritance did not control the single vs. paired spikelets, and the genetic control of the trait was more complex. Mangelsdorf (1947) showed the molecular marker loci on chromosomes 1L, 2S, 3L, and 4S link with single vs. paired spikelets. Subsequently, with the development of molecular marker technology, a RFLP marker was used for mapping QTL for single vs. paired spikelets. A total of 18 QTL for single vs. paired spikelets were identi ed, being located on all 10 chromosomes. The phenotypic variation explained by individual QTL ranged from 4.9-46.1% (Doebley et al. 1990;Stec 1991, 1993;Lauter and Doebley 2002b). Two major QTL located on chromosomes 1L (QTL-1L) and 3L (QTL-3L) were detected in three different genetic backgrounds and the effect of these two QTL was con rmed by analyzing of near-isogenic lines (NILs) (Doebley et al. 1995a).
Epistasis is used to denote the interaction between alleles from different loci and plays an important role in quantitative genetics analysis (Mackay 2014;Vazquez et al. 2015;Shang et al. 2016). Recently, epistasis was reliably detected to have an effect on complex trait variation, including disease resistance (Roncallo et al. 2012;Singh. et al, 2013;Vazquez et al. 2015) and yield-related traits Soyk et al. 2017;Sundaram et al. 2018). Based on several experiments in which the total phenotypic variance of epistasis ranged from 16-81%, Carlborg and Haley argued that epistasis should be accounted for in order to understand the genetic basis of complex traits (Carlborg and Haley 2004). In previous studies, the signi cant epistasis between loci umc107 (QTL-1L) and umc92 (QTL-3L) for single vs. paired spikelets was detected in the maize × teosinte (Zea mays ssp. mexicana) F 2 population (Doebley and Stec 1991). In the maize × teosinte (Zea mays ssp. parviglumis) F 2 population, the proportion of paired spikelets increased from 0-60% due to the combined effect of these two QTL (Doebley et al. 1995a). Then, the near-isogenic lines (NILs) of QTL-1L and QTL-3L were constructed in the maize and teosinte genetic backgrounds, respectively. The mean trait values of maize allele of QTL-1L and QTL-3L in the teosinte genetic background were 0.8% and 1.9%, respectively; the value was 7.3% when these two QTL were combined, which was 4.6% higher than the sum of two individual QTL. In addition, compared with the combined effect of these two QTL in the F 2 population (60%), about 52.7% reduction was observed in the teosinte genetic background, suggesting that the higher-order epistatic interactions may be responsible for PEDS (Doebley et al. 1995a).
In this study, a two-ranked (four rows) primitive maize was chosen to develop an F 2 population by crossing it with two-ranked, two-rows teosinte (Z. mays ssp. mexicana). The oligonucleotides pool assay (OPA) with 3072 well-distributed and high-quality SNPs was used for genotyping F 2 plants. QTL mapping for PEDS was performed, identifying one minor QTL explaining 1.93% of the phenotypic variance, one major QTL explaining 23.79% of the phenotypic variance and 50 epistatic QTL explaining 94.40% of the total phenotypic variance. The major QTL was validated in three environments. These results suggested the importance of not only major and minor QTL, but the genomewide epistasis as well, in in uencing the variation of PEDS.

Plant materials
An F 2 population was developed through a cross between waxy maize inbred line SICAU1212 and teosinte MT1 (Z. mays ssp. mexicana).
SICAU1212 ear has two ranks with two spikelets per cupule (paired spikelets), whereas MT1 ear has two ranks with one spikelet per cupule (single spikelet), thus excluding the effect of number of ranks, which is essential for investigation of single vs. paired spikelets. SICAU1212 was derived from a waxy maize landrace Silunuo by 10 consecutive generations of self-pollination, and had multiple primordial traits of maize, including small stature, rami cation, stooling, narrow leaves, and small spike. The history of Silunuo cultivation in Yunnan Province of China can be traced back to at least 1890 (Zeng and Pu 1981;Tian et al. 2009).

Field trials and trait measurements
The PEDS was de ned by Doebley et al. (1990) as the percentage of cupules lacking the pediculate spikelet, and we followed it. In order to easily and accurately investigate the phenotype, the single and paired spikelets on the basal-most secondary lateral in orescence were counted as soon as the silks of plant were visible, avoiding the interference of oral organs with developmental failure. The phenotypic data were cubic root transformed to reduce skewness and kurtosis (Doebley et al. 1990). Basic statistical analysis of PEDS was performed using SPSS19.0 software (http://www.spss.com).
Genotyping and construction of genetic linkage map Genomic DNA was extracted from young leaves of F 2 population (162 randomly selected plants from the 409-plant population, and 43 plants with PEDS >20% selected from the rest of population), their F 1 and parents in the 11YA environment using the cetyl trimethyl ammonium bromide (CTAB) method (Saghai-Maroof et al. 1984). The oligonucleotides pool assay (OPA) consisting of 3072 welldistributed and high-quality SNPs was used to genotype the plants mentioned above using the procedure described by Fan et al (2006) and Hou et al (2015). The SNP markers showed polymorphism between the parents SICAU1212 and MT1 and heterozygotes allele in the 162 F 2 plants, excluding obviously segregation distortion markers in SICAU1212 or MT1 allele less than 21, or heterozygotes allele less than 41.
The chi-square analysis was conducted for the segregation ratio of the remaining SNP markers. The genetic linkage map was constructed by using MAPMAKER/EXP 3.0 with LOD threshold >3.0, and the genetic distances were calculated by the Kosambi mapping function (Lander et al. 1987;Kosambi 2016). QTL analysis QTL for PEDS were identi ed using two approaches: one was QTL mapping in 162 randomly selected plants, and another was selective genotyping in 86 plants with extreme traits. Combined with the SNP-based linkage genetic map and phenotyping of 162 F 2 plants, the main QTL and epistatic effects for PEDS were identi ed in the QTL IciMapping 4.1 (http://www.isbreeding.net/) using the inclusive composite interval mapping of additive (ICIM-ADD) and two epistatic QTL (ICIM-EPI) modules . Before QTL mapping, the missing phenotypes were removed; to detect additive QTL, walking speed was set at 1.0 cM, probability in the stepwise regression was 0.001, and threshold LOD scores were calculated using 1000 permutations with a type 1 error of 0.05. Epistatic QTL were identi ed by using walking speed of 5 cM, probability of 0.0001 in stepwise regression, and threshold LOD of 5.0.
A subset of 86 plants, namely 43 plants with PEDS = 0 chosen from the 162 F 2 plants and 43 plants with PEDS >20% selected from 409 plants representing the most extreme phenotypes, were used for selective genotyping to identify QTL for PEDS. The detection of QTL for PEDS was performed in the QTL IciMapping 4.1 using the selective genotyping mapping (SGM) and the inclusive composite interval mapping of additive (ICIM-ADD) modules. The parameters of threshold LOD were the same as that for additive QTL detection described above. QTL designations followed the standard nomenclature of McCouch et al. (1997). For example, for QTL qPEDS3.1, q represents a QTL, peds is the abbreviation for "Percentage of cupules lacking the pedicellate spikelet", 3 indicates that the QTL was located on the chromosome 3, and 1 is the serial number of that QTL; the QTL had the same serial number when mapped within the same marker interval or when sharing a common marker.

Validation of QTL regions
The two QTL on chromosomes 1L and 3S identi ed in the 11YA environment were validated in three environments. First, a total of 36 InDel markers  located in the QTL regions were chosen to detect polymorphism between SICAU1212 and MT1. The procedure of PCR ampli cation was conducted as described in our previous reports (CHEN et al. 2017), and denatured ampli ed products were separated on 6% polyacrylamide gels and visualized by silver staining (Sanguinetti et al. 1994). Then, the co-dominant and non-segregation distortion markers were utilized to construct the local genetic linkage map. With the phenotypic data, QTL mapping was performed following the procedures described above.

Phenotypic variation in the trait
The ear of maize SICAU1212 is about 12 cm, presenting two ranks with 21-24 cupules per rank and two spikelets per cupule (paired spikelets); in contrast, the ear of teosinte MT1 is about 7 cm, presenting two ranks with six cupules per rank and one spikelet per cupule ( Fig. 1). The differences in PEDS between SICAU1212 and MT1 were distinct, namely PEDS = 0% in SICAU1212 versus PEDS = 100% in MT1; the PEDS of F 1 plants equaled 0%, indicating that maize was completely dominant to teosinte (Table 1). The mean values of PEDS were 7.63%, 8.77%, 8.32%, and 13.35% in the F 2 populations, with a range from 0-100%, in the four environments, respectively. The ratio of individuals with PEDS > 0% in the F 2 population ranged from 1/6 (17WJ) to 1/4 (17JH).In addition, at least one plant with all single spikelets (PEDS = 100%) was observed in each environment. However, the phenotypes of F 2 populations did not follow a normal distribution, which was instead signi cantly skewed toward the maize phenotype (Fig. 2). with the ratio of 27.31%. To exclude the obviously segregation distortion markers, 586 SNPs were used to construct the genetic linkage map. The linkage map covered all 10 maize chromosomes, and the number of markers ranged from 36 on chromosome 10 to 81 on chromosome 8. It spanned a total length of 1758.29 cM, with an average marker interval of 3.00 cM; the biggest linkage distance between two adjacent markers was 32.12 cM, located on chromosome 3 between PZE-103118170 and PZE-103150482 (Fig. 3).

Main And Epistatic Effect Qtl For Peds
Taking advantage of the SNP-based linkage map and precise phenotypic data, two QTL for PEDS were identi ed (Table 3, Fig. 3) in the F 2 population with 162 randomly selected plants in the 11YA environment; they were located on chromosomes 1L (qPEDS1.1) and 3S (qPEDS3.1), were anked by markers PZE-101196838-PUT-163a-76013247-3735 and PZE-103018221-SYN28119, and accounted for 1.93% and 23.79% of the phenotypic variance, respectively. In addition, the SNP marker with the highest LOD score of 4.485 was located at SYN28119, as identi ed using the SGM module in the population of 86 extreme phenotypes. Moreover, the major QTL (qPEDS3.1) was also detected in the same population using the ICIM-ADD module; it was located between markers PZE-103018221 and SYN28119, and explained 64.01% of the phenotypic variance. The additive effect of all QTL was negative, indicating that the alleles that increased the expression of PEDS were from MT1. The qPEDS1.1 showed the additive effect due to the dominance/additive ratio being 0.17, whereas qPEDS3.1 showed the dominance effect due to the dominance/additive ratio being 0.84 .  Fig. 4). The total phenotypic variance explained (PVE) by these 50 epistatic QTL was 94.40%, and the phenotypic variance explained by each epistatic QTL ranged from 0.68-4.12%. However, these 50 epistatic QTL were divided into two categories based on the additive-by-additive effect . One category, including 29 epistatic QTL, had a negative effect that decreased the expression of PEDS, and accounted for 55.59% of the phenotypic variance; another category consisting of 21 epistatic QTL had a positive effect that increased the expression of PEDS, and accounted for 38.81% of the phenotypic variance. Among them, ve epistatic QTL (EPqpeds-6 to EPqpeds-10) were located in the region of qPEDS1.1, and nine epistatic QTL (EPqpeds-1, EPqpeds-13, and EPqpeds-17 to EPqpeds-23) were in the region of qPEDS3.1.

Validation Of The Main Qtl
To validate the main QTL (qPEDS1.1 and qPEDS3.1), 11 polymorphic InDel markers were selected from 36 InDel markers, and were utilized to construct the local genetic linkage map in the qPEDS1.1 and qPEDS3.1 regions. No signi cant QTL was detected in the qPEDS1.1 region, whereas one major QTL (named qPEDS3.1) was identi ed in the qPEDS3.1 region in the three environments (Table 3, Fig. 5). The QTL was anked by markers chr3-14843386 and chr3-19759100, accounting for 31.05% and 38.94% of phenotypic variance in the 17WJ and 17YA environments, respectively. In addition, the QTL located between chr3-11534053 and chr3-14843386 explained 23.16% of phenotypic variance in the 17YA environment. All QTLs had a negative additive effect, indicating that alleles from MT1 contributed predominantly to the phenotype.

Discussion
The advantage of four-row waxy maize landrace with relatively primitive morphology in studying PEDS In order to dissect the genetic basis of single vs. paired spikelets (PEDS), the two-ranked (four rows) maize SICAU1212 (derived from fourrow waxy maize landrace Silunuo) was crossed to the two-ranked (two rows) teosinte MT1 to develop an F 2 population. The SNP arrary was used for genotyping of the F 2 plants, resulted in one major, one minor and 50 epistatic QTL for PEDS were identi ed. There are several advantages in using SICAU1212 as a parent in the study of the origin of paired spikelets during maize domestication. Firstly, the ears of both SICAU1212 and teosinte MT1 had two ranks, avoiding the confusion potentially arising from the segregation of ranks (rows) in a multi-ranked maize-teosinte F 2 population; notably, the fact all plants of the F 2 population of SICAU1212 × MT1 had only two ranks meant the only difference was paired spikelets in SICAU1212 ear versus single spikelets in MT1, making it easy to investigate the PEDS phenotype, as shown in Fig. 1. For instance, the SICAU1212 ear has an average ear length of 12 cm and approximately 80 kernels in four rows, whereas the earliest maize ear was 6 cm in length with 28 kernels in four to eight rows. According to the records, the maize landrace Silunuo (SICAU1212) has been cultivated in Yunnan Province, China, since 1890 (Zeng and Pu 1981). In addition, the cluster analysis of SICAU1212, 368 maize inbred lines and eight teosintes resulted in classifying SICAU1212 and eight teosintes into the same sub-group ( Figure S1, unpublished data). In other words, SICAU1212 is a relatively primitive maize in both phenotypic and genetic relationship analyses. The other advantage was that the differences in the genetic background between SICAU1212 and teosinte was small, which reduced the effect of complex genetic background on QTL mapping; the QTL for PEDS identi ed in this study were probably related to the transformation of spikelets from single to paired during maize domestication rather than due to maize improvement.
The Qtl (Qpeds3.1) Was A Novel Major Qtl In this study, two QTL for PEDS were identi ed, one major stable QTL (qPEDS3.1) located on the short arm of chromosome 3 (3S) and one minor QTL (qPEDS1.1) located on the long arm of chromosome 1 (1L). In previous studies, the highest frequency of QTL for PEDS was located on chromosomes 1L and 3L (Doebley et al. 1990;Stec 1991, 1993). The QTL qPEDS1.1 was consistent with the QTL-1L detected by Doebley et al., but the phenotypic variance (1.93%) of qPEDS1.1 was much lower than that of QTL-1L (~ 19.5%); the reason may be that the mapping populations and the method of QTL mapping were different. The QTL qPEDS3.1 was mapped on chromosome 3S, accounting for 23.16-38.94% of the phenotypic variance in the four environments tested in this study, whereas the QTL-3L located on chromosome 3L explained 13.0-46.1% of the phenotypic variance as reported by Doebley et al (Doebley and Stec 1993;Doebley et al. 1990;Doebley and Stec 1991). Similarly, one QTL explaining 5.8% of the phenotypic variance and located on chromosome 3S was also identi ed by Lauter and Doebley (2002a); that QTL was approximately 4 Mb away from the QTL qPEDS3.1 based on the physical position of markers, and further veri cation is needed to ascertain whether they are the same. Notably, the phenotypic variance associated with qPEDS3.1 was 64.01% (detected using selective genotyping), which was 2-fold higher than that of QTL detected in the 162 randomly selected F 2 plants.
One reason may be that the selective genotyping eliminated the effect of the high segregation distortion in the maize-teosinte population.
Another reason may be that the obstacles in QTL mapping arising from epistatic effects can be overcome using selective genotyping analysis ).
Using the high-density linkage map detected more segregation distortion regions potentially related to maize domestication Different types of molecular markers have been used for constructing the linkage maps in the maize × teosinte F 2 or BC populations, including RFLP, AFLP, SSR and SNP, with the number of markers ranging from 58 (RFLP) to 338 (AFLP) (Mano et al. 2005;Briggs et al. 2007;Jun-qing et al. 2011;Wang et al. 2012 (Mano et al. 2005) and a W22 × parviglumis F 2 population (1474.9 cM) (Briggs et al. 2007). In addition, the SNP array was also utilized to genotype the F 2:3 population derived from a maize × maize (R08 × Ye478) cross, and the total length of the linkage map was 2007.9 cM, which was longer than in this study (Hou et al. 2015). The linkage map constructed in this study was a relatively high-density linkage map, which was helpful in detection of QTL and the segregation distortion regions.

Epistasis Contributed Signi cantly To The Variation Of Peds
In this study, 50 epistatic QTL were identi ed, explaining 94.40% of the total phenotypic variance, which was much larger than that of the main effect QTL (~ 25.72%), suggesting that epistasis played an important role in the variation of PEDS as also described in previous studies (Doebley and Stec 1991;Doebley et al. 1995a). Doebley and Stec (1991) detected signi cant epistasis between loci umc107 and umc92 in the F 2 population derived from maize Sin2 × teosinte Doebley643, which increased the expression of PEDS; and epistasis was con rmed in the later study (Doebley et al. 1995a). In contrast, no signi cant epistasis was identi ed in the F 2 population derived from maize Nay15 × teosinte Iltis&Cochrane 81 Doebley et al. 1995a). We also did not detect epistasis between these two loci, but the locus umc107 was close to the marker PZE-101196838 that interacted with loci on chromosomes 2, 6, 7, and 10 (EPqpeds-6 to EPqpeds-9). In addition, the QTL anked by bnl5.46-umc42a and bnl8.32-umc151 detected by Doebley et al (Doebley et al. 1990;Doebley and Stec 1991) overlapped with the epistatic QTL EPqpeds-27 and EPqpeds-44 based on the physical position, respectively. However, no signi cant epistasis was identi ed in these two QTL regions in previous studies. The possible reason was that epistasis was affected by the genetic background (Doebley and Stec 1991). Epistasis for PEDS was detected in the teosinte genetic background and signi cantly increased the expression of PEDS, whereas PEDS were invariant in the maize genetic background (Doebley et al. 1995a;Lukens and Doebley 1999).
The main force of crop domestication probably came from mutations, recombination and genetic drift. In conditions where epistasis is common; mutation may be neutral or bene cial in one genetic background but deleterious in other genetic backgrounds (De Visser et al. 2011;Lehner 2011;Breen et al. 2012). In other words, epistasis could determine a consequence of mutation effect, resulting in different phenotypes of the trait (Wagner et al. 1997;Rice 1998;Nijhout 2002;Azevedo et al. 2006). Epistasis is regarded as a byproduct of the evolutionary process (Hurst 2000). In contrast, some researchers think it can determine the outcome of many evolutionary processes (Kondrashov 1988;Azevedo et al. 2006). The single spikelet in teosinte ear was compelely tranformated into paired spikelets in maize ear during maize domestication. The paired spikelets are dominant to single spikelet and are invariant in modern maize after stabilizing selection by breeders, suggesting that the process of selection led to the genetic robustness of paired spikelets as the bene cial trait (Waddington 1953;Azevedo et al. 2006). In the present study, 50 epistatic QTL for PEDS were identi ed, suggesting that epistasis may be a remnant of evolution that played an important role in selection for the genetic robustness of paired spikelets; therefore, PEDS have had a complex evolutionary process.