Xiaowei Li, Min Wang, Renyu Zhang, Hui Fang, Xiuyi Fu, Xiaohong Yang, Jiansheng Li

State Key Laboratory of Plant Physiology and Biochemistry, National Maize Improvement Center of China, MOA Key Laboratory of Maize Biology, Beijing Key Laboratory of Crop Genetic Improvement, China Agricultural University, Beijing 100193, China

Keywords:Embryo size Kernel size Quantitative trait loci Genetic relationship

ABSTRACT The embryo in maize has a critical role in controlling kernel nutrition components and grain yield.We measured five embryo weight and size traits, six kernel weight and size traits, and five embryo-tokernel ratio traits in a nested association mapping (NAM) population of 611 recombinant inbred lines(RILs)derived from four inbred lines including the high-oil,giant-embryo line BY815 as the common parent.Using three statistical methods, we identified 5-22 quantitative trait loci (QTL) for each trait,explaining 4.7%-46.7% of the phenotypic variation.The genetic architecture of maize embryo size and its related traits appeared to be dominated by multiple small-effect loci with little epistasis, and the genetic control underlying embryo size appeared to be distinct from that underlying kernel size.A trait-QTL association network included 205 nodes and 439 edges and revealed 28 key loci associated with at least three traits.Cloned maize genes including ZmUrb2, Emp12 and ZmBAM1d, maize orthologs of known rice genes that control seed size including BG1, XIAO and GS9, and 11 maize orthologs of Arabidopsis EMBRYO-DEFECTIVE (EMB) genes were identified as underlying these key loci.Further, the phenotypic and genetic relationships between embryo size and kernel size were evaluated, and genetic patterns for identified loci that control embryo size and its related traits were proposed.Our findings reveal distinct genetic architectures for embryo size, kernel size, and embryo-to-kernel ratio traits and establish a foundation for the improvement of embryo-size-mediated kernel nutrition and grain yield.

1.Introduction

Maize has an irreplaceable role in ensuring world food and nutrition security[1].The maize kernel is composed of the embryo,developing from the zygote; the endosperm, developing from the fertilized central cell;and the pericarp,developing from the maternal integuments.In contrast to the starch-rich endosperm, the embryo accounts for ~10% of kernel dry weight and includes~85% of the kernel oil, 10%-20% of the kernel protein, and ~50%of the kernel vitamin E[2-4].As oil contains 2.25 times more calories than starch on a weight basis[2],greater energy yield could be obtained from maize by increasing embryo size.However,a reduction in seed size has been associated with an increase in the proportion of the embryo in previous high-oil maize selection experiments [5,6], suggesting complex relationships between embryo size, kernel size, and kernel nutrition components.An understanding of the genetic architecture responsible for maize embryo size and its genetic contribution to kernel size and kernel nutrition variation would accelerate the improvement of maize.

Substantial efforts [7-14] have been made to characterize and investigate the genetic control underlying various maize kernel mutants includingdefective kernel(dek),embryo specific(emb),andempty pericarp(emp).To date, about 50 genes controlling maize embryo or kernel development have been cloned using mutants.Among them, >30 genes encoding pentatricopeptide repeat (PPR) proteins are involved in regulation of the RNA metabolism of mitochondrial or chloroplast genes [15], indicating the important role of these genes as essential regulators in controlling maize embryo or kernel development.However, there have been few attempts to elucidate the genetic architecture of embryo size based on natural variation.Thus, the genetic basis of the natural variation of embryo size in maize is poorly understood [16,17].

Linkage mapping and association mapping are the two main methods for investigating the genetic architecture of complex quantitative traits.The first method uses recombination events that have occurred during the construction of the mapping population.The clear genetic structure of the linkage population is conducive to controlling the number of false positives, and this design can effectively reveal low-frequency variations.Association mapping, in contrast, affords greater mapping resolution based on the use of historical recombination events that accumulate in natural populations,but the complex population structure often yields false positives or false negatives.Association mapping also has little power to detect low-frequency alleles and small-effect QTL[18,19].Nested association mapping simultaneously exploits the advantages of both linkage and association mapping [20] and has been successfully applied [21] to investigate the genetic basis of complex quantitative traits in maize.

In this study, we combined the high statistical power of the maize NAM panel of 611 RILs with high-density SNPs (singlenucleotide polymorphisms) to identify QTL and underlying candidate genes responsible for natural variation in maize embryo size and related traits and to study the genetic relationship between embryo size and kernel size in maize with the aim of clarifying the role of embryo size in controlling kernel nutrition components and grain yield.

2.Materials and methods

2.1.Materials, phenotype measurements, and phenotypic analysis

The common parent BY815 (a high-oil/giant-embryo inbred line) was crossed with the maize inbred lines DE3, K22, and KUI3(Fig.S1) and the progeny were advanced by selfing for six generations to yield 207, 207, and 197 RILs, respectively (Fig.S1).The resulting NAM population of 611 lines was planted in three locations in China:Yunnan province in 2012,Beijing in 2013,and Hainan province in 2013.Across all of the trials, a randomized complete block design was used.All lines were planted in unreplicated single-row plots consisting of 2.5-m rows with 11 plants per row and 0.67 m between rows.The 11 plants in each row were self-pollinated, and 300 kernels were sampled in bulk for each row, with equal numbers of kernels being collected from the middle parts of five ears harvested from different plants.Twenty kernels were selected randomly from each bulk sample for phenotyping and dried for 60 h at 45 °C.

Five embryo weight and size traits were measured:hundredembryo weight (HEMW), embryo volume (EV), embryo width(EW), embryo length (EL), and embryo thickness (ET).Six kernel weight and size traits included hundred-kernel weight (HKW),hundred-endosperm weight (HEDW), kernel volume (KV), kernel width (KW), kernel length (KL), and kernel thickness (KT).Five embryo-to-kernel ratio traits included embryo-to-endosperm weight ratio (EER), embryo-to-kernel volume ratio (EVR),embryo-to-kernel width ratio (EWR), embryo-to-kernel length ratio (ELR), and embryo-to-kernel thickness ratio (ETR).These traits are described in Table S1.The KL and KW of single kernels were measured with a semi-automatic imaging system (National Engineering Research Center for Information Technology in Agriculture,Beijing,China),and the KT of single kernels was measured with an electronic vernier caliper to 0.01 mm precision.The 20-kernel volume was estimated by the ethanol displacement method using a burette with 0.1 mL precision and converted to single KV.After these measurements were made, the kernels were soaked in deionized water for 24 h at 65°C followed by manual dissection into embryo and endosperm with a scalpel.The dissected embryos and endosperms were again dried for 60 h at 45 °C and then were weighed using an electronic scale to 0.001 g precision.The 20-embryo weight and 20-endosperm weight were converted to HEMW and HEDW by multiplying by 5, and HKW was calculated as the sum of HEMW and HEDW.The EL,EW,and ET of single kernels were measured using the same methods as the kernel-size traits.The 20-embryo volume was estimated by the ethanol displacement method using a burette with 0.05 mL precision and converted to single EV.The five embryo-to-kernel ratio traits were calculated from the embryo and kernel trait values.

Variance decomposition and broad-sense heritability estimation [22] were applied to the phenotypic data across the three environments using theaovfunction in R (version 3.5.1) [23].Because the NAM population was phenotyped in all three environments without replication, the following model was fitted to obtain unbiased estimates of each trait value in each line:yij=μ+ei+fj+εij,whereyis the observed phenotype for thejthfamily in theithenvironment,μ is the grand mean(the total of all the data values from all three environments divided by the total sample size),eiis the random effect of theithenvironment,fjis the random effect of thejthfamily,and εijis an error term that follows an independently and identically distributed N (0,) distribution.Broad-sense heritability was calculated asH2= [(MSF-MSε)/e]/(MSF/e), whereMSFis the mean square between families,MSε is the mean square error, andeis the number of environments.The best linear unbiased prediction (BLUP) value for each line was calculated across all environments using a mixed linear model with both genotype and environment treated as random effects in the R packagelme4[24],and the trait BLUP values were used for subsequent analyses.Owing to the nonzero intercorrelations among predictors,it is difficult to determine the relative weights of predictors in multiple regression.Trait-trait phenotypic contribution analysis was performed in R according to a previous method proposed by Johnson[25] that can solve the correlated-predictor difficulties, producing results very similar to those produced by more complex methods.Principal component analysis (PCA) of 16 traits was performed usingprcompfunction in R.

2.2.Genotyping

The NAM population was genotyped with the MaizeSNP50 BeadChip [26], which contains 56,110 SNPs.Among the RIL population, respectively 13,729 SNPs, 13,603 SNPs, and 14,024 SNPs were polymorphic within the BY815 × DE3, BY815 × K22, and BY815 × KUI3 populations.For the three populations, highdensity genetic maps were constructed as described previously[27], containing respectively 2382, 2263, and 2683 genetic bins(intervals within which no recombination events occur).A joint linkage map for the 611 RILs was constructed using the CMP module in QTL IciMapping(version 4.1)software[28]and was 1771 cM in length with 2630 bins.Genotype data from the 611 RILs with 2630 bins were used for the genome-wide association study(GWAS).

2.3.Separate-linkage mapping (SLM), joint-linkage mapping (JLM),and GWAS

Within each RIL population, SLM was performed using the CIM model [29] in Windows QTL Cartographer 2.5 [30].The program parameter settings were as follows:-Model 6 -Control marker 5 -Window size 10 cM-Regression method 2.Within each population,500 permutations were performed for each trait to determine the logarithm of odds (LOD) threshold for the QTL significance test,and the resulting LOD threshold ranged from 2.7 to 3.3(α=0.05).For simplicity,LOD=3.0 was used as the global threshold and the QTL support interval was determined by the width of the LOD peak minus 1.5-LOD[31].QTL detected for different traits were integrated into a ‘‘consensus QTL” if their support intervals overlapped, as described previously [32].

JLM was performed using the NAM module in QTL IciMapping(version 3.3).The mapping method was ICIM-ADD [33] and the step was set to 1.0 cM.The QTL detection threshold was manually set to LOD = 3.0, to be consistent with that used in the SLM.The QTL support interval was determined as described above for the SLM.Epistasis analysis of JLM QTL was performed as described previously [34,35].A two-way analysis of variance was used to test pairwise additive-by-additive epistatic interactions for all identified QTL for each trait (P< 0.05) from the JLM, and the proportion of variance explained by epistasis were estimated by comparing the residual variance of the full model containing all single-locus effects and two-locus interaction effects with that of the reduced model excluding two-locus interaction effects.

GWAS was performed using the stepwise regression and resampling methods reported previously [36] with minor modification.For each of 100 iterations,80%of the RILs were randomly selected with replacement from the population, and the analysis was repeated.The proportion of times that each SNP was included in the multiple-SNP model was defined as its bootstrap posterior probability (BPP) score.Finally, the SNPs retained in the model were subjected to linkage disequilibrium analysis, and the other SNPs in the block where the peak SNP(with most significant associationP)was located were removed usingr2=0.7 as the threshold.The peak SNP in the block was defined as the trait-associated SNP.The physical position of the trait-associated SNP was extended to 2 Mb (determined by the average bin size of the three separate linkage maps) on both sides, and the resulting region was defined as the support interval of that site.QTL identified by GWAS were used to calculate the inter-trait genetic contribution values.The genetic contribution of trait x to trait y(expressed as trait y-trait x)was defined as the proportion of reduced PVE(phenotypic variance explained) for trait y after treating trait x-associated QTL as covariates.

In this study, provided that the physical intervals of different QTL overlapped (whether partially or completely), they were defined as overlapping loci.QTL detected by the three methods were compared for each trait, and the numbers of overlapping and distinct QTL were collected for each trait.The numbers of overlapping and distinct QTL detected over all of the traits were then summed to yield the overall proportion of overlapping loci.

2.4.QTL effect estimation

The effects of all QTL identified by GWAS were normalized as described previously [37] with modification, as follows:scaled effect = GWAS effect/SD × heritability, where the scaled effect is the normalized effect, the GWAS effect is the original QTL effect detected by GWAS,SD is the standard deviation of the BLUP values for each trait in the NAM population,and heritability is the broadsense heritability of each trait in the NAM population.

2.5.Construction of trait-QTL association networks

Trait-QTL association networks were constructed using Cytoscape(version 3.7.2)[38]as described previously[35,39]with modification.Mapping results from GWAS were used to construct trait-QTL association networks, as the uniform support intervals of QTL detected by GWAS were more convenient for identifying ‘‘consensus QTL” than QTL from JLM.A total of 222 QTL detected by the GWAS were integrated into 81 consensus QTL based on the overlaps of QTL support intervals.First, all 16 traits and their corresponding QTL were assigned as nodes.Links between individual traits and QTL and links within consensus QTL were assigned as edges.The network was first processed by applying the‘‘Preferred Layout”option,and the distribution of nodes in the plane was then adjusted according to the degree of association to avoid overlaps.Finally,the sizes of the nodes corresponding to QTL were adjusted to correspond to the strength of their association with the corresponding traits.

2.6.Method of nominating candidate genes

First, all the genes contained in the QTL support interval were extracted using the physical position of the QTL (B73_RefGen_v3)[40].The functions of these genes were determined according to protein information for the genes or annotation information from their orthologs in rice orArabidopsis[41].Gene with the most relevant functional annotation information in the QTL support interval was nominated as the candidate gene.

3.Results

3.1.Phenotypic variation in 16 embryo size and related traits

The 16 traits analyzed in the NAM population showed abundant variation (Fig.1; Table S2).The values of RILs in the NAM population for embryo weight and size traits generally fell between the low parent, DE3, and the high parent, BY815, whereas kernel weight and size traits showed high transgressive segregation(Fig.1).This finding indicated that the genetic controlling for individual embryo traits might differ from those of kernel traits.Broadsense heritability ranged from 0.40 (KV) to 0.87 (EWR), and the heritabilities of embryo and kernel width traits (EW, 0.86; KW,0.73; EWR, 0.87) were greater than that of embryo and kernel length traits (EL, 0.70; KL, 0.73; ELR, 0.55) and embryo and kernel thickness traits (ET, 0.47; KT, 0.42; ETR, 0.48).

The correlation between EL and KL (r= 0.86,P< 2.20 × 10-16)was much stronger than those between EW and KW (r= 0.38,P< 2.20 × 10-16) or between ET and KT (r= 0.58,P< 2.20 × 10-16), indicating that the strength of the constraints of the kernel on the embryo differs among the three dimensions.EER was much more strongly correlated with HEMW (r= 0.71,P< 2.20 × 10-16)than with HKW (r= 0.14,P= 4.83 × 10-4), indicating that most of the embryo-to-kernel weight ratio variation might be from embryo weight instead of kernel weight (Fig.S2).

In PCA of 16 traits, PC1 separated the embryo traits from the embryo-to-kernel ratio traits, PC2 separated embryo traits and embryo-to-kernel ratio traits from kernel traits,and PC3 separated embryo length and kernel length from the other traits (Fig.S3),indicating that the 16 measured traits comprehensively covered the morphological variations of the embryo and kernel.The PCA results revealed the relationships among the 16 traits at the phenotypic level:the variation in embryo traits differed from that in kernel traits, the variation in embryo-to-kernel ratio traits differed from that in embryo traits, and the variation in EL was similar to that in KL.

3.2.Genetic architecture of 16 embryo size and related traits

SLM revealed 202 QTL (Table S3) for the 16 measured traits,including respectively 66, 68, and 68 QTL for embryo, kernel, and embryo-to-kernel ratio traits.The phenotypic variation explained by each QTL ranged from 3.4% to 19.1%, with a mean of 7.7%.The 202 QTL identified by SLM were merged into 47 consensus QTL.With JLM (Fig.2), 210 QTL (Table S4) were identified, including 64,78,and 68 QTL for embryo,kernel,and embryo-to-kernel ratio traits, respectively.The 210 QTL from JLM were merged into 78 consensus QTL.With the GWAS (Fig.2), 222 QTL (Table S5) were identified, including 55, 84, and 83 QTL for embryo, kernel, and embryo-to-kernel ratio traits, respectively.The 222 QTL from the GWAS were merged into 81 consensus QTL.

Fig.1.Phenotypic variation of 16 embryo size and related traits in the NAM population.Solid triangles represent parental lines.Embryo traits are indicated with red labels,kernel traits with blue labels, and embryo-to-kernel ratio traits with green labels.HEMW, hundred-embryo weight; EV, embryo volume; EW, embryo width; EL, embryo length; ET, embryo thickness; HKW, hundred-kernel weight; HEDW, hundred-endosperm weight; KV, kernel volume; KW, kernel width; KL, kernel length; KT, kernel thickness; EER, embryo-to-endosperm weight ratio; EVR, embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-to-kernel length ratio; ETR,embryo-to-kernel thickness ratio.

The loci detected by the three methods overlapped extensively,indicating the robustness of the QTL analyses(Fig.S4).The ratio of QTL overlapped between SLM and JLM was 71.1%, the ratio of QTL overlapped between SLM and GWAS was 61.8%, and the ratio of QTL overlapped between JLM and GWAS was 72.9%.Of the QTL identified with SLM, 68.8% were also detected by JLM and/or GWAS.Similarly, 68.6% of JLM QTL and 70.7% of GWAS QTL were identified by the other two methods.These results demonstrated the reliability of the identified QTL and indicated that the integrated use of these three complementary models could improve our ability to elucidate the genetic architecture of embryo size and its related traits.

The 5-22 QTL associated with each trait by JLM explained 4.7%-46.7% of the phenotypic variation, and negligible epistasis was detected across the three families for these 16 measured traits(Fig.3).The 210 QTL identified by JLM were grouped into classes based on their significance level:seven major-effect QTL were detected within each of the three separate families (LOD values of these seven QTL within each separate family and across all three families were >3); 20 and 105 moderate-effect QTL were detected within two separate families and one family, respectively; 78 minor-effect QTL were detected only by joint-linkage analysis using all 611 RILs across all three families simultaneously.

Fig.2.QTL mapping results for 16 embryo size and related traits by JLM (upper) and GWAS (lower) in the NAM population.In the upper panel, colored rectangles indicate joint-linkage QTL regions,and the color density indicates the LOD value of the identified QTL.In the lower panel,triangles indicate QTL identified by the GWAS.Red triangles indicate that the allele in the common parental line BY815 was associated with an increase in the trait,and blue triangles indicate that the BY815 allele was associated with a decrease in the trait.Embryo traits are indicated in red,kernel traits in blue,and embryo-to-kernel ratio traits in green.HEMW,hundred-embryo weight;EV,embryo volume;EW,embryo width;EL,embryo length;ET,embryo thickness;HKW,hundred-kernel weight;HEDW,hundred-endosperm weight;KV,kernel volume;KW,kernel width;KL,kernel length; KT, kernel thickness; EER, embryo-to-endosperm weight ratio; EVR, embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-tokernel length ratio; ETR, embryo-to-kernel thickness ratio.

The 210 QTL identified by JLM could also be grouped by their effect directions:for 134 QTL, the BY815 alleles showed the same effect direction(95 positive and 39 negative)across the three families;for the other 76 QTL,the BY815 alleles showed distinct effect directions within different families.For above 76 QTL, BY815 alleles showed positive or negative effects at the same locus in different separate families,indicating the importance of allele series for variations in embryo size and its related traits (Table 1).In total,78.6% of the QTL identified by JLM explained <5% of phenotypic variation, indicating many minor-effect QTL underlying natural variation in embryo size and its related traits.The detection of a large number of minor-effect QTL was expected, as the NAM approach should increase the detection efficiency of minor-effect QTL by integrating recombination information from multiple crosses [33].

Table 1 QTL identified for 16 embryo size and related traits by JLM.

Based on the JLM results, the 64 QTL for five embryo traits, 68 QTL for five embryo-to-kernel ratio traits,and 78 QTL for six kernel traits were merged into respectively 40,48,and 41 consensus QTL.The ratio of QTL overlapped between embryo traits and kernel traits was 59.3%, the ratio of QTL overlapped between embryo traits and embryo-to-kernel ratio traits was 53.4%, and the ratio of QTL overlapped between kernel traits and embryo-to-kernel ratio traits was 44.5%.The three-way overlap for the above 40,48, and 41 consensus QTL was 31.4% (Fig.S5).The mean effect of EW loci was greater than that of EL loci and ET loci, and the same trend was observed for loci controlling KW, KL, and KT (Fig.4).Thus, although there were close relationships among the 16 embryo and kernel traits and evidence of pleiotropic loci, the genetic architecture underlying each of these traits was distinct.This is similar to the control of male and female inflorescence traits in maize, in which inflorescence traits had larger QTL effects than flowering and leaf traits [37].

Fig.3.Effects of QTL identified for 16 embryo size and related traits in the NAM population.(a) Effect size (represented by PVE) and the origin of the increasing alleles for individual identified QTL.Bars above and below zero indicate that the increasing alleles came from BY815 and other parental lines,respectively.Orange bars indicate major QTL that were detected within each of the three separate families(LOD values of these QTL within each separate family and across three families were>3);red and blue bars indicate moderate QTL that were detected within two separate families and one family, respectively.Green bars indicate minor QTL that were detected by joint linkage analysis only when all 611 RILs across the three families were analyzed simultaneously.(b) Total PVE contributed by individual (blue bars) and epistatic (red bars) QTL for each trait.The green line indicates broad-sense heritability.HEMW, hundred-embryo weight; EV, embryo volume; EW, embryo width; EL, embryo length; ET, embryo thickness;HKW,hundred-kernel weight;HEDW,hundred-endosperm weight;KV,kernel volume;KW,kernel width;KL,kernel length;KT,kernel thickness;EER,embryo-toendosperm weight ratio; EVR, embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-to-kernel length ratio; ETR, embryo-to-kernel thickness ratio.

3.3.Distinct effect patterns of 28 key loci as determined by a trait-QTL association network

Based on the trait-QTL associations and QTL-QTL colocalizations from the GWAS, a trait-QTL association network was constructed, in which 205 nodes and 439 edges were included, and all traits had at least one QTL that was associated with another trait (Fig.5).In agreement with the phenotypic correlation pattern of the traits (Fig.S2), QTL controlling correlated phenotypes, such as EL and KL, tended to co-localize and cluster within more closely connected sub-networks.As expected, traits from the same class, such as kernel weight and size traits, including HKW, KV, KL, and KW, tended to be closely located in the network and showed a high level of connectivity, indicating genetic co-regulation.There were 15 loci, 8 loci, 5 loci, 2 loci, and 1 locus involved in at least four, five, six, seven, and eight traits,respectively.One noteworthy example is theZmBG1locus, which was associated with 11 traits and located in the center of the network.

Fig.4.Comparison of the effect size of QTL identified by GWAS for 16 embryo size and related traits.Background colors indicate different traits,black dots represent detected QTL,orange horizontal lines represent the median effect of the QTL for each trait,and the red dashed lines encircle traits for which a trend of decreasing QTL effect sizes was noted.HEMW,hundred-embryo weight;EV,embryo volume;EW,embryo width;EL,embryo length;ET,embryo thickness;HKW,hundred-kernel weight;HEDW,hundredendosperm weight; KV,kernel volume; KW, kernel width; KL, kernel length; KT, kernel thickness; EER, embryo-to-endosperm weight ratio; EVR, embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-to-kernel length ratio; ETR, embryo-to-kernel thickness ratio.

For the 28 key loci that were associated with at least three traits,four loci displayed a‘‘uE-uK-uR”(upregulated in association with the embryo, kernel, and embryo-to-kernel ratio traits) effect pattern, four a ‘‘uE-uK” pattern, eight a ‘‘uE-uR” pattern, three a‘‘uR”pattern,five a‘‘dK”(downregulated in association with kernel traits)pattern,three a‘‘dK-uR”pattern,and one a‘‘uE-uK-dR”pattern.For the four loci that displayed a ‘‘uE-uK-uR” pattern, the BY815 alleles increased embryo size and the embryo-to-kernel size ratio without reducing kernel size.These consistent effect patterns for both embryo and kernel regulation indicated the potential of simultaneous improvement of maize kernel nutrition and grain yield.In the Beijing High-Oil (BHO) population, kernel oil content reached 15.5%after 17 cycles of selection,with only a slight reduction in kernel weight[42].The common parent line BY815 that was used to construct our NAM population was selected from the BHO population.By pyramiding favorable alleles from multiple loci under comprehensive breeding schemes, the goal of cultivating high-nutrition and high-yield maize varieties might be achieved.

Cloned maize genes includingZmUrb2,Emp12,ZmBAM1d,Dek15,andDek605[43-47]; maize orthologs of known rice genes that control kernel size includingBG1,XIAO,GS9,andIPA1[48-51]; and 11 maize orthologs ofArabidopsis EMBgenes [52] were identified as underlying the 28 key loci (Table S6).A QTL cluster was identified at ~22 Mb of chromosome 1 that controlled seven traits including HEMW, ET, EER, and KT (Fig.6a).The most significant SNP associated with kernel thickness in the GWAS (BPP of 0.48)was located upstream ofZmUrb2, which was recently cloned by mutant analysis [43], and the most significant SNP associated with ET (BPP of 0.5) was 115.8 kb fromZmUrb2.LikeZmUrb2, the candidate geneZmBG1was coincident with a QTL cluster controlling multiple traits (Fig.6b).BG1affects the kernel size of rice and was cloned by mutant analysis [48].In this study,ZmBG1was associated with 11 traits including HEMW and HKW,implying that the function of this gene is conserved among grass crops.The SNP associated with KW (BPP of 0.18) was located upstream ofZmBG1and was 375 bp from the gene.ForZmUrb2andZmBG1,the alleles from BY815 showed positive effects on all associated traits, that is, these two loci showed consistent effect patterns for both embryo and kernel regulation, the ‘‘uE-uK-uR” effect pattern(Table S6).

The locus at ~120 Mb of chromosome 9 was associated with multiple traits including EW, ET, EV, HEMW, KW, KT, KV, and HKW (Fig.6c), and the effects of this locus on the traits were consistent.Thus, this locus may control the variation of HEMW and HKW via regulation of EW and ET.However, none of the three QTL mapping methods detected strong associations between this locus and embryo-to-kernel ratio traits except for a minor EER QTL.A PPR gene was identified as underlying this locus, and this gene is the homolog of the embryo-defective geneEMB975inArabidopsis[52].The locus identified near 20 Mb on chromosome 4 affected embryo size and embryo-to-kernel ratio traits (Fig.6d).This locus conferred consistent effects on embryo traits (EW, EV,HEMW) and embryo-to-kernel ratio traits (EWR, EVR, EER),although it showed no effects on kernel traits.We speculated that this locus regulates embryo size specifically,without affecting kernel size.The maize homolog ofArabidopsisembryo defective geneEMB1395was identified as underlying this locus, and this gene encodes anS-adenosyl homocysteine hydrolase [52].These two novel loci aroundEMB975andEMB1395showed ‘‘uE-uK” and‘‘uE-uR” effect patterns, respectively (Table S6).

3.4.Genetic relationships between embryo and kernel-size traits

Fig.5.Trait-QTL association network for 16 embryo size and related traits.Large and smaller dots indicate traits and their corresponding QTL,respectively.A trait and a QTL are connected by a solid black line in the network if the QTL was associated with that trait.QTL are connected by dashed magenta lines in the network if two QTL co-localized at the same locus.The colors of dots indicate different traits:red for embryo traits, blue for kernel traits, and green for embryo-to-kernel ratio traits.Information for the presented candidate genes is given in Table S6.HEMW, hundred-embryo weight; EV, embryo volume; EW, embryo width; EL, embryo length; ET, embryo thickness; HKW,hundred-kernel weight; HEDW, hundred-endosperm weight; KV, kernel volume; KW, kernel width; KL, kernel length; KT, kernel thickness; EER, embryo-to-endosperm weight ratio; EVR, embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-to-kernel length ratio; ETR, embryo-to-kernel thickness ratio.

High positive correlations among five embryo traits and five embryo-to-kernel ratio traits were observed (Fig.S2), with the exception of the correlation between EL and EWR (r= -0.16,P=8.57×10-5).This negative correlation was expected,as variations in length and width are likely mutually restrictive.High positive correlations were also found among all six kernel traits,except for the correlation between kernel length and thickness (r=-0.19,P=3.97×10-6).Interestingly,there was a marked positive correlation between EW and KW (r= 0.38,P< 2.20 × 10-16),whereas EWR was negatively correlated with KW (r= -0.28,P=1.01 × 10-11), indicating that a narrower kernel tended to have a larger embryo-to-kernel width ratio.Similar relationships among EL,KL,and ELR were observed.EER,the embryo-endosperm weight ratio,is determined by both the embryo and endosperm of an individual kernel.The correlation between EER and embryo weight was strong(r=0.71,P<2.20×10-16),whereas almost no correlation was found between EER and endosperm weight (r= 0.019,P= 0.66), indicating that EER variation in our panel was almost independent of variation in endosperm weight at the population level.

Fig.6.Associations around known maize genes and novel loci.(a)Associations identified at the known locus ZmUrb2.(b)Associations identified at the novel locus ZmBG1,the maize homolog of the known rice gene BG1.(c) Associations identified at the novel locus EMB975.(d) Associations identified at the novel locus EMB1395.The upper plot in each panel shows the association signals from the linkage analysis and the lower plot shows the association signals from the GWAS.Dashed lines indicate the physical position of each candidate gene.HEMW,hundred-embryo weight;EV,embryo volume;EW,embryo width;EL,embryo length;ET,embryo thickness;HKW,hundred-kernel weight; HEDW, hundred-endosperm weight; KV, kernel volume; KW, kernel width; KL, kernel length; KT, kernel thickness; EER, embryo-to-endosperm weight ratio; EVR,embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-to-kernel length ratio; ETR, embryo-to-kernel thickness ratio.

The genetic contributions from embryo to kernel variations(indicated by green dots in Fig.7) were <50%, except for the KLEL pair.EL explained 77.7% of the variation in KL at the genetic level, ET explained 48.1% of the variation in KT, and EW explained 31.1% of the variation in KW.These results suggest that EL and KL share a greater genetic basis than ET and KT or EW and KW, indicating that the variation in EL might be strictly constrained by KL,whereas EW and ET might be less constrained by KW and KT and thus have a greater possible variation.EW explained >80% of HEMW variation at the genetic level, much more than ET (61.9%)or EL (37.5%).For kernel weight variation, kernel width made the greatest contribution, whereas the contribution from the embryo weight was <30%.

4.Discussion

4.1.Distinct genetic architectures for embryo size and kernel size

In this study, the genetic architecture associated with five embryo weight and size traits, six kernel weight and size traits,and five embryo-to-kernel ratio traits was investigated using a NAM population of 611 RILs derived from crosses between a high-oil/giant embryo inbred line and three regular inbred lines.Shared and unique loci were identified for embryo size and kernel size (Fig.S5).More than half of the loci (59.3%) identified for embryo size and kernel size were shared.Finding both shared and unique genetic control with respect to the natural variation of embryo size and kernel size was expected and is consistent with the genetic and functional evidence found in maize kernel mutant analyses that used different categories of mutants, includingdekmutations, which affect both endosperm and embryo, andembmutations,which result in more or less normal endosperm formation and embryo-specific effects [53].Different effect sizes for the QTL underlying embryo size traits were observed:the mean effect of embryo width loci was greater than that of embryo length loci and embryo thickness loci.The effect size trend was also observed for loci controlling kernel width, length, and thickness (Fig.4).

4.2.Candidate genes underlying key loci revealed by the trait-QTL association network

In this study,information from trait-QTL associations and QTLQTL co-localization was integrated to construct a trait-QTL association network that includes 16 embryosize and related traits and 222 QTL identified by GWAS (Fig.5).In this network, 28 key loci associated with at least three associations were revealed,including eight loci associated with at least five traits.Underlying these key loci,cloned maize genes that control kernel development includingZmUrb2,Emp12,ZmBAM1d,Dek15,andDek605[43-47]were identified.For example, the protein encoded byZmUrb2is involved in the biogenesis of ribosomes and is required for kernel development.Mutation of this gene resulted in a dramatic decrease in kernel thickness as well as in embryo size[43].Among the kernel-size traits, only KT was associated withZmUrb2, indicating the robustness and accuracy of the QTL analysis.We also identified maize homologs of cloned rice genes that affect seed development as underlying the key loci, includingBG1,XIAO,GS9,andIPA1[48-51].One noteworthy example is theZmBG1locus,which was associated with 11 traits, including embryo weight and kernel weight,and is located in the center of the network.BG1affects the kernel size of rice by regulating auxin transport, and this gene has been cloned by mutant analysis [48].We also identified 11 novel loci that were maize homologs ofArabidopsis EMBgenes (Fig.5;Table S6).These 11 genes inArabidopsiscontrol seed development,and their functions have been confirmed by knockout using the TDNA insertion method [52].

Compared with previous studies, which focused only on embryo size [16,17] or kernel size [54,55], the advantage of our approach is the simultaneous analysis of the genetic architectures of embryo size and kernel size.The 16 embryo size and related traits were grouped into three categories:embryo, kernel, and embryo-to-kernel ratio traits.For the 28 key loci associated with at least three traits in the trait-QTL association network,five were associated with traits from three categories, 15 with traits from two categories, and the remaining eight loci with traits from only one category.For example,ZmUrb2andZmBG1were associated with seven traits (including four from the embryo trait group,two from the embryo-to-kernel ratio trait group, and one from the kernel trait category) and 11 traits (including three from the embryo trait group, three from the embryo-to-kernel ratio trait group, and five from the kernel trait group), respectively (Fig.5).ZmUrb2andZmBG1, which showed consistent effect patterns for both embryo and kernel regulation, might be exploited in future breeding programs.We also identified three loci that showed the‘‘dK-uR” effect pattern,ZmIPA1,EMB2296, andEMB2369.For these three loci, the alleles from BY815 increased the embryo-to-kernel ratio,whereas they reduced the kernel size.The key loci associated with multiple traits,especially those that resulted in the repulsion effect pattern, may correspond to pleiotropic regulation or closely linked genes.The key loci identified by the trait-QTL association network, with their specific effect patterns, offer greater application potential in breeding programs than the single-trait QTL reported previously [16,17].

4.3.Variation in embryo size and improvement of nutrition and grain yield

In previous high-oil maize breeding programs, artificial selection for higher oil content led to an increase in the proportion of the embryo in the maize kernel [5,6,42].In this NAM population,EW explained >80% of HEMW variation at the genetic level, much more than ET (61.9%) or EL (37.5%), suggesting that embryo width makes the greatest contribution to embryo weight variation at the genetic level (indicated by red dots in Fig.7).Embryo weight contributed~50%of the variation in EER.As EER is a direct measure of the proportion of embryo tissue in the kernel,it should be possible to regulate EER effectively by modifying EW loci to further increase embryo size.The contribution of embryo weight(HEMW)to kernel weight (HKW) variation at the genetic level is small (27.1%), suggesting that variation in kernel weight is due mainly to the starch-rich endosperm tissue and in turn that the kernel allows us space for increasing the embryo-to-kernel ratio by pyramiding the favorable alleles at loci that showed the‘‘uE-uR”effect pattern to further increase kernel oil and other nutritional components.

Negative correlations between kernel oil content and starch concentration and kernel size were observed in Illinois High-Oil and Alexho Synthetic populations [5,6].This phenomenon is consistent with the negative correlations observed between EWR and KW and between ELR and KL in our NAM population,and could be explained by loci showing the repulsive effect pattern‘‘dK-uR”.It should thus be possible to increase kernel oil content moderately by increasing the oil concentration in the embryo without changing embryo size to avoid the negative effects of smaller endosperm and seed size.A successful example of this approach is the application of the favorable allele ofDGAT1-2,which was used to improve the hybrid Zhengdan 958[56,57].In this hybrid,kernel oil content was increased from 3.8% to 4.5% without affecting the yield.

Fig.7.Phenotypic and genetic relationships between embryo size and related traits.Traits were considered in pairs,and the variation explained by the second trait in each pair with respect to the first trait was calculated at the phenotypic(y axis)and genetic(x axis)level,respectively.The colors of dots indicate different trait combinations,and the dot size corresponds to the total phenotypic variation explained by all the QTL identified in the GWAS for the first trait in each pair.The level of the confidence interval used for linear model fitting was 95% and is indicated by gray shading.HEMW, hundred-embryo weight; EV, embryo volume; EW, embryo width; EL, embryo length; ET,embryo thickness; HKW, hundred-kernel weight; HEDW, hundred-endosperm weight; KV, kernel volume; KW, kernel width; KL, kernel length; KT, kernel thickness; EER,embryo-to-endosperm weight ratio; EVR, embryo-to-kernel volume ratio; EWR, embryo-to-kernel width ratio; ELR, embryo-to-kernel length ratio; ETR, embryo-to-kernel thickness ratio.

Maize starch production accounts for 91% of global starch production,and the demand for maize starch production is increasing for both food and non-food applications such as the bioethanol industry[58,59].In contrast to breeding efforts for high-oil maize,high-starch maize cultivars with higher yield potential could be developed by reducing embryo size to increase the endosperm by exploiting the unfavorable alleles at the loci that showed the ‘‘uEuR” effect pattern.For example, theZmGE2locus on chromosome 1, which affects natural variation of EER and kernel oil content[17], was associated with multiple embryo traits including EW,EWR, and EER in the present study (Fig.S6).This locus was also associated with kernel traits including KW, HEDW, and HKW.Although the allele from BY815 increased the size of the embryo,it also decreased endosperm and kernel weight.Interestingly,two other genes,ZmGS3andZmDof3, were identified surroundingZmGE2.ZmGS3was reported[60]to affect kernel size based on candidate gene association analysis, andZmDof3was reported [61] to regulate the development of the endosperm in maize kernels based on mutant analysis.Thus, the effects from this locus on chromosome 1 were likely caused by multiple genes, and the alleles that reduce embryo size and increase endosperm and kernel weight might be exploited to develop high-starch maize.

Besides identifying candidate genes for 28 key loci that control the natural variation in maize embryo size and its related traits,our study revealed the genetic relationship between embryo size and kernel size and suggested the role of embryo size in controlling kernel nutrition components and grain yield.These findings shed light on the genetic basis of embryo size and its related traits in maize and provide loci for potential use by maize breeders to manipulate embryo size in developing ideal cultivars.

CRediT authorship contribution statement

Xiaohong Yang and Jiansheng Lidesigned the research plan;Xiaowei Li,Hui Fang,and Xiuyi Fuperformed the research work;Xiaowei Li,Min Wang,and Renyu Zhanganalyzed the data;Xiaowei Li and Jiansheng Liwrote the manuscript.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgment

The research was supported by the National Natural Science Foundation of China (31421005).

Appendix A.Supplementary data

Supplementary data for this article can be found online at https://doi.org/10.1016/j.cj.2021.03.007.