Mapping of the QTLs governing grain micronutrients and thousand kernel weight in wheat (Triticum aestivum L.) using high density SNP markers

Biofortification is gaining importance globally to improve human nutrition through enhancing the micronutrient content, such as vitamin A, iron, and zinc, in staple food crops. The present study aims to identify the chromosomal regions governing the grain iron concentration (GFeC), grain zinc concentration (GZnC), and thousand kernel weight (TKW) using recombinant inbred lines (RILs) in wheat, developed from a cross between HD3086 and HI1500. The experiment was conducted in four different production conditions at Delhi viz., control, drought, heat, and combined heat and drought stress and at Indore under drought stress. Grain iron and zinc content increased under heat and combined stress conditions, while thousand kernel weight decreased. Medium to high heritability with a moderate correlation between grain iron and zinc was observed. Out of 4,106 polymorphic markers between the parents, 3,407 SNP markers were used for linkage map construction which spanned over a length of 14791.18 cm. QTL analysis identified a total of 32 chromosomal regions governing the traits under study, which includes 9, 11, and 12 QTLs for GFeC, GZnC, and TKW, respectively. A QTL hotspot was identified on chromosome 4B which is associated with grain iron, grain zinc, and thousand kernel weight explaining the phenotypic variance of 29.28, 10.98, and 17.53%, respectively. Similarly, common loci were identified on chromosomes 4B and 4D for grain iron, zinc, and thousand kernel weight. In silico analysis of these chromosomal regions identified putative candidate genes that code for proteins such as Inositol 1,3,4-trisphosphate 5/6-kinase, P-loop containing nucleoside triphosphate hydrolase, Pleckstrin homology (PH) domains, Serine-threonine/tyrosine-protein kinase and F-box-like domain superfamily proteins which play role in many important biochemical or physiological process. The identified markers linked to QTLs can be used in MAS once successfully validated.


Introduction
Hidden hunger, often referred as micronutrient deficiency, occurs due to inadequate intake of minerals like iron, zinc, iodine, and vitamin A. The main cause of hidden hunger is the consumption of energy-rich but nutrient-poor diet. According to reports, 2 billion people suffer from micronutrient deficiency worldwide (1). The devastating effects of malnutrition are seen in underdeveloped and developing countries with low-income levels and reduced dietary diversification (2). Essential micronutrients are required in minute quantities for a healthy life, otherwise, deficiencies will lead to morbidity and mortality (3).
Among essential micronutrients, iron is the most important element needed for oxygen transportation, synthesis and repair of DNA (4). The deficiency of iron causes anemia around 32.9% of people suffer from anemia worldwide (5). The recommended average daily intake of iron is 27 mg per day for pregnant women, 18 mg for women, 9 mg for lactating mother, and 8 mg for men (6). Zinc is a trace mineral that acts as a catalyst for many enzymes, restore impaired energy metabolism, and helps to regulate gene expression (7). Zinc deficiency leads to alopecia, poor growth, and sexual dysfunction (8). It's reported that 17.3% of the world's population is affected by zinc deficiency (9). The average recommended intake of zinc is 8 mg for women, 11 mg for men and pregnant women, 12 mg for lactating mother (6).
Wheat is a staple food crop and the second-largest grain in both acreage and production worldwide. It's the second largest consumed food grain and accounts for 20% of protein and calories (10). The wheat cereal-based diets fail to provide the necessary amount of minerals like iron and zinc (11). Hidden hunger is more pronounced in the areas with cereals as sole source of food supply, low dietary diversification and unavailability of biofortified varieties (12). Biofortification is the most ideal as it is cost-effective and sustainable way to improve grain micronutrient concentration (13). Therefore, attention should be given to develop biofortified varieties using breeding and molecular techniques (14).
Enrichment of grain micronutrients through traditional breeding approach is difficult because of polygenic inheritance and interaction with environment. With the advent of novel molecular techniques breeders have put a step forward in utilizing markers for assessment of genetic diversity, germplasm characterization, identification of QTLs and their utilization in practical plant breeding. Identification of loci governing GFeC and GZnC with high PVE and the markers linked to them will be very useful in molecular breeding. Several QTLs have been identified previously from the populations derived from crosses between tetraploid and hexaploid wheat varieties (15)(16)(17)(18). Few biofortified wheat varieties released so far in India include, WB2 and HPBW 01 (high Zn and Fe), PBW1Zn (high Zn content), Pusa Tejas and Pusa Ujala (high protein content along with Zn) (19). However, only a few QTLs are identified in wheat with high phenotypic variance and are hardly being used in Indian molecular breeding for grain nutrient improvement. Hence, the present study aims to identify the novel and stable QTLs for grain zinc concentration, grain iron concentration and thousand kernel weight using mapping population derived from a cross between HD3086 and HI1500. The study was conducted across the environments viz, Delhi and Indore, in control and stress conditions (drought, heat, and combined drought and heat).

Plant material and environment
The mapping population consists of 166 recombinant inbred lines (RILs) derived from a cross between HD3086 and HI1500. HD3086 is an Indian high-yielding hexaploid wheat variety suitable for timely sown, irrigated condition developed at IARI, New Delhi (20,21). HI1500 is a popular variety recommended for cultivation under restricted irrigation in central zone harboring many important traits related to drought and heat tolerance (21,22). The variety HI1500 has higher GFeC, GZnC, and TKW (44.50 mg/kg GFeC, 53.50 mg/kg GZnC, and 35.88 g TKW) in all the conditions as compared to HD3086 (38.60 mg/kg GFeC, 38.60 mg/kg GZnC, and 35.58 g TKW). The RILs population along with parents were evaluated under four conditions namely timely sown irrigation (TSIR) taken as control, timely sown restricted irrigation (TSRI), late sown irrigation (LSIR), and late sown restricted irrigation (LSRI) conditions at Delhi, and under restricted irrigation condition at Indore. Two irrigations were given in TSRI, one at germination and other at 21 days after sowing; six irrigations were provided during cropping period in irrigated condition (TSIR). Late sown trials (LSIR and LSRI) were planted in second fortnight of December to expose plants to heat stress, and under LSRI condition irrigation was withheld to expose plants to both terminal heat and drought stress. The genotypes were evaluated in an alpha-lattice design with two replications. Each genotype was sown in 3 rows of 1 m each with 22.5 cm distance between rows and 10 cm distance between plants. Uniform agronomic practises were practiced for proper establishment of crop stand. The details of the sowing conditions and locations are presented in Table 1.

Phenotyping for grain micronutrients and thousand kernel weight
From each plot, 20 random spikes were harvested and spikes from each plot were threshed separately. While cleaning, care was taken to prevent metal and dust contamination. The grain iron concentration (GFeC) and grain zinc concentration (GZnC) were measured using Energy Dispersive X-ray Fluorescence (ED-XRF) machine (model X-Supreme 8000 M/s Oxford Inc., USA). The thousand kernel weight (TKW) was recorded by counting 1,000 grains manually and weighted with an electronic balance.

Phenotypic data analysis
Analysis of variance was done using PBTools v1.4 software (23). Heritability and correlations among traits were calculated using the MetaRv6.0 (Multi Environment Trial Analysis with R) software (24). From phenotypic data best linear unbiased predictors (BLUPs) were calculated for the individual conditions and combined across all production conditions and environments for further QTL mapping.

DNA extraction and genotyping
DNA was extracted from 21 days old seedling using the CTAB method (25). Genomic DNA quality was determined using 0.8% agarose gel electrophoresis with λ DNA as the standard and quantified using nanodrop. The 35K SNP Axiom breeders' array was used for genotyping of parents and the RILs population.

Linkage map construction
A total of 4,106 SNP markers were polymorphic between the two parents. Among these, the redundant markers were removed by binning and markers deviating from mendelian segregation were deleted. Finally, a set of 3,407 non-redundant SNP markers spanning all over the chromosomes were used for the linkage map construction using IciMapping v4.2 software (26). Kosambi mapping function was used to calculate map distances between markers. Marker grouping was done using a rec value of 0.37 and the linear order of markers was determined. Ordering is done using function K-optimality 3-optMAP with NN initials of 10. Rippling was carried out using recombination with a window size of 5 cm. The final generated linkage map is used for QTL mapping. QTL analysis and identification of candidate genes QTL mapping was performed using IciMapping 4.2 software (26) with inclusive composite interval mapping (ICIM-ADD) model (27). BLUP values were calculated in the individual environment and pooled over environments along with genotypic data used for QTL analysis. "Mean replacement" was used to address missing phenotypic data. The walking speed of 1.0 cm, with P = 0.001 was used in stepwise regression. The LOD score of 3.0 along with 1,000 permutations was chosen for the declaration of the QTL. Identified QTLs were named following standard nomenclature available in the catalog of wheat gene symbols (28). The candidate genes (CGs) were identified based on the positions of flanking markers of the QTL. BLAST search was done to identify putative candidate genes in the physical location of markers against IWGSC wheat (Triticum aestivum L.) reference genome embedded in the Ensembl Plants database. 1
The mean TKW of the RILs was highest in the timely sown condition compared to drought, heat and combined stress, whereas the GFeC and GZnC were low under control condition as compared to combined stress. The mean TKW in control was 38.28 g, which decreased to 37.36 g under drought, 29.99 g under heat, and 25.49 g under combined stress condition, while TKW in restricted irrigated condition at Indore increased to 41.88 g ( Table 3)     Depiction of the distribution of GFeC, GZnC, and GZnC under control, drought, heat, and combined stress conditions through violin plots.
The mean decrease in TKW and increase in GFeC and GZnC was found significant (p-value of 0.05) using one tailed and two tailed z-test ( Table 3). Analysis of variance showed significant variation among the RILs for the traits GFeC, GZnC, and TKW in all the conditions viz., control (TSIR), drought (TSRI), heat (LSIR), and combined stress (LSRI) in Delhi and drought stress (TSRI) in Indore.
The coefficient of variation was high under drought and combined stress conditions compared to the control but heritability was found to be lower under drought and combined stress conditions for all the traits. Broad sense heritability of the traits studied was medium to high except for the GZnC which was low under combined stress condition. The range, heritability, coefficient of variation and mean decrease/increase are given in Table 3. Significant association among GFeC and GZnC under control (0.51 * * * ), drought (0.57 * * * ) and heat (0.50 * * * ) at Delhi location and GFeC and GZnC under drought stress (0.41 * * * ) at Indore (at p-value of 0.001) was observed. Traits viz., GFeC and GZnC did not show significant association with TKW except positive significant (0.20 * * ) association between TKW and GZnC (at p-value of 0.01) under the drought ( Table 4).

Linkage map
A total of 4,106 polymorphic markers were identified between HD3086 and HI1500, out of which 3,407 non-redundant markers uniformly distributed over all the 21 chromosomes were used for the linkage map construction. The linkage map spanned a genetic length of 14791.18 cm, ranging from 2.03 cm/marker in 1B to 9.41 cm/marker in 4D chromosome with an average marker density of 4.34 cm/marker. The highest number of markers mapped on genome B with a total of 1386 followed by A and D with 1,046 and 975 markers, respectively. Chromosome 2B had the highest number of markers 239 while chromosome 4B had only 61 markers ( Table 5).

QTL mapping
A total of 32 QTLs were mapped on 17 different chromosomes for the traits GFeC, GZnC, and TKW across the environment. Out of 32, nine QTLs were mapped for GFeC, eleven for GZnC, and twelve for TKW. Chromosome 4B carried the highest number of QTLs i.e., 5, followed by chromosome 2D which carried 4 QTLs, while chromosomes 1D, 2A, 3B, 4D, 5D, 6A, 7A, and 7B carried 2 QTLs each, and the remaining chromosomes 3A, 4A, 5A, 5B, 6B, 6D, and 7D carried only one QTL. A list of QTLs identified along with flanking markers, LOD score, PVE (%), and additive effects are given in Table 6. In silico analysis of QTL regions identified a few important candidate genes which are having various roles in different pathways of growth and development were given in Table 7.

Discussion
Micronutrient elements such as grain iron, grain zinc, are governed by many genes whose expression is influenced by the external environment. In the present study, to identify QTLs expressed in different stress conditions the experiment was carried out in drought, heat, and combined stress conditions taking timely sown environment as control. The treatment-wise highest mean for GFeC 41.51 mg/kg, GZnC 59.33 mg/kg, and TKW 41.88 g was observed in drought (Indore), heat, and drought (Indore) conditions, respectively, as observed in previous studies (29,30). The GFeC and GZnC were increased in late sown conditions but the TKW decreased because the high temperature during grain filling caused forced maturity, reduced photosynthetic assimilation and starch content leading to the development of undernourished shriveled seeds (31). Although the grain micronutrients concentration was highest in late sown condition, the total micronutrient yield per area was highest in timely sown condition as the grain yield was highest in timely sown condition similar kind of results were also observed in a study by Velu et al. (30). GCV for GFeC ranged from 2.33 to 10.00, GCV for GZnC ranged from 12.26 to 29.92 and GCV for TKW ranged from 7.18 to 29.18 in different treatment conditions. The lower the value of GCV, the greater the influence of the environment on the expression of the particular trait whereas higher value of GCV, indicates the variation in the population mainly attributes to the genetic makeup of the individual. The CV for GFeC was ranged from 5.66 to 6.62, GZnC ranged from 1.63 to 48.90, TKW ranged from 6.97 to 15.44, similar findings were observed in previous studies by 29, 32, 33. A wide variation was observed for GFeC, GZnC, and TKW (Table 3) (17,34,35). Heritability is important selection parameter that aids plant breeders in determining the characters for which selection is to be performed (36). The heritability of the traits studied was found to be medium to high indicating a predominance of genotypic variance which can be exploited by selection in crop improvement, (17,37). GFeC is significantly correlated with GZnC in all the conditions, indicating they are positively associated with each other and the selection for improvement of one character will simultaneously bring improvement in another. The results were similar to previous studies by 38, 39. In the present study a total of 32 QTLs were identified out of which 9, 11, and 12 for GFeC, GZnC, and TKW, respectively (Figure 2). Most of the QTLs for GZnC and TKW were derived from parent HI1500 while QTLs for GFeC were derived from parent HD3086. Interestingly, one pleiotropic QTL QGFe.iari-4B/QGZnC.iari-4B.2/QTKW.iari-4B.1 has been identified at position 335 cm which explained the phenotypic variance of 29.28% for GFeC, 10.98% for GZnC, and 17.53% for TKW. Two Stable QTLs QGFeC.iari-1D.1 and QTL QGFeC.iari-4D.1 were identified for GFeC, and one for TKW QTKW.iari-4D.1 and. were also reported in earlier studies (40)(41)(42). QTLs for grain iron content are mapped most on the D genome while B genome carried the highest number of QTLs for grain zinc content and A genome carried highest number of QTLs for TKW.
The stable QTLs for GFeC was identified on chromosome 1D and 4D as also reported by previous workers on 4D (32,37), 1D Genetic linkage map and QTL positions identified on A, B, and D genomes of RILs derived from the cross HD3086/HI1500. Red color indicates QTLs for GFeC; blue color indicates QTLs for GZnC; green color indicates QTLs for TKW; purple color indicates QTL governing GFe, GZn and TKW; yellow color indicates QTL governing both GFe and TKW. (18). The putative candidate gene TraesCS4D02G03950 in the QTL QGFeC.iari-4D.1 region coded Pleckstrin homology (PH) domains; these are protein modules made up of 100-120 amino acids and has the potential to bind phosphoinositides (43). Phosphoinositides are known to regulate intracellular membrane trafficking by providing intrinsic membrane signals (44) and also help in regulation of ion channel function (45). Another stable QTL QGFeC.iari-1D.1 harbors a putative candidate gene TraesCS1D02G240700 which codes for Regulator of chromosome condensation 1/beta-lactamase-inhibitor protein II. Three major QGFeC.iari-4D.1, QGFeC.iari-4B, and QGFeC.iari-1D.1 for GFeC were also reported by previous workers on 4D (32,37), 1D (18), and 4B (37) in different mapping populations. In silico analysis of chromosomal regions harboring QTL QGFeC.iari-4B found a putative candidate gene TraesCS1D02G240700 near right marker AX-95215762 which codes for P-loop containing nucleoside triphosphate hydrolase, and belong to the special class of metallochaperones (46). Metallochaperones are specific class of molecular chaperones that mediate the intracellular transport of metal ions to metalloproteins and metalloenzymes through proteinprotein interactions (47). These proteins are Fe 2+ chaperones and important multifunctional adaptors that function in the nuclear and cytosolic compartmentalization, storage and export of iron in the form of ferritin and ferroportin (48). Another putative candidate gene TraesCS4B02G056800 codes for Inositol 1,3,4-trisphosphate 5/6-kinase, it acts as a chelator of metal ions such as iron and zinc.

Conclusion
The number of QTLs identified for GFeC, GZnC and TKW are 9, 11, and 12, respectively. One pleiotropic QTL was identified on chromosome 4B which is associated with grain iron, grain zinc concentration and thousand kernel weight. The common locus was identified for Zn content and TKW on chromosome 4B, while common locus for Fe and TKW was found on 4D. The same genes might be functioning in the accumulation of Fe, Zn, and TKW which need to be studied further in detail at molecular and biochemical levels. Genomic regions on chromosomes 1D and 4D were associated with Fe content, 2A and 6A associated with TKW, while on chromosome 4B were associated with Zn content, therefore fine mapping of the such regions may be rewarding. Further in silico analysis of QTL regions identified putative candidate genes which are directly or indirectly related to these traits. So, the identified QTLs can be used in practical plant breeding to develop biofortified varieties after the successful validation of these identified markers.

Data availability statement
The original contributions presented in this study are included in the article/Supplementary material, further inquiries can be directed to the corresponding authors.