Four maize inbred lines (X987F, AJ5001, AJ7001, and AJ9010) from Maize Yufeng Biotechnology LLC were collected to construct three DH populations with X987F as the common parent, forming a NAM population. The plants were cultivated using a randomized complete block design across three different locations in 2022: Beijing (BJ, 40°08′N, 116°10′E), Neimeng (NM, 40°31′N, 107°05′E), and Liaoning (LN, 40°82’N, 123°56’E). Each line was grown in single-row plots, 250 cm in length, with 11 plants per row and 60 cm spacing between rows, under natural field conditions. The 11 plants in each row were self-pollinated, and 300 kernels were bulk-sampled from five moderate-sized ears, with an equal number of kernels collected from each ear. These bulk kernels were then used for phenotyping.

Phenotypic data were collected on three phenotypic quality traits of maize, including starch content, oil content, and lysine content. Specifically, a near-infrared reflectance (NIR) spectrometer (DA 7250, Perten Instruments Inc., Sweden) was utilized to assess the quality traits in maize kernels. Each sample underwent scanning three times to obtain an average value. The phenotypic variation of the quality traits was analyzed using R software version 4.0.1 with the “AOV” function (ANOVA). The ANOVA model employed is expressed as y = μ + α_g + β_e + ϵ, where α_g is the effect of the g^th line, β_e is the effect of the e^th environment, and ϵ is the error. All of the effects were considered to be random. These variance components were used to calculate the broad-sense heritability as h² = σ _g ² /(σ _g ² + σ _ϵ ²/e), where σ _g ² is the genetic variance, σ _ϵ ² is the residual error, and e is the number of environments. To eliminate the influence of environmental effects, the best linear unbiased predictor (BLUP) value for each line was calculated using a linear mixed model that considered both genotype and environment as random effects in the R function “LME4,” and the trait BLUP values were used for subsequent analyses.

The genotype of the three DH populations with their parents was obtained by utilizing the GenoBaits Maize 1K marker panel (Mol Breeding Biotechnology Co., Ltd., Shijiazhuang, China). A total of 4,589 SNP markers were identified based on genotyping using a target sequencing platform. The minor allele frequency (MAF) and missing rate were estimated for each population, with SNPs exhibiting MAF < 0.1 or a missing rate > 0.6 being filtered out (Purcell et al., 2007). Following quality control measures, the polymorphic SNPs between the two parental lines were utilized to construct genetic linkage maps through the functions est.rf and est.map from the R/qtl package, employing the Kosambi mapping method (Broman et al., 2003). A joint linkage map for the three DH populations was developed using the CMP module in QTL IciMapping software (version 4.2).

Using the genetic linkage maps derived from three DH populations, SLM was conducted utilizing composite interval mapping as previously described (Zeng, 1994; Hu et al., 2021), which was implemented in Windows QTL Cartographer 2.5, according to the methodology of Wang et al. (2010) for each DH population. Model 6 of the Zmapqtl procedure (Basten et al., 1997) in the Composite Interval Mapping module was employed to identify QTLs across the genome by surveying intervals of 1.0 cM between markers within a window size of 10 cM. Forward–backward stepwise regression with five controlling markers was used to account for background contributions from adjacent markers. To determine the threshold logarithm of odds (LOD) value for putative QTLs, 1,000 permutations were performed for each trait in each DH population, and the resulting LOD score threshold ranged from 2.76 to 3.27 (α = 0.05). For simplification purposes, a LOD score threshold of 3.0 was adopted globally. The confidence interval concerning the position of the QTL was estimated using the 1.5-LOD support interval method as articulated in Liu’s study (Liu et al., 2017). Additionally, we utilized the R function “LM” to calculate total PVE by significant individual QTLs.

The NAM module in QTL IciMapping (version 4.2) was used to conduct the nested association mapping. The mapping approach employed was ICIM-ADD (Li et al., 2011a), with a step size set at 1.0 cM. In the same way, 1,000 permutations were performed for each trait. The resulting LOD score threshold ranged from 4.84 to 5.36 (α = 0.05), and we adopted a LOD score of 5.0 as the global threshold.

All genes located within the QTL support interval were extracted based on their physical position using the Zea mays L. B73 reference genome version 4.0. The functional annotation of genes and information on starch, oil, and lysine metabolism pathways were obtained from the maizeGDB database (http://www.maizegdb.org/). The expression levels of candidate genes in different tissues are derived from public transcriptome data (Yi et al., 2019).

Four inbred lines exhibiting starch content ranging from 70.99% to 74.22%, oil content between 4.01% and 5.04%, and lysine content varying from 0.27% to 0.35% were collected for the development of three DH populations. The three DH populations, namely, DHPop1, DHPop2, and DHPop3, comprised 163, 219, and 320 lines, respectively ( Table 1 ). The distributions of all traits within the DH populations showed continuous normality without significant skewness ( Figures 1A–C ). Analysis of variance (ANOVA) revealed that genotype variance exceeded environmental variance across almost all populations ( Table 1 ), indicating that phenotypic variations were predominantly controlled by genetic factors. In addition, starch content, oil content, and lysine content exhibited average broad-sense heritabilities of 78.40%, 71.69%, and 74.89%, respectively ( Table 1 ). Therefore, the abundant phenotypic variation along with high heritability observed in maize kernel quality traits across the three DH populations provides a solid genetic foundation for identifying new QTLs in maize.

The three DH populations were genotyped by using the GenoBaits Maize 1K marker panel including 4,589 SNP markers. After quality control, DHPop1, DHPop2, and DHPop3 contained 1,266, 1,210, and 1,249 high-quality SNP markers, respectively, and covered all ten maize chromosomes. Based on the reference parental polymorphic loci, three linkage maps were independently constructed, with genetic map lengths of 728, 734, and 749 cM ( Supplementary Figure S1 ). The average genetic distance between adjacent markers was 0.58, 0.61, and 0.62 cM in each DH population, respectively. The combined linkage map for the three DH populations using the CMP module in QTL IciMapping, included 1,993 SNP markers spanning a total genetic distance of 1,343 cM.

We performed SLM analysis in each DH population with the composite interval mapping method (Zeng, 1994). In total, 37 QTLs for quality traits were detected, including 14 QTLs for starch content, 13 QTLs for oil content, and 10 QTLs for lysine content ( Supplementary Table S1 , Supplementary Figure S2 ). The 1.5-LOD supporting QTL interval averaged 46.26 cM, with a range of 1.7–168.2 cM. The total PVE by all identified QTLs in a population ranged from 20.67% to 53.51% for starch content, 28.71% to 33.98% for oil content, and 7.53% to 37.13% for lysine content ( Supplementary Table S1 ). The PVE for the three traits showed far less than broad-sense heritability ( Figure 2A ), which suggested that some minor QTLs for these traits cannot be detected in bi-parent populations. Of these QTLs, 75.68% had the PVE <10% ( Supplementary Table S1 ). For starch content, the PVE for each QTL ranged from 3.97% (qSC14 in DHPop3) to 17.24% (qSC04 in DHPop1), and only 28.57% (4/14) of the QTLs had large effects with PVE ≥ 10% ( Figure 2A , Supplementary Table S1 ). For oil content, the PVE for each QTL ranged from 3.16% (qOiL13 in DHPop3) to 14.27% (qOiL07 in DHPop2), and only 23.08% (3/13) of the QTLs had large effects with PVE ≥ 10% ( Figure 2A , Supplementary Table S1 ). For lysine content, the PVE for each QTL ranged from 5.14% (qLys05 in DHPop3) to 12.06% (qLys08 in DHPop3), and only 20.00% (2/10) of the QTLs had large effects with PVE ≥ 10% ( Figure 2A , Supplementary Table S1 ). For each DH population, both parents in each DH population contained the alleles that increased starch content, oil content, and lysine content, respectively ( Figures 2B–D , Supplementary Table S1 ). In addition, the QTL co-localization analysis among these three populations showed partial overlaps among more than two populations. Moreover, 10 QTLs uniquely detected in a given population underscored the genetic diversity of the founders of the DH populations ( Figure 2E ).

Subsequently, we identified a total of 17 QTLs using the NAM analysis, including seven QTLs for starch content, four QTLs for oil content, and six QTLs for lysine content ( Table 2 , Supplementary Table S2 ). The total PVE averaged 5.42%, with a range of 2.60%–9.72% in all populations for the three traits. Compared with the SLM results, the QTL interval was expectedly small, with 82.35% (14/17) of the QTL intervals falling within 50 Mb. The QTL co-localization analysis showed that almost all QTL from NAM were overlapping intervals with these from SLM, which indicates that the QTLs mined from the three DH populations established are representative ( Figure 2E ).

The genes located within the QTL support interval were extracted, resulting in the identification of 20,367 genes ( Supplementary Table S3 ). According to the annotation in the MaizeGDB database (www.maizegdb.org), genes with the most relevant functional annotation information was nominated as the candidate genes. Phosphoglucomutase2 protein-coding gene pgm2 was the most abundant annotations gene for starch content positioned in 11,195,010–11,201,017 interval of chromosome 5 ( Figure 3A , Supplementary Table S3 ). Expression pattern analysis showed that pgm2 was expressed in all tissues of maize, and the expression levels were relatively the same except in pollen ( Figure 3B ). There were two annotations genes for oil content. Acyl carrier protein coding gene ACP was positioned in 252129992–252138538 interval of chromosome 1 and expressed in all tissues ( Figures 3C, D ; Supplementary Table S3 ). Acyl-coenzyme A oxidase coding gene ACOX was positioned in 7209473–7214489 interval of chromosome 4 ( Figure 3E , Supplementary Table S3 ). ACOX was also expressed in almost all tissues of maize except in pollen ( Figure 3F ). In addition, the expression levels in the embryo and endosperm vary with different developmental stages ( Figure 3F ). For lysine content, glyceraldehyde-3-phosphate dehydrogenase 4 coding gene gpc4 and Acyl-activating enzyme-like protein coding gene o7 were positioned in 186104976–186109690 interval of chromosome 5 and 149382212–149384169 interval of chromosome 10, respectively ( Figure 3G, I ; Supplementary Table S3 ). gpc4 was expressed in all tissues of maize, and the expression level of o7 was significantly lower in the shoot apical meristem (SAM), ear, cob, and embryo tissues compared to other tissues ( Figures 3H, J ).

SNP markers represent the most prevalent variations within genomes, and their application in plant breeding has significantly enhanced the precision of QTL mapping and genetic analysis (Bhattramakki et al., 2002; Mammadov et al., 2012; Kaur et al., 2021). By conditioning linked markers during testing, the sensitivity of the test statistic to the positions of individual QTLs is heightened, thus improving the accuracy of QTL mapping (Zeng, 1994). With advancements in sequencing technology, an increasing array of molecular markers has been utilized for QTL mapping, leading to substantial improvements in its accuracy (Schnable et al., 2009; Chia et al., 2012; Bukowski et al., 2018; Fang et al., 2021). In this study, the average physical intervals of the QTLs were 42.11 Mb; 29.79% (14/47) spanned physical intervals of <10 Mb, and 68.09% (32/47) spanned <50 Mb. These findings indicate a considerable enhancement in resolution due to both the large number of employed markers and the suitable population type used in this research. The resolution is likely on the order of 2–3 cM, since pairs of markers that are farther apart rarely exhibit significant levels of linkage disequilibrium.

The maize kernel quality traits including starch content, oil content, and lysine content are complex quantitative traits, and the genetic basis are still unresolved issues at present (Holding et al., 2008; Huang et al., 2021). A set of parameters that elucidates the genetic component underlying trait variation within or among populations will influence the genetic architecture of quantitative traits. Mapping populations represent one of the critical parameters capable of capturing genetic diversity and possessing the power to detect QTLs with small effects (Odell et al., 2022). DH populations have the advantages of retaining homozygosity and genetically identical replicates (Bordes et al., 2006; Foiada et al., 2015; Yan et al., 2017; Chaikam et al., 2019) and nested association mapping simultaneously exploits the advantages of both linkage and association mapping and has been successfully applied to investigate the genetic basis of complex quantitative traits in maize (Yu et al., 2008; Gage et al., 2020). In this study, we developed a NAM population formed by three DH populations for detection of QTLs and identification of candidate genes. The results of the phenotypic and genetic detection showed that there was a wide phenotypic variation extent in starch content (64.64%–76.83%), oil content (3.60-5.96%), and lysine content (0.17-0.46%) and high broad-sense heritability in the DH populations. It illustrated that these traits were mainly controlled by genetic factors. Totally, 47 QTLs were found in the three DH populations and distributed on chromosomes 1–10. Among the identified QTLs, three QTLs (qSC02 in DHPop1 and qSC09 in DHPop3 for starch content, qLys06 in DHPop3 for lysine content) spanned a 104.93 Mb physical interval on chromosome 4. Four QTLs (qSC07 in DHPop2 and qSC11 in DHPop3 for starch content, qOiL08 in DHPop3 and qOiL12 in DHPop3 spanned 86.51 Mb physical interval on chromosome 6. These results indicated that there is cross-over in the regulatory network of starch content, oil content, and lysine content in maize kernel. The average PVE per QTL was 7.89% using SLM in each DH population, whereas it was 5.42% using NAM method across three DH populations. Similar to previous studies, our findings align with the quantitative nature of these quality traits, which are acknowledged to be governed by a multitude of genes and QTLs with minor effects.

For starch content, pgm2 (Zm00001d013428) located in qSC10 was the most abundant annotations among the annotated genes ( Figure 3A ), which encoded phosphoglucomutase 2. Phosphoglucomutase (PGM) is a phosphoenzyme (EC 5.4.2.2) and catalyzes an important trafficking point in carbohydrate metabolism in cells of prokaryotic and eukaryotic organisms (Manjunath et al., 1998). In maize, cytosolic isozymes of PGM are encoded by pgm1 and pgm2 (Goodman et al., 1980; Stuber and Goodman, 1983), and the existing levels of PGM are sufficient to maintain the flux of Glc-1-phosphate into glycolysis under O₂ deprivation (Manjunath et al., 1998).

There were three most abundant annotations genes for oil content ( Figure 3A ). ACP (Zm00001d033149) in qOiL10 interval was annotated to be one of the important acyl carrier protein encoding gene. ACP is linked to the malonyl group via a transacylation reaction when acetyl-CoA is carboxylated to malonyl-CoA and participates in de novo fatty acid synthesis in plants and prokaryotes (Kalinger and Rowland, 2023). In qOiL03 interval, there were two abundant annotations genes both related to lipid synthesis metabolic pathways. ACOX (Zm00001d048890) was identified as acyl-coenzyme oxidase encoding gene that involved in the fatty acid β-oxidation pathway in plant peroxisomes (Li et al., 2024; Ruan et al., 2024). ACOX family members may specifically recognize distinct chain lengths of fatty acids and catalyze the first step of peroxisomal fatty acid β-oxidation to converse fatty acyl-CoAs to trans-2-enoyl CoAs (Chu et al., 1994). This step is thought to be critical the rate of carbon flux in the β-oxidation pathway (Pinfield-Wells et al., 2005).

Zm00001d017121 in qLys08 is one of the small multi-gene family, which encoded a glyceraldehyde-3-phosphate dehydrogenase. In maize cytoplasm, glyceraldehyde-3-phosphate dehydrogenase protein is synthesized in roots during anoxic conditions and is known to be one of the “anaerobic polypeptides” (Manjunath and Sachs, 1997). In plants, glyceraldehyde-3-phosphate dehydrogenase (GAPDH), a key enzyme in the glycolytic pathway, reversibly converts the glyceraldehyde-3-phosphate to 1,3-bisphosphoglycerate by coupling with the reduction in NAD1 to NADH (Sirover, 1997) and plays an important role in several cellular processes, including plant hormone signaling, plant development, and transcriptional regulation (Kim et al., 2024). Zm00001d026649 in qLys04 was annotated as opaque endosperm7 (o7) gene. o7 is one of the three most important high-lysine mutants in maize, alongside opaque2 (o2) and floury2 (fl2) (Misra et al., 1972; Miclaus et al., 2011). The o7 mutant endosperm was characterized as having significantly more lysine due to a general reduction in zein levels and an increase in other lysine-rich proteins like the albumins and globulins (Misra et al., 1972; Azevedo et al., 2004).