QTL for Maize Midparent Heterosis in the Heterotic Pattern American Dent × European Flint under Corn Borer Pressure

Despite the importance of heterosis and the efforts to comprehend this phenomenon, its molecular bases are still unknown. In this study, we intended to detect Quantitative trait loci (QTL) for mid-parent heterosis under infestation with the Mediterranean corn borer (MCB, Sesamia nonagrioides Lef.) using a North Carolina design III approach with a RIL population derived from a European flint inbred (EP42) × American dent inbred (A637) cross. QTL for heterosis of kernel yield have been positioned in regions corresponding to previously identified QTL for the same trait in different backgrounds. These results reinforce the high congruency of genes controlling heterosis across populations, even when populations have been developed from different heterotic patterns. A high percentage of genetic variation for mid-parent heterosis (Z2) for kernel yield could not be explained. Furthermore, genomic regions involved in heterosis for yield and plant height were not found despite the high genetic correlation between Z2 transformations for kernel yield and plant height. The moderate power in detecting QTL for mid-parent heterosis suggests that many genes with low augmented dominance effects contribute to the genetic architecture of mid-parent heterosis; dominance and additive-additive epistatic effects could also contribute to heterosis. However, results from this and previous studies suggest that the region 8.03–8.05 deserves special attention in future works in order to fine map loci involved in mid-parent heterosis for yield.


INTRODUCTION
Heterosis is defined as the greater biomass, speed of development or fertility in heterozygotes, polyploids or hybrids compared to their parents (Schnable and Springer, 2013). Despite the importance of this phenomenom, its molecular basis is still unknown. In maize, mid-parent heterosis for yield, defined as the superiority of the hybrid over the mean of the parental inbreds, is particularly important because maize hybrids are extensively grown around the world. Most studies on the basis of maize mid-parent heterosis have been specifically focused in a unique heterotic pattern, Stiff Stalk Synthetic × Lancaster Sure Crop (Stuber et al., 1992;Frascaroli et al., 2007;Franco Garcia et al., 2008;Springer et al., 2009;Schon et al., 2010;Barber et al., 2012;Ding et al., 2014;Paschold et al., 2014), and conclusions from those studies cannot be generalized to other hybrids. Therefore, identifying quantitative trait loci (QTL) involved in heterosis in different genetic backgrounds is recommended as it has been rarely done (Tang et al., 2010;Lariepe et al., 2012;Guo et al., 2014). The heterotic pattern American dent × European Flint is particularly interesting because is the preferred pattern for hybrid development in Northern and Central Europe and Japan (Messmer et al., 1992;Moreno-Gonzalez et al., 2000;Technow et al., 2014).
The North Carolina design III and the immortilized F 2 population design have been used to study the QTL involved in maize heterosis (Stuber et al., 1992;Cockerham and Zeng, 1996;Frascaroli et al., 2007;Tang et al., 2010;Guo et al., 2014). To detect QTL for mid-parent heterosis with an immortalized F 2 population, a population of random crosses among RILs along with the RILs per se should be tested. Therefore, although the North Carolina design III with RILs does not allow to separate the contribution of dominance and epistasis to heterosis (Melchinger et al., 2007); it has been more extensively used because it does not require to test the RILs per se, but the backcrosses of those RILs to both parents (Frascaroli et al., 2007;Schon et al., 2010;He et al., 2012;Lariepe et al., 2012). The original North Carolina III mating design was developed by Comstock and Robinson (1948) to estimate the average level of dominance. This mating design has primarily been used in F 2 populations to determine the effect of linkages on the estimates of additive and dominance variances, and the average level of dominance using progenies developed by backcrossing individual S 0 plants of the F 2 population to both parents (Hallauer and Miranda, 1981). In the present study, we attempted to detect QTL for mid-parent heterosis using a North Carolina design III approach with a RIL population derived from a cross between a European flint (EP42) and an American dent (A637) inbreds. The evaluations were conducted under infestation with the Mediterranean corn borer (MCB, Sesamia nonagrioides Lef.) to induce stalk damage, a major limiting factor for maize yield in the Mediterranean region. The cross EP42 × A637 presented high specific combining ability for yield under MCB infestation  and the backcrosses of RILs derived from this cross are ideal materials to study quantitative genetics of mid-parent heterosis for yield under infestation with MCB. Genes involved in heterosis under high pest pressure could be clearly different from those implicated in heterosis at optimal conditions because experimental data suggests that a suppression of defense response gene activities are important for generating the hybrid vigor phenotype due to the possible competition for resource distribution between plant defense and plant growth (Groszmann et al., 2015).

Plant Material and Genotypic Data
A population of 146 RILs derived from the cross of the European flint inbred line EP42 and the American dent inbred line A637 was developed (Samayoa et al., 2014)  successfully backcrossed to both parents for mapping QTL for heterosis under infestation with MCB eggs. The mechanisms of resistance can be classified into three groups: antixenosis, antibiosis and tolerance (Painter, 1951). Inbreds EP42 and A637 were both classified as susceptible under MCB attack attending to their antibiosis levels, but A637 was tolerant (Butron et al., 1998) and showed good general combining ability for yield ; while the cross EP42 × A637 presented significant and favorable specific combining abilities for yield and yield loss under MCB infestation . The linkage map used has been previously published by Samayoa et al. (2014).

Phenotypic Data
The backcrosses of the 136 RILs to both parents were evaluated at Pontevedra (42 • 24 ′ N, 8 • 38 ′ W, and 20 m above sea level), Spain, in 2010 and 2011. Backcrosses of RILs were arranged in 17 sets comprised each one of 8 RILs backcrossed to both parents, backcrosses of the same RIL were kept in adjacent plots. Sets and RILs within sets were randomized in each replication (two per trial). The trials were hand planted and each experimental plot consisted of one row spaced 0.8 m apart from the other row with 13 two-kernel hills spaced 0.18 m apart. Plots were overplanted and thinned, obtaining a final density of ∼70,000 plant ha −1 . The evaluations were performed under artificial infestation with eggs of MCB obtained at the Misión Biologica de Galicia by rearing MCB as described by Eizaguirre and Albajes (1992) and Khan and Saxena (1997). Five plants of each plot were infested with ∼40 MCB eggs placed between the stem and the sheath of a basal leaf. Data collected were: days to silking computed as the number of days from planting to 50% of plants in the plot showing silks; ear and plant heights on five representative plants as the length (in cm) from the ground to the main ear and from the ground to the top of the plant, respectively; stem lodging defined as the percentage of plants in the plot with the stem broken below the main ear; kernel resistance to MCB attack on the ears of the infested plants collected at harvest according to a subjective visual resistance scale from 1 to 9 (1 indicates completely damaged and 9 indicates no damage); tunnel length as the average stem tunnel length made by borers on the five infested plants; kernel yield of the plot expressed in Mg ha −1 at 140 g H 2 O kg −1 (infested and non-infested plants were considered); and kernel moisture at harvest as g of water in 100 g of kernels.

Statistical Analysis
The individual analysis of each trial was made with the SAS mixed model procedure (PROC MIXED) (SAS Institute Inc, 2011) considering replications as random effects and sets and RILs within sets as fixed effects. The combined analysis across years was made considering years and replications as random effects. A Best linear unbiased estimate (BLUE) was obtained to estimate mean phenotypic values of each RIL backcrossed to P 1 (EP42) and P 2 (A637) in each year and across years.
To express the mid-parent heterosis as the sum of individual QTL effects, Melchinger et al. (2007) defined a new type of heterotic genetic effect d * i that includes the dominance effect of QTL i minus half the sum of its additive × additive epistatic interactions and was named augmented dominance effect. The same authors defined another genetic effect a * i (named augmented additive effect) that includes the additive effect for QTL i minus half the sum of dominance × additive epistatic interactions with all other QTL and corresponds exactly to the net contribution of QTL i to parental differences. They also demonstrate that contrasts Z 1 and Z 2 in a modification of design III of Comstock and Robinson (Comstock and Robinson, 1952) using RILs instead of a F 2 account for augmented additive and dominance effects, respectively. The linear transformations Z 1 and Z 2 are defined as follow: Z 1 = (H 1 + H 2 )/2 (trait mean across each pair of RIL backcrosses), Z 2 = (H 1 -H 2 )/2 (half the trait difference between each pair of RIL backcrosses), H 1 and H 2 being the phenotypic observations on backcrosses of each RIL with the parental inbreds P 1 and P 2 , respectively.
As a genome scan with Z 2 localizes genomic regions affecting mid-parent heterosis and with Z 1 identifies QTL that contribute to parental differences, QTL analyses were performed with values of Z 1 and Z 2 for each RIL and trait across years. In addition, a joint QTL analysis of Z 1 and Z 2 was implemented. Genetic variances and heritabilities (ĥ 2 ) across environments were estimated for Z 1 and Z 2 for each trait on a family-mean basis as described previously by Holland et al. (2003). The augmented degree of dominance (D * ) was estimated as the square root of the division between genetic variances for Z 2 and Z 1 . The genetic correlation coefficients among Z 1 tranformations for the different trais, among Z 2 tranformations and between Z 1 and Z 2 for each specific trait were computed following Holland (2006).
The QTL analyses were performed with the software package PlabMQTL (Utz, 2012). Composite interval mapping approach (CIM) was conducted for the QTL detection and to estimate QTL effects. According to a previously executed permutation test with 1000 random reshuffles, a LOD threshold of 2.5 was chosen to declare significant a putative QTL assuming an experimentwise error <30%. The phenotypic variance explained by the QTL model was estimated by the adjusted coefficient of determination (R 2 adj ) which accounts for the total proportion of the phenotypic variance explained by all detected QTL in the final fit. Then, the phenotypic variance explained by an individual putative QTL i was calculated as: Where n = total of QTL in the final fit and part R 2 i = partial coefficient of determination, estimated for the ith QTL detected (Zhu et al., 2004). The proportion of the genotypic variance (p) explained by all detected QTL in final fit was estimated from the ratiop = R 2 adj /ĥ 2 and the proportion of the genetic variance explained by each individual QTLi was estimated asp i = R 2 i /ĥ 2 . A 5-fold cross validation (CV) approach was employed for evaluating QTL mapping results (Utz et al., 2000). The data set (DS) was split at random into 5 data subsets; four subsets were combined to form de estimation set (ES) and the remaining subset formed the test set (TS). 1000 runs of CV were performed in order to determine the QTL frequency and bias of QTL effect estimate at the position of a QTL detected in the original data set (Melchinger et al., 2004). The proportion of the bias in the effect estimation of each individual QTL due to genotypic samplings was calculated as the difference between the average estimates in ES and TS divided by the estimate in ES.

RESULTS AND DISCUSSION
Genetic variances for Z 1 were significantly different from zero for all traits except for stalk lodging, tunnel length and kernel resistance; while genetic variances for Z 2 were significantly different from zero for yield, yield-related traits, such as plant and ear heights, days to silking and stalk lodging ( Table 1). The augmented degree of dominance was significantly higher than 1 for kernel yield, and significantly lower than 1 for days to silking, kernel moisture, and plant an ear heights. These results emphasize the great importance of heterosis for kernel yield and suggest that overdominance could play an important role as it has been reported in similar studies, although it should be taken into account that the augmented degree of dominance is biased by QTL linkage and epistasis and, consequently, does not exactly reflect gene action at individual loci (Melchinger et al., 2007;Schon et al., 2010;Lariepe et al., 2012). Heterosis was also important for yield-related traits, such as plant and ear heights and days to silking, but, as the augmented degree of dominance was significantly <1, dominance would be more important than over-dominance effects for these traits. In accordance to Lariepe et al. (2012), the higher importance of heterosis for days to silking in this cross compared to that in crosses between Reid and Lancaster inbreds could be due to the earliness of materials used and could emphasize the importance of studying heterosis in different backgrounds to get a better picture of the genetic architecture of the trait. Heritabilities for Z 1 of agronomical traits, such as days to silking, kernel moisture, plant and ear heights and kernel yield were low to moderate; while heritabilities for Z 2 were higher for plant and ear heights and kernel yield, slightly lower for days to silking, and considerably lower for kernel moisture than heritability estimates for the corresponding Z 1 transformations ( Table 1). Heritabilities did not differ from zero for Z 1 and Z 2 for related-resistance traits, although in a previous work heritability for tunnel length of the RILs per se was low but significantly different from zero (Samayoa et al., 2014). The discrepancy between the heritabilities computed for tunnel length of the RILs per se and for Z 1 in this study could be consequence of opposite effects for additive × dominance epistasis and additivity. It could also be due to a higher experimental error in the current experiment based on the number of genotypes tested (backcrosses of RILs to both parents) that doubled the number of RILs per se. As additive × additive epistasis effects were not significant for tunnel length in a previous study (Butron et al., 2009), the absence of heritability for Z 2 for tunnel length would indicate that dominance effects are not important for resistance to stem tunneling by MCB as it has been already suggested Cartea et al., 1999).
As the environment × genotype interaction was not significant for Z 1 and Z 2 for yield and resistance-related traits, QTL analysis was not independently done for each environment (Data not shown). No QTL for Z 1 for kernel yield were found, while two were found for Z 1 for days to silking that explained the 32% of the genetic variation, and one for Z 1 for kernel moisture, tunnel length, plant and ear heights and kernel resistance that explained the 17, 72, 20, 21, and 18% of genetic variation, respectively ( Table 2). The QTL identified for Z 1 for tunnel length at 9.03-9.04 co-localized with QTL for Z 1 for both plant and ear heights. There was not a complete agreement between QTL detected for Z 1 and QTL detected for the per se traits in a previous study (Samayoa et al., 2014); although evaluations a Definitions of transformations Z 1 and Z 2 can be found in the materials and methods' section. b Effects estimated with the complete data set (DS), the estimation set (ES), the validation set (TS), and the proportion of the bias in the effect estimation of each individual QTL due to genotypic samplings calculated as the difference between the average estimates in ES and TS divided by the estimate in ES; R 2 i , the adjusted proportion of phenotypic variance explained by QTL i ; p i , the adjusted proportion of the genetic variance explained by QTL i ; Freq, percentage of cross-validation runs in which the QTL has been detected.
Frontiers in Plant Science | www.frontiersin.org of the RILs per se and backcrossed to the parental inbreds were performed in adjacent trials. Higher experimental errors expected in the current experiment compared to those in which the RILs per se were evaluated could be responsible for detecting fewer QTL for Z 1 . In addition, additive and additive × additive epistasis effects were confounded in the study by Samayoa et al. (2014); while additive and additive × dominance epistasis effects are undistinguished in the current study and both contribute to the augmented additive (a * ) effect. Therefore, epistasis effects in this particular cross could be important enough to counteract additive effects. Previous studies have (Yan et al., 2006;Guo et al., 2014) already pointed out that digenic interactions at the twolocus level might play an important role in the genetic basis of maize heterosis. In any case, QTL for Z 1 for tunnel length, and plant and ear heights positioned at 9.03-9.04, QTL for Z 1 for kernel resistance positioned at 5.05-5.06, and QTL for Z 1 for days to silking positioned at 8.05-8.06 were placed in the same bin or in adjacent bins to those where QTL for these traits per se were previously found (Samayoa et al., 2014) corroborating the presence of genes with additive effects for these traits in those regions.
Quantitative trait loci (QTL) for mid-parent heterosis (Z 2 transformation) were found for the agronomical traits stem lodging, kernel moisture, kernel yield, plant height, and ear height and explained the 42, 100, 39, 11, and 48% of genetic variance, respectively ( Table 2). QTLs for mid-parent heterosis for kernel yield were located at bins 1.02-1.03, 304-3.05, 6.06, and 8.05, were stable across years, and explained the 8.5, 6.2, 8.9, and 15.3 of the genetic variance, respectively (Tables 2, 3). The accuracy of effect estimation was low for QTL at 1.02-1.03 and 3.04-3.05, moderate for the QTL at 6.06 and high for the QTL at 8.05 because the estimations of bias proportions derived from cross-validation runs were 0.56, 0.58, 0.30, and 0.08, respectively. In addition, the level of occurrence of these QTLs in crossvalidations was moderate and high for QTL at 6.06 and 8.05, respectively. When the Akaike information criterion (AIC), a less conservative criterion, is used instead of the modified Bayesian information criterion (BIC) for fitting a model that includes all QTLs for heterosis for kernel yield, the percentage of genetic variance explained increases to 52% and additional QTL were found at 2.04-2.06, 3.03, and 5.03-5.04 (Data not shown). Schon et al. (2010) also detected QTLs for the Z 2 transformation for kernel yield in all these regions (same or adjacent bins) in at least one of three populations derived from Reid × Lancaster inbred crosses. However, Lariepe et al. (2012) detected nine QTL for Z 2 for kernel yield that explained the 53% of phenotypic variance in crosses among an European flint and two dent inbreds, but only three of the genomic regions detected (the centromeric regions of chromosomes 3 and 5 and the long arm of chromosome 8) bear QTL for Z 2 in the current study. Lastly, Tang et al. (2010) located QTL for mid-parent heterosis for yield in bins 1.01, 1.10, and 8.03 using an immortalized F 2 population derived from the cross between a Chinese dent inbred and an exotic flint inbred (Tang et al., 2010;Guo et al., 2014). Then,  QTL for kernel yield heterosis under high insect pressure were located in genomic regions known for including loci involved in heterosis at optimal conditions, but these results cannot serve to declare that, in general, competition between hybrid vigor and defense mechanisms does not exist because, in the current study, genetic variability for resistance-related traits was low. The number of QTLs detected for yield heterosis were lower than those detected by Schon et al. (2010) in dent crosses, similar to those detected by Lariepe et al. (2012) in early populations derived from European flint × dent and dent × dent crosses and larger than those detected by Tang et al. (2010) in a dent × flint cross. Assuming similar population sizes, the higher power to detect QTL for mid-parent heterosis for yield in the Reid × Lancaster heterotic pattern than in other heterotic patterns evaluated would relate to the larger magnitude of heterosis (absolute value) in the Reid × Lancaster pattern compared to the others. Results from this and previous studies suggest that the region 8.03-8.05 deserves special attention in future works in order to fine map loci involved in mid-parent heterosis for yield. This genomic region harbors QTL with important additive augmented (a * ) effects for days to silking and dominance augmented (d * ) effects for ear height and kernel yield in this study. Previously, in the same region, significant QTL for Z 2 for yield were identified in several populations (Schon et al., 2010;Tang et al., 2010;Lariepe et al., 2012) and is one of the genomic regions with high density of QTL peaks for yield based on a MetaQTL analysis involving results from 44 manuscripts (Martinez et al., 2016). In addition, Samayoa et al. (2014) also detected important QTL for days to silking and plant and ear heights in the same region. This region presents high homology with the region of chromosome 6 [these regions correspond to duplicated regions of an ancient tetraploid (Schnable et al., 2011)] where a QTL for yield heterosis has been located in the current study and in that performed by Schon et al. (2010).
No overlapping QTL for Z 1 and Z 2 for kernel yield were found although genetic correlation between both parameters was moderate (0.50) suggesting that genomic regions containing QTL with augmented additive and dominance effects for yield went undetected due to the low magnitude of genetic effects (Table 4). In a similar way, although the correlation coefficients between Z 1 for yield and Z 1 for ear height and between Z 2 for yield and Z 2 for plant and ear heights were high and significant (Table 4), only two QTL for Z 2 for kernel yield (at 8.05 and 5.03-5.04) co-localized with QTL for Z 2 for ear height. These facts along with the lack of power for detecting QTL for Z 1 for kernel yield and Z 2 for days to silking, although genetic variability for these transformations significantly differed from zero (Table 1), suggest that many genes with low augmented additive and/or dominance effects can contribute to kernel yield and mid-parent heterosis for yield, at least in this particular cross. In addition, the augmented dominance effects of some loci that contribute to heterosis could be consequence of additiveadditive epistatic effects. Lariepe et al. (2012) showed that when a QTL for Z 2 for yield overlap with a QTL for yield with no significant dominance effect, the QTL probably presents significant additive-by-additive interaction effects. Therefore, the QTL detected for Z 2 for yield at 5.03-5.04 that was already reported as important for yield per se (Samayoa et al., 2014) could rather contribute to Z 2 due to additive-by-additive interactions than to dominance effects. The significant contribution of additive-by-additive interactions to heterosis in this particular cross was already reported in a previous study (Butron et al., 2009).
As conclusions, the low to moderate power of the experiment to detect QTL for mid-parent heterosis suggests that many genes with low augmented dominance effects can contribute to the genetic architecture of mid-parent heterosis. Furthermore, dominance and additive-additive epistatic effects could also contribute to heterosis. The genomic regions detected for midparent heterosis for kernel yield in this particular cross (European flint × American dent inbreds) had been detected as significant for heterosis in other heterotic patterns reinforcing the idea of high congruency of heterosis genes across populations. Results from this and previous studies suggest that the region 8.03-8.05 deserves special attention in future works in order to fine map loci involved in mid-parent heterosis for yield.

AUTHOR CONTRIBUTIONS
RM and AB designed the work; LS did the collection and statistical analysis of data; and AB drafted the manuscript that was critically revised by LS and RM.

FUNDING
This work was supported by the National Plan for Research and Development of Spain (projects AGL2012-33415 and AGL2015-67313-C2-1-R, both funded in part by the European Regional Development Fund).

ACKNOWLEDGMENTS
LS acknowledges a contract JAE-Predoc from the Spanish Council for Scientific Research (CSIC). We thank Jeffrey Dunne for his careful language review and comments.