Genetic dissection of yield-related traits and mid-parent heterosis for those traits in maize (Zea mays L.)

Background Utilization of heterosis in maize could be critical in maize breeding for boosting grain yield. However, the genetic architecture of heterosis is not fully understood. To dissect the genetic basis of yield-related traits and heterosis in maize, 301 recombinant inbred lines derived from 08 to 641 × YE478 and 298 hybrids from the immortalized F2 (IF2) population were used to map quantitative trait loci (QTLs) for nine yield-related traits and mid-parent heterosis. Results We observed 156 QTLs, 28 pairs of loci with epistatic interaction, and 10 significant QTL × environment interactions in the inbred and hybrid mapping populations. The high heterosis in F1 and IF2 populations for kernel weight per ear (KWPE), ear weight per ear (EWPE), and kernel number per row (KNPR) matched the high percentages of QTLs (over 50%) for those traits exhibiting overdominance, whereas a notable predominance of loci with dominance effects (more than 70%) was observed for traits that show low heterosis such as cob weight per ear (CWPE), rate of kernel production (RKP), ear length (EL), ear diameter (ED), cob diameter, and row number (RN). The environmentally stable QTL qRKP3–2 was identified across two mapping populations, while qKWPE9, affecting the trait mean and the mid-parent heterosis (MPH) level, explained over 18% of phenotypic variations. Nine QTLs, qEWPE9–1, qEWPE10–1, qCWPE6, qEL8, qED2–2, qRN10–1, qKWPE9, qKWPE10–1, and qRKP4–3, accounted for over 10% of phenotypic variation. In addition, QTL mapping identified 95 QTLs that were gathered together and integrated into 33 QTL clusters on 10 chromosomes. Conclusions The results revealed that (1) the inheritance of yield-related traits and MPH in the heterotic pattern improved Reid (PA) × Tem-tropic I (PB) is trait-dependent; (2) a large proportion of loci showed dominance effects, whereas overdominance also contributed to MPH for KNPR, EWPE, and KWPE; (3) marker-assisted selection for markers at genomic regions 1.09–1.11, 2.04, 3.08–3.09, and 10.04–10.05 contributed to hybrid performance per se and heterosis and were repeatedly reported in previous studies using different heterotic patterns is recommended. Electronic supplementary material The online version of this article (10.1186/s12870-019-2009-2) contains supplementary material, which is available to authorized users.


Background
Utilization of heterosis in maize is of great importance for boosting grain yield [1][2][3][4]. Hybrids in maize accounted for 65% of total maize cultivation by the late twentieth century and had contributed to a quadrupling of annual maize production [1]. Numerous breeders have been interested in heterosis for many years; however, comprehension of the genetics implicated in heterosis for grain yield still remains elusive. Therefore, investigating and assessing genetic mechanisms of maize heterosis for grain yield would provide a paved route for understanding those phenomena and help to optimize breeding for grain yield in different heterotic groups.
Through earlier inbred selection and experimental breeding [5,6], maize breeders realized the importance of germplasm resources and heterosis. Previous studies [2,[7][8][9][10][11][12][13] used molecular markers, phenotypic identification, and pedigree information to divide different germplasm into different heterotic groups. The heterotic pattern Stiff Stalk Synthetic (SS) × Non-Stiff Stalk (NSS) has been widely used for maize breeding [12,[14][15][16]. Li and Wang [17], based on former studies, proposed five heterotic groups, namely Tangsipingtou (TSPT), Lancaster, Lancaster-like, the improved Reid (PA), and Tem-tropic I (PB). Extensive research has been focused on the TSPT × Reid heterotic pattern because it is one of the predominant heterotic patterns in northern China and has been used for the production of elite hybrids such as Zhengdan958 [18]. However, this heterotic pattern is rarely used in southwestern China, where alternative, scarcely investigated, heterotic patterns such as PB × PA are more suitable [17]. Therefore, the current study investigated the genetics of heterosis for yield-related traits in a PB × PA hybrid and the genetics of the yield-related traits per se.
Heterosis is defined as the superiority of F 1 hybrids in their performance over their parents [19,20] and has been attributed to three different genetic hypotheses: dominance [21,22], overdominance [19,23], and/or epistasis [24,25]. Specific loci for heterosis reported in previous studies corroborated each one of these hypotheses [26][27][28][29][30][31][32][33][34][35][36][37][38]. Those loci were identified using different designs such as the North Carolina design III with F 3 families or recombinant inbred lines (RILs), immortalized F 2 (IF 2) populations, or chromosome segment substitution lines (CSSLs). Those studies suggested different degrees of dominance for grain yield and inconsistency of heterotic loci across different heterotic patterns. For example, Frascaroli et al. [29] used the RILs from a cross between two elite inbred lines, B73 and H99, and three testcross populations obtained by crossing 142 RILs to each parent and F 1 found a significant correlation between heterozygosity level and phenotypic performance. Schön et al. [32] demonstrated high congruency of QTLs for heterosis of grain yield in three populations derived from the crosses of the heterotic pattern Iowa Stiff Stalk × Lancaster Sure Crop and proposed that diverse alleles were fixed in each heterotic group. Tang et al. [31] identified 13 heterotic loci for yield-related traits, and Guo et al. [34] reported that the genetic basis of heterosis for grain yield of the elite maize hybrid Yuyu22 relies on the cumulative effects of dominance, overdominance, and epistasis. Li et al. [37] detected 38 heterotic loci for ear-related traits and suggested a genotype-dependent combination of heterotic loci in different hybrids. Notably, Wang et al. [35,38] identified 169 heterotic loci associated with grain yield and yield components in two CSSL test populations. Therefore, the results pertaining to genetics of grain yield and heterosis using crosses from well-studied heterotic patterns, such as Reid × Lancaster and TSPT × Reid, cannot be extrapolated to other heterotic patterns such PA × PB.
In this study, in order to dissect genetic architecture of grain yield-related traits in maize and the contribution of different genetic effects to heterosis for grain yield-related traits, we evaluated 301 RILs from 08 to 641 × YE478 in four environments and an IF 2 population produced from those RILs under three environments. The elite line YE478, belonging to the PA heterotic group [17], is one of the hybrid foundation parents in China. The other parent, 08-641 (SAU08-641), which is derived from PB (or Tem-tropic I heterotic group) germplasm, is a typical elite inbred line widely used in southwest China since 1998. PA × PB is one of the predominant heterotic patterns in southwestern China. The main objectives of this research were (1) to detect and estimate QTLs with additive and additive-additive effects on yield-related traits using an RIL population (2) to identify QTLs with additive, dominant, and/or epistatic effect on yield-related traits using a hybrid population derived from the RIL population (IF 2 population) (3) to identify QTLs for mid-parent heterosis (MPH) using both mapping populations and (4) to identify the most suitable genomic regions associated with grain yield and MPH to be managed by marker-assisted selection.

Results
The performance of grain yield-related traits and heterosis across environments The parental line Ye478 had significantly higher mean RKP and kernel number per row (KNPR) compared with those of 08-641 (Additional file 1: Figure S1). The RIL and IF 2 populations showed considerable variations for yield-related traits (Fig. 1). The means of the F 1 were significantly higher than the means of any of the parental lines for ear weight per ear (EWPE), cob weight (CWPE), ear length (EL), ear diameter (ED), kernel weight per ear (KWPE), KNPR, and RKP. MPH of the hybrid Ye478 × 08-641 was < 16% for ED and RKP, 43.5% for EL, and > 70% for EWPE, KWPE, CWPE, and KNPR. Average MPH values across hybrids of the IF 2 population were similar to MPH values for F 1 (Additional file 4: Table S2). These results suggested a strong heterosis level for EWPE, KWPE, KNPR, and CWPE. The RIL and IF 2 populations thus seem suitable for dissection of genetic effects involved in yield-related traits and heterosis for these traits.
In the combined analyses of variance of the RIL and IF 2 populations across all environments, the sources of variation for genotype, environment, and genotype × environment interaction were highly significant (p < 0.01) for most yield-related traits ( Table 1). The broad-sense  The hierarchical cluster analysis based on the standardized data from the RILs classified the yield-related traits into two unrooted groups, with traits within each group being, in general, moderately to highly correlated based on Pearson's correlation coefficients (Fig. 2b).

Identification of QTLs
In the current study, 156 QTLs were found for averaged yield-related traits and MPH for those traits across all environments ( Table 2; Fig. 3). Ninety-six QTLs detected across environments were also identified in at least one environment via single environment analysis (Additional file 6: Table S4). These QTLs were distributed on all chromosomes, although two QTL-concentrated regions were identified at bin 3.08-3.09 and 10.04-10.05. Each QTL explained 0.38 to 33.99% of phenotypic variance (PVE), and over half of QTLs accounted each for less than 5% of PVE (Fig. 2c). The QTL number and the total PVE greatly differed among traits and populations ( Table 2; Fig. 2d, e). Nine QTLs were detected simultaneously in the RIL and IF 2 populations, while 13 QTLs detected in the IF 2 population and in the MPH dataset co-localized. Only two QTLs found in the RIL population were also positioned in the same regions when QTL analyses were performed for MPH. These results suggested that genetics involved in hybrid performance per se could be, in part, related to that implicated in MPH and it is in agreement with the moderate to high correlation coefficients between hybrid performance per se in the IF 2 population and the corresponding MPH levels reported for all traits, except for RKP (Additional file 7: Table S5).        The contribution of dominance and overdominance to the heterosis phenomena was assessed by the (|D/A|) ratio for each locus [34,39] using single-marker analysis (SMA) and inclusive composite interval mapping (ICIM) (Fig. 4) Among traits, we also found a strong association between the proportion of QTLs with overdominance and the ratio of MPH/ mid-parent for F 1 (r 2 = 0.951 *** ) or the MPH level across hybrids from the IF 2 population (r 2 = 0.966 *** ). We found, in general, significant but low correlations  Table S6).

Analysis of EPI and QEI
We estimated the total PVE by the epistatic loci in the whole-genome search using ICIM-EPI (Additional file 2: Figure S2). In the IF 2 population , epistatic loci had higher PVE than single-locus QTLs for CWPE (50% vs 41%), ED (24% vs 20%), CD (56% vs 32%), and KNPR (53% vs 37%). In the MPH dataset, the total PVE by epistatic loci was greater than by single-locus QTLs for all traits except for EWPE (20% vs 49%), CWPE (23% vs 28%), KWPE (43% vs 52%), and RKP (16% vs 43%), indicating that EPI had a greater impact on these traits than single-locus QTLs, regardless of the statistical significance. The current findings also suggest that the cumulative effects resulting from combining the epistatic loci and single-locus QTLs explained a large proportion of PVE. Meanwhile, based on mixed-model-based composite interval mapping (MCIM) with the full QTL model, 28 pairs of marker intervals showed significant epistatic interaction effects; 13 additive by additive (AA) interactions were detected in the RIL population, 10 epistatic interactions in the IF 2 population, and five digenic interactions in the MPH dataset (Additional file 9: Table S7). We observed 23 pairs of significant AA interactions, three for additive by dominance (AD) interactions, seven for dominance by additive (DA) interactions, and six for dominance by dominance (DD) interactions. Each marker interval interacted with one to five loci. In addition, a total of nine loci related to four traits exhibited significant additive × environment interactions in the RIL and IF 2 populations, seven in the RILs, and two in the IF 2 population; only one locus, qKWPE2-2, for the MPH level for KWPE showed a significant dominance × environment interaction in the IF 2 population (Additional file 10: Table S8). Many of these marker intervals involved in epistasis and QEI also contained QTLs with main effects on yield-related traits or MPH. For instance, the interval PZE-101184757/SYN25826 was involved in three significant AA interactions and contains qRKP1-1.

QTL co-localization across traits
Overlapping QTLs for different traits were assigned to the same QTL cluster. In total, 95 QTLs were grouped into 33 QTL clusters distributed across all chromosomes (Table 3). Each QTL cluster contained QTLs for two to six traits. The cluster QC3-3 harboring six QTLs affected the means of EL, RN, ED, KNPR, RKP, and KWPE, as well as the MPH for RKP. QC10-3, within a 5 cM region on chromosome 10, had a pleiotropic effect on EWPE, KNPR, CWPE, and the MPH level of ED and CD. As expected from the moderate to high correlation between hybrid performance per se and MPH, some loci associated with yield-related traits and MPH gathered together. Nearly 60% of QTL clusters were associated with two traits. Some of these clusters were specific for a determined population or for MPH. For example, clusters QC1-5 (bin 1.09-1.11) affected mid-parent heterosis for EWPE and KWPE. Contrarily, many other clusters comprised QTLs for yield-related traits and heterosis for those traits. Cluster QC1-4 (bin 1.08) affected KWPE and MPH for RKP. Genomic region 3.09 affected KNPR, RKP, and the MPH level of KWPE and EWPE. Interestingly, some of these clusters co-located in the interval of loci showing epistasis and QTL × environment interactions. For instance, QC7-1 (bin 7.02) had a pleiotropic effect on RKP, EWPE, KNPR, KWPE, and the MPH level for KWPE. This locus interacted with other loci exhibiting epistatic effects for RN in the RILs and for KWPE in the IF 2 population. This hotspot also interacted with E1 (2014JH) for RKP in the RILs. These results further confirmed the complex genetic architecture of yield-related traits.

Discussion
Heterotic loci with dominance and epistasis and the complex genetic network of "crosstalk" contribute to grain yield Previous studies [29,31,34,40,41] proposed that better-parent heterosis was yield-related loci displaying dominance, overdominance, and/or epistasis. In the present study, over half of QTLs exhibited overdominance for KWPE, EWPE, and KNPR, which is similar to previous results reported by Frascaroli et al. [29] for overdominant loci for grain yield and kernel number per plant and by Guo et al. [34] for single markers and QTLs for grain yield. Such a high average degree of overdominance for grain yield and kernel number per plant could be, in part, explained by linkage in the repulsion of loci of the dominance type in mapping populations [29]. In accordance with this, greater number of genome markers with overdominance for EL and ED was previously observed in the F 2:3 population, from which the current RIL population has been derived [42], than in the present mapping populations. The degree of dominance for the traits was well in accordance with the level of F 1 superiority over mid-parent performance and the   type of inheritance of different traits varied among those traits, which was consistent with other studies [29,33,34,37,43]. High correlation coefficients were observed between the proportion of QTLs displaying overdominance and the MPH levels across different traits. For example, the F 1 did not differ from the inbred parent lines 08-641 and Ye478 for RN and CD, and accordingly QTLs, such as qRN1-2 and qRN7-2 displaying similar additive effects, detected for those traits showed low or no significant dominance effects and were mainly detected in RIL and IF 2 populations [34]. The QTLs for EWPE, KNPR, and KWPE in IF 2 population and the MPH dataset displayed high dominance effects. In general, similar dominance effects were found in the IF 2 population and in the MPH dataset, but epistatic effects could interfere in the estimation of dominance effects and cause discrepancies between results obtained in different populations [44]. Discrepancies between additive effects estimated in the RIL and IF 2 populations were even more conspicuous. For example, additive effects estimations for qCWPE6 were different when obtained from RIL or IF 2 evaluations. These findings demonstrated that the loci displaying dominance contributed to yield-related traits and heterosis for those traits, although overdominance had an outstanding contribution to heterosis for EWPE, KWPE, and KNPR, which presented the highest values for MPH.
The level of whole marker loci with heterozygosity and phenotypic variations (hybrid performance per se or MPH) were weakly correlated, whereas Frascaroli et al. [29] and Larièpe et al. [33] demonstrated moderate-tohigh correlations between them. Such low correlations indicated that the overall heterozygosity contributed little to hybrid performance per se or heterosis [45,46]. Huang et al. [46] demonstrated that phenotypic variations for yield-related traits are highly correlated with the number of accumulated superior gene alleles. In contrast, numerous epistatic effects and QTL × environment interactions for yield-related traits have been addressed in the literature [29,31,34,37,[47][48][49][50]. The total PVE values by single-locus QTLs and epistasis together approach the high heritability for these yield-related traits [34], indicating the cumulative effects of different genetic effects [34,41]. The QTL model, which was used to estimate the loci showing additive, dominance, and EPI effects, identified 28 pairs of significant loci with epistatic effects and 10 loci with QTL × environment interaction in the two mapping populations. Of those, only five pairs of epistatic loci were observed for MPH. In addition, many of these loci showing interactions did not co-map with putative QTLs. Similar features were also reported by Yang et al. [51]. EPI and QEI are also involved in the complex genetic web for grain yield via "cross-talk" among loci [45,47,48,52]. Collectively, the cumulative effects of various genetic effects discussed above seem to explain the genetic mechanism of grain yield and betterparent performance [34,41].
QTL hotspots regulate the architecture of yield-related traits and heterosis Over 60% of QTLs gathered together and were integrated into 33 clusters distributed on 10 chromosomes. We found a good consistency for QTL clusters between our study and former studies (Table 3). Thirteen pleiotropic regions in the present study were consistently mapped to the clusters in all chromosomes, except in chromosome 4, that were defined by Chen et al. [53] in their meta-analysis of yield-related QTLs. Notably, many of these clusters related to MPH co-located with QTLs found in previous studies [29,31,34,35] for better parent performance (Table 3). Cluster QC1-2 associated with EL, CD, and CWPE overlapped with the hotspot region for ear weight, grain weight, and KNPR between 38 and 42 Mb on chromosome 1, which was highlighted by Guo et al. [49], and with a genomic region containing QTLs for heterosis for grain yield, ear weight, kernel weight, and ear length reported by Wang et al. [35]. The heterosis locus QC1-5 (1.09-1.11) for EWPE and KWPE co-located with QTLs for heterosis for EWPE and KRP presented by Wang et al. [35]; with a heterosis locus for grain yield at approximately 271 Mbp found by Tang et al. [31] using the heterotic pattern Mixed × Flint and by Larièpe et al. [33] using the pattern European flint × Dent; and with QTLs for 100-kernel weight and grain yield in an IF 2 population reported by Guo et al. [34]. The genomic region 2.04 (QC2-2) controlling EWPE and the MPH level for EWPE and KWPE co-located with a QTL involved in heterosis for EL, RN, KWPE, KRP, EWPE, and grain yield by Wang et al. [35], and it was nearby the region implicated in heterosis for grain yield in the Reid × Lancaster pattern reported by Frascaroli et al. [29] and the meta-QTL for grain yield, ear-related traits, and kernel-related traits highlighted by Chen et al. [53]. Cluster QC3-3, being particularly interesting in the current study because of its association with heterosis for several traits, overlapped with a QTL hotspot for heterosis for EL, EWPE, RN, KRP, KWPE, and grain yield located between 209 and 213 Mbp on chromosome 3 reported by Wang et al. [35]. Genomic region 3.09 (QC3-4), related to KNPR and RKP and the MPH level for KWPE and EWPE overlapped with the environmentally stable region for heterosis for KNPR, ear weight, and kernel weight between 218 and 228 Mbp on chromosome 3 reported by Li et al. [37] in the heterotic pattern Reid × TSPT. Another hotspot region, QC10-3, related to five traits co-mapped with a heterotic locus that affects grain yield and kernel number presented by Frascaroli et al. [29]. In addition, genomic regions 2.04, 7.02, and 9.03, associated with heterosis for KWPE, contain QTLs for grain yield heterosis in the Reid × Lancaster heterotic pattern published by Schön et al. [32]. The cluster QC6-3 for RKP, KWPE, and the MPH level for EL overlapped with the heterotic loci for grain yield found by Samayoa et al. [36] in the population derived from an American dent × European flint cross. These results suggested the consistency of hotspot regions across different heterotic patterns, which may be used for genetic improvement of grain yield in the future.
In contrast, hybrid performance per se for yield-related traits and the MPH level for those traits were moderately to highly correlated and 13 QTLs found in the IF 2 population co-localized with QTLs for MPH. These findings indicated the likely common genetic basis between grain yield and its MPH level. The results also partly supported that heterosis loci for grain yield are not independent from loci for grain yield, thus agreeing with Guo et al. [34] and refuting the conclusions by Tang et al. [31] who proposed independence between genomics involved in performance per se and heterosis.

Genetic basis for grain yield and its components in different populations from the same parents
Our study revealed that QTLs found in different populations derived from the same cross are not always consistent with previous studies [34,[54][55][56]. Among the 156 detected QTLs for yield-related traits and heterosis, only one QTL, qRN10-2, was previously found across F 2:3 families derived from the cross 08-641 × YE478 [60]. In addition, only nine QTLs with |D/A| < 1 were simultaneously detected in the RILs and the IF 2 population. Lack of congruency for QTLs found in populations derived from the same cross but with different homozygous levels could be due to the different genetic effects acting in each population, but also to the different environments used for evaluations and biases caused by population sampling [29].

Conclusions
Heterosis for yield-related traits in maize is one of the main issues for maize breeders, and investigation and assessment of heterotic loci involved in specific heterotic patterns are of particular interest. The inheritance of yield-related traits and MPH varied among different traits. A large proportion of the loci with dominance effect had a greater effect on most traits compared with loci with other genetic effects, whereas overdominance also contributed greatly to MPH for KNPR, EWPE, and KWPE in the current heterotic pattern. As QTL hotspots at 1.09-1.11, 2.04, 3.08-3.09, and 10.04-10.05 encompass genomic regions where other authors have found QTLs for heterosis across different heterotic patterns, markers in those regions could be used in marker-assisted selection programs for increased yield.

Plant material and field experiments
The mapping populations used in this study consisted of 301 RILs derived from a cross between the maize inbred lines 08-641 and Ye478. The elite inbred 08-641 comes from the Maize Research Institute of Sichuan Agricultural University, and Ye478 was provided by the Project of National Major Basic Dairy Research "973" Plan, which is developed by the Laizhou Academy of Agricultural Sciences [17]. The parental line 08-641 was crossed with Ye478 in the winter of 2010 to produce a set of RILs (F 8 ) by using the single seed descent method. sub-tropical humid monsoon climate), in Sichuan province, April 2016 (CZ2016). An IF 2 population was developed with 320 RILs following the procedure of Hua et al. [57] and Tang et al. [31] The 320 RILs were randomly divided into two groups of 160 RILs. Then, the lines of both groups each RIL was only involved in a single cross, yielding 160 F 1 crosses. We repeated the procedure mentioned above twice and constructed 320 crosses that constituted the IF 2 population. After genotyping the RILs, those with heterozygosity above 20% and their respective crosses were removed. After filtering the RILs, the IF 2 population consisted of 298 crosses that were evaluated under three environments (2015JH, 2016JH, and 2016CZ). In each environment, each trial was performed following a randomized complete block design with two repetitions. Plant density was 57,000 plants per ha. Each plot consisted of a single row. Each row with 14 plants was 3 m long, with 0.80 m space between rows. The trials with two mapping populations were adjacent and occupied a uniform and square parcel. Field management was in accordance with local practices.

Phenotypic measurements and analysis
In each plot, 10 well-pollinated ears were randomly chosen from all ears harvested. The nine yield-related traits measured (Additional file 3: Table S1) were: EL (cm), ED (mm), CD (mm), RN (count), KNPR, EWPE (g), CWPE (g), KWPE (g), and RKP. Means for each replication were used for further analyses. The distribution of the phenotypic traits and the Pearson's correlation, using PROC CORR from the statistical software package SPSS 17.0 (SPSS, Inc., Chicago, IL, USA), were conducted with the means of the phenotypic traits in the 301 RILs and 298 hybrids of the IF 2 population across all environments. A hierarchical cluster analysis ("hclust") was carried out for phenotypic traits, based on the standardized data from the RILs [58,59]. Combined analyses of variance of the RILs and the IF 2 population were computed for each trait using the GLM procedure in SPSS 17.0 (SPSS, Inc., Chicago, IL, USA) with genotype as a fixed effect and replication and environment as random effects. Broad-sense heritability (h 2 ) for the RILs and IF 2 population was computed on an entry mean basis as described by Hallauer and Miranda [60]: where σ 2 g represents the genetic variance, σ 2 ge represents the genotype × environment interaction variance, σ 2 represents the error variance, b represents the number of replications, and n is the number of environments. The 90% confidence interval of heritability (h 2 ) was determined according to Knapp et al. [61]. Each phenotypic trait of the 301 RILs and each trait of the IF 2 population were analyzed following two mixed models fitted by restricted maximum likelihood: where Y mk is the phenotypic value of genotype m in replication k; Y mik indicates the phenotypic value of genotype m in environment i and replication k; μ is the overall mean of the RIL population or the IF 2 population; G m is the random effect of genotype m; GE mi is the random effect of the interaction between genotype m and environment i; E i is the random effect of environment i; R k is the random effect of replication k; and ε mk and ε mik denote the random experimental error. The best linear unbiased predictor values (BLUPs) with Eqs.
(2) and (3) were used for single and combined environment analysis, respectively. The adjusted mean BLUP values of each RIL across four environments and each hybrid across three environments were used for QTL mapping. The analyses were conducted with the R program for statistical computing [59,62]. Besides, MPH was estimated as [40]: where F12 is the genotypic value of the F12 individual from the IF2 population, and (P1 + P2)/2 is the average of BLUP values of the corresponding parents (P1 and P2 from the RIL population) estimated from the RIL evaluation. The genotypic dataset for MPH depends on dominance and epistatic effects, but epistatic effects were not included and an additive-dominance model was assumed for the IF2 population [31,34].

Molecular linkage map construction
Following the modified CTAB protocol [63], DNA was isolated from 7-day-old seedling leaves of the 301 RILs and the parents grown in shade (10 plants per RIL as a bulk). The oligonucleotide pool assay consisted of 3072 well-distributed, high-quality single nucleotide polymorphisms (SNPs) from all 10 maize chromosomes that were selected from 56,110 SNPs in 513 maize inbred lines developed by the National Maize Improvement Center of China using Illumina GoldenGate technology. The protocol for genotyping SNPs using an Illumina BeadStation 500 G (Illumina, San Diego, CA, USA) was conducted as described by Fan et al. [64]. Marker data from the RILs were filtered for heterozygous data points (< 20%), missing data points (< 20%), and segregation distortion (in accordance with the expected Mendelian segregation ratio of 1:1). A total of 683 SNPs selected for their uniform distribution throughout all 10 maize chromosomes were used to construct the linkage map with a total genetic length of 1786.1 cM and an average interval distance of 2.61 cM. The genotypes of each cross from the IF 2 population were deduced from the marker genotypes of their RIL parents, with QTL mapping for the IF 2 population and the MPH dataset performed using the molecular linkage map of the RIL population [31,57]. The genetic map for the RILs and IF 2 population was developed using MapDisto 1.7.5 [65]. The Kosambi mapping function was used for converting recombination frequencies to genetic distances [66]. The average proportion of homozygous loci for YE478 variants, homozygous loci for 08-641 variants, and heterozygous loci in each individual of the IF 2 population was 25.2, 24.7, and 50.1%, respectively, according to the expected ratio in an F 2 population.

QTL mapping
QTL mapping for each trait in the RIL and IF 2 populations and for MPH using both populations was performed using QTL ICIMapping software [67] with ICIM [68,69]. For QTL detection on ICIM-ADD, the largest p-value for entering variables in stepwise regression of residuals on marker variables, PIN was set to 0.001. The largest p-value for removing variables is assumed to be two times the PIN. The threshold logarithm of the odds (LOD) score to declare significant QTL was calculated using the 1000 permutation test at Type I error of 0.05 with the step size of 1 cM was 2.5 for the RILs, IF 2 population, and the MPH dataset since very minor differences were observed in LOD score among them. Under the additive-dominance model assumed, only QTLs with dominant effect were detected for MPH [34]. We also performed single marker analysis (SMA) to detect QTLs with the same empirical threshold LOD mentioned previously, and the degree of dominance was calculated as a ratio of dominance to additive effects (|D/A|) for each QTL [34]. Loci with |D/A| values greater than 1.26 displayed overdominance; otherwise, the QTL was treated as dominance type [39]. The QTLs were nominated as follows: a lowercase letter "q" standing for QTL, the trait abbreviation, one number standing for the chromosome and another for the physical order or the QTL within the chromosome. The QTLs for different traits detected in overlapping or adjacent intervals were regarded as a QTL cluster. The graphical presentation of linkage maps and QTLs were generated with the R program [59]. To scan for significant marker interval interactions, epistasis interaction (EPI, digenic interactions; AA, AD, DA, and DD) was identified with the step in scanning of 5 cM and PIN = 0.0001 on ICIM-EPI. An empirical threshold LOD for digenic epistasis interaction was set at 5. Meanwhile, the QTL Network program 2.1 based on MCIM [70][71][72] was also used to identify QTL with AA effects in the RILs; AA, AD, DA, and DD effects in the IF 2 population using combined data across all environments and the MPH dataset using the average data across environments. QTL × environment interaction effects (QEI) in the RIL and IF 2 population were estimated for yield-related traits via combined analysis across all environments. The testing window and filtration window size was set at 10 cM with a walk speed set to 2 cM for one-dimensional genome scans. The F-statistic and the critical F-value were estimated with the help of 1000-permutations test [73] with p < 0.05 as experiment-wise significance level for candidate interval selection. The final full QTL model for each trait incorporated significant additive, dominance, and epistatic effects, as well as their interactions with environments. Bonferroni correction was used to compute the comparison-wise significance threshold assuming an experiment-wise error < 0.05 [74]. We also calculated the percentage of heterozygous loci for each RIL analyzed the correlation relationships between phenotypic performance and heterozygosity, and conducted regression analyses of phenotypic performance on heterozygosity [59].