Genome-Wide Association Mapping Revealed SNP Alleles Associated with Spike Traits in Wheat

Abstract

Wheat (Triticum aestivum L.) is one of the most important crops in the world. Four spike-related traits, namely, spike weight (SW), spike length (SL), the total number of spikelets per spike (TSNS), total kernels per spike (TKNS), and thousand-kernel weight (TKW), were evaluated in 270 F_3:6 Nebraska winter wheat lines in two environments (Lincoln and North Platte, NE, USA). All genotypes in both locations exhibited high genetic variation for all yield traits. High positive correlations were observed among all yield-related traits in each location separately. No or low correlation in yield-related traits was observed between the two environments. The broad-sense heritability estimates were 72.6, 72.3, 71.2, 72.3, and 56.1% for SW, SL, TSNS, TKNS, and TKW, respectively. A genome-wide association study (GWAS) was used to identify SNPs associated with yield traits. In the Lincoln environment, 44 markers were found to be significantly associated with spike-related traits (SW, SL, TSNS, TKNS, and TKW), while 41 were detected in North Platte. Due to the strong significant genotype x environment, no common SNP markers were found between the two locations. Gene annotation of the significant markers revealed candidate genes encoded for important proteins that are associated directly or indirectly with yield traits. Such high genetic variation among genotypes is very useful for selection to improve yield traits in each location separately.

1. Introduction

Bread wheat (Triticum aestivum L.) is one of the most important crops worldwide. Wheat has the fourth-highest global production of all grains and the second-highest net production value of any crop. Global wheat production in the marketing year 2019/2020 was over 765 million metric tons. In comparison to the previous marketing year, this was an improvement of over 30 million tonnes [1]. Due to the diverse genetic architecture and low heritability of wheat grain yield (GY), improving GY is one of the most difficult goals in wheat breeding, and improving GY is usually the most important breeding goal [2]. Wheat yield is considered to have three main components: spikes per unit area, kernel number per spike (KNS), and thousand kernel wheat (TKW). Other traits that affect grain yield include spike length (SL), spike weight (SW), and total spikelet number per spike (TSNS). All of these yield-related traits are a little more heritable (h2) than GY, making them simpler to select in small plots early in breeding programmers [3,4,5]. Only a few yield-related trait results have been used in the selection of wheat lines in ongoing wheat breeding programs. This is due to the great effort and higher labor costs necessary to collect these types of data.

Marker-assisted selection (MAS) is thought to be a crucial strategy for helping to break through the GY bottleneck and increase wheat GY potential [2]. The number and diversity of available alleles and closely related molecular markers determine the potential applications of MAS [6,7]. A better understanding of yield-related traits combined with marker-assisted selection is thought to be the most effective way to increase GY potential [7]. Approximately 65 genes in wheat have been cloned to date, with 40 of them linked to GY and its components [1,7,8,9,10,11]. About 150 functional markers have been converted to kompetitive allele-specific PCR (KASP) formats, which are useful for low-cost, high-throughput genotyping for many of the cloned genes [10].

QTL mapping has been widely utilized to unravel the genetic architecture of agronomically important traits, as most are likely controlled by a large number of genes [12,13,14,15,16,17]. However, there are significant disadvantages to QTL mapping, including poor resolution and a restricted number of alleles that may be examined in each study. By producing high-resolution linkage maps, high-throughput genotyping methods that provide a large number of single nucleotide polymorphism (SNP) data have dramatically enhanced the resolution of QTL mapping [7,9,18,19,20]. GWAS data analysis is focused on linkage disequilibrium (LD), which has a far higher precision than conventional QTL mapping using biparental populations for capturing insights into the genetic architecture of complex traits [7,21,22]. Sallam et al. [22] contrast QTL mapping using biparental populations, GWAS, which makes use of readily available germplasm and avoids the time spent developing segregating populations. There have been multiple publications on quantitative trait loci (QTL) mapping and genome-wide association studies (GWAS) of GY and its related traits, such as [9,11,13,14,23,24,25,26,27,28]. We now believe it has had little impact on breeding programs, and additional research is needed to clarify the function of yield-related traits in the development of new cultivars.

In large genome species, novel genotyping approaches based on next-generation sequencing (NGS) have recently been developed as a viable tool for producing high-density genome-wide markers at a cheap cost per sample [29,30]. Poland et al. [31] developed a GBS methodology that employs two restriction enzymes (PstI and MspI) to significantly reduce genome complexity and provide more uniform libraries for sequencing than the single enzyme protocol [30]. The GBS approach offers several benefits, including low cost, ease of sample handling, and fewer purification processes [29,31,32].

A single nucleotide polymorphism (SNP) can induce genetic variation in crop traits, but it is more likely that many SNPs within a haplotype block are responsible [33]. As a result, SNP-GWAS and haplotype-GWAS are both used to classify genes that regulate complex traits. SNP-GWAS is often used in crop genetic experiments, while Haplotype-GWAS is often used in cross-pollinated crops to identify heterozygous chromosome fragments [34,35]. The plant material in this study has successfully been applied for the GWAS to identify the significant SNP markers associated with grain yield [18] and stem rust resistance [18,36]. The present study aimed to study the genetic variation among GY and its related traits in two different locations, and identify alleles associated with increased GY-related traits for future application in MAS.

2. Materials and Methods

2.1. Plant Materials

In this study, 270 winter wheat genotypes were selected from 2017 F_3:6 nurseries (the preliminary yield trial, referred to as the Nebraska Duplicate Nursery; DUP2017,). These genotypes are the progeny of selections from segregating populations derived from 800 to 1000 crosses primarily among genotypes adapted to the American Great Plains, with an emphasis on genotypes adapted to Nebraska [18,29]. Pedigree information for F_3:6 Nebraska Duplicate Nursery winter wheat is presented in (Supplementary Table S1).

2.2. Experimental Design and Layout

The 270 genotypes investigated in this study were cultivated in two environments in Nebraska: Lincoln (latitude 40.8136° N, longitude 96.7026° W) and North Platte (latitude 41.1403° N, longitude 100.7601° W). The experimental design was an augmented incomplete block design with a single replication and ten incomplete blocks in each location. Each incomplete block included three check cultivars (Goodstreak, Freeman, and Camelot) and 27 experimental lines (a total of 30 lines per incomplete block). As such, the experiment had 300 plots, 270 experimental lines, and the three check cultivars were ten times to provide an estimate of error and spatial variability. Each plot consisted of five 3 m long rows and 0.23 m between rows in each location. The seedling rate was 54 kg ha⁻¹. At the end of the growing season, the GY of each genotype was harvested by a Wintersteiger Classic combine harvester (Wintersteiger Inc., Ankeny, IA, USA), bagged, and stored at room temperature before weighing. The GY was converted to kg ha⁻¹.

2.3. Yield Related Traits

Six spikes from six different single plants were randomly selected from each plot, and the following traits were measured: spike weight (SW, g; measured after physiological maturity by measuring the weight of each individual spike), spike length (SL, cm; measured from the base of the spike to the tip, excluding awns if any), total spikelet number per spike (TSNS; measured by counting the number of spikelets per spike), the total number of kernels per spike (TNKS; was measured by threshing each spike and counting the number of kernels), and 1000-kernel weight (TKW, g; was measured by weighing 1000 randomly selected kernels from each plot). Eltaher et al. [18] previously evaluated grain yield (kg ha⁻¹) for all genotypes in nine locations, including North Platte and Lincoln. We used the grain yield data scores in these two locations with the other data scores in this study to analyze the correlations between spike traits and the grain yield in depth.

2.4. Statistical Analysis of the Studied Yield Components

For each studied trait, data from all the tested two environments were combined and analyzed using the lmerTest R package. The analysis of variance model was:

Y = Check + Environment + Iblock (Environment) + Genotype + G × E + Error

(1)

All factors except checks were fitted as random effects in this model, whereas the check was fitted as a fixed effect. The variance component was also used to estimate broad sense heritability using the following formula:

H2 = Var(G)/[(Var(G) + Var(G × E)/E + Var(E)/(E × R)]

(2)

The R software package “maten” was used to calculate Pearson’s correlations among all GY-related traits based on the genotype performance of each experimental genotype in both environments.

2.5. DNA Extraction, Genotyping-By-Sequencing and SNP Calling

Following the manufacturer’s instructions, DNA was extracted and purified from the 2–3 wheat leaves of two-week old seedlings using BioSprint 96 DNA Plant Kits (Qiagen, Valencia, CA, USA). The GBS method was performed as described by Poland et al. [31]. The SNPs were called using Tassel v5.2.40 software [37]. SNPs were called using the reference genome v1.0 of the “Chinese Spring” genome assembly from the International Wheat Genome Sequencing Consortium (IWGSC). To increase genome coverage and read depth, the raw sequences of the 270 genotypes in this investigation, as well as 6791 other genotypes previously genotyped in our method, were combined and analyzed together [38,39]. The GBS methods identified 206,620 SNPs, which were filtered according to the following criteria: minor allele frequency (MAF > 0.05), maximum missing sites per SNP <20% and maximum missing sites per genotype <20% [38,40]. To avoid miscalculation of allele effects, heterozygous loci were treated as missing [41,42]. As a result of these filters, 28,568 SNPs remained and were used for GWAS in this study.

2.6. Genome-Wide Association Study (GWAS) for the Studied Yield Components

Genome-wide association mapping analysis between all the studied yield components and the filtered SNPs was carried out using TASSEL v5.0 software [37]. The mixed linear model (MLM) with population structure and kinship coefficients was used. The threshold for the p-value (1.98 × 10⁻⁵) was calculated based on the number of markers (P = 1/n, where n is the total number of SNPs used) according to the method reported by Eltaher et al. [18]. For multiple comparison adjustments, the marker-trait associations (MTAs) were tested against Bonferroni corrections (BC) at a significance level of 5%. For all significant MTAs, the percentage of explained phenotypic variation (R²), major and allelic effects were reported. Manhattan plots for yield-related traits were visualized using http://www.bioinformatics.com.cn/en (accessed on 10 June 2022), an online platform for data analysis and visualization. Linkage disequilibrium (r²) was estimated using TASSEL 5.0 between each pair of significant SNPs located on the same chromosome.

2.7. Candidate Genes and Gene Annotation for Yield Component Traits

All significant SNP markers were detected by GWAS in the Ensemble Plants genomic database (http://plants.ensembl.org/Triticum_aestivum/Info/Index (accessed on 10 June 2022) to see if they are located within gene models using the reference genome assembly (IWGSC Ref Seq v2.0). The physical position of each significant SNP was used to find the gene model associated with the target significant SNP. The expression of the candidate genes was tested through the Wheat eFP browser at http://bar.utoronto.ca/efp_wheat/cgi-bin/efpWeb.cgi (accessed on 10 June 2022).

3. Results

3.1. Analysis of Variance for the Yield Components Traits

The analysis of variance results for the 2017 growing season for SW, SL, TSNS, TKNS, and TKW (

Table 1. Analysis of variance (Mean of Square and Properties of significance for F value) and Broad sense heritability for yield-related traits in Lincoln and North Platte.

Source	SW		SL		TNSN		TKNS		TKW
	Mean of Square	Pr (>F)	Mean of Square	Pr (>F)	Mean of Square	Pr (>F)	Mean of Square	Pr (>F)	Mean of Square	Pr (>F)
Genotypes	0.5599	2.2 × 10−15 ***	3.079	2.2 × 10−15 ***	6.55	2.2 × 10−15 ***	118.86	2.2 × 10−15 ***	16.819	6.74 × 10−5 ***
Environment	18.0789	0.0001648 ***	76.152	0.000273 ***	671.28	1.09 × 10−5 ***	788.51	0.00824 ***	39.589	4.57 × 10−5 ***
Iblock (Environment)	0.4371	0.1797082	1.889	0.239888	9.87	0.07659	96.45	0.20626	26.308	0.372
Genotypes × Environment	0.583	<2.2 × 10−15 ***	2.986	2.2 × 10−15 ***	6.94	2.2 × 10−15 ***	144.16	2.2 × 10−15 ***	16.836	6.46 × 10−5 ***
Residuals	0.22		1.14		2.8		55		13.145
Heritability	72.6		72.3		71.2		72.3		56.1

*** Significant at p-value < 0.001.

) revealed significant differences among genotypes, between environments Lincoln and North Platte, and the GEI were significant for all five traits: SW, SL, TSNS, TNKS, and TKW. The broad-sense heritability estimates were 72.6, 72.3, 71.2, 72.3, and 56.1% for SW, SL, TSNS, TKNS, and TKW, respectively.

3.2. Phenotypic Analysis for the Yield Related Traits

The phenotypic distribution was analyzed and presented using boxplots and histograms for all yield-related traits in both environments, Lincoln, and North Platte (Figure 1). The individual genotypic values for GY-related traits in Lincoln and North Platte can be found in (Supplementary Table S2). Continuous and wide-ranging distributions were detected in all yield-related traits under investigation, as one would expect for QTLs. The continuous distributions indicated that the characters are likely polygenic in nature, quantitatively inherited, and measured with some variability. In Lincoln, the highest spike weight of 4.02 g was found in NE17464, and the lowest spike weight of 1.65 g was found for NE17550 with an average of 2.50 g. In North Platte, the highest spike weight of 2.81 g was found in NE17532, and the lowest spike weight of 1.36 g was found in NE17487. For spike length, the longest spike, 12.35 cm, was found in NE17464, and the shortest spike, 6.78 cm, was found in NE17617, with an average of 9.10 cm in the Lincoln, while in the North Platte, the longest spike of 9.88 cm was found in NE17532, and the shortest spike length of 6.81 cm was found in NE17660, with a mean of 8.17 cm. For TSNS, the maximum number of spikelets of 20.33 was found in NE17566, and the minimum number of spikelets of 13.67 was found in NE17563, with an average of 16.87 in the Lincoln. In North Platte, the highest number of spikelets, 17.67, was found in NE17598, and the lowest number of spikelets, 12.66, was found in NE17660, with a mean of 14.97. For TKNS, in Lincoln, the maximum number of TKNS of 58.67 kernels was found in NE17431, and the minimum number of TKNS of 30.17 kernels was found in NE17425, with an average of 42.87 kernels. In North Platte, the highest number of kernels per spike of 52.33 kernels was detected in NE17665, and the lowest number of kernels per spike of 29.66 kernels was noticed in NE17487, with a mean of 36.66 kernels. In Lincoln, the highest value of TKW, 44.87 g, was weighted in the genotype NE17438, and the lowest value of thousand kernel weight, 25.91 g, was detected in genotype NE17609 with a mean of 34.29 g. In North Platte, the greatest value of the thousand kernels weight, 44.17 g, was found in the genotype NE17404, and the smallest value of the thousand kernels weight, 25.46 g, was found in genotype NE17623, with a mean of 32.77 g. Overall, the box plot analysis revealed that all genotypes had higher yield attributes in Lincoln than in North Platte.

3.3. Correlation Coefficients for Yield-Related Traits

Eltaher et al. [18] previously investigated the grain yield of 270 F_3:6 Nebraska winter wheat genotypes in eight Nebraska and one Kansas environment. In each environment, for each trait (Figure 2), there was no significant correlation in the yield-related traits between the two environments. In each environment, highly positive and significant correlations were found among spike-related traits. In Lincoln, the highest positive significant correlation was found between GY and TKNS (r = 0.85 **), while the correlation between SL and TSNS was the highest with r of 0.82 ** in North Platte. GY was significantly correlated with all traits in Lincoln, while, in North Platte, it was significantly correlated with SL, TSNS, and TKNS. Notably, TKW was not significantly correlated with any trait in Lincoln, while it was negatively correlated with TKNS, TSNS, and SL in North Platte.

3.4. Genome Wide Association Studies for Yield-Related Traits

The analysis of GWAS revealed 44 significant SNPs associated with yield-related traits in Lincoln and 41 in North Platte (Figure 3a). SW, TKNS, and TSNS had higher QTL in Lincoln than in North Platte, while SL and TKW had a higher number of QTL in North Platte than in Lincoln (Figure 3b). There were no common markers between the two environments (Lincoln and North Platte), as shown in Figure 3c.

In Lincoln, 44 markers were identified as being significantly associated with yield-related traits SW, SL, TSNS, TKNS, and TKW (Supplementary Table S3). Additionally, the summarized GWAS analysis for yield-related traits is presented in

Table 2. Significant SNP loci, phenotypic variation (R²) and allele effect identified for yield-related traits in Lincoln.

Traits	Number of SNPs	Log 10 p Value		Range of R2		Range of Allele Effect
		Min	Max	Min R2	Max R2	Min Allele Effect	Max Allele Effect
SW	8	3.63 × 10−7	3.27 × 10−5	5.75	8.87	0.19 A	0.42 A
SL	5	1.16 × 10−7	1.17 × 10−5	6.07	9.01	0.43 T	0.62 T
TSNS	10	3.63 × 10−7	1.01 × 10−5	6.03	8.69	0.64 T	1.44 A
TKNS	15	3.63 × 10−7	3.27 × 10−5	5.70	8.96	2.90 G	6.71 A
TKW	6	1.1 × 10−6	1.77 × 10−5	5.92	8.74	2.40 C	2.89 C

. Manhattan plots for yield-related traits are presented in Figure 4. Significant markers were found on nine different chromosomes, including 2A, 3B, 5A, 5D, 6A, 6D, 7A, 7B, and 7D. SL had the fewest significant SNPs for yield-related traits, while TSNS had the most significant SNPs. The −log 10 p value varied from 1.161 × 10⁻⁷ to 3.27 × 10⁻⁵ for markers S7B_165529101 and S5A_47458032, respectively. The phenotypic variation R² varied from 5.69 to 9.01%. The three significant markers S5A_46628103, S7B_607427421, and S5A_380823821 were found to be significantly associated with SW. The lowest allele effect was accounted for in marker S6D_469537865, which was significantly associated with TNKS. The maximal allele effect (A) of 6.71 was observed.

In North Platte, 41 markers were found to be significantly associated with yield-related traits SW, SL, TSNS, TKNS, and TKW, as described in (Supplementary Table S3). Additionally, the summarized GWAS analysis for yield-related traits is presented in (

Table 3. Significant SNP loci, phenotypic variation (R²) and allele effect identified for yield-related traits in North Platte.

Traits	Number of SNPs	Log 10 p Value		Range of R2		Range of Allele Effect
		Min	Max	Min R2	Max R2	Min Allele Effect	Max Allele Effect
SW	5	3.72 × 10−7	2.02 × 10−5	6.91	8.26	0.40 C	0.65 C
SL	10	3.91 × 10−7	3.06 × 10−5	5.90	8.36	0.12 T	0.25 A
TSNS	9	3.72 × 10−7	2.67 × 10−5	5.45	8.03	0.46 G	1.55 C
TKNS	8	3.91 × 10−7	2.11 × 10−5	7.00	13.36	1.85 A	3.06 C
TKW	9	1.58 × 10−7	1.88 × 10−5	9.54	12.03	1.44 G	2.84 T

). Manhattan plots for yield-related traits are presented in Figure 5. All these significant markers were distributed on nine different chromosomes: 1A, 2A, 3B, 4A, 5D, 6A, 6B, 7A, and 7B. The lowest number of significant SNPs, 5 SNPs for yield-related traits, was detected in SL, while the highest number of SNPs, 15, was observed in TKNS. The −log 10 p value ranged from 3.72 × 10⁻⁷ to 3.0618 × 10⁻⁵ for markers S5D_61792984 and S3B_60737182, respectively. The R² ranged from 5.45 to 13.36% for marker S5D_62479367 and marker S3B_62315382, respectively. The minimum allele effect (T) of 0.12 was observed in four markers, S3B_60737182, S6A_19036296, S7A_610993044, and S3B_64172577, which were significantly associated with SL. However, the maximum allele effect (C) of 3.06 was found in both markers S5D_72377429 and S5D_61792984, which were significantly associated with TNKS.

3.5. Common Markers Associated with Yield-Related Traits

Table 4. Common marker associated with yield-related traits using the mixed linear model (MLM) at the significance level of 5% bonferroni correction (BC 5%) using 11991 SNPs in Lincoln.

Traits	Marker ID	Chromosome	Position	Log 10 p Value	R2	Target Allele	Allele Effect
SW	S5A_380823821	5A	380823821	1.17 × 10−5	6.06	T	0.43
SL				4.8 × 10−6	5.9		0.19
TNKS				1.06 × 10−5	6.32		3.06
TSNS				1.01 × 10−5	6.03		0.64
SW	S5A_46628103	5A	46628103	4.68 × 10−7	5.82	A	0.19
TNKS				3.27 × 10−5	6.47		3.00
SW	S5D_548379143	5D	548379143	1.05 × 10−6	9.11	T	0.50
SL				3.63 × 10−7	8.87		0.23
TNKS				3.63 × 10−7	8.96		3.53
TSNS				3.63 × 10−7	8.69		0.74
TNKS	S6D_469537865	6D	469537865	4.8 × 10−6	6.53	A	6.71
TSNS				1.06 × 10−5	6.56		1.44
TNKS	S7B_164151731	7B	164151731	1.06 × 10−5	5.87	C	3.46
TSNS				1.29 × 10−5	6.51		0.78
SW	S7B_165529101	7B	165529101	1.16 × 10−7	6.19	C	0.49
SL				1.06 × 10−5	5.78		0.22
TNKS				7.8 × 10−6	6.50		3.52
TSNS				4.68 × 10−6	7.12		0.79
TNKS	S7B_181032630	7B	181032630	1.29 × 10−5	6.17	A	3.50
TSNS				4.8 × 10−6	6.84		0.79
SW	S7B_329792071	7B	329792071	5.5 × 10−6	6.91	A	0.54
SL				1.29 × 10−5	6.56		0.24
TNKS				4.68 × 10−6	6.87		3.80
TSNS				3.77 × 10−6	7.79		0.87
SW	S7B_607427421	7B	607427421	3.27 × 10−6	5.70	G	0.19
TNKS				1.29 × 10−6	5.70		2.90
SW	S7D_485517060	7D	485517060	4.83 × 10−6	6.67	T	0.62
SL				3.77 × 10−6	6.61		0.28
TNKS				3.77 × 10−6	7.01		4.43
TSNS				7.8 × 10−6	6.67		0.92

describes ten common markers found to be significantly associated with yield-related traits SW, SL, TSNS, and TKNS in Lincoln. Chromosome 3B had four comment markers associated with more than one trait. Five common markers, S5A_380823821, S5D_548379143, S7B_165529101, S7B_329792071, and S7D_485517060, were found to be significantly associated with all yield-related traits, except TKW. Two markers were on chromosome 7B, and one marker was on each of chromosomes 5A, 5D, and 7D. The R² ranged from 5.78% in marker S7B_165529101 with trait SL, to 9.11% in marker S5D_548379143 with trait SW. Two markers, S5A_46628103 and S7B_607427421, were found to be significantly associated with SW and TNKS. Three markers, S6D_469537865, S7B_164151731, and S7B_181032630, were found to be significantly associated with TNKS and TSNS.

Ten markers were also found to be significantly associated with yield-related traits SW, SL, TSNS, and TKNS in North Platte (

Table 5. Common marker associated with yield-related traits using the mixed linear model (MLM) at the significance level of 5% bonferroni correction (BC 5%) using 11991 SNPs in North Platte.

Trait	Marker ID	Chromosome	Position	Log 10 p Value	R2	Target Allele	Allele Effect
SL	S3B_60737182	3B	60737182	3.06 × 10−5	7.00	A	0.12
TNKS				2.11 × 10−5	7.00		1.85
SW	S3B_62315382	3B	62315382	3.91 × 10−7	8.36	C	0.15
SL				1.96 × 10−6	8.26		0.40
TNKS				3.91 × 10−7	13.36		1.85
TSNS				3.3 × 10−6	8.03		1.55
SW	S3B_62315407	3B	62315407	2.5 × 10−7	8.36	T	0.15
SL				2.02 × 10−6	8.26		0.40
TNKS				2.5 × 10−6	9.23		2.23
TSNS				1.96 × 10−6	8.03		0.55
SL	S3B_64172577	3B	64172577	2.11 × 10−5	6.01	T	0.12
TNKS				1.61 × 10−5	8.63		1.90
TSNS				2.02 × 10−5	6.31		0.48
SW	S5D_61792984	5D	61792984	1.09 × 10−5	6.17	C	0.21
SL				3.72 × 10−7	7.96		0.65
TNKS				1.09 × 10−5	9.09		3.06
TSNS				3.72 × 10−7	6.85		0.83
SW	S5D_62479367	5D	62479367	2.02 × 10−5	6.91	G	0.58
TSNS				2.52 × 10−5	5.45		0.71
SW	S5D_72377429	5D	72377429	1.61 × 10−5	5.96	C	0.21
SL				3.3 × 10−6	7.77		0.64
TNKS				1.61 × 10−5	8.63		3.06
TSNS				2.02 × 10−5	6.18		0.80
SL	S6B_668517613	6B	668517613	1.61 × 10−5	5.90	A	0.20
TNKS				2.11 × 10−5	7.56		3.05
TSNS				2.09 × 10−5	5.75		0.77
SL	S7A_610993044	7A	610993044	2.11 × 10−5	6.00	G	0.12
TSNS				2.67 × 10−5	5.87		0.46
SL	S7B_729441244	7B	729441244	2.11 × 10−5	6.62	A	0.25
TSNS				1.87 × 10−5	6.60		0.94

). Chromosome 7B had a set of six markers associated with more than one trait. Four markers, S3B_62315382, S3B_62315407, S5D_61792984, and S5D_72377429, were found to be significantly associated with all yield-related traits, except TKW. Two markers, S3B_64172577 and S6B_668517613, were significantly associated with SL, TNKS, and TSNS. Three markers, S3B_60737182, S7A_61099304,4 and S7B_729441244, were associated with SL and TNKS. The marker S5D_62479367 was associated with SW and TSNS.

By considering the SNP markers associated with (GY), which were described in Eltaher et al. [18], four markers were found in common between GY and yield-related traits in Lincoln (

Table 6. Repetitive SNPs and their chromosome revealed by GWAS identified in GY described by Eltaher et al. [18] and yield-related traits of the investigated environments (Lincoln and North Platte).

Grain Yield	Polygenic Marker	Chromosome	Yield-Related Traits
Lincoln	S5D_548379143	5D	(SL, SW, TNKS and TSNS)
	S6D_469537865	6D	(SW, TNKS and TSNS)
	S7B_329792071	7B	(SL, SW, TNKS and TSNS)
	S7D_485517060	7D	(SL, SW, TNKS and TSNS)
North Platt	S3B_60737182	3B	(SL and TNKS)
	S3B_62315382	3B	(SL, SW, TNKS and TSNS)
	S3B_62315407	3B	(SL, SW, TNKS and TSNS)
	S3B_64172577	3B	(SL TNKS and TSNS)
	S5D_61792984	5D	(SL, SW, TNKS and TSNS)
	S5D_62479367	5D	(SW and TSNS)
	S5D_72377429	5D	(SL, SW, TNKS and TSNS)
	S7B_729441244	7B	(SL and TSNS)

). Four markers were significantly associated with GY and spike-related traits in Lincoln, while eight markers were common between GY and spike-related traits in North Platte. Notably, no shared markers were found on the A genome and most of the shared markers between GY and spike-related traits were located on the D and B genomes in both locations, with six markers each.

3.6. Gene Annotation for Yield-Related Traits

The candidate genes associated with the significant SNPs detected by GWAS in both locations are presented in Supplementary Table S4. The result of gene annotation revealed nineteen candidate genes.

In Lincoln, ten candidate genes were detected on different chromosomes. Out of the ten candidate genes, six were related to protein-coding with unknown function or non-translating coding sequences (CDS). The important common SNP marker S7D_485517060 (located on chromosome 7D), which was found to be significantly associated with all yield-related traits except TKW, was also found to be associated with grain yield in Lincoln by Eltaher et al. [18]. This SNP marker was annotated to TraesCS7D02G375100, and this gene had an effective value in the spike traits and was turned on in the spike (Figure 6). This gene translated to CDP-choline: 1,2-diacylglycerol cholinephosphotransfer, which plays a great role in the accumulation and deposition of triacylglycerols in the starchy endosperm of wheat grain, especially in the aleurone layers.

A set of 14 gene models were identified in North Platte. Of the 14 SNPs, nine genes were related to protein-coding with unknown function or non-translating CDs. The S3B_64172577 marker, which is associated with SL, TSNS, and TKNS, was found to be TraesCS3B02G095300.1, which encodes protein kinase-like domain superfamily.

4. Discussion

4.1. Genetic Variation for Yield-Related Traits

Analysis of genotypic responses across diverse production environments is valuable to classify suitable genotypes and test environments for enhanced breeding and cultivar development [18,19,43,44,45,46]. The analysis of variance in the present study found highly significant (p ≤ 0.001) differences among genotypes, test environments, and GEI effects. The variable genotypic response across test environments was recommended by a significant GEI in this study. Hence, identifying appropriate wheat genotypes that perform well in test environments is critical. Furthermore, environmental effects contributed much more to total variability than genotype and GEI effects, as indicated by the largest sum of squares for the traits studied. Moreover, extremely significant variation in genotypes, environments, and GEI reflected the differences in genotypes within a single environment, as well as between two environments. This result suggests that there is a high level of genetic variation, allowing plant breeders to use the full potential of genetic and environmental variation while supporting the selection process between genotypes [18,29]. The highest genotype for yield-related traits differed by environment, and the phenotypic correlation among environments varied due to the highly significant GEI interaction. As a result, the breeding program for increasing high grain yield may fluctuate depending on the environment [18,44]. All traits had high heritability estimates, except TKW, which had moderate heritability, indicating that the selection for these traits will be fruitful for the breeding program. Such high genetic variation among genotypes for each trait, in addition to high heritability, made the identification of genetic variants using GWAS feasible.

The significant GEI interaction hindered finding the most suitable genotypes for both environments, and specific breeding programs should be performed according to the environment. This notion is supported by the lack of correlation between the two environments for each trait. Significant correlations were only found among the traits within each environment separately. The same finding was observed in Eltaher et al., [18,44]. This highly significant GEI can be explained by the differences in the precipitation, snow cover, and temperature in each environment during the growing season. The high GEI is frequently interpreted as indicating that a site-specific breeding program for increasing grain yield is necessary for optimal improvement in each target environment [19,44].

4.2. Genome Wide Association Mapping for Yield-Related Traits

The most interesting SNP markers linked with yield-related traits SW, SL, TNKS, TSNS, and TKW were found on chromosomes 1A, 2A, 3A, 5A, 5D, 7B, and 7D in both environments in the present study. Environments have a major impact on grain yield and its components, and it is hard to select high-yielding lines in small plots early in a breeding program. Environments, on the other hand, have a significantly smaller impact on yield components, and some more stable QTL for these traits have been discovered, which is consistent with previous results [4,18,19,47,48,49,50,51,52,53,54].

All markers detected by GWAS in Lincoln had an R² of 10%, indicating that these markers had minor effects on spike-related traits in Lincoln, while a set of nine markers in North Platte had a major effect on spike-related traits (R² > 10%).

Markers with pleotropic effects, which were associated with more than one trait, were detected in each environment, separately. The non-shared markers for yield and its component traits between the two environments were due to the strong GEI. Many previous studies have found a strong relationship between grain yield and its component traits [4,51,52,55]. However, several markers for different traits were detected in the same or neighboring positions (3A, 5A, 5D, 6A, 7A, and 7B) as those identified in previous studies [56,57,58]. These results revealed the pleiotropism of markers for the GY and related traits, which may be due to the complex and often compensatory relationships among these traits [52,56,57,58,59,60].

Fifteen makers were found to be significantly associated with TKW. These SNPs were located on chromosomes 2A, 3B, and 7A in Lincoln and on 1A, 2A, 4A, 6A, and 7B in North Platte. Previous studies found QTLs associated with TKW traits were mainly distributed on chromosomes 2A, 3B, and 7A [51,61,62,63,64,65]. Significant markers for TKW were recorded on chromosome 7B, which was previously reported to have significant associations with TKW, using the markers Xwmc606, Xgwm537, wPt1715 and wPt2449 in a collection of tetraploid durum (T. durum L.) wheat genotypes [52,66].

Remarkably, significant markers were found to be associated with more than one trait, indicating that these markers have pleotropic effects. However, these markers have pleotropic effects in the specific environment due to the highly significant effects of G × E. These markers can be converted to Kompetitive-specific, allele-specific PCR for further validation in different genetic backgrounds before using them in marker-assisted selection.

4.3. Gene Annotation for Yield-Related Traits

The candidate gene TraesCS2A02G477600 was translated into the Bifunctional inhibitor/plant lipid transfer protein/seed storage helical domain superfamily for TKW. These protein families are large. The proteins are annotated as bifunctional inhibitors or cereal seed allergens that belong to the seed storage helical domain. These proteins can be found at high concentrations in the seeds of both mono-and dicotyledonous plants and are an important component of the normal human diet [67,68]. Seed storage proteins are proteins that accumulate significantly in the growing seed and serve as nitrogen, carbon, and sulfur storage reserves. These proteins are rapidly activated during seed germination and serve as the seedlings’ primary source of nutrients. After learning how important this gene is for preserving protein inside seeds, it is no surprise that it is also associated with TKW. As a result, it should be stated that selecting for this gene should increase protein storage within the grain of wheat, resulting in grains with high-quality features in terms of protein quantity and possibly quality.

The common candidate gene TraesCS6D02G398100.1 encodes the nucleic acid-binding, Oligonucleotide/Oligosaccharide-Binding (OB) Fold Proteins, Replication Factor A protein-like and winged helix DNA-binding domain superfamily. Within the Oligonucleotide/Oligosaccharide-Binding (OB) Fold Proteins, the nucleic acid-binding superfamily is the largest, and proteins with this motif are involved in practically every single-stranded DNA or RNA (ssDNA/ssRNA) that is present or needs to be manipulated [69]. DNA replication, recombination, repair, and telomere homeostasis are just a few of the biological processes that OB-fold proteins have been shown to have a role in. The common candidate gene, TraesCS7B02G135400.1, is translated into the BURP protein domain. The BURP domain is a C-terminal protein domain with four common members: BNM2, USP, RD22, and PG1. Plant-specific BURP domain-containing proteins have only been discovered so far, implying that their functions may be plant-specific [70]. BURP domain-containing proteins have been found in many species, such as rice (Oryza sativa L.), soybean (Glycine max (L.) Merr.), maize (Zea mays L.), and sorghum (Sorghum bicolor (L.) Moench) [71,72,73]. Plant BURP domain-containing proteins play a vital role in plant metabolism and development, although their expression patterns are variable, and several of their functions are unclear [74,75]. For example, VfUSP is an abundant non-storage seed protein from field beans (Vicia faba L.) with unknown activities that is expressed during zygotic embryogenesis and in vitro embryogenesis [76,77]. The two common candidate genes, TraesCS3B02G095300.1 and TraesCS3B02G095300.2, translate into the protein kinase like-domain superfamily. They play a role in diurnal and circadian regulation; cell cycle regulation; developmental processes; vesicle transport and channel activity modulation, and cellular metabolic regulation [78,79,80,81].

The common gene TraesCS7D02G375100 is associated with spike-related traits and GY in Lincoln, which is translated to CDP-choline: 1,2-diacylglycerol cholinephospho transfer. This gene plays a role in the storage of lipids inside the wheat grain [82,83]. Triacylglycerols (TAGs) are the major storage lipids in seeds, although they are only minor components in cultivated cereals (with the exception of oats, Avena sativa L.). They are abundant in the aleurone layer and scutellum of the wheat embryo, accounting for 60–80% of total lipids in these tissues [83]. On the other hand, TAGs make up around a third of the total lipids in the starchy endosperm tissue that white flour is produced from [83,84]. Their synthesis and deposition are poorly understood. Although lipid droplets have been observed in starchy endosperm cells, it is also possible that certain lipids (including TAGs) are transferred from the aleurone and embryo to the flour during milling [83,85]. This gene was found to be highly expressed in spike and grain, which agrees with the strong association between this gene and spike-related traits detected by GWAS in our study.

In conclusion, the highly significant G × E interaction found between the two locations hindered selecting the highest yielding genotypes for breeding programs, as well as sharing markers for marker-assisted selection in both locations. Therefore, it is highly recommendable that each location should have its own specific breeding program supported by promising markers to accelerate breeding programs for wheat improvement.