Genome-Wide Association Study in Bread Wheat Identifies Genomic Regions Associated with Grain Yield and Quality under Contrasting Water Availability

Abstract

The objectives of this study were to identify genetic loci in the bread wheat genome that would influence yield stability and quality under water stress, and to identify accessions that can be recommended for cultivation in dry and hot regions. We performed a genome-wide association study (GWAS) using a panel of 232 wheat accessions spanning diverse ecogeographic regions. Plants were evaluated in the Israeli Northern Negev, under two environments: water-limited (D; 250 mm) and well-watered (W; 450 mm) conditions; they were genotyped with ~71,500 SNPs derived from exome capture sequencing. Of the 14 phenotypic traits evaluated, 12 had significantly lower values under D compared to W conditions, while the values for two traits were higher under D. High heritability (H² = 0.5-0.9) was observed for grain yield, spike weight, number of grains per spike, peduncle length, and plant height. Days to heading and grain yield could be partitioned based on accession origins. GWAS identified 154 marker-trait associations (MTAs) for yield and quality-related traits, 82 under D and 72 under W, and identified potential candidate genes. We identified 24 accessions showing high and/or stable yields under D conditions that can be recommended for cultivation in regions under the threat of global climate change.

Keywords: GWAS; Triticum aestivum L.; grain protein content (GPC); grain yield; marker trait association (MTA); quantitative trait loci (QTL); water-limited.

Figure 1

Comparison of averages of phenotypic traits in the wheat panel of 232 accessions, measured under well-watered (W) and water-limited (D) conditions. Key of traits and scales measured units: grain yield (GY, g/m²), grain protein content (GPC, %), spike weight (SW, g), grain weight per spike (GWpS, g), grain number per spike (GNpS, #), number of spikes per plot (NSpP, #), thousand-grain weight (TGW, g), the efficiency of grain filling (GFE, #), leaf area (LA, cm²), specific leaf weight (SLW, mg/cm²), plant height including spike (PH, cm), peduncle length (PedL, cm), spike length (SL, cm), days-to-heading (DH, #).

Figure 2

(a) Surface temperature map of the experimental field consisting of well–watered (W), and water-limited (D) treatments; (b) Experimental field, where the border between two water regimes is marked with red flag; (c) Distribution of crop water stress index (CWSI) within the two treatments.

Figure 3

Principal component analysis (PCA), and biplot vectors of 14 traits measured in the 232 wheat accessions in two water regimes: (a) well-watered (W) and (b) water-limited (D). Key: grain yield (GY), grain protein content (GPC), spike weight (SW), grain weight per spike (GWpS), grain number per spike (GNpS), number of spikes per plot (NSpP), thousand-grain weight (TGW), the efficiency of grain filling (GFE), leaf area (LA), specific leaf weight (SLW), plant height including spike (PH), peduncle length (PedL), spike length (SL), days to heading (DH).

Figure 4

Geographic origins of the best performing accessions. Accessions with high GY in D (red); in D + W (purple); in W (green). The origins of the 232 accessions are marked as small brown circles.

Figure 5

Distribution and density of SNPs mapped to the wheat genome: (a) Distribution of SNPs on 7 chromosomal groups, and the number of mapped SNPs on each wheat genome A, B, and D; (b) Density plot of SNPs mapped to the different wheat chromosomes.

Figure 6

Population structure and kinship matrix of the wheat diversity panel based on 71,571 SNP markers: (a) STRUCTURE analysis and plot of delta K (1 to 10), and the presence of a peak at K = 4 hint at four subgroups. The 1st cluster (red), 2nd cluster (green), 3rd cluster (blue), and 4th cluster (yellow); (b) Heat map of the identity-by-descent based on genomic relationship matrices (GRM).

Figure 7

Multilocus model for comparing GY in the wheat panel: (a) Quantile-quantile (QQ) plot of five GWAS models (GLM, MLM, MLMM, BLINK, and FarmCPU) showing the expected versus observed -log10 (p-value) of each SNP marker (shown as dots). The red line is the expected distribution under the null hypothesis; (b) Minor allele frequency (MAF) plot –log₁₀(p). The p-values for each SNP marker were plotted against their MAF for the FarmCPU model.

Figure 8

The distribution of MTAs for the different traits were revealed under water-limited (D) and well-watered (W) conditions: (a) MTAs in 21 wheat chromosomes; the colors represent an increasing number of MTAs for each trait (e.g., 1 MTA (red), 2 MTAs (green); 3 MTAs (dark green)); (b) distribution of MTA distribution under D and W conditions, in A, B, and D genomes of wheat.

Figure 9

Manhattan Plots across all 21 wheat chromosomes, for nine traits: GY—Grain yield; SW—spike weight; GWpS—grain weight per spike; NSpP—number of spikes per plot; PH—plant height; PedL—peduncle length; GNpS—grain number per spike; DH—days to heading; GPC—grain protein content. SW, GNpS, GWpS, NSpP, PH, PedL, DH (FDR < 0.05), and GPC (FDR < 0.1). The (a) well-watered (W), and (b) water-limited (D) conditions. The FarmCPU model was used for GY, SW, GNpS, GWpS, NSpP, PH, PedL, DH, and the GLM model for GPC. The two horizontal lines indicate 10⁻⁴ and 10⁻⁶ thresholds of significance. The x-axis located 21 chromosomes in wheat. The y-axis located p-values (−log transformed).

Figure 9

Manhattan Plots across all 21 wheat chromosomes, for nine traits: GY—Grain yield; SW—spike weight; GWpS—grain weight per spike; NSpP—number of spikes per plot; PH—plant height; PedL—peduncle length; GNpS—grain number per spike; DH—days to heading; GPC—grain protein content. SW, GNpS, GWpS, NSpP, PH, PedL, DH (FDR < 0.05), and GPC (FDR < 0.1). The (a) well-watered (W), and (b) water-limited (D) conditions. The FarmCPU model was used for GY, SW, GNpS, GWpS, NSpP, PH, PedL, DH, and the GLM model for GPC. The two horizontal lines indicate 10⁻⁴ and 10⁻⁶ thresholds of significance. The x-axis located 21 chromosomes in wheat. The y-axis located p-values (−log transformed).