QTL mapping of agronomic traits in wheat using the UK Avalon ×  Cadenza reference mapping population grown in Kazakhstan

Avalon × Cadenza mapping population

The original mapping population (MP) Avalon × Cadenza (A × C) consisted of 201 samples including 100 winter type lines and 101 spring type lines The MP was produced from a cross between widely grown British wheat cultivars Avalon (winter wheat) and Cadenza (spring wheat). The A × C DH population was developed as part of the Wheat Genetic Improvements Network (WCIN) (http://www.wgin.org.uk/) (Allen et al., 2011). The parental cultivars differ in their photoperiod sensitivity alleles by: Ppd-A1, Ppd-D1, Ppd-B1, and vernalization genes Vrn-A1, Vrn-B1, Vrn-D1 (Avalon) and Vrn-A1a (Cadenza). They also differ in terms of reduced height genes as Avalon carries the allele Rht-D1b, while Cadenza carries the wild type allele Rht-D1a (Farré et al., 2016). In this study, only spring-type DH lines and Cadenza were subjected for the analysis along with local standards.

Evaluation of the MP for variation in agronomic traits

The studied traits were formally divided into two groups: plant adaptation-related traits and yield components. The plant adaptation traits included the heading time (HT, days), seed maturation time (SMT, days), and plant height (PH, cm). The yield components, including the spike length (SL, cm), productive tillering (PT, pcs), number of kernels per spike (NKS, pcs), thousand kernel weight (TKW, g), and number of kernels per meter (NKM, pcs), were calculated as PT × NKS, yield per square meter (YM2, g). These A × C spring lines were evaluated in three regions of Kazakhstan, at the Karabalyk Agricultural Experimental Station (North Kazakhstan), the Karaganda Institute of Agriculture (Central Kazakhstan), and the Kazakh Rice Research Institute (South Kazakhstan) over three years, 2013–2015 (Fig. S1). In Northern and Central Kazakhstan, DH lines were grown in non-irrigated plots, while in Southern Kazakhstan, plants were grown in an irrigated field. DH lines and Cadenza were planted in three replications at each location in randomized 1 m² plots. In addition, local standards “Karabalykskaya 90”, “Karagandinskaya 22”, and “Kazakhstanskaya 4”, were planted in Northern, Central, and Southern Kazakhstan, respectively. The distance between rows was 15 cm, and the distance between plants in a row was 5 cm, respectively (Dospekhov, 1985). The climate conditions recorded during the trials were shown in Table 1, and more extended climate information for the last eight years was provided in the Raw data file.

Linkage mapping and QTL analysis

The genetic map for A × C was developed by Griffiths et al. (2009) and Griffiths et al. (2012), and details of this map are available from the WGIN website (WGIN, 2008). The genetic map was previously reported to consist of 3647 polymorphic DNA markers, including 22 amplified fragment length polymorphisms (AFLPs), 16 COS (conserved orthologous sequences), 88 DArT (diversity array technology), 3325 SNPs (single-nucleotide polymorphisms), 153 SSRs (simple sequence repeats), 3 STSs (sequence-tagged sites), and 12 genes (Table S1). The total map length was 3246.9 centiMorgan (cM), with an average chromosome length of 154.6 cM; range: 16.80 cM (chromosome 6D) to 264.8 cM (chromosome 5B) (Griffiths et al., 2009; Griffiths et al., 2012). QTL identification was conducted using the composite interval mapping (CIM) methods of Windows QTL Cartographer v2.5 software (Wang, Basten & Zeng, 2012). A logarithm of the odds ratio (LOD) threshold of 3.0 was used to determine a significant QTL. MapChart v2.32 software was used to draw the genetic map (Voorrips, 2002). The correlation analysis was calculated using the Rstudio software (R Studio Team, 2020). The GGE (genotype plus genotype-by-environment interaction) effects were analyzed using GenStat software (International, 2019).

Evaluation of agronomic traits of the A ×C population in three regions of Kazakhstan

The duration of HT differed sharply across the three regions based on the analysis of means over three years (Table 2). The earliest HT was registered in the Karaganda Institute of Agriculture (KA) region (42.3 ± 1.11 days), followed by Karabalyk Agricultural Experimental Station (KB) (49.4 ± 1.32 days), and then Kazakh Rice Research Institute (KO) (56.6 ± 3.09 days). The duration of the mean SMT showed a different trend for the three regions, with the earliest seed maturation observed in KO (22.5 ± 2.18 days), followed by KA (47.2 ± 1.16 days), and then KB (49.0 ± 1.35 days). The mean PH ranged from 47.1 ± 5.6 cm in KA to 58.1 ± 5.9 cm in KO (Table 2).

As the KB station represents Northern Kazakhstan, where wheat is grown on over 80% of the total sowing area in Kazakhstan, it was essential to compare the MP to the local standard (check cultivar) “Karabalykskaya 90”. The comparison showed that the mean performance for the HT, SMT, and PH in the A × C lines was less optimal than for the local standard. Notably, the average HT was 9.1 days, and the average SMT was 2.7 days longer in comparison with Karabalykskaya 90. This pattern was also observed at KA (Central Kazakhstan), but was reversed at the KO station (Southern Kazakhstan) as in the latter case; the SMT was shorter than in the check cultivar for south Kazakhstan, “Kazakhstanskaya 4” (Table 2). Under the irrigated conditions of Southern Kazakhstan (KO), the mean yield of the A × C lines was comparable with the check cultivar, and Cadenza showed even better productivity in comparison to Kazakhstanskaya 4 (Table 2).

An analysis of the means for YM2 revealed that nine DH lines exceeded the YM2 of the local standard cultivar, Karabalykskaya 90 (556.7 g/m²) in Northern Kazakhstan. Similar calculations performed for the trials in the Central and Southern regions suggested that 22 and 26 DH lines, respectively, had heavier yields than the corresponding local standard cultivars. Two particular lines, A × C52 and A × C55, demonstrated higher productivity than the check cultivars in all three regions. The averaged YM2 over three years in non-irrigated sites of KB and KO were significantly correlated (P < 0.01), and the averages in both locations were not correlated with the irrigated sites in KO (P < 0.81).

Pearson’s correlation index showed that in Northern Kazakhstan, the yield was not correlated with HT, SMT, and PH (Fig. 1A). However, TKW, which is one of the important agronomic traits, was negatively correlated with PH (P < 0.05), suggesting that plant height is favorable for wheat productivity in this region. A negative correlation of TKW with both HT and SMT and positively correlated with PH was recorded at the KA station, and a negative correlation was revealed between YM2 and SMT at the KO station (Figs. 1B, 1C).

A GGE biplot into YM2 divided the three regions for the four mega-environments. PC1 (25.81%) effectively separated KO2014 and KO2015 from the KB and KA sites, and PC2 (14.78%) separated KA2014 from the remaining environments (Fig. 2).

Two particular lines, A × C52 and A × C55, demonstrated adaptability to both non-irrigated and irrigated sites and showed higher productivity than the check cultivars in all three regions located in a biplot between two mega-environments (Fig. 2).

Identification of QTLs for agronomic traits in three regions of Kazakhstan

The QTL analysis in three studied regions led to the identification of 232 QTLs for nine agronomic traits. The number of QTLs per trait varied from 17 for PH to 40 for NKS (Table 3). Only 83 out of the 232 QTLs found were statistically significantly associated in two and more environments, suggesting that only 1/3 of associations were stable in the three regions. Among the nine traits, the number of identified QTLs varied from two for SMT to 12 each for both PH and SL. The largest LOD score 25.3 for the traits was recorded for PT on chromosome 2D for the Central and Southern regions. The numbers of QTLs among the three genomes A, B, and D were 31, 21, and 27, respectively, suggesting that the A and D genomes were the main locations of the stable associations. The number of stable QTLs identified for the group of adaptation-related traits was eighteen, and the number of QTLs for the group of traits for yield components was seventy-seven (Tables 3 and 4).

The number of identified QTLs found in the data from the three different regions varied significantly, and most MTAs were found in the Northern Kazakhstan data (50 QTLs), followed by the Southern (45 QTLs) and then Central sites (37 QTLs), of Kazakhstan (Table 3). Despite that more QTLs were identified in the Northern station data, the number of associations for the PT and NKS was nearly twice as high as that found in the Central and Southern regions (Table 3).

QTL mapping for traits related to plant adaptation in the Avalon × Cadenza DH population

A total of 18 QTLs was identified for plant adaptation-related traits, 12 of them detected in the PH trait (Table 4). The majority of QTLs were detected at the irrigated KO site (eight QTLs), while at the non-irrigated KB and KA sites, five QTLs were recorded at each location. Field trials in the three sites led to the identification of only four common QTLs (1D, 2D, 5B, and 6A), and one of them, QHt-AxC.ippb-1D, was identified over the years in KB sites. Another QTL for SMT, QHt-AxC.ippb-2D, was identified at the KO region, and Avalon was the donor of the increasing alleles. For SMT, we detected only two QTLs (1A), which were identified in the KB region, and in both cases, Cadenza was the donor of the increasing alleles (Table 4).

Unlike the HT and SMT analyses, where only a few QTLs were identified, twelve QTLs were genetically mapped for PH (Table 4). The R ² for the PH ranged from 10% (for QPh-AxC.ipbb-3A) to 48% (QPh-AxC.ipbb-4D) (Table 4), where the latter QTL was mapped in the vicinity of the Rht1 gene. QPh-AxC.ipbb-4D had the highest LOD score (21.5) compared to the other PH-associated QTLs; however, this QTL was significant only at the KO irrigated site (Table 4, Fig. 3). At the Northern KB site, five QTLs on chromosomes 2A, 2B, and 5A were identified. For those five QTLs, four alleles for increasing height were from Avalon, and only the QTL on 5A (135.2–183.7 cM) had the increasing height allele from Cadenza (Table 4).

QTL mapping for yield components in the Avalon × Cadenza DH population

A total of 65 stable QTLs were identified for six traits directly related to grain yield. The number of stable QTLs per trait is ranged from 9 in TKW to 12 in SL and NKM. Three QTLs for NKS (Avalon), TKW (Cadenza/Avalon), and YM2 (Cadenza) were mapped on chromosome 4D in the vicinity of the RhtMrkD1 gene associated with reduced plant height. The largest number of QTLs were identified at the KB site (39 QTLs) followed by KO (33 QTLs) and KA (31 QTLs) (Table 3). Twelve QTLs were identified for SL, and their R ² ranged from 11% (QSl-AxC.ipbb-6A) to 53% (QSl-AxC.ipbb-2D.1). The locus QSl-AxC.ipbb-2D.1 was detected at both the KB and KO sites. However, the largest QTL effect for SL was QSl-AxC.ipbb-5A.2 (0.85) , with Avalon being the donor of the favorable allele (Table 4, Table S2, and Fig. S2).

Eleven QTLs were identified for PT, where the R² values were ranged from 17% (for QPt-AxC.ipbb-1A) to 66% (for QPt-AxC.ipbb-2D). Three QTLs (QPt-AxC.ipbb-3D.1, QPt-AxC.ipbb-3D.2, and QPt-AxC.ipbb-3D.3) were genetically mapped to different regions of chromosome 3D by using trial data from the KA and KB sites (Table 4). Of the eleven QTLs for NKS, the largest QTL effects were due to QNks-AxC.ipbb-6A (3.83) and QNks-AxC.ipbb-7B (2.35), and in both cases, Avalon was the donor of these alleles associated with an increasing effect. The largest number of QTLs for NKS was identified in the data from KB (nine QTLs), where three QTLs were located in each of the A and D genomes, and four QTLs were in the B genome (Table 4). The 101 studied DH lines were separated into groups with high (from 6 to 9), middle (from 4 to 5), and low (from 1 to 3) numbers of positive QTLs (favorable alleles of significant SNP), and groups was represented by 19, 41, and 41 accessions, respectively. The unpaired t-test for DH lines harvested in Northern Kazakhstan suggested that the YM2 of the group with the high number of positive QTLs for NKS was significantly superior in comparison to the middle (P < 0.05) and low (P < 0.01) groups (Fig. 4). A similar outcome was recorded in Central Kazakhstan, where only three QTLs for NKS were identified, and samples with three positive QTLs (n = 8 samples) were having significantly higher YM2 (P < 0.05) in comparison to the group of DH lines with none or one positive QTL (n = 63 samples). A different result was recorded in Southern Kazakhstan, as the groups with more positive QTLs for NKS showed no statistical advantages in averaged YM2 over the groups with less positive QTLs (Fig. 4, Table S3).

Of the twelve QTLs for NKM, the largest QTL effects were observed for QNkm-AxC.ipbb-1D (−8.03) and QNkm-AxC.ipbb-3A (−6.03), and in both cases, Cadenza was the donor of these alleles associated with an increasing effect. The R² values for nine QTLs for TKW ranged from 9% (for QTkw-AxC.ipbb-7D) to 26% (for QTkw-AxC.ipbb-4D). The QTL with the largest LOD score (9.0) for TKW was QTkw-AxC.ipbb-4D. However, the most important QTL with the largest effect (1.84 g, Avalon) was QTkw-AxC.ipbb-5B, where the allele associated with an increasing effect came from Avalon (Table 5). Finally, eleven QTLs were identified for YM2, including QYM2-AxC.ipbb-5A.1, which was the locus with the largest LOD score (11.9) and QTL effect (615.5) at the KB/KA sites (Table 4, Table S2 , Fig. S2).

Yield assessment of the A × C DH population in three contrasting regions of Kazakhstan during the period 2013–2015

The field performance of the studied population significantly depended on geographic locations and key environmental parameters, including mean temperature, average rainfall, day length, soil quality, and etc. (Table 1). Therefore, these factors, particularly temperature and amount of precipitation in key stages of plant growth, may lead to different plant performances of the same collection of samples in different wheat-growing regions (Two-factorial ANOVA in Raw meteorological data file). In the present study, the correlation analysis showed negative influence of late heading time on major yield components-YM2 and TKW in Central and South regions, but not in the North (Fig. 1). Traditionally, the requirements for the early development of wheat in the Northern region were negated by the fact that local breeders were mostly focusing on grain quality parameters and, therefore, they targeted lines with an early flowering time (Kamran, Iqbal & Spaner, 2014; Tshikunde et al., 2019) and late SMT. However, the analysis of meteorological data revealed that heavy rains in early September were occurring more often in Northern Kazakhstan than in previous decades, which might result in a change in the breeding goals toward an early SMT as well. Hence, the negative correlation of yield-related traits with HT and SMT observed in this study, although not significant, is additional evidence of the necessity to adjust the local breeding priorities in northern parts of Kazakhstan.

The analysis of the averaged YM2 revealed 9, 22, and 26 DH lines that exceeded the YM2 of the local standard cultivars in the Northern, Central, and Southern regions, respectively. In addition, 86, 24, and 3 DH lines exceeded the YM2 of the Cadenza (parent) in the KA, KB and KO, respectively. Two particular lines, A × C52 and A × C55, demonstrated adaptability to both non-irrigated and irrigated sites and showed higher productivity than the local standards in all three of the studied regions (Fig. 2). The application of a GGE biplot analysis suggested some more insights into the assignment of particular DH lines for their possible usage in breeding projects at the three different regions using YM2. This result is particularly important as the correlation test suggested that average YM2 in the non-irrigated sites KB and KA were highly correlated (P < 0.01), and the yield in both locations was not correlated with the irrigated site in KO (P < 0.811).

Comparative analysis with associations revealed in previously published reports

The QTL analyses in the three regions led to the identification of 83 stable QTLs that were significant for nine agronomic traits in two and more environments (Table 4). Notably, the least number of associations was identified for HT and SMT, which is an indication of a narrow range of heading times in the population tested under these new environments. On the other hand, the determination of only a few QTLs associated with HT and SMT suggests that the majority of those identified for yield-related traits were not associated with the pleiotropic effects of major genes.

A comparison of the mapped QTLs analyzed in this study with those from other previous studies indicated that 40 QTLs matched known associations (Table 5). Twelve associations matched the results from studies of the SQ1 × CS DH mapping population (Quarrie et al., 2005; Abugalieva, 2007), where five associations with PT, and two with SL, NKS, and TKW, were identified through studies in Southeastern Kazakhstan. Another nine associations were identical to the genetic positions of QTLs identified with the analyses of six traits using GWAS based on the assessment of common wheat in three different regions of Kazakhstan (Turuspekov et al., 2017). Notably, five of those nine associations were also genetically mapped in other GWAS conducted around the world.

The literature survey demonstrated that 16 out of the 84 QTLs identified in our study had also been detected in previous QTL mapping studies for PH, NKM, TKW, and YM2 traits using the A ×C population in Europe (Griffiths et al., 2012; Ma et al., 2015). The majority of those matches were found for PH (nine QTLs), followed by TKW (three QTLs), NKS, and YM2 with two QTLs for each trait (Ma et al., 2015) (Table 5).

Assessment of presumed novel QTLs based on the field trials of the A ×C DH population

The identification of 43 novel putative QTLs identified in this work underlines the importance of collaborative efforts as the A ×C was developed as a reference DH population within the UK Wheat Genetic Improvement Network (http://www.wgin.org.uk). These results are additional evidence of the importance of extensive germplasm exchange. On the other hand, the identification of new highly significant MTAs underlies the significance of field trials under diverse environmental niches, particularly in those countries where cultivation plays an enormous role in global food security. Hence, the combination of these two factors may lead to the discovery of new important MTAs controlling both plant adaptation-related traits and yield-related traits.

For instance, CIM allowed for the identification of seven novel putative QTLs for PH, including four associations revealed in the Northern region of Kazakhstan (Table 4). One of those QTLs, QPh-AxC.ipbb-5A.1, possibly affects both NKS and TKW in the Northern region as their mapping intervals on chromosome 5A co-localized (Table 4). Similar findings were found from studies in the Southern region as QPh-AxC.ipbb-3D and QPh-AxC.ipbb-4D share locations with QTLs for NKS and TKW, respectively (Table  4). In the search for novel QTLs in TKW, three out of five QTLs were revealed in the Northern region (Table 4). Notably, QTkw-AxC.ipbb-1D (Cadenza) had a matching QTL position with the association for YM2 (QYM2-AxC.ipbb-1D, (Avalon)) in Northern Kazakhstan, and QTkw-AxC.ipbb-3B (Cadenza) and QTkw-AxC.ipbb-6A (Avalon) matched corresponding QTLs for YM2 in Southern Kazakhstan. As seen from Fig. 1, NKS was a highly significant trait for yield performance in Northern Kazakhstan; therefore, it was important to assess whether the identified QTLs for NKS have contributed to the average YM2 over three years. Hence, the DH mapping population was partitioned into groups with the high, middle, and low number of QTLs that carry SNPs with favorable alleles. The evaluation of average YM2 in those groups has demonstrated that having more positive QTLs is highly advantageous over lines with a low number of positive QTLs for plant performances (Fig. 4). Therefore, the study is another confirmation that the accumulation of favorable QTLs is a promising approach in wheat breeding conducted in specific environments (Würschum, Leiser & Langer, 2018; Tshikunde et al., 2019). Still, the results suggest that the higher number of favorable QTLs for NKS does not always seem significant for increased yield, as it was notable for KO (South Kazakhstan) site (Fig. 4). Evidently, despite the benefit of having more positive QTLs for NKS in non-irrigated KB and KA sites, the irrigated KO condition have some masking effect on this advantage for yield performance.