QTL mapping for quality traits using a high-density genetic map of wheat

Protein- and starch-related quality traits, which are quantitatively inherited and significantly influenced by the environment, are critical determinants of the end-use quality of wheat. We constructed a high-density genetic map containing 10,739 loci (5,399 unique loci) using a set of 184 recombinant inbred lines (RILs) derived from a cross of ‘Tainong 18 × Linmai 6’ (TL-RILs). In this study, a quantitative trait loci (QTLs) analysis was used to examine the genetic control of grain protein content, sedimentation value, farinograph parameters, falling number and the performance of the starch pasting properties using TL-RILs grown in a field for three years. A total of 106 QTLs for 13 quality traits were detected, distributed on the 21 chromosomes. Of these, 38 and 68 QTLs for protein- and starch-related traits, respectively, were detected in three environments and their average values (AV). Twenty-six relatively high-frequency QTLs (RHF-QTLs) that were detected in more than two environments. Twelve stable QTL clusters containing at least one RHF-QTL were detected and classified into three types: detected only for protein-related traits (type I), detected only for starch-related traits (type II), and detected for both protein- and starch-related traits (type III). A total of 339 markers flanked with 11 QTL clusters (all except C6), were found to be highly homologous with 282 high confidence (HC) and 57 low confidence (LC) candidate genes based on IWGSC RefSeq v 1.0. These stable QTLs and RHF-QTLs, especially those grouped into clusters, are credible and should be given priority for QTL fine-mapping and identification of candidate genes with which to explain the molecular mechanisms of quality development and inform marker-assisted breeding in the future.


Introduction
Common wheat (Triticum aestivum L., 2n = 6x = 42, AABBDD) is one of the most widely grown crops worldwide, with an annual output of more than 620 million tonnes. As a staple food crop for 30% of the human population, approximately 20% of the protein and calories consumed by humans come from wheat (FAO, http://www.fao.org/faostat/). Food products, such as bread, cake, noodles and steamed bread, made from processed wheat grain are a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 consumed daily across the globe. Improving processing and end use quality is an important objective in wheat breeding programs [1,2].
In wheat, protein and starch functionality are the critical determinants of end-use quality. Protein-and starch-related quality traits include grain protein content, dough rheological properties, sedimentation value, falling number (FN), and starch pasting characteristics. Grain protein content can affect the nutritional value and baking properties of wheat, which also influence the dough rheological properties and are critical marketing characteristics [3,4]. Rheological properties govern the performance of wheat flour dough during mechanical treatment [5,6]. Farinograph parameters, including water absorption (WA), development time (DT), and stability time (ST), are used to assess these rheological properties, providing insight into baking performance [7]. The sedimentation value is measured to evaluate the gluten strength, which is correlated with dough rheological properties [8]. The falling number is an important quality trait that determines the effect of -amylase on wheat grain starch [9][10][11]. Starch pasting properties are highly related to the quality of steamed bread, noodles, and other products [12].
Quality traits are quantitatively inherited and significantly influenced by the environment [12]. Quantitative trait locus (QTL) analysis based on high-density genetic maps provides an effective approach to dissect complicated traits into component loci and study the relative effects of the loci on a specific trait [13]. In recent decades, multiple QTL mapping studies have been performed to understand the genetic architecture underlying end-use quality, including grain protein content [2; 8; 12; 14-20], dough rheological properties [2; 12; 19; 21; 22], sedimentation value [4; 23-26], falling number [11; 27; 28] and starch pasting properties [24; 29-32]. Some stable or relatively high frequency QTLs (RHF-QTL) and important QTL clusters in the same genomic regions have been detected. For example, McCartney et al. [23] identified 99 QTLs over 18 chromosomes for 41 quality traits. One stable locus, which included twenty of those QTLs, was mapped on chromosome 4D. Another 10 QTLs, which were mainly associated with traits for mixograph and farinograph performance, baking performance, and starch functionality, were mapped near a QTL for time to maturity on chromosome 7D. Mann et al. [2] detected two relatively high frequency QTLs for grain protein content on chromosomes 3A and 7A at three of five sites, and two dough rheology QTLs were highly consistent across the 5 sites, with major effects associated with the Glu-B1 and Glu-D1 loci identified. Carter et al. [8] identified two QTL clusters located on chromosomes 3B and 4D and associated with traits for milling quality and starch functionality. Boehm Jr. et al. [18] identified four stable genomic regions on chromosomes 1BL, 2DS, 7AS and 7BS, with numerous colocalizing QTLs influencing multiple end-use quality traits, such as grain protein content, water absorption, and flour yield.
Single nucleotide polymorphisms (SNPs), the discovery of which was largely based on sequence information, are abundant, codominant and uniform throughout the genome [33]. Recently, SNP genotyping arrays, such as the 90K iSelect and 660K SNP array, have been used to construct high-density maps and detect QTLs in wheat [31; 34-39]. High-density genetic maps are useful in detecting QTLs related to important agronomical traits, identifying and localizing candidate genes, and map-based cloning of genes of interest and in marker-assisted selection (MAS) and predictive molecular breeding [33; 40; 41].
In this study, we performed QTL analysis of 13 quality traits measured with the seeds sample harvested from the three growing season and their mean values using a high-density genetic map containing 10,739 loci (5,399 unique loci) covering all 21 chromosomes [42]. The aims of the present study were to identify QTLs for quality traits using the RIL population and find relatively stable QTL and QTL clusters, which may be used in QTL cloning and wheat breeding programs.

Plant materials and field trials
The RIL population used in the study was derived by single-seed descent (SSD) methods from a cross of "Tainong18 × Linmai6" (TL-RIL, F 8 in 2013). Tainong18 (TN18) is a cultivated variety developed by our research group (Approved Number: Lunongshen 2008056), which was released in 2008 and is planted on approximately 30 thousand hectares per year in the Huang-huai winter wheat region, China. TN18 possesses several salient features, such as dwarfness and excellent lodging resistance, high grain yield and fine quality. TN18 became a core parent with the development of the cultivars Shannong 29 (Approved Number: Lunongshen 2016002, Guoshenmai 2016024) and Shannong 30 (Approved Number: Guoshenmai 20170019) and over 20 strains under regional trials. The male parent Linmai6 (LM6) is an elite wheat line developed by the Linyi Academy of Agricultural Science, the mother of which is a sister line of the famous cultivar "Jimai 22" (Approved Number: Lunongshen 2006050, Guoshenmai 2006018). A total of 184 lines of the RILs, which were randomly selected from the original 305 lines, were used to construct the genetic map and investigate the phenotypes.
The field trials were conducted at the experimental farm of the Shandong Agricultural University in Tai'an for three growing seasons with two replications: 2011-2012 (E1), 2012-2013 (E2), and 2013-2014 (E3). The seeds of the TL-RILs and their parents were sown in early October every year and harvested in early June the next year. Each plot consisted of 3 rows that were 1.2 m long and spaced 25 cm apart. Fifty seeds were planted in each row. The experimental fields had loamy soil, and the grain yield was approximately 9 000 kg�ha −1 . During the growing season, cultural management were carried out according to the local cultivation practices and the plants were not affected by insect pests, disease, or lodging, and the rainfalls in Tai'an was 190-220 mm.

Quality trait measurements
The seed samples obtained from the harvested populations were stored at room temperature for approximately three months and then milled using a Bühler experimental mill (Bühler mill, Bühler-Miag Company, Braunschweig, Germany). The whole flour extracted from the wheat grain with natural drying, was obtained with a Perten 3100 experimental mill (Perten Instruments AB, Huddinge, Sweden). Meanwhile the moisture content of whole flour was measured by direct drying method according to GB/T 5009.3-2016, which was used to adjust the GPC values based on dry matter. The flour extraction rate was approximately 70%, and the flours were stored at 4˚C until use. The whole flour was used to test the grain protein content (GPC), and the flour was used to test the sedimentation value (SV), farinograph parameters, falling number (FN) and starch pasting properties. The GPC was measured by the Kjeldahl method using an NC analyzer (KDY-9820, Tongrun Ltd., China). The SV was determined with a sedimentation volume instrument (BAU-A type) [30]. The farinograph parameters, including development time (DT), stability time (ST), and water absorption (WA), were determined by a farinograph (Brabender GmbH and Co KG). The FN was determined using a falling number instrument (FN1500 type, Perton Co., Sweden). The starch pasting properties, including peak viscosity (PV), trough viscosity (TV), final viscosity (FV), breakdown (BD), setback (SB), peak time (PTi) and pasting temperature (PTe), were measured with a Rapid ViscoAnalyser (RVA-Super 3, Newport Scientific, Australia), following the manufacturer's instructions.

Data analysis
The analyses of variance (ANOVA) were calculated using SPSS 17.0 software (SPSS Inc., Chicago, IL, USA) and Excel 2010. Pearson correlations coefficients (r) between different traits under each environment were calculated and plotted in R version 3.2.2 (https://www.r-project. org/). The frequency distribution of different traits was also plotted using R software. Heritability was calculated using a model in which the three environments were regarded as three replications, and the genotype × treatment interaction was used as the error term. The broadsense heritability (h B 2 ) was estimated according to the following formula: h B 2 = σ g 2 / (σ g 2 +σ e 2 ), where σ g 2 is the genotypic variance and σ e 2 is the total error variance.

The high-density genetic map
The TL-RILs and their parents were scanned with a total of 156 thousand markers, including the wheat PstI (TaqI) V3 and wheat GBS1.0 DArT arrays (Diversity Arrays Technology, https://www.diversityarrays.com/), the iSelect 90K SNP array [34], and SSR (http://wheat.pw. usda.gov) and EST-SSR markers [43][44][45]. The anchor information for the SSR markers was obtained from the GrainGene website (http://wheat.pw.usda.gov), and the anchor information for the DArT was provided by DArT arrays canning results (https://www.diversityarrays.com/ ) and SNP arrays was provided by Wang et al [34]. A linkage analysis was performed in two steps using JoinMap 1 4.0 [46]. The map construction was performed using the Kosambi mapping function with the following JoinMap parameter settings: Rec = 0.35, LOD = 4.0 and Jump = 5. A regression mapping algorithm was used to order the markers in the maps and to calculate the genetic distances between the markers. The linkage groups were named according to the wheat chromosome nomenclature, followed by a number. Among the 184 RILs, the markers with identical genotypes were defined as cosegregated and considered unique loci, which were named using randomly selected marker names. Finally, 10 739 markers were mapped to the genetic map, including 5548 DArTs, 5085 SNPs, and 106 SSR and EST-SSRs (S1 Table). Eliminating the colocated markers, the whole map consisted of 5399 unique loci, including 3788 DArTs, 1506 SNPs, and 105 SSRs and EST-SSRs. The genetic map spanned a total map length of 3394.47 cM and included 43 linkage groups covering all 21 chromosomes, with an average chromosome length of 161.64 cM and a marker density of 0.63 cM/marker (S1 Table).

QTL mapping
Both the phenotypic values obtained from the three environments (E1, E2, and E3) and their average values (AV) were used for QTL mapping analyses. Windows QTL Cartographer 2.5 software [47] was used to perform the QTL mapping. Composite-interval mapping (CIM) was selected to search for QTL of each trait separately. The parameter setup ''model 6 standard analysis" was used with a walk speed of 0.5 cM, ''forward and backward" regression for the selection of the markers to control for the genetic background, up to five control markers, and a blocked window size of 10 cM to exclude closely linked control markers at the testing site. The threshold for declaring the presence of a significant QTL was defined by 1,000 permutations at p�0.05 [48] and the minimum LOD score was chosen. We defined the relatively high frequency QTL (RHF-QTL), which was detected in more than two environments. A QTL cluster was defined as two or more traits with significant QTLs (at least one RHF-QTL contained) having overlapping confidence intervals.

Identification of candidate genes related to QTLs
To determine whether these flanking chromosomal region markers are associated with candidate genes for quality traits, all the marker sequences from the QTL cluster confidence interval were aligned with the wheat reference genome (IWGSC RefSeq v1.0). In addition, according to the results of the alignment of the genetic linkages and physical maps based on IWGSC RefSeq v1.0, we predicted all the candidate genes within the confidence intervals associated with the 12 QTL clusters for quality traits. BLAST analyses were carried out on Chinese Spring IWGSC RefSeq v1.0 using tools available at the URGI (https://wheat-urgi.versailles.inra.fr/ Tools/BLAST-Public). Gene annotations were obtained from https://wheat-urgi.versailles. inra.fr/Seq-Repository/Annotations.

Phenotypic variation and correlation analysis
The parents of the RIL population, TN18 and LM6, showed obvious differences for most of the investigated traits (Table 1). TN18 had higher SV, DT, ST than LM6, but lower GPC, WA, and starch pasting properties except SB than LM6. For the RIL population, there was a wide range of variation with from 1.35% to 47.48% (Table 1). Transgressive segregation was observed for almost all of the quality traits under E1, E2, E3 and AV. All 13 investigated traits in each environment (including AV) exhibited a continuous distribution (S1 Fig). The h B 2 for the investigated traits ranged from 29.23 (Pte) to 70.05% (SV) ( Table 1). The results of the ANOVA showed that the variance of the genotype and the environmental effects of the 13 investigated traits were significant at p � 0.001 (Table 2). Most pearson correlation coefficients (r) among protein-related traits (GPC, SV, WA, DT and ST) were significant, except for the r between ST and GPC, between WA and ST (Fig 1). Likewise, most pearson correlation coefficients (r) among the starch-related traits (FN, PV, TV, FV, BD, SB, PTi and PTe) were significant, except for the r between BD and TV, between Pte and FN/PV/TV/FV/BD/SB/PTi (Fig 1). In addition, the r values between the protein-and starch-related traits were relatively less significant; for example, the GPC was significantly negatively correlated with the FV and SB traits, the SV was significantly correlated with the PTi traits, the WA was significantly negatively correlated with the PV, BD, FV and SB traits, the DT was significantly correlated with the FN, TV, FV and PTi traits, and the ST was significantly correlated with the TV, FV, SB and PTi traits (Fig 1).

QTL analysis
According to the results of 1,000 permutations at p�0.05, the minimum LOD score for each trait-environment is listed in S2 Table. A total of 106 QTLs for 13 quality traits were detected, distributed on the 21 chromosomes (S3 Table; Fig 2). The QTLs explained 5.32-35.09% of the phenotypic variation. The highest LOD value for a single QTL was 30.2.

QTL mapping compared with that of a previous study
In the present study, a total of 107 QTLs for GPC, SV, FN, three farinograph parameters and seven starch pasting properties were detected (S3 Table). The stable loci, including 26 RHF-QTLs, were mainly located on chromosomes 1B, 4B and 6A for GPC; 1A, 1D, 4B, 6D for SV; 1A and 1D for DT; 1A, 1D and 6A for ST; 1B, 4B and 6A for WA; 1D for FN; 1A for FV; 5D, 6B, 7A and 7D for BD; and 1A and 6D for Pti. By comparing our QTL mapping results with those of previously mapped QTLs, we found that some QTLs were detected on the same or adjacent marker regions as those in previous results (Table 4). For example, the RHF-QTLs for SV were detected on chromosome 1D in all the environments tightly linked with the SNP marker BS00018250_51, for which there was also a detected QTL for DT reported by Jin et al.  [11], and for starch pasting properties on 1A, 1D, 7A, and 7D [30; 31]. However, due to the very large differences in the markers used to create the genetic maps, the background of the populations and the environments of the trials, the majority of the QTLs on the same chromosome were mapped to new marker regions in the present study.

The stable QTLs for quality traits and predication of candidate genes
In wheat, a large number of stable QTLs, especially those distributed into clusters, have been mapped in the same genomic regions [2; 8; 17; 18; 23; 31]; these regions are credible and should be the priority for fine-scale mapping and identification of candidate genes with which to elaborate on the molecular mechanisms of quality development. The QTL clusters that include stable QTLs for different traits in several environments are usually detected at the same marker interval or adjacent loci on the same chromosome, which may be due to the linkage of several genes or to the pleiotropic effects of a common gene. Therefore, these stable QTL regions have a very large significance for fine-scale gene mapping, cloning and molecular marker-assisted breeding. Moreover, with the completion of the IWGSC RefSeq v1.0 assembly, it is possible to dissect quantitative traits genetically and implement modern breeding strategies for future wheat improvement by resolving the inherent complexity of gene families related to important agronomic traits [49]. In the present study, twenty-seven RHF-QTLs were identified in more than two environments (Table 3). Twenty-one of the 26 RHF-QTLs were grouped into 12 QTL clusters (C1-C12), mapped on eight chromosomes, including 1A, 1B, 1D, 4B, 5D, 6A, 6D and 7A (Table 4, Fig 2), and classified into three types: I, II and III. By aligning with the wheat reference genome (IWGSC RefSeq v1.0), a total of 339 markers on our genetic linkage map, flanked with 11 QTL clusters (all except C6), were found to be highly homologous, with 282 high confidence (HC) and 57 low confidence (LC) candidate genes (S4 Table). In addition, according to the results of the alignment of the genetic linkages and physical maps (based on IWGSC RefSeq v1.0) (S5 Table), we predicted all the candidate genes within the confidence intervals associated with the 12 QTL clusters for quality traits (Table 5). Several important QTL clusters (e.g., C2, C3, and C4) were linked with known quality genes, and other clusters (e.g., C5, C6, and C10) could be linked with new candidate genes, which are discussed as follows. Cluster C2, located on chromosome 1A between two flanking markers (BS00089894_51~wmc312), represents a genetic distance of 6.07 cM. This corresponds to an approximate physical interval of 16.2 Mb in a region between 497.52 and 513.75 Mb, which contains 185 HC genes and 264 LC genes based on IWGSC RefSeq v 1.0 (S5 Table,  Table 5). Of the candidate genes, TraesCS1A01G317500, TraesCS1A01G466100LC, TraesC-S1A01G466400LC and TraesCS1A01G466500LC, which encode the high molecular weight glutenin subunit (HMW-GS), were detected in the confidence interval of cluster C2. Cluster C3, located on chromosome 1B between two flanking markers (D-996292~D-3024514), represents a genetic distance of 6.44 cM. This corresponds to an approximately physical interval of 32.0 Mb in a region between 532.69 and 564.65 Mb, which contains 271 HC and 452 LC genes (S5 Table, Table 5). Three HMW-GS genes, including TraesCS1B01G330000, TraesCS1B01G569900LC, and TraesCS1B01G570600LC, were detected in the confidence interval of cluster C3. Cluster C4, located on chromosome 1D between two flanking markers (swes1148~Ex_c68201_1182), represents a genetic distance of 10.97 cM. This corresponds to an approximately physical interval of 28.1 Mb in a region between 407.81 and 435.94 Mb, which contains 373 HC and 444 LC genes (S5 Table, Table 5). A total of three Table 5. QTL clusters for more than two traits in more than two environments. candidate genes related to the functions of the HMW-GS were predicted on chromosome 1D, including TraesCS1D01G317300, TraesCS1D01G435600LC, and TraesCS1D01G435700LC. In the present study, the SNP markers BS00022768_51 and BS00082503_51, which are tightly linked with QSv-1A-1, were further identified in the coding region of the TraesCS1D01G317300 and TraesCS1D01G435700LC genes, respectively. According to the BLASTN program from the National Center for Biotechnology Information nucleotide sequence database (http://blast. ncbi.nlm.nih.gov), the candidate genes TraesCS1A01G317500, TraesCS1B01G330000, and TraesCS1D01G317300 were consistent with the globulin loci reported by Gu et al. [50]; TraesC-S1A01G466100LC and TraesCS1A01G466500LC were consistent with the 1Ax and 1Ay null genotpyes of the Glu-A1 locus [50,51], respectively; TraesCS1B01G570600LC was consistent with the y9 subunit of the Glu-B1 locus [52]; and TraesCS1D01G435600LC and TraesCS1D01G435700LC were consistent with the x5 and y10 subunits of the Glu-D1 locus [53,54]. Numerous studies have identified associations between the HMW-GS and the quality of wheat [2; 16; 17; 23; 25; 26; 55; 56]. Furthermore, many previous studies have reported that QTLs for some quality traits (such as sedimentation value, development time, and stability time) were detected on the Glu-1 loci, for example Glu-A1 [30,31], Glu-B1 [17], and Glu-D1 [12,24,30,31]. In the present study, clusters C2 and C4 were determined to be colocations of stable QTLs for the SV and farinograph parameter (ST and DT) traits in more than three environments. Cluster C3 on chromosome 1B was determined to be a colocation of stable QTLs for the GPC traits under E3/AV environments. This result suggests the reliability of QTL mapping using our high-density genetic map and helps us to facilitate further fine-scale mapping and discovery of candidate genes with other QTL clusters.
In addition, we found that clusters C5 and C6, located on chromosome 4B, had colocalized stable QTLs for both the protein-and starch-related traits (e.g., GPC, SV, farinograph parameters and starch pasting properties) in one or more environments. For these two clusters, there was a genetic distance of 2.32/4.92 cM (for cluster C5/C6) corresponding to two approximate physical intervals of 1.  Table, Table 5). These results may provide a molecular foundation for the annotation of QTL and QTL clusters. Moreover, we found another cluster C10, located on chromosome 6A between two flanking markers (D-1074651~D-3943958), representing a genetic distance of 9.04 cM. This corresponds to two approximate physical intervals of 60.0/1.4 Mb region between 38.72 and 98.75 Mb and 594.53 and 595.93 Mb, which contain 541/31 HC genes and 841/27 LC genes based on IWGSC RefSeq v1.0 (S5 Table, Table 5). Cluster C16 was detected based on the colocation of stable QTLs for the GPC, WA, SV, ST traits in more than two environments. However, there were no known quality genes, such as α-and β-gliadins (Gli-A2, Gli-B2 and Gli-D2), detected based on the IWGSC_RefSeq_Annotations_v1.0 (https://urgi. versailles.inra.fr/). The α-and β-gliadins (Gli-A2, Gli-B2 and Gli-D2) are located on the short arms of group 6 chromosomes and it has been speculated that there may be new genes for protein-related metabolism between the physical intervals of cluster C10. Therefore, these regions may facilitate further fine-scale mapping and the discovery of candidate genes. Otherwise, SNP or DArT markers closely linked to quality traits and identified by QTL mapping can be successfully converted into assays for marker-assisted breeding in the future.

Conclusions
In this study, we used a high-density genetic map containing 10,739 loci (5,399 unique loci) of TL-RILs population to identify QTLs associated with grain protein content, sedimentation value, farinograph parameters, falling number and the performance of the starch pasting properties in wheat. A total of 106 QTLs for 13 quality traits were detected. Of these, 38 and 68 QTLs for protein-and starch-related traits were detected in three environments and the average values (AV), respectively. Twenty-six relatively high-frequency QTLs (RHF-QTLs) that were detected in more than two environments. Twelve stable QTL clusters containing at least one RHF-QTL were detected and classified into three types: type I, type II, and type III. A total of 339 markers flanked with 11 QTL clusters (all except C6), were found to be highly homologous with 282 high confidence (HC) and 57 low confidence (LC) candidate genes based on IWGSC RefSeq v1.0. Our findings might be of great usefulness for marker-assisted breeding, and could provide detailed information for the QTL fine-mapping and candidate gene discovery.