Abstract Brace roots are the primary organs for water and nutrient absorption, and play an important role in lodging resistance. Dissecting the genetic basis of brace root traits will facilitate breeding new varieties with lodging resistance and high yield. In present study, genome-wide association study (GWAS) and genomic selection (GS) for brace root penetrometer resistance (PR), root number (RN), and tier number (TN) were conducted in a multi-parent doubled haploid (DH) population. Although the three brace root traits were significantly influenced by genotype-by-environment interactions, relatively high heritabilities were observed for RN (79.25%), TN (81.00%), and PR (74.28%). At the threshold of P-value < 7.42 × 10^− 6, 22 SNPs distributed on all 10 chromosomes, except for chromosome 7, and 50 candidate genes were identified using the FarmCPU model. Zm00001d028733, Zm00001d002350, Zm00001d002349, and Zm00001d043694 may be involved in brace root growth. The prediction accuracies of RN, TN, and PR estimated by five-fold cross-validation method with genome-wide SNPs were 0.39, 0.36, and 0.51, respectively, which can be greatly improved by using 100–300 significant SNPs. This study provides new insights into the genetic architecture of brace root traits and genomic selection for breeding lodging resistance maize varieties. Keywords: Maize, Brace root, Lodging resistance, Genome-wide association study, Genomic selection Subject terms: Genetics, Molecular biology, Plant sciences Introduction Lodging poses a significant threat to maize grain yield, ranking as one of the key factors leading to a decrease in photosynthesis, nutrient transport, and grain quality^[38]1. Globally, the annual yield loss caused by lodging is approximately 5–20%^[39]2. Maize brace roots, which penetrate the soil and slant towards the ground, play a crucial role in plant anchorage and lodging resistance^[40]3,[41]4. Dissecting the genetic basis of brace roots contributes to breeding new varieties with lodging resistance and high yield. Previous studies have reported numerous QTL (QTL, Quantitative Trait Locus) related to brace root traits. Ku et al. (2012) conducted linkage mapping for total brace root tier number (TBRTN) and effective brace root tier number (EBRTN)^[42]5. Ten QTL on chromosomes 1, 3, 5, 6, 7, 9, and 10 were identified in a RIL (Recombinant Inbred Lines) population and 15 QTL on chromosomes 1, 2, 5, 6, 7, 8, 9, and 10 were identified in an IF[2] population. A B[2]F[1] population developed from the cross between Yi17 and Yi16 were used for explore QTL relate to seven maize brace root traits^[43]6. A total of 21 QTL were detected for root dry weight, brace root dry weight, brace-root tier number, brace-root surface area, brace-root average diameter, brace-root volume, and brace root total length. The genomic regions at bins 3.05 and 8.04-8.05 were rich in QTL for brace root traits. Zhang et al. (2018) identified three QTL for the number of tiers, three QTL for the total number of roots, and two QTL for the radius of branched roots in 207 recombinant inbred lines of BY815/K22^[44]7. Although, linkage mapping is a powerful tool for QTL analysis, the mapping resolution is low and fine mapping is required. Genome-wide association study (GWAS), as another potent tool for genetic analysis with high mapping resolution, has been widely used in maize^[45]8,[46]9. Wang et al. (2021) detected 126 SNPs (SNP, Single Nucleotide Polymorphism) and 113 candidate genes for crown-root angle, crown-root diameter, and crown-root number using GLM, MLM, FarmCPU, FASTmrMLM, and FASTmrEMMA models^[47]10. A GWAS of bracing root angle (BRA) and diameter (BRD) were conducted by both MLM and BLINK models in 508 inbred lines of maize. Ten SNPs and 27 candidate genes were identified^[48]6. These identified loci can improve the breeding efficiency of new varieties with lodging resistance by marker assisted selection (MAS). Genomic selection (GS), using whole genome markers for selection, is an upgrade of MAS. It was first proposed by Meuwissen et al. in 2001^[49]11. RRBLUP is one of the most used models for GS in maize. A best linear unbiased prediction (BLUP) model was established in a training population with genotype and phenotype data to estimate the effect of each marker, and to predict the breeding value of individuals in the prediction population based on genotypes^[50]12. GS has the potential to greatly accelerate genetic gains per unit time and reduce the costs of plant breeding programs involving complex traits^[51]13. However, no research has been reported on GS of maize brace root traits. Brace root penetrometer resistance (PR), root number (RN), and tier number (TN) are important components for brace root. However, the genetic basis is still unclear. In this study, PR, RN, and TN of 379 multi-parent doubled haploid (DH) lines was evaluated in four environments for GWAS and GS. The objectives of our study were as follows: (1) to identify SNPs and candidate genes associated with PR, RN, and TN through GWAS; (2) to explore the potential of GS in breeding resistance lodging resistant varieties; and (3) to assess the impact of various factors on prediction accuracy. Results Genetic diversity and heritability of brace roots traits A total of 379 DH lines from 21 different maize hybrids were evaluated for PR, RN, and TN of brace root in four environments. PR, RN, and TN showed wide phenotypic variation across environments (Fig. [52]1). PR ranged from 19.99 N in QT to 38.66 N in DF, with a BLUP value of 30.41 N across environments. RN ranged from 14.79 in LD to 18.40 in SG, with a BLUP value of 16.52 across environments. TN ranged from 1.36 in SG to 1.63 in LD, with a BLUP value of 1.43 across environments. Fig. 1. [53]Fig. 1 [54]Open in a new tab Boxplot of PR、RN and TN in the four environments. (a) PR, brace root penetrometer resistance; (b) RN, number of brace roots; (C) TN, brace root tier number. *** represents significant at P ≤ 0.001. Analysis of variance (ANOVA) showed that genetic variance, and genotype × environment (G × E) interaction variance were significant (p < 0.01) across environments. This indicating that PR, RN and TN were not only controlled by genetic factors, but also greatly influenced by genotype × environment interactions. The broad-sense heritabilities (H^2) of PR, RN and TN were 74.28%, 79.25% and 81.00%, respectively, indicating that the phenotypic data was reliable for GWAS and GS (Table [55]1). Table 1. Descriptive statistics and heritability of brace root traits. Trait^1 Mean Env Min Max SD^2 CV Variance components^3 H^2(%) Inline graphic Inline graphic Inline graphic PR 42.90 SG 37.04 47.93 1.91 0.04 186.79*** 225.56*** 66.08*** 74.28 46.29 DF 38.66 53.87 2.49 0.05 38.08 LD 32.67 44.31 2.08 0.05 21.18 QT 19.99 23.08 0.6 0.03 37.08 BLUP 30.41 41.99 1.81 0.05 RN 18.40 SG 12.09 26.08 2.08 0.11 23.80*** 18.92*** 6.40*** 79.25 17.01 DF 11.08 23.06 2.1 0.12 14.79 LD 10.86 19.31 1.53 0.1 15.92 QT 13.49 18.52 0.87 0.05 16.52 BLUP 14.35 20.25 0.85 0.05 TN 1.36 SG 0.74 2.1 0.26 0.19 1.34*** 1.01*** 0.25*** 81.00 1.37 DF 0.4 2.34 0.34 0.25 1.63 LD 1.01 2.15 0.2 0.12 1.38 QT 0.78 2.23 0.23 0.17 1.43 BLUP 0.92 2.03 0.22 0.15 [56]Open in a new tab PR brace root penetrometer resistance, RN number of brace roots, TN brace root tier number, SD standard deviation. Inline graphic , genotypic variance, Inline graphic genotypic × environment interaction variance, Inline graphic error variance. **Significant at P < 0.01. ***Significant at P < 0.001. Correlation of brace root traits with other agronomic traits The results of correlation analysis are shown in Fig. [57]2. For the three brace root traits, TN was significant positive correlate with PR and RN. The Pearson correlation coefficient between RN was -0.008, which was not significant. PR was positive correlate with PH, LW, and HKW. RN and TN were significant positive correlate with all the five agronomic traits. All the three brace root traits were positive correlated with HKW, indicating that brace roots have important impacts on yield. Fig. 2. [58]Fig. 2 [59]Open in a new tab Schematic diagram and correlation analysis for brace root traits in maize. (a) Schematic diagram of maize brace rooting field; (b) Correlation between brace root traits and other agronomic traits. PR, Brace root penetrometer resistance; RN, number of brace roots; TN, number of brace root layers; PH, plant height; EH, ear height; LW, ear leaf length; LL, ear leaf width; HKW, hundred kernel weight; * represents significant at P ≤ 0.05, ** represents significant at P ≤ 0.01, *** represents significant at P ≤ 0.001. Genome-wide association analysis of brace root traits related to lodging resistance Based on the BLUP values across environments, GWAS was conducted using the FarmCPU method with a threshold of P-value < 7.42 × 10^−6. The quantile-quantile (q-q) plots showed that population structure was well controlled using the FarmCPU method (Fig. [60]3). Fig. 3. [61]Fig. 3 [62]Open in a new tab Manhattan and QQ plots of brace root traits in maize. The dashed line is the significnt threshold at P-value = 7.42 × 10^−6. Eight significant SNPs located at bins 1.09, 3.07, 4.05, 4.09, 6.01, 6.07, and 9.02 were identified for PR. The P-value of significant SNPs ranged from 8.96×10^−12 to 6.59 × 10^−6, the MAF ranged from 0.08 to 0.47, and the additive effect ranged from − 0.55 to 0.49. The most significant SNP 4_213271033 was located on chromosome 4 with a MAF of 0.47, and a SNP effect of 0.50. Four significant SNPs located at bins 5.04, 6.06, 9.05, and 10.04 were detected for RN. The P-value of these SNPs ranged from 7.59 × 10^−8 to 2.69 × 10^−6, MAF ranged from 0.15 to 0.50, and additive effect ranged from − 0.23 to 0.25. The most significant SNP 5_145812474 was located on chromosome 5 with a MAF of 0.42 and a SNP effect of 0.23. Ten significant SNPs located at bin 1.03, 2.02, 2.05, 3.02, 3.05, 4.07, 4.09, 6.06, 8.03, and 10.06 were detected for TN. The P-value of the SNPs ranged from 3.63 × 10^−11 to 2.69 × 10^−6, and MAF ranged from 0.15 to 0.50, and additive effect ranged from − 0.06 to 0.07. The most significant SNP 10_141744253 was located on chromosome 10 with a MAF of 0.35 and a SNP effect of 0.06. Haplotypes analysis Haplotype analysis of the loci that have been localised for maize aerial root puncture resistance strength showed that the SNP marker (6_169021369) located on chromosome 6 has one Block. Haplotype analysis of SNP loci that have been localised for maize aerial root number showed that the SNP marker (10_117250342) located on chromosome 10 has one Block (Fig. [63]4). The frequency of haplotypes corresponding to SNPs on chromosome 6 was 84.43% in the test material, containing 320 self-crosses; The frequency of haplotypes corresponding to SNPs on chromosome 10 was 78.10% in the test material, containing 296 self-crosses. Differential analysis of the self-crosses containing SNP markers showed that for 6_169021369, the PR of maize plants with the CC allele was significantly greater than that of maize plants with the TT allele; For 10_117250342, the RN of maize plants with the CC allele was significantly greater than the PR of maize plants with the AA allele (Fig. [64]5). Fig. 4. [65]Fig. 4 [66]Open in a new tab Block analysis of the located SNPs. Fig. 5. [67]Fig. 5 [68]Open in a new tab Differential analysis of SNP marker-containing inbred lines. Candidate genes associated with significant SNPs A total of 50 candidate genes were identified within 10 kb upstream and downstream of the significant SNPs. Thirty-three candidate genes were annotated with known functions, including 15 candidate genes for PR, 4 for RN, and 14 for TN (Table [69]2). GO (GO, Gene Ontology) analysis showed that these candidate genes were mainly involved in metabolic and synthetic processes of biological systems. Genes affecting cellular components were mainly concentrated in intracellular structures and various organelles. Genes affecting molecular functions were mainly related to catalytic activity and binding (Fig. [70]6A). KEGG (KEGG, Kyoto Encyclopedia of Genes and Genomes) analysis identified 16 pathways including purine metabolism, pyrimidine and porphyrin metabolism, metabolic pathways, secondary metabolite biosynthesis, starch and sucrose metabolism as well as biosynthesis of triterpenoids and sterols. Secondary metabolite biosynthesis and metabolic pathways appeared to be particularly important and may be related to brace root traits (Fig. [71]6B). Table 2. Significant SNPs and candidate genes associated with the brace root. Trait SNP^1 P-value MAF^2 SNP effect^3 Candidate gene Gene annotation PR 1_253175968 5.97E-06 0.47 0.29 Zm00001d033181 Fructokinase-like 2, chloroplastic Zm00001d033180 Cytochrome P450 superfamily protein PR 3_207562844 6.51E-06 0.37 0.31 Zm00001d043694 Putative magnesium transporter NIPA9 Zm00001d043693 Casein kinase 1-like protein 2 PR 4_139080812 4.55E-07 0.08 -0.55 Zm00001d051053 Uncharacterized LOC100193062 PR 4_213271033 8.96E-12 0.47 0.50 Zm00001d053091 Adagio protein 3 Zm00001d053090 14-3-3-like protein PR 6_167242430 6.60E-06 0.17 0.36 Zm00001d038932 Histone chaperone ASF1B Zm00001d038933 Ribosomal protein L7/L12-like Zm00001d038931 Syntaxin-71 Zm00001d038934 Unknown PR 6_169021369 1.02E-08 0.35 -0.41 Zm00001d039038 Eukaryotic translation initiation factor 3 subunit D Zm00001d039037 Mitogen-activated protein kinase 7 Zm00001d039039 Cold acclimation protein COR413-TM1 Zm00001d039040 Chlorophyll a-b binding protein, chloroplastic PR 6_67355900 8.52E-07 0.41 -0.34 Zm00001d036031 Serine/threonine-protein phosphatase 2A 65 kDa regulatory subunit A PR 9_12697915 1.08E-06 0.25 0.37 Zm00001D045104 Ubiquitin carboxyl-terminal hydrolase 23 Zm00001d045105 Unknown RN 5_145812474 7.59E-08 0.42 0.23 Zm00001d016131 Transcription factor GTE4 Zm00001d016130 GDSL esterase/lipase At5g45910 RN 6_161969932 2.69E-06 0.17 -0.23 Zm00001d038660 Uncharacterized LOC100284420 Zm00001d038661 Signal recognition particle 54 kDa protein 2 RN 9_140011269 1.54E-06 0.15 0.25 Zm00001d047726 Unknown Zm00001d047725 Unknown Zm00001d047724 Unknown Zm00001d047723 Unknown RN 10_117250342 5.78E-07 0.50 0.18 Zm00001d025402 Peroxidase 72-like TN 1_44557007 1.87E-06 0.30 0.05 Zm00001d028733 SNF1-related protein kinase Zm00001d028732 Unknown TN 2_10688903 1.05E-07 0.15 -0.06 Zm00001d002350 Ent-kaurene synthase B Zm00001d002349 Ent-kaur-16-ene synthase, chloroplastic TN 2_133457516 5.08E-07 0.31 0.05 Zm00001d004716 Probable NAD kinase 2, chloroplastic TN 3_173524912 1.18E-06 0.19 0.05 Zm00001d042610 Methyltransferase PMT26 Zm00001d000647 Unknown TN 3_4459301 3.31E-10 0.26 0.07 Zm00001d039439 Protein BRASSINAZOLE-RESISTANT 1 Zm00001d039440 Unknown TN 4_178021073 9.60E-08 0.49 -0.04 Zm00001d052059 RING-H2 finger protein ATL13 Zm00001d052060 Probable alpha,alpha-trehalose-phosphate synthase [UDP-forming] 9 TN 4_230820977 2.69E-06 0.38 -0.05 Zm00001d053416 Unknown Zm00001d027011 Unknown Zm00001d053417 Unknown TN 6_159738234 2.21E-08 0.34 0.05 Zm00001d038547 Glutamate-1-semialdehyde 21-aminomutase 2 chloroplastic Zm00001d038548 Putative disease resistance protein RGA3 Zm00001d038549 Cysteine-rich receptor-like protein kinase 10 Zm00001d001460 Unknown TN 8_102989335 6.92E-07 0.37 0.05 Zm00001d010174 COV1-like protein TN 10_141744253 3.63E-11 0.35 0.06 Zm00001d026242 Uncharacterized LOC100280313 Zm00001d026243 Multiple organellar RNA editing factor 2 chloroplastic Zm00001d026244 Putative putrescine-binding domain family protein Zm00001d026245 Trihelix transcription factor GT-3a [72]Open in a new tab SNP chromosome_position, MAF minor allele frequency. Positive values indicate that the major allele increase the phenotype value, and the negative values indicate that the minor allele reduce the phenotype value. Fig. 6. [73]Fig. 6 [74]Open in a new tab Enrichment (a) and Pathway (b) analysis of candidate genes for brace rooting-related traits. Further tissue-specific analysis of the expression of these candidate genes at each period of B73 (Fig. [75]4) showed that 22 genes were expressed at all periods. Among them, ZM00001D053090, ZM00001d038931, ZM00001D043693, ZM00001D016130, ZM00001d036031, ZM00001d039038, ZM00001d038932ZM00001d039037. ZM00001D039439, ZM00001D025402 and other genes were highly expressed in in roots (Fig. [76]7). Fig. 7. [77]Fig. 7 [78]Open in a new tab Expression analysis of candidate genes at different stages of B73. Effects of training population size, marker density, and significant markers on prediction accuracy The r[MG] evaluated from 5-fold cross validation was 0.51, 0.39 and 0.36 for PR, RN, and TN, respectively. The prediction accuracies evaluated with different TPS using all the SNPs are shown in Fig. [79]8a. As TPS increased from 10 to 90%, the r[MG] was 0.42, 0.47, 0.48, 0.48, 0.50, 0.51, 0.51, 0.51 and 0.52 for PR, 0.19,0.26, 0.30, 0.33, 0.34, 0.36, 0.38, 0.39 and 0.39 for RN, and 0.18, 0.23, 0.27, 0.30, 0.32, 0.33, 0.35, 0.36 and 0.36 for TN. The results indicated that the r[MG] increased with the increase of TPS. Fig. 8. [80]Fig. 8 [81]Open in a new tab Genomic prediction accuracy for brace root traits when the number of SNPs ranged from 10 to 10,000, the training population size ranged from 10 to 90% of the total population size, and the number of significant SNPs varied from 1 to 500. (a) GS with different marker density; (b) GS with different training population sizes; (c) GS with different number of significant SNPs. The prediction accuracies evaluated with different MD using 5-fold cross validation are shown in Fig. [82]8b. As MD increased from 10 to 10,000 (10, 30, 50, 100, 300, 500, 1,000, 3,000, 5,000, and 10,000), the r[MG] was 0.40, 0.40, 0.42, 0.45, 0.48, 0.49, 0.51, 0.51, 0.51 and 0.52 for PR, 0.14, 0.21, 0.23, 0.27, 0.32, 0.34, 0.36, 0.38, 0.39 and 0.40 for RN, and 0.11, 0.17, 0.20, 0.24, 0.29, 0.32, 0.33, 0.34, 0.34 and 0.35 for TN. The results indicated that the r[MG] increased with the increase of MD. The prediction accuracies evaluated from different number of significant SNPs using 5-fold cross validation are shown in Fig. [83]8c. As the number of significant SNPs increased from 1 to 500 (1, 3, 5, 10, 30, 50, 100, 300, 500), the r[MG] was 0.07, 0.33, 0.36, 0.46, 0.66, 0.73, 0.75, 0.71 and 0.70 for PR, 0.19, 0.28, 0.40, 0.57, 0.71, 0.74, 0.75, 0.73 and 0.71 for RN, and 0.26, 0.39, 0.49, 0.61, 0.73, 0.74, 0.75, 0.76 and 0.71 for TN. The highest r[MG] of PR, RN and TN were obtained, when the number of significant SNPs was 100, 100, and 300, respectively. Discussion In this study, a multi-parent DH population was used to explore the genetic basis of brace root traits in multi-environment trials. ANOVA analysis revealed significant differences in genetic variance and genotype × environment variance for all the three brace root traits, indicating that the BLUP values across environments were more reliable for genetic analysis. The broad-sense heritability for PR (H^2 = 74.28%), RN (H^2 = 79.25%) and TN (H^2 = 81.00%) was moderate, significant differences were found in both genetic variance and genotype × environment variance, which was consistent with previous results^[84]6,[85]7. In the study of Zhang et al. (2018), the heritability of brace-root tier number, brace root number, and radius of the brace roots was 0.80, 0.69, and 0.76, respectively. In the study of Sun et al. (2022), the broad-sense heritability of bracing root angle and diameter was 71% and 52%, respectively the results of heritability indicating that the brace root-related traits were mainly affected by genetic factors. TN was significant positive correlate with PR and RN, suggesting that these traits may coordinate and influence each other. Significant positive correlation was observed between brace root traits and agronomic traits, indicating that brace roots have a significant impact on plant growth and yield. These findings are consistent with previous research, which have shown that root systems can influence stem structure formation^[86]20, leaf area, leaf water content^[87]21, leaf blade hierarchy^[88]22, as well as nitrogen balance^[89]23. Brace roots affect the growth and development of aboveground plant parts by regulating the absorption and transport of water and nutrients. These findings provide a basis for future root-focused breeding strategies. GWAS based on a population with broad genetic diversity can effectively identify the genetic loci of target traits with high mapping resolution. It has been widely used for genetic dissection of complex quantitative traits. In this study, GWAS identified 22 SNPs associated with brace root traits distributed on chromosomes 1, 2, 3, 4, 5, 6, 8, 9, and 10, including eight SNPs for PR, four SNPs for RN, and ten SNPs for TN. These results indicated that brace root traits are complex traits that controlled by multiple loci. It is consistently with previous studies. Numerous QTL linked to brace root traits have been discovered in maize. Some significant SNPs detected in present study were located at the same region of QTL reported before. For PR, the significant SNP 6_67355900 located at bin6.01 was mapped in the mapping interval of qRPR6 (bin6.01), a QTL related to rind penetrometer resistance (RPR) of stalk in a high-oil maize population^[90]24. The significant SNP 3_207562844 located at bin3.07 was mapped in the same region of qEAR3-1, a QTL for effective brace root tier number detected in a RIL population^[91]5. For RN, significant SNP 9_140011269 located at bin9.05 was in the same region of qRPR9 for rind penetrometer resistance^[92]24, qTAR9-1 for total brace root tier number, qEAR9-1 for effective brace root tier number (Ku et al., 2012), and qRN9-1 for root number per tier^[93]7. SNP 10_117250342 at bin 10.04 was mapped to the QTL interval of qInD10 for internode diameter, qFreW10 for fresh weight of internode, and qDryW10 for dry weight of internode^[94]24. SNP 5_145812474 located at bin5.04 coincided with qTAR5 for total brace root tier number and qEAR5 for effective brace root tier number reported by Ku et al. (2012). For TN, SNP 8_102989335 located at bin8.03 was in the same region as qRPR8、qInD8、qFreW8, and qDryW8-1 reported by Hu et al. (2012). These common loci detected in different studies are stable QTL and deserve further analysis. Some loci showed pleiotropism effect, consistent with the correlation analysis. Increasing the frequency of favorable alleles in breeding germplasm can improve lodging resistance and increase yield. Candidate gene analysis can better understand the genetic architecture of brace root traits. In present study, 50 candidate genes were found and 33 of them were annotated with known function. According to gene annotation, some functional genes were involved in the growth and development of brace roots. The candidate gene Zm00001d039037 encodes a mitogen-activated protein kinase (MAPK). Mitogen-activated protein (MAP) kinases are key regulators of cellular signal transduction systems. The MAPK pathway is a signaling cascade that is differentially regulated by growth factors, mitogens, hormones, and stress, mediating cell growth, differentiation, and survival. MAPK activity is regulated through a three-tiered cascade reaction consisting of MAPK, MAPK kinase (MAPKK, also known as MEK), and MAPK kinase kinase (MAPKKK, or MEKK). Substrates of MAPK include other kinases and transcription factors^[95]25. The plant MAPK pathway encompasses a large number of distinct components of the MAPK cascade response, and studies of the functions of these components have demonstrated that MAPK plays a critical role in signaling and developmental processes in response to various stresses and most plant hormones^[96]26. The genome of Arabidopsis (Arabidopsis thaliana) encodes approximately 80 potential MAPKKKs (including 21 from the MEKK subfamily), 10 MAPK kinases (MKKs), and 20 MAPKs (MPKs) 27]. Arabidopsis MKK20 is involved in osmotic stress response by regulating MPK6^[97]28. Recent studies have shown that insertion mutants of MKKK20 exhibit insensitivity to abscisic acid (ABA) in root growth and stomatal responses, indicating that MKKK20 regulates ABA response by acting upstream of MKK5 and MPK6^[98]29. Analysis of the transcriptome data from ‘B73’ revealed that this gene is highly expressed in maize roots (Fig. [99]8). Therefore, it is hypothesized that the Zm00001d039037 candidate gene may influence root development by regulating cell division. The candidate gene Zm00001d025402 encodes a peroxidase, which is found in bacteria, fungi, plants, and animals, and can be classified as a member of a superfamily consisting of three major classes. Class III encompasses secreted plant peroxidases, which serve a variety of tissue-specific functions, including the removal of hydrogen peroxide from chloroplasts and cytoplasmic lysates, the oxidation of toxic compounds, cell wall biosynthesis, defense responses to injury, indole-3-acetic acid (IAA) catabolism, and ethylene biosynthesis, among others^[100]30. In rice, increased activity of cell wall peroxidase has been shown to inhibit root growth^[101]31. In Arabidopsis roots, catalase plays a role in the low potassium signal transduction pathway and is involved in the production of reactive oxygen species (ROS) under conditions of potassium deprivation, which is critical for many signal transduction pathways^[102]32. Transcriptome analysis of ‘B73’ revealed that this gene is highly expressed in maize roots (Fig. [103]8). Therefore, it is hypothesized that the Zm00001d025402 candidate gene may indirectly influence root growth by regulating catalase levels. Over the past decade, GS has been widely used in plant and animal molecular breeding. Predictive accuracy is an important indicator for evaluating the application of GS in breeding, which is influenced by several factors, including population size, population structure, MD, significant markers, genetic variance, statistical models, TPS, and the genetic relationship between training and breeding populations^[104]33–[105]37. The effect of MD, TPS, and significant markers on prediction accuracy was evaluated in this study. As the MD and TPS increased, the prediction accuracy first increases and then stabilizes. When the TPS reached 50-60% or the MD reached 3000, the r[MG] tends to stabilize. The similar trend has been reported by Ren et al. (2020) for common rust resistance, Guo et al. (2020) for kernel Zinc concentration, Cao et al. (2017) for tar spot complex resistance, and Cui et al. (2020) for hust traits in Maize^[106]38–[107]40. As the cost of genotyping decreases, numerous loci related to different traits have been identified in maize. Integrating significant SNPs of target traits into GS can improve the prediction accuracy^[108]41. Liu et al. (2020) conducted GS with significantly associated SNPs and genome-wide SNPs for Fusarium ear rot resistance in tropical maize germplasm. GS with significantly associated SNPs showed higher prediction accuracy than that with the genome-wide SNPs. In present study, as the number of significant marker increased from 1 to 500, the prediction accuracy first increased and then decreased. The highest prediction accuracy was obtained with 100-300 significant markers. The prediction accuracy was improved from 0.51 to 0.75 for PR, from 0.39 to 0.75 for RN, and from 0.36 to 0.76 for TN than using all the SNPs. Our results provide a new perspective on GS of brace root traits. Materials and methods Experimental materials and field trial design The 379 DH lines derived from 21 hybrids were used for GWAS and GS, including XY1466, XY1366, XY1266, XY1148, XY1140, XY047, XY1225, XY335, XY1224, J1652, DK653, C6361, Lidan771, Lidan638, Lidan618, DK516, Lihe869, DK159, ZD958, Lidan295 and JK968. The DH population were planted in Sangong Town experimental station (SG, 87° 12′57″ E, 43° 56′ 54″ N) and Dafeng Town experimental station (DF, 86° 34′ 49″ E, 44° 10′ 47″ N), Changji, Xinjiang, in May 2022; Ledong experimental station (LD, 108° 57′ 14″ E, 18° 27′ 14″ N), Hainan, in November 2022, and Qitai County experimental station (QT, 89° 44′ 19″ E, 44° 5′ 77″ N), Xinjiang, in May 2023. A randomized com-plete block design was used. The DH population was grown in a single row plot, with two replicates per site. The row length was 2.5 m, plant spacing was 0.25 m, and row spacing was 0.6 m. Standardized field management was followed. To eliminate boundary effects, the test area was set up with four protective rows. In each plot, 3 plants with intact and undamaged roots were measured for penetrometer resistance (PR), root number (RN), and tier number (TN) of brace root 30 days after silking. Three brace roots were sampled from each plant and measured for PR using a stalk strength penetrometer HJ03-YYD-1B (Beijing Beixin Keji Analytical Instrument Co., Ltd, Beijing, China). RN is the tier number of brace roots striking deep into the soilthe total number of brace roots. TN is the tier num-ber of brace roots initiated from above-ground nodes of the stem (Fig. [109]2a). Statistical analysis of phenotype data The statistical analysis of phenotypic data was conducted using lme4 package in R software ([110]www.R-project.org). The best linear unbiased prediction values (BLUP) of maize root related traits were calculated using a mixed linear model as follows: graphic file with name M7.gif where Inline graphic is the overall average of target traits, Inline graphic is the genetic effect of the Inline graphic th line, Inline graphic is the effect of the Inline graphic th environment, Inline graphic is the interaction effect between Inline graphic th genetic and lth environmental, Inline graphic is the kth repetitive effect in the lth environment, and Inline graphic is the residual error. The broad-sense heritability (H^2) on an entry-mean basis was evaluated as follows: graphic file with name M17.gif where Inline graphic is the genetic variance, Inline graphic is the genotype-environment interaction variance, Inline graphic is the residual variance, n is the number of environments, and r is the number of replications^[111]14. To explore the relationship of brace root traits with other agronomic traits, a correlation analysis between PR, RN, TN, plant height, ear height, ear leaf length, ear leaf width, and hundred kernel weight was conducted using R. Genotyping and genome wide association analysis The DH population was genotyped using the 48K liquid-phase hybridization probe capture technique at China golden marker biotech Co. (Beijing, China). Trimmomatic-0.36 software was used for quality control^[112]15. BWA 0.7.17 software was employed to align clean reads to the Maize Zm-B73-REFERENCE-GRAMENE-4.0^[113]16. The recommended basic filtering parameters of GATK were applied for filter, and 1,322,279 SNPs were obtained. Vcftools (Variant Call Format tools) was used to remove markers with missing rate over 20%, and minimum allele frequency below 0.05^[114]17. A total of 134,785 SNPs were left for further analysis. The number of SNPs on each chromosome ranged from 21,818 on chromosome 1 to 9,411 on chromosome 10. GWAS was performed using the Fixed and Stochastic Model Cyclic Probability Unified (FarmCPU) method in GAPIT R package. The first three principal components and kinship matrix were evaluated for false-positive control. The significant threshold of P-value was 7.42 × 10^−6 calculated by 1/n (n = 134,785). Enrichment and pathway analysis of candidate genes Candidate genes were retrieved within 10 kb upstream and downstream of the significant SNPs and annotated on maize GDB ([115]http://www.maizegdb.org/) and NCBI ([116]https://www.ncbi.nlm.nih.gov/). GO enrichment analysis was performed using AgriGO v2.0 ([117]http://systemsbiology.cau.edu.cn/agriGOv2/). KEGG pathway enrichment analysis was performed using KOBAS version 3.0 ([118]http://kobas.cbi.pku.edu.cn/kobas3)^[119]18. Genomic prediction analysis Genomic prediction analysis of PR, RN, and TN was conducted using the ridge regression best linear unbiased prediction (RRBLUP) model in the rrBLUP package^[120]19. A five-fold cross-validation method was used, with 80% of the randomly selected population set as the training population and the other 20% set as the prediction population. GS was conducted 100 times and the mean prediction accuracy was evaluated. The impact of marker density (MD), training population size (TPS), and significant markers (SM) on prediction accuracy (r[MG]) was assessed. To estimate the effect of TPS on r[MG], GS was conducted using all the markers with the TPS varied from 10% to 90% in increments of 10% of the total population size. To estimate the effect of MD on r[MG], the number of SNPs ranging from 10 to 10,000 (10, 30, 50, 100, 300, 500, 1,000, 3,000, 5,000, and 10,000) was used to estimate the prediction accuracy of the DH population. To investigate the effect of SM on r[MG], SNPs were ranked by P-value estimated by GWAS and 1, 3, 5, 10, 30, 50, 100, 300, and 500 SNPs with the smallest P-value were selected for GS. The prediction accuracies were obtained by 100 replicates of five-fold cross-validation. Author contributions Conceptualization, S.L., J.R., XX. and P.W.; Data curation, S.L., Z.F., J.J., Y.Z. and Y.M.; Funding acquisition, JR. and P.W.; Investigation, S.L. and Z.F.; Methodology, S.L., J., X.X., J.J. and P.W.; Project administration, J.R. and P.W.; Resources, J.R. and P.W.; Visualization, S.L., X.X., Z.F., Y.Z. and Y.M.; Writing – original draft, S.L.; Writing – review & editing, J.R. and P.W. Funding This research was funded by National Natural Science Foundation of China (32060484, 32360491, U2003304); Tianshan innovation team (2022D14017); Xinjiang Uygur Autonomous Region Natural Science Foundation key project (2022D01D34); Tianshan Yingcai (2022TSYCJU0003); Xinjiang Uygur Autonomous Region Major Science and Technology Special Projects (2022A02003-4); XinJiang Agriculture Research System (XJARS-02). Data availability The datasets generated during and/or analysed during the current study are available from the corresponding author on rea-sonable request. Declarations Competing interests The authors declare no competing interests. Footnotes Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Shaohang Lin and Xiaoming Xu contributed equally to this work and share the first authorship. Contributor Information Jiaojiao Ren, Email: renjiaojiao789@sina.com. Penghao Wu, Email: craie788@126.com. References