QTL mapping for micronutrients concentration and yield component traits in a hexaploid wheat mapping population

Bread wheat is a major staple cereal provides more than 20% of dietary energy and protein supply to global population. However, with increasing population growth, the problem of nutritional deficiencies is increasingly affecting the health of resource people with predominantly cereal-based diet. Therefore, the development of wheat genotypes with micronutrient-dense grains along with high-yield potential is one of the major priorities of wheat biofortification program at CIMMYT. We conducted a QTL mapping study using a recombinant inbred line (RIL) population derived from a cross between a Chinese parental line with highGZnC and a Mexican commercial bread wheat cultivar Roelfs F2007 to identify QTLs that could potentially be integrated in mineral nutrient concentrations and agronomic-related traits breeding. We evaluated 200 RIL lines for mineral nutrient concentrations and agronomic-related traits over two years. A total of 60 QTLs were detected, of which 10 QTLs for GZnC, 9 for GFeC, 5 for GPC and 36 for agronomic-related traits. Moreover, a total of 55 promising candidate genes were identified from the list of associated markers for GFeC and GZnC using the recently annotated wheat genome sequence. We identified the promising genomic regions with high mineral nutrient concentrations and acceptable yield potential, which are good resource for further use in wheat biofortification breeding programs.


Introduction
Poor quality diets dominated by food staples, are often deficient in minerals and vitamins, causing mineral malnutrition problem or 'hidden hunger', which affects more than 2 billion people globally. Especially iron (Fe) and zinc (Zn) deficiencies cauing serious health problems in pregnant women and childrens below 5 years of an developing regions. Anemia, whose largest cause is Fe deficiency, afflicts 24.8% of the global population and more than 65% of the preschool-age children in Africa and South-east Asia (HO, 2008). Besides, Fe deficiency is the 13th most significant global health risk factor in 2010, accounting for 48,225,000 disability-adjusted life years (Murray and Lopez, 2013). Meanwhile, inadequate Zn intake threatens estimated 17.3% of the global population, which mostly distributes in Asian and African developing countries (Wessells and Brown, 2012). In addition, Zn deficiency is responsible for the annual mortality of 433,000 children under the age of five (WHO, 2009).
One of the measures to alleviate Fe and Zn deficiencies is to increase the GFeC and/or GZnC within edible parts of agricultural crops, which is named asbiofortification. Among various means to tackle Fe/Zn deficiencies (Stein, 2010), breeding high-GFeC/GZnC crops through genetic biofortification has emerged as a most cost-effective and sustainable approach. Although being a poor source of micronutrients such as Fe and Zn, cereals are dominant in the diet of most resource-scarce population in developing regions (Graham et al., 2001). Wheat (Triticum aestivum) is the cereal with the third largest global production after maize and rice (FAO, 2014), contributing to 17% of the calories in the human diet in the developing world (Braun et al., 2010). Hence, biofortification of wheat by means of breeding is a promising strategy to ameliorate Fe and Zn deficiencies in the developing countries, and therefore increasing the GFeC/GZnC has been an important breeding target in the recent wheat breeding programs. Modern wheat varieties have limited variation for grain Zn and Fe, hence a large scale screening of wheat genetic resources at the germplasm bank of the International Maize and Wheat Improvement Center (CIMMYT) showed large variation for Zn and Fe amongst the wild relatives and landraces (Ortiz-Monasterio et al., 2007;Velu et al. 2011Velu et al. , 2019. The available genetic variation from different species and landraces were utilized in the biofortification breeding program to develop nutrient-enriched wheat germplasm with competitive yield potential and stress tolerance (Velu et al., 2014).
While GFeC and GZnC biofortification is an important objective in wheat breeding programs, conventional traits such as grain yield and GPC should not be compromised in exchange of mineral concentrations as the farmers income primarily depends on the grain yield. In addition, wheat is a major source of protein accounting for 20% of human protein intake in the developing regions (Braun et al., 2010). As been reviewed by Cakmak et al. (2000), there is generally no or only a slight negative correlation between the yield and GFeC/GZnC. However, recent research by Velu et al. (2016) involving advanced lines under heat and drought stresses showed a significant negative correlation between the GZnC and grain yield per unit area, indicating that breeding for high GZnC would not optimize the yield under such stress conditions. Nonetheless, it has been reported that protein concentration of the wheat grain decreases under elevated air CO 2 concentration of 550 μmol/mol (Fernando et al., 2012). By the end of the 21 st century, it is predicted that the global temperature will rise by 1.1-3.1°C and atmospheric CO 2 concentration will reach above 550 ppm under intermediate scenarios (IPCC, 2014). Moreover, it is expected that heat waves will occur with a higher frequency and longer duration (IPCC, 2014). Given these climatic variations predicted in the coming decades and how wheat grain yield and GPC would be affected by such changes, it is paramount to improve grain yield potential and stress tolerance with additional levels of GFeC and GZnC in the biofortification breeding. Therefore, identifying the quantitative trait loci (QTL) that regulate the accumulation of high mineral nutrient concentration in the wheat grain and yield related agronomic traits would allow breeders to develop biofortified cultivars more efficiently by using closely linked molecular markers to screen and select the most favourable genotypes.
In this study, a recombinant inbred line (RIL) population derived from a cross between a Mexican commercial bread wheat cultivar (Roelfs F2007) and a highGZnC cultivar with Chinese origin was developed. The main objective of the study was to: (a) explore the genotypic potential of a biofortified wheat cultivar with Chinese progenitor in its background for GFeC, GZnC and GPC and grain yield; (b) identify genetic loci associated with high GFeC, GZnC, GPC and yield component traits to provide evidence for breeding biofortified wheat varieties with high mineral nutrient concentration and TKW; and (c) determine the relationship between mineral nutrients concentration and other associated agronomic traits.

Plant materials
An F6 recombinant inbred lines (RIL) population between a Mexican commercial bread wheat cultivar Roelfs F2007 as female and a Chinese parental line whose lineage is Hong Hua Mai/. ../Blouk #1 consisting of 200 recombinant lines were grown at Norman E. Borlaug Research Station, Ciudad Obregon, Sonora, Mexico along with their parental lines. Hong Hua Mai, which is one of the Chinese landrace previously reported to have high GFeC/GZnC.

Phenotyping and analysis of phenotypic data
The RIL population was evaluated forGZnC, GFeC, GPC concentrations and other agronomic traits in Ciudad Obregon, Mexico (Norman E. Borlaug Experiment Station) for two successive crop seasons during the 2016-2017 (Y16-17) and 2017-2018 (Y17-18). Each entry was grown on a double-row plot of 1 m long and 0.8 m wide in a bedplanting system; a randomized complete block design with two replicates was utilized to conduct the evaluations. Diseases and pests were controlled chemically, whereas weeds were controlled manually and chemically. Plant materials were harvested after physiological maturity when grains were totally dry in the field. Grain samples of about 20 g for each entry were carefully cleaned to discard broken grains and foreign material, and used for micronutrient and protein analysis. Mineral concentrations in the grains (GZnC and GFeC) were measured with a 'bench-top', non-destructive, energy-dispersive X-ray fluorescence spectrometry (EDXRF) instrument (model X-Supreme 8000, Oxford Instruments plc, Abingdon, UK) that was standardized for high-throughput screening of GZnC and GFeC (unit: mg/kg) content in whole grain wheat (Paltridge et al., 2012). About three laboratory commercial checks (PBW343, Borlaug100 F2014 and Kachu #1) were used as in-house quality control checks. These checks were chosen to ensure wide range of GZnC and GFeC in wheat samples. The XRF machine at CIMMYT, Ciudad Obregon, Mexico showed within the acceptable levels of accuracy in a proficiency study conducted by Flinders University, Adelaide by comparing XRF data with the Inductively Coupled Plasma (ICP) results and confirmed that XRF data is reproducible and as accurate as ICP results. GPC was assessed using nearinfrared reflectance (NIR Systems 6500, Foss Denmark) according to the official method AACC 39-70A.
Thousand-kernel weight (TKW) were measured with SeedCount digital imaging system (model SC5000, Next Instruments Pty Ltd, New South Wales, Australia) that was standardized to measure TKW (g 10 −3 kernels). The Seed-Count system can rapidly and accurately analyze wheat grain samples, determining the grain number and grain physical characteristics based on software and flatbed scanner technology. The data on days to maturity (MD) were measured in Y16-17 (E1) and Y17-18 (E2) environments by counting the number of days from germination to physiological maturity when more than 50% of spikes were ripe and had turned yellow. The data on days to heading (HD) were measured in Y16-17 (E1) and Y17-18 (E2) environments by counting the number of days from germination to physiological heading. Plant height (PH) was measured from the ground to the tip of the spike excluding awns at the late grain-filling stage.
Analyses of variance (ANOVA) and pearson correlation coefficients for traits in each year were performed using the SPSS version 22.0 for Windows (SPSS Inc., Chicago, IL, USA).

Genotyping, linkage mapping and QTL analysis
DNA was extracted according to standard procedures (Dreisigacker et al., 2013) and genotyping was done with DArT sequence at CIMMYT, SAGA platform. The genotypic information consisted of an initial number of 45021 markers codified as 0 and 1 for homozygote and 2 for heterozygote. Parental non-polymorphic markers were discarded. Thus, a total of 3556 DArT markers were used for the QTL analyses.
The linkage map was assembled from the genotypic data using QTL IciMapping v3.2 software (http://www.isbreeding.net), applying a LOD threshold of 3.0 between adjacent markers (Li et al., 2007). QTLs were identified with the inclusive composite interval mapping (ICIM) algorithm for additive gene effects implemented in QTL IciMapping v.3.2 software. The QTL expressed in each environment were defined, as well as set of QTL which were stable across all the environments. For both procedures, the walking step was set to 1 cM and a relaxed LOD threshold of 2.5 was applied to call significance. QTL nomenclature was standard (http://wheat.pw.usda.gov/ggpages/wgc/98/Intro.htm).

Phenotypic evaluations
Mean values of traits for the parents 'Roelfs F2007' and Chinese parental line, and for the recombinant inbred lines (RIL) over environments are shown in Table 1. Large differences between the parents were observed for GZnC and plant height traits, whereas small differences between the two parents were found for GFeC, GPC, TKW, HD and MD. Mineral nutrient concentrations and agronomic traits in the RIL population showed wide range of variation. These results indicated that RIL population showed transgressive segregation, suggesting polygenic inheritance of the traits. In all environments, the mean value of GZnC and GFeC was higher than that of both parents and three commercial checks (Table 1). Among them, the highest mean GZnC (89.10 mg/kg) was recorded in E2 environment, whereas the lowest mean for GZnC (39.60 mg/kg) was recorded in E1 environment and the highest GFeC (63.10 mg/kg) was also recorded in E1. In addition, the mean value of GPC was higher than that of both parents and PBW343, whereas the mean value of plant height was lower than that of both parents. For the other agronomic trait, there was no significant difference between both parents and RILs mean value in all environments (Table 1).
The correlation coefficients between two mineral nutrient traits, GPC and four agronomic traits were calculated (Table 2). For the mineral nutrient traits, the strongest positive correlation (r = 0.517; P < 0.01) was observed between GZnC and GFeC. Moreover, slight negative correlations were found between PH and TKW (r = -0.032), and between PH and GZnC (r = -0.009), but were not significant. HD showed a significant positive correlation with MD and PH, but a negative correlation with TKW. There was a significant positive correlation between MD and PH. Significant negative correlation were also found between MD and TKW. Both GZnC and TKW showed a significant positive correlation with GPC. However, both GPC and TKW showed no significant correlation with GFeC. (Table 2). The distribution of traits across two years was continuous (Fig. 1). The highest frequency distribution region of each trait was different, such as the GZnC was about 60 mg/kg, GFeC was about 43 mg/kg, GPC was around 13-14%, HD was around 84 days, MD was 125 day, PH was 99 cm and TKW was 50 g.

Construction of genetic map
The linkage map was constructed using 3556 informative markers. A total of 1198 marker mapped to A genome chromosomes, 1955 to B genome chromosomes and 403 to D genome chromosomes. Among them, the 1B chromosome has the most markers (520), and the 3D chromosome has the least markers (41) ( Table S1). The full map covered a genetic length of 5747 cM with a mean inter-marker distance of 1.6 cM.

QTL analysis
QTLs detected by ICIM were shown in Table 3, and their map positions provided in Fig. 2. A total of 60 QTLs were detected in this study, the highest number of QTLs was found in the B genome, with 30 QTLs (50%); 17 (28%) and 13 (22%) QTLs were found in genomes A, B and D, respectively (Table 3). Of these, ranging from 4 to 14 QTLs for each trait. They were mapped to 16 chromosomal locations with 1-3 QTLs per cluster, the largest gene cluster was located between 100008884|F|0 and 1237294|F|0 on the 6D chromosome (Fig. 2).

QTL analysis of GFeC
Nine QTLs on chromosomes 1A, 2A, 3B, 3D, 4B, 5A and 6B were detected for GFeC. The phenotypic variation explained by these QTLs ranged from 2.10 to 14.56%. For three of nine QTLs, the 'Chinese parental line' allele increased GFeC. The largest portion of the total phenotypic variation (R 2 = 14.56%) was explained by QGFe.co-3B.1 with an additive effect increase from the 'Chinese parental line' allele (Table 3).

QTL analysis of GPC
Five QTLs were associated with GPC. These five QTLs associated with a PVE of between 3.16 and 10.84% (Table 3). For the QTLs located on chromosomes 2A, 2B, and 4A, the Roelfs F2007 allele increased GPC, whereas for the QGpc.co-2B.1 on 2B, the 'Chinese parental line' allele increased GPC. QGpc.co-2A was the stably expressed QTL in E2, and also when performance was averaged across the all environments. The PVE for QGpc.co-2A was 5.76% across all environments. However, no QTL was detected in E1 environment.

QTL analysis of TKW
Fourteen QTLs for TKW were identified in all, associated with a PVE of between 3.24 and 13.56% (Table 3). QTKW.co-2B was the most stably expressed QTL, followed by QTKW.co-6D.2. Both of these QTL were detectable in each of the two environments, and when performance was averaged across the two environments. The PVE for QTKW.co-2B was 13.22% across all environments, while the overall PVE for QTKW.co-6D.2 was 7.01%. For the two QTLs, the Roelfs F2007 allele increased TKW. QTKW.co-5A was identified in E2 and across all environments. The QTKW.co-1A on chromosome 1A was detected in only one environment (E1), but the largest portion of the total phenotypic variation (R 2 = 13.56%) was explained by QTKW.co-1A. In addition, other QTLs were located on chromosomes 2A, 3B, and 5B.

QTL analysis of PH
Four QTLs on chromosomes 2A and 3B were detected for PH (Table 3). The phenotypic variation explained by these QTLs ranged from 5.12to 9.94%. QPH.co-2A mapping to chromosome 2A (flanked by the SNP marker 1010905 and 1109890) was detected in E1 environments and its PVE was 5.12% (Table 3). Other three QTLs mapping to chromosome 3B was detected in E1 and E2 environments, respectively. For QPH.co-2A, QPH.co-3B.1 and QPH.co-3B.2, the 'Chinese parental line' allele played a role in reducing PH.

QTL analysis of HD
A total of six QTLs were associated with HD. These six QTLs associated with a PVE of between 5.21 and 7.99% (Table 3). QHD.co-6D was the most stably expressed QTL in all environments and explained 7.97, 6.40 and 6.63% of the phenotypic variance, respectively. For this QTL, the 'Chinese parental line' allele delayed HD. Moreover, the other three QTLs were detected on 1A, 1B and 3B. Finally, a QTL on chromosome 1A was detected in only E2 environment, but the largest portion of the total phenotypic variation (R 2 = 7.99%) was explained by QHD.co-1A. For QHDco-1A the Roelfs F2007 allele delayed heading date HD.

QTL analysis of MD
Twelve QTLs on chromosomes 2B, 3B, 4A, 6B, 6D, 7B and 7D were detected for MD (Table 3). The phenotypic variation explained by these QTLs ranged from 2.66 to 12.28%. QMD.co-6D was the most stably expressed QTL in all environments and explained 4.24, 5.94 and 7.40% of the phenotypic variance, respectively. For this QTL, the 'Chinese parental line' allele delayed MD. The QMD.co-3B.2 on 3B explained 12.28 and 5.75% of the phenotypic variation in E2 and E3, respectively. For QMD.co-6B and QMD.co-4A, the ROELFS F2007 allele showed a delayed effect on MD.

Pleiotropic and multigenic effects revealed by QTL mapping
Some of the genomic regions identified were found to be associated with more than one trait based on multi-trait QTL mapping. Of these, two genomic regions on 3B and 6D were simultaneously associated with HD and MD. This is consistent with the significant positive correlation between HD and MD (Fig. 2). In particular, the genomic region on 6D was not only significantly associated with HD and MD, but associated with TKW. The pleiotropic effects were supported by correlations between agronomic traits.
The GZnC QTL on chromosome3D was found at the same location as the QTL for GFeC. Moreover, genomic region on chromosome 2B was simultaneously associated with GZnC and GPC, which was closer to the genomic region associated with TKW (Fig. 2). Multigenic effects were also observed in this study, where each of traits GZnC, GFeC, GPC, PH, MD, HD and TKW was significantly associated with multi-trait markers (Table 3, Fig. 2).

Candidate genes that may be linked to traits
The sequence of the flanking marker of each QTL were entered in the International Wheat Genome Sequencing Consortium database (IWGSC; http://www.wheatgenome.org/) and Ensambl Plants database of the bread wheat genome sequence (http://plants.ensembl.org/index. html). The results showed that a total of 55 QTLs appeared to be located in regions were genes coded for diverse proteins, including mainly: .cytochrome P450, Leucine-rich repeat receptor-like protein kinase, Protein kinase family protein, Acid phosphatase/vanadium-dependent haloperoxidase-related protein, protein embryonic flower 1, etc (Table S2).

Discussion
The QTL mapping analysis is a useful tool to identify genomic regions associated with mineral nutrient concentrations and yield component traits. In addition, the validated QTLs and associated markers can be integrated into the elite germplasm through marker assisted breeding.

QTLs for mineral nutrient concentrations traits
GZnC, GFeC and GPC are quantitatively inherited traits, as shown by their continuous distribution across the RIL population (Fig. 1). Our results indicated that the RILs also showed some transgressive segregation in both directions suggesting that both parents carried a few different genes with alleles contributing to increased GZnC, GFeC and GPC, and this is in accordance with the previous results presented by others (Xu et al., 2012). A large number of previous QTL mapping studies have also mapped QTL for GZnC, GFeC and GPC in 21 chromosomes of wheat and wheat related species (Peleg et al., 2009;Krishnappa et al., 2017). In our study the largest PVE (14.56%) was displayed by QGFe.co-3B.1 on chromosome 3B, followed by QGZn.co-5A with PVE 14.22% on chromosome 5A. In addition, some co-location of QTLs occurred for the three traits on chromosomes 2B (517 cM) and 3D (125 cM) (Fig. 2). Among them, QTL was detected for the GZnC and GPC on chromosomes 2B and another QTL was detected for GZnC and GFeC on chromosomes 3D. On chromosome 2B, the favourable parent allele was contributed by 'Chinese parental line' for GPC and GZnC. On chromosome 3D, the favourable allele was from Roelfs F2007 for the GZnC and GFeC. These genetic regions may have pleiotropic effects. Meanwhile, significant and positive correlations among GZnC, GFeC and GPC were found in this study, which was consistent with the studies in wheat (Xu et al., 2012;Zhao et al., 2009). All the results indicated common physiological and/or genetic basis for GZnC, GFeC and GPC in wheat, suggesting that these three traits could be improved simultaneously (Zhao et al., 2009).
The famous Gpc-B1 gene on chromosome 6BS conferring stable expression of GPC, GZnC and GFeC from wild emmer has been cloned (Uauy et al., 2006a). In this study, the highly associated markers were used to predict putative candidate genes for QTLs associated with mineral nutrienttraits. For the two traits, the results showed that on chromsomes 1A, 1B, 2A, 2B, 3A, 3B, 3D, 4A, 4B, 5A, 6B and 7A there are genes encoding for the cytochrome P450, acid phosphatase, tyrosine decarboxylase, zinc-binding in reverse transcriptase, GRF zinc finger protein, auxin response factor, heavy metal-associated domain containing protein, protein kinase and so on (Table S2) Among them, the cytochrome P450 is reported to be related to Zn and Fe homeostasis, and frequently expressed under high Zn conditions in Arabidopsis (van de Mortel et al., 2006). The heavy metal-associated domain containing protein, zinc-binding in reverse transcriptase, GRF zinc finger protein and acid phosphatase may play major roles in accumulating additional Zn and Fe to wheat grain. The QGZn.co-5A on chromosome 5A displayed a large PVE 14.22% across two years, the BLAST results showed that there are genes encoding for the FAR1 protein and zinc-binding in reverse transcriptase, which are related with zinc ion binding (Hudson et al., 2003). The QGpc.co-4A displayed its location in a region where genes code for the protein kinase and, which are reported to be involved in catalyze phosphorylation processes in which some protein structures are Zn related (Scheeff and Bourne, 2005). Furthermore, for the QGFe.co-3D and QGZn.co-3D on the same location, the BLAST results showed that on chromosome 3D there are genes encoding for the retrotransposon protein, which is related with zinc-binding in reverse transcriptase. The QGZn.co-2BQ and Gpc.co-2B.2 displayed its location in a same region on chromosome 2B where genes code for the purple acid phosphatase, NBS-LRR disease resistance protein-like protein and auxin response factor.

QTLs for agronomic-related traits
Four QTLs for PH were mapped on chromosomes 2A and 3B (Fig. 2). These QTLs mapped in the similar locations of chromosomes 2A and 3B were also reported by Cui et al. (2011). Our BLAST results showed that candidate genes related to PH are involved in cysteine-rich receptor kinase (CRK), hydrolases superfamily protein and NBS-LRR resistance-like proteinetc (Table S2). Some related studies have demonstrated that CRK5 affects plant growth and development in Arabidopsis, which further proves that this gene is closely related to PH (Burdiak et al., 2015).
Of four QTLs detected for HD, one QTL on chromosomes 6D was coincident with QTLs for MD and TKW, and another QTL on chromosome 3B was located in the same position as the QTL for MD (Fig. 2). This is consistent with our results that HD was significantly correlated with MD. Candidate genes associated with the multiple QTLs on chromosomes 3B is involved in protein embryonic flower 1 (EMF1) and auxin response factor 11 (Table S2). The former researches have indicated that the plant-specific protein embryonic flower1 (EMF1) functions in maintaining the repression of the flower homeotic gene AGAMOUS (AG) during vegetative development in Arabidopsis thaliana by acting in concert with the EMF2 complex (Calonje et al., 2008). Wang et al. (2014) showed that plant-specific embryonic flower1 (EMF1) is involved in promoting Arabidopsis flowering under long daylight. Another multiple QTLs on chromosomes 6D displayed its location in a region where genes code for the protein kinase family protein, endonuclease/exonuclease/phosphatase and UDP-glucosyl transferase (Table  S2). Previous related research showed that the expression abundance of UDP-glucosyl transferase in the starch synthesis pathway is different during the endosperm filling period (Yamakawa et al., 2007). Furthermore, putative protein phosphatase with kelch-like repeat domain is related to grain length, and usually thousand kernel weight is positively correlated with grain length, suggesting that phosphatase is also related to TKW . The stable QTKW.co-2B was located to a QTL for thousand kernel weight on chromosome 2B also identified by Ramya et al. (2010) in a bread wheat population.

Conclusions
In conclusion, few QTLs for mineral nutrient concentrations and yield componenttraits were co-located, indicating that there is a great opportunity for simultaneous selection which was evident from the significant positive correlation between GPC and GZnC, GZnC and GFeC. Although there was a significant positive correlation between GPC and TKW, but no co-location QTL was found in this study. Meanwhile, the mineral nutrient concentrations of the RIL population were significantly improved compared to the parents, and these results suggested that the rich genetic diversity for mineral nutrient concentrations in landraces provides novel alleles for genetically enhancing wheat grain mineral nutrient concentrations. The promising QTLs identified for GZnC, GZnC and GPC offers further possibilities for selective biofortification breeding of wheat with enhanced GZnC and GFeC. The major and intermediate effect QTLs from diverse sources were effectively utilized in Biofortification breeding program by pyramiding QTL regions through marker assisted breeding.

Conflicts of interest
Authors declare there are no conflicts of interest.