Genome‐Wide Association Mapping of Leaf Rust Response in a Durum Wheat Worldwide Germplasm Collection

Thirteen durum wheat accessions showed resistance to all Puccinia triticina races tested GWAS revealed 88 SNPs (37 loci) associated with leaf rust response in durum wheat Associations were identified on all chromosomes except 1B and 7B GWAS revealed 14 previously uncharacterized loci for leaf rust resistance

D urum wheat (2n = 4x = 28) is an important cereal crop grown in many parts of the world, especially in the Mediterranean basin where ~50% of global production and 75% of the growing area are located (Elias and Manthey, 2005). The world durum wheat production was estimated at 33.5 million metric tons in 2014.
Durum wheat is an important crop that is concentrated in localized areas, often in developing countries, where it represents a large portion of total wheat planted as well as a major staple food used for pasta, couscous, and flat bread. Moreover, because of its adaptability to arid climate conditions, marginal soils, and relatively low water requirements, improvement of durum wheat production should be an agricultural and economic priority to ensure food security in these regions. Despite the broad adaptability of durum wheat, its production is often limited by different fungal diseases including rusts, septoria leaf blotch, fusarium head blight, and root rot (Nachit 2000;Nsarellah et al., 2000;Singh et al., 2005).
Leaf rust is a significant disease affecting wheat production worldwide. Durum wheat has been traditionally considered more resistant to Pt than bread wheat (T. aestivum L.; 2n = 6x = 42) in many regions worldwide. However, races of the leaf rust pathogen, virulent to widely grown durum cultivars in several production areas, are increasingly impacting durum production (Singh et al., 2004;Goyeau et al., 2006;Huerta-Espino et al., 2009). A Pt race, BBG/BN, in northwestern Mexico that appeared in 2001, with virulence to Lr72, overcame the resistance of widely adapted durum wheat cultivars from the CIMMYT breeding program that had been effective in Mexico for >25 yr (Singh et al., 2004;. Increased susceptibility of durum wheat to leaf rust has also been observed in the Mediterranean basin (Martinez et al., 2005;Goyeau et al., 2012), the Middle East (Ordoñez and Kolmer, 2007a), and Chile (Singh et al., 2004). In the United States, a Pt race, designated BBBQD using the North American Pt nomenclature system (Long and Kolmer, 1989) and possessing a similar virulence pattern to previously identified Mexican races, was collected from California durum fields in 2009. In 2013, an isolate of the same race was found in Kansas. The occurrence of this race in these regions increases the likelihood of its spread to the northern Great Plains (Kolmer, 2015a) and most importantly, to North Dakota where 58% of the total US durum wheat is produced (USDA-National Agricultural Statistics Service, 2015).
A number of studies have determined that the Pt populations predominant on bread wheat differ from those found on durum wheat. The Pt isolates from bread wheat are often avirulent on durum wheat (Singh, 1991;Huerta-Espino and Roelfs, 1992;Ordoñez and Kolmer, 2007a). In contrast, Pt isolates collected from durum wheat are avirulent to all but the most susceptible bread wheat cultivars (J.A. Kolmer, unpublished data, 2013). Huerta-Espino and Roelfs (1992) determined that the majority of the Pt collections from durum wheat were race BBB, indicating that they were avirulent to all resistance genes carried by the first three sets of international leaf rust differentials (Long and Kolmer, 1989). Isolates collected from durum wheat in several countries share highly similar virulence patterns on 'Thatcher' near-isogenic lines, suggesting a common origin (Ordoñez and Kolmer, 2007a;Ordoñez and Kolmer, 2007b). However, some Pt isolates collected from durum wheat in Ethiopia have a distinct virulence phenotype with avirulence to commonly susceptible bread wheat hosts including Thatcher. These isolates also grouped separately for simple-sequence repeat (SSR) genotypes (Kolmer and Acevedo, 2016) and DNA sequence compared with other isolates collected on durum and bread wheat in Europe, South America, and Mexico (Liu et al., 2014).
In addition, the adult plant resistance slow rusting gene Lr46 has been reported in CIMMYT durum wheat (Herrera-Foessel et al., 2011). However, the extensive use of single-race-specific resistance genes caused rapid selection of new virulent Pt races only a few years after their deployment in commercial cultivars. For instance, in 2008, a Pt race was detected in Mexico that carried virulence on Lr27+Lr31 (Huerta-Espino et al., 2009). Similarly, a variant of the old Mexican race BBG/BN, designated as BBB/ BN_Lr61Vir , has virulence to Lr61. Both Lr14a and LrCamayo are still effective against the current Mexican races, however virulence on Lr14a has been reported in France (Goyeau et al., 2006), Spain (M. Acevedo, unpublished data, 2015), Tunisia (Gharbi et al., 2013), and Morocco (M. Acevedo, unpublished data, 2015).
The mapping of plant disease resistance genes and quantitative trait loci (QTL) has been traditionally performed through biparental mapping populations. Recently, as a result of the progress in high-throughput genotyping technologies and the improvement of statistical programs, GWAS or AM has been used as an alternative approach to biparental mapping (Yu et al., 2006;Zhu et al., 2008). Association mapping is based on revealing correlations between phenotype and genotype in a germplasm collection (Zondervan and Lon 2004). The GWAS takes advantage of the linkage disequilibrium (LD) between alleles to identify molecular markers significantly associated with a trait response. The GWAS uses the recombination events that occur throughout the evolutionary history of a germplasm. This allows the breakup of the LD blocks within the genome and results in a faster decay of the LD in the AM panels than in recombinant inbred lines and doublehaploid populations, in which only the allelic diversity that segregates between the parents can be assessed. Therefore, GWAS can identify associated loci with the trait response at a much higher mapping resolution than biparental mapping (Rafalski, 2002;Nordborg and Weigel, 2008;Zhu et al., 2008;Neumann et al., 2011).
In durum wheat, AM has been used to dissect the genetic basis of important agronomic traits including grain yield, yellow pigment, root architecture, plant height, and drought and salinity tolerance (Reimer et al., 2008;Canè et al., 2014, Hu et al., 2015Maccaferri et al., 2010a;Turki et al., 2015). Moreover, GWAS has been used to identify markers associated with disease resistance to Fusarium head blight, leaf rust, and stem rust (Ghavami et al., 2011;Maccaferri et al., 2010b;Letta et al., 2014). However, the previous AM analysis for leaf rust response only included elite germplasm (cultivars and advanced lines) that was genotyped using only 225 SSR markers. Wheat landrace germplasm collections may carry new genes for resistance to diseases, including leaf rust, stem rust, and stripe rust, since the use of landraces in the modern breeding programs has been relatively rare (Reif et al., 2005;Bonman et al., 2007;Newton et al., 2010;Bux et al., 2012;Gurung et al., 2014). Recently, a high-density, SNP-based consensus map of tetraploid wheat was developed by Maccaferri et al. (2015), which will increase the effectiveness of GWAS and QTL meta-analyses.
The current study describes the first durum wheat leaf rust GWAS using a highly diverse germplasm panel comprised of a worldwide collection of landraces, cultivars, and breeding lines genotyped using the Illumina iSelect 9K wheat array. This study provides relevant insight into the genetic basis underlying resistance in durum wheat to North American Pt races.

Genetic Material
A durum wheat collection of 496 accessions from the USDA-Agricultural Research Service (ARS)-National Small Grain Collection (NSGC) was evaluated for leaf rust resistance in this study. These accessions were originally sourced from 67 countries including accessions from Africa (96), Asia (172), Australia (7), North America (42), South America (34), Europe (140), and unknown origin (5). The collection includes 234 landraces, 55 cultivars, 128 cultivated lines, 77 breeding lines, and 2 of unknown accession type. Seeds used in this project were obtained from single plant selections increased in a nursery at the USDA-ARS-NSGC and Potato Germplasm Research Unit, Aberdeen, ID (Supplemental Table S1).

Phenotyping at Seedling Stage
The disease screenings at the seedling stage of plant development were conducted at three locations: (i) North Dakota (ND) Agricultural Experiment Station Greenhouse Complex, Fargo, ND, USA; (ii) the USDA-ARS Cereal Disease Laboratory (CDL) in Saint Paul, MN, USA; and (iii) the CIMMYT, Mexico (El Batán, State of Mexico).
The durum wheat collection was evaluated for response to Pt races BBBQD (California isolate), BBBQJ (Mexican isolate, BBG/BP), BBBDB (USA, Race 1), MBDSD (North Dakota isolate), MCDSS (North Dakota isolate), and a mixture of bread-wheat-type races predominant in the United States (MHDSB, MFPSB, MLDSB, TBBGG, TFBJQ, and TFBGQ). The virulence-avirulence profile of the rust races was based on infection types on seedlings of Thatcher wheat differentials that are near isogenic for single-resistance genes following the race nomenclature of Long and Kolmer (1989) (Table 1). The isolates BBBQD, BBBQJ, MBDSD, and MCDSS were collected from durum wheat. Races MBDSD and MCDSS are known to be virulent on bread wheat as well. All other isolates were collected from bread wheat germplasm.
The experiment performed at the CDL with races BBBDB, BBBQD, and the race mixture (MHDSB, MFPSB, MLDSB, TBBGJ, TFBJQ, and TFBGQ) were conducted in augmented design with one replicate. Five to seven seeds of each durum wheat accession were planted in 50-cell trays containing vermiculite. The susceptible checks Little Club and RL6089 were included in the experiment several times. Urediniospore increase, inoculation, incubation, and greenhouse conditions were as previously described by Kolmer and Hughes (2013). The screening using race BBBQJ was conducted at CIM-MYT in augmented design and the experiment was repeated twice. The seedlings were grown as previously described by Loladze et al. (2014). Seven to 10 seedlings from each accession were grown at 22C in the greenhouse in pots with four accessions planted per pot with a soil mix consisting of one part peat moss, one part sand, and one part black soil. Pots were fertilized twice with a urea fertilizer (5 g per 10 L of water) 4 to 5 d after planting and 2 to 3 d after inoculation. Seedlings were inoculated with a suspension of urediniospores and light mineral oil (Soltrol 170; 1 mg of urediniospores per 1 mL of oil) using a hydrocarbon propellant pressure pack. The oil was allowed to evaporate from the leaves for at least 30 min before incubating in dark at 20C for 18 h. Following incubation, plants were transferred to the greenhouse at 23°C with 16 h photoperiod.
Leaf rust infection types (ITs) were assessed 12 d postinoculation on the second leaves of seedlings using a 0-to-4 scale (Long andKolmer, 1989, McIntosh et al., 1995) where IT 0 = no visible symptom, ; = hypersensitive flecks, 1 = small uredinia with necrosis, 2 = small-to medium-size uredinia surrounded by chlorosis, 3 = medium-size uredinia with no chlorosis or necrosis, and 4 = large uredinia with no necrosis or chlorosis. Larger and smaller uredinia than expected for each IT were designated with + and −, respectively. Symbols C and N are used to indicate more than usual degrees of chlorosis and necrosis, respectively. Heterogeneous infection type evenly distributed over the leaf surface was designated as X (mesothetic reaction). Accessions with ITs of 0 to 2+ and X were considered resistant, while 3 and 4 scores were considered susceptible (McIntosh et al., 1995;Long and Kolmer, 1989).
To account for multiple ITs in a single plant at seedling stage, the 0-to-4 scale for leaf rust was converted to a linearized 0-to-9 scale using the weighted mean of the most and the least predominant IT on the same leaf surface . Ratings of 0 to 6 were classified as resistant IT and 7 to 9 were considered as susceptible IT.
Analysis of variance (ANOVA) on the linearized 0-to-9 ratings was performed using SAS software 9.3 (SAS Institute, 2011) to test for the homogeneity of variances for experiments and replicates. The homogeneity of variance test resulted in combining IT data from the two North Dakota experiments on race BBBQD and the mean was used for the GWAS, while the CDL experiment for the same race was analyzed as a separate experiment. For all other experiments with two or three replicates on the rest of the races, the mean across replicates for each accession was used.

Evaluations of Adult Plants
The accessions were field tested in 2012 at St. Paul (MN_ StP_F) and Crookston (MN_Cr_F) in Minnesota and at CIMMYT Centro Experimental de Norman E. Borlaug station in Ciudad Obregón, Sonora (MX_Ob_F) and at CIMMYT Headquarters, El Batán experimental station (MX_EB_F) in 2013.
In MN_StP_F and MN_Cr_F, 50 to 60 seeds of each accession were planted in a 3-m row and 30 cm apart. Rows of the susceptible wheat cultivars, Thatcher, Little Club, and Morocco were planted perpendicular to the durum wheat entries. The spreader rows were inoculated with a mixture of six Pt races (MHDSB, MFPSB, MLDSB, TBBGJ, TFGJQ, and TFBGQ) from bread wheat that were present in the Great Plains region. The common spring wheat cultivars Thatcher, Tom, Verde, and Knudson were planted every 100 entries as checks. Field inoculation was conducted as described by Kolmer (2015b).
The germplasm at all locations was evaluated when at least 70 to 80% of the flag leaf area of the susceptible checks was covered by uredinia. The accessions were evaluated using the modified Cobb scale (Peterson et al., 1948) ~2 mo after planting. The scoring was based on both the disease severity (the percentage of tissue infected) and the plant response to infection. Plant response was recorded as resistant (R), moderately resistant (MR), moderately susceptible (MS), and susceptible (S) reactions (McIntosh et al., 1995). In some cases, infection responses were a combination between any two categories on the same leaf with most predominant infection response first followed by the least predominant one. For instance, MSMR referred to overlapping of MS and MR categories, where MS was observed more frequently than MR response. All accessions with R, RMR, MR, and MRMS infection responses were considered resistant. In addition, accessions with MSMR and MS infection response but severity lower than 20% were also considered resistant.
For AM analysis, disease severity and infection response data were combined in a single value as the coefficient of infection, which was calculated by mul-

Single-Nucleotide Polymorphism Marker Genotyping and Analysis
The durum wheat collection were genotyped through the Triticeae Coordinated Agricultural Project using the Illumina iSelect 9K wheat array (Cavanagh et al., 2013) at the USDA-ARS genotyping laboratory in Fargo, ND. A total of 5490 high-quality polymorphic SNPs were originally selected. Marker data are available through the USDA-NIFA-funded T-CAP project (http://www. triticeaecap.org) and only 0.7% of SNP data was missing. However, only 3569 SNP markers (i.e., 1.37 marker cM −1 ), which were in common with those included in the tetraploid wheat consensus map of ~2600 cM (Maccaferri et al., 2015), were used in the GWAS (Supplemental Table S2). Markers with minor allele frequency (MAF) of <5.0% were eliminated to reduce the chance of detecting false positives. In addition, markers that had >10% missing data were discarded from further analysis. The genetic position of the SNP markers was estimated based on the tetraploid wheat consensus map developed from Illumina iSelect 90K wheat array (Maccaferri et al., 2015). Redundant accessions sharing exactly the same SNP genotypes were subsequently excluded, which resulted in eliminating 64 accessions from the analysis.
Linkage disequilibrium for all pairwise comparisons between intrachromosomal SNPs was computed and the genome-wide LD decay was estimated using JMP Genomics 6.1 software (SAS Institute, 2012). The LD was computed as the squared correlation coefficient (R 2 ) for each of the marker pairs. The genome-wide LD decay was estimated by plotting LD estimate (R 2 ) from all 14 chromosomes against the corresponding pairwise genetic distances (cM). Smoothing spline Fit (lambda = 114,551.3) was applied to LD decay.

Association Analysis
The population structure (Q matrix) (Price et al., 2006) was assessed through principal component (PC) analysis. The familial relatedness was estimated using an identityby-state matrix (K matrix) (Zhao et al., 2007). Both K and Q matrices were generated using JMP Genomics 6.1. Four regression models were used to analyze marker-trait association using JMP Genomics 6.1. They included (i) naïve model that did not account for kinship and population structure, (ii) kinship, (iii) kinship plus population structure (first two PCs that collectively explained 21.8%), and (iv) kinship plus population structure (first three PCs that collectively explained 26.6%). The K and the Q matrices were included in the regression equation to ensure that only genetically significant associations were detected from the GWAS and were not spurious associations resulting from population structure or familial relatedness. The general formula of the mixed linear model used for the GWAS follows the regression equation y = Xb + Qu + Iµ + e, where y is a vector of phenotypic values, X is a vector of SNP marker genotypes, b is a vector of fixed effects as a result of the genotype, Q is matrix of principle component vectors estimating population structure, u = vector of fixed effects resulting from population structure, I is an identity matrix, µ is a vector of random effects that estimates the probability of coancestry between genotypes, and e is a vector of residuals. The variances of µ and e effects are based on these assumptions; Var(u) = 2KV g and Var(e) = V R , where K is the kinship matrix deduced from genotypes based on the proportion of shared allele values, V g is the genetic variance, and V R is the residual variance (Yu et al., 2006;Zegeye et al., 2014). Each SNP marker was then fitted into the regression equation to generate a P-value. Marker-trait associations were considered significant at a P ≤ 0.001. For each regression model, the SNP markers were ranked from smallest to largest P-values. The best model for each leaf rust race-field location was chosen based on the mean squared difference (MSD) between observed and expected P-values  because of the uniform distribution of random marker P-values (Yu et al., 2006). The MSD was calculated using the following formula: where n is the number of markers, i is the rank number that is from 1 to n, and pi is the probability of the ith-ranked P-value. Significant markers associated with response to leaf rust were selected only from the model with the lowest MSD value. The P-values of the selected model were later adjusted by calculating the corresponding positive false discovery rate (pFDR) (Benjamini and Yekutieli, 2001) using JMP genomics 6.1. Marker-trait associations were finally considered significant at a pFDR ≤0.1. Furthermore, we performed a stepwise regression using JMP genomics 6.1 on the significant SNPs of each trait. This allowed determining the minimum number of SNPs independently associated with leaf rust response (Gurung et al., 2014). The selected SNPs from the stepwise regression explain similar phenotypic variation as that described by all the significant SNPs considered together for each trait .

Phenotypic Data Analysis
A total of 496 durum wheat accessions were evaluated for response to different Pt races at seedling stage in the greenhouse and at adult plant stage in the field. Most of the accessions were phenotyped in the experiments conducted in the United States. In Mexico, 364 were evaluated with race BBBQJ at the seedline stage in the greenhouse, while 371 and 383 accessions were screened at adult-plant stage in MX_EB_F and MX_Ob_F, respectively (Table 2; Supplemental Table S1). Phenotypic data for seedling stage were homogeneous based on Levene's test (Levene, 1960) on ITs for each Pt race except for BBBQD. Therefore, phenotypic data of replicates and experiments were pooled for each race and arithmetic means were calculated and used for AM. The North Dakota and CDL experiments for race BBBQD (BBBQD (ND) and BBBQD (CDL), respectively) were analyzed separately in the GWAS because of the nonhomogeneity of variance (Supplemental Table  S1). The discrepancies between experiments with race BBBQD could be due to different experimental conditions at the CDL and North Dakota sites and differing interpretation of ITs.
A high percentage of resistance was observed among the 496 accessions when evaluated at adult plant stage in both Minnesota trials that were inoculated with bread wheat specific races. At MN_Cr_F, 85.90% (371 accessions) of the durum germplasm was resistant, while 89.60% (406 accessions) of the collection screened at MN_StP_F was classified as resistant. On the contrary, there was much lower proportion of resistant accessions at both Mexico locations, where races virulent to durum wheat were used for inoculation and present as natural inoculum in the field. The percentage of resistant germplasm was estimated at 53.79% (206 accessions) and 39.89% (148 accessions) in MX_Ob_F and MX_EB_F, respectively (Table 2).

Population Structure, Kinship Analysis, and Regression Model Selection for Association Mapping
Population structure was inferred using principal component analysis (PCA). The PCA shows that two, three, and 10 PCs explain a cumulative 21.2, 26.1, and 41.5% of the genotype variation, respectively. The first three PCs clustered the collection into three subpopulations (Fig. 2). Subpopulation 1 (in red color) contains the largest number of accessions (376), which are mainly from Europe (121 accessions), Asia (123 accessions), North America (45 accessions), and South America (27 accessions). Subpopulation 2 (blue color) was the smallest with only 27 accessions mainly from Africa (7 accessions), Asia (7 accessions), and Europe (8 accessions). Subpopulation 3 (green color) includes 29 accessions mostly from Africa (21 accessions, 18 from Ethiopia and three from Eritrea). Thus, clustering of these accessions was not based on the geographic location except for Subpopulation 3. In Subpopulation 1, 59% were cultivars and breeding lines, while in Subpopulations 2 and 3, 63.0 and 76.0% were landraces, respectively.
The familial relatedness was estimated using an identity-by-state matrix (K matrix). Kinship between accessions was calculated and a heat map of the markerbased K matrix is illustrated in Supplemental Fig. S1. The durum collection had intermediate familial relationships as some hotspots with related lines were observed on the heat map. Accounting for the population structure and familial relationship between individuals in the AM analysis reduces the number of false-positive associated markers. Based on MSD values of the four regression models tested, no single model fits best for all traits. For instance, the kinship model was used for MX_ Ob_F, MX_EB_F, MN_Cr_F, and MN_StP_F, while 2PCs+kinship was the model used for BBBQD (ND) and BBBQD (CDL). Model 3PCs+Kinship fits best for MBDSD, MCDSS, BBBDB, and race mixture (Table 4).

Marker-Trait Associations for Leaf Rust Response
Association mapping based on the infection type to different Pt isolates at seedling stage in the greenhouse and at adult stage in the field revealed 88 significant SNP markers. The leaf rust response-SNP associations were distributed across all chromosomes except for 1B and 7B. Twenty markers were associated with leaf rust response at seedling stage, while 68 markers were associated with leaf rust response at adult plant stage. Of these 88 markers, 33 were located on chromosomes 2A and 2B. Based on the pairwise LD, the 88 SNPs represented 37 different loci (Supplemental Table S3, S4). Evaluation of the 88 SNPs in all resistant accessions allowed verification of the association of the SNP markers with disease response and identification of alleles associated with the resistance to Pt (Supplemental Table S4).
At seedling stage, 14 SNPs on chromosomes 2A, 2B, 3B, 4A, and 5A were associated with response to BBBQD (ND). These markers represent six different loci. Three of the 14 markers fit into a stepwise regression, accounting for 11.9% of the phenotypic variation. The experiment at the CDL with the same race (BBBQD) revealed a total of nine associated SNPs located on 2A, 2B, 5A, and 6A, representing six loci. Four of the nine markers fit into a stepwise regression and accounted for 11.9% of the phenotypic variation (Table 5; Fig. 3, 4). The experiments with BBBQD at the North Dakota and CDL sites shared five associated SNPs corresponding to three loci located on 2A, 2B, and 5A (Table 6).
Two markers, IWA7547 and IWA757, on chromosomes 2A and 5B, respectively, were associated with response to Pt races MCDSS and MBDSD and to the Pt race mixture at seedling stage. The marker IWA7547 was the only one generated from the stepwise regression model in races MCDSS, MBDSD, and Pt race mixture.
Only one significant marker, IWA757 (located on chromosome 5B), was found associated with response to race BBBDB. This marker was considered significantly associated despite a pFDR of 0.2 (pFDR  0.1 cutoff used for significant associated markers across all experiments) because it was also associated with response at seedling stage to races MBDSD, MCDSS, BBBQD (ND), BBBQD (CDL), and the race mixture.
Three markers on chromosomes 4B and 6A were associated with response to BBBQJ representing two loci on chromosomes 4B and 6A. Two markers fit into a stepwise regression and accounted for 7.8% of the phenotypic variation (Table 5; Fig. 3, 4).
At the adult plant stage in the field experiments, 12 markers on chromosomes 1A, 2A, 3A, 4A, and 6A were significantly associated with response to Pt races in MN_StP_F, representing six loci. Furthermore, six of the 12 markers fit into a stepwise regression and accounted for 52.7% of the phenotypic variation.
Fourteen markers on chromosomes 2A, 3A, 3B, 4A, and 5A were significantly associated with response to Pt races in MN_Cr_F. These SNPs represent six loci. Five of the 14 significant markers fit into a stepwise regression and accounted for 27.1% of the phenotypic variation. Six markers, representing two loci, on chromosomes 2A (Lr.locus-2A4) and 3A (Lr.locus-3A2), respectively, were associated with response to Pt races at both locations in Minnesota. Similarly, 25 markers were significantly associated with response to Pt races in the MX_EB_F trial and were located on chromosomes 1A, 2A, 3A, 3B, 5A, 6B, and 7A. Twelve out of these 25 SNPs represent different loci. Eight of the 25 significant markers fit into a regression model and together accounted for 48.0% of the phenotypic variation.
Two markers, IWA1387 and IWA7151, on chromosome 1A and a marker, IWA757, on chromosome 6B were detected in response to Pt races across the Mexican locations. Moreover, markers IWA3512 on 3A and IWA4618 on 4B were associated with resistance or susceptibility to the races in field trials in Minnesota and Mexico (Table  6). A total of 19 SNPs, representing 11 loci, were significantly associated with leaf rust response in two or more tests (Table 6). Marker IWA7547 located on chromosome 2AL at 137.1-cM position based on the consensus map   # pFDR, positive false discovery rate. † † Markers in bold were maintained after stepwise regression.   (Maccaferri et al., 2015) was associated with response to all the Pt races tested at seedling stage except BBBQJ.

Discussion
In this study, we identified durum wheat accessions carrying race-specific and broad-spectrum resistance to leaf rust. The GWAS showed 88 significant SNPs associated with leaf rust response. These markers represent 37 loci distributed across all chromosomes except 1B and 7B. Interestingly, chromosomes 1B and 7B are known to carry adult-plant leaf rust resistance genes in bread wheat including Lr46 (1BL) and Lr68 (7BL) (Singh et al., 1998;Herrera-Foessel et al., 2012).
Some of the identified loci in this study appeared to be located in genomic regions of previously characterized Lr genes and QTL in earlier AM studies and biparental populations in durum and bread wheat. However, several other loci are previously unknown to be associated with leaf rust response. These findings will be valuable to both durum and bread wheat breeding programs.
Difference in frequency of resistance to races collected from common and durum wheat were evident in the germplasm collection. The majority of accessions were resistant to races collected from bread wheat, while most were susceptible to races from durum wheat, which agrees with previous studies (Singh, 1991;Huerta-Espino and Roelfs, 1992;Ordoñez and Kolmer, 2007a). Unexpectedly, low levels of resistance to the isolates collected from durum wheat fields in North Dakota in 2012 (race MBDSD and MCDSS) were observed. Isolates of these races are also commonly found on bread wheat in North Dakota and other regions in the United States (Kolmer and Hughes, 2014). Screening of 60 cultivars and breeding lines from the North Dakota durum wheat breeding program with the MCDSS revealed that 51.7% of accessions were susceptible at seedling stage (data not shown). This suggests that Pt races currently present in North Dakota with a combined virulence on both common and durum wheat may become a threat to durum wheat production in the United States, especially since race MBDSD has been common in recent US surveys (Kolmer and Hughes, 2014).
Of all evaluated genotypes, 13 accessions were resistant across all experiments. The SNP genotype showed that none of these accessions were duplicates. However, Fig. 4. (continued on next three pages) Chromosome locations of significant single-nucleotide polymorphism (SNP)-leaf rust associations in this study relative to mapped known Lr genes. Marker order and locations (left side of chromosome bars) are as reported by Maccaferri et al. (2015). The association mapping results are reported for Puccinia triticina races used at seedling stage in the greenhouse and at artificially inoculated field trials performed in the United States and Mexico. Markers in red are the significant SNP-leaf rust association observed in this study. † Simple-sequence repeat (SSR) marker associated with leaf rust response in durum wheat in Maccaferri et al. study (2010b). The genetic locations of Lr genes are indicated with arrows. The Lr genes assigned to chromosomes but not yet mapped are in a box at the bottom of each chromosome. Centromere position is indicated with a black circle on the chromosome bar. Not all SSR in the tetraploid consensus map (Maccaferri et al. 2015) are presented in this figure, a more saturated genetic map is presented in Supplemental Table S3. genetic mapping and allelism tests are required to verify whether these accessions have different or the same resistance genes. Accession data in the NSGC records indicates that seven of these resistant accessions were originally collected from Mediterranean countries and the Fertile Crescent, while four accessions were from Ethiopia. This is not surprising, as these countries belongs to the center of origin and diversity of durum wheat (Demissie and Habtemariam, 1991;Tesemma and Belay, 1991;Salamini et al., 2002). The coevolution of the leaf rust pathogen and wheat in these areas is believed to exert selection pressure, resulting in fitness advantages and maintenance of disease resistance diversity in wheat (Newton et al., 2010).
Further screening at seedling stage with races collected from durum wheat from Argentina (Arg. 9.3: BBBQD), France (FRA4.3: BBBQD), Ethiopia (E11D2-1: MCDSB, E114-1: BBBQD, and E125-1: EEEEE, avirulent on Thatcher), Arizona (09AZ103A:BBBQB), Mexico (LCJ/ BN and BBB/BN_Lr61vir), and Italy (PSB7: FGBQ) and at adult stage in Ethiopian and Moroccan fields in 2014 showed that eight out of the 13 previously mentioned accessions were resistant across all experiments (unpublished data). The remaining five accessions were resistant to most races. PI 278379 was susceptible only to Ethiopian isolates (E125-1 and E114-1), PI 519832 was susceptible to a French isolate (FRA 4.3) and an Ethiopian isolate (E114-1), PI 387263 was susceptible to only one Ethiopian isolate (E125-1), PI 324928 was susceptible to only Mexican races (LCJ/BN and BBB/BN_Lr61vir), and CItr 14623 was susceptible to only an Ethiopian isolate (E125-1) (Supplemental Table S5). Thus, these accessions In this study, the identified loci were not associated with response across all experiments and Pt races tested in North America. For instance, 26 of the discovered loci in this project were race or assay specific, while other loci appear to be more stable. This suggests that breeders may need to pyramid multiple loci to achieve broad-spectrum and stable leaf rust resistance.

Comparison to Previously Mapped Lr Genes
Chromosomes 3A, 5A, and 6A have not been previously shown to carry known Lr genes in either bread or durum wheat (McIntosh et al., 2014). However, some Lr genes from the ancestors of wheat were reported on chromosomes 3A and 6A. This includes Lr63 (from T. monococcum L.) and Lr66 (Aegilops speltoides Tausch), both on 3A, and Lr64 [from T. dicoccoides (Korn. ex Asch. & Graebn.) Schweinf.] on 6A (Kolmer et al., 2010;Marais et al., 2010, Kolmer, 2008. The Lr genes introgressed from alien species were not considered in this comparison because of the low chance of their presence in the current can be used to diversify leaf rust resistance in different breeding programs globally. Out of the 37 loci associated with leaf rust resistance or susceptibility in durum wheat germplasm in this study, 11 were associated with response to at least two Pt races and experiments. Loci associated with both seedling and adult-plant response were not observed at the cutoff used, P ≤ 0.001 and pFDR ≤ 0.1. This probably is due to the difference in the races used for screening of the durum accessions at seedling stage in the greenhouse and at adult stage in the field in addition to the presence of natural inoculum in field trials. However, increasing the pFDR to 0.22 to 0.29 showed SNPs (136.9-137.1 cM) on Lr.Locus-2A.4 associated with response to races in MN_Cr_F, MN_StP_F, and MX_EB_F. The same locus was also significantly associated with response to US races (BBBQD, MBDSD, MCDSS, and race mixture) tested at seedling stage. Similarly, marker IWA4187 (176.2 cM) on 7A was associated (at pFDR = 0.3) with response to MBDSD and races in MN_StP_F. durum wheat collection used. In addition, distantly located known Lr genes from the associated loci in this study were regarded as different.
Two loci on chromosome 2BL, Lr.locus-2B2 (95.2 cM) and Lr.locus-2B3 (100.9 cM), associated with response to leaf rust in MX_Ob_F were in close proximity to the adult plant resistance gene Lr35. However, this gene has been introgressed from A. speltoides (Seyfarth et al., 1999). On chromosome 4AS, the gene Lr30 was located close to the identified Lr.locus-4A1 (59.9 cM); however, the mapping information of Lr30 was not sufficient to make inferences.
On chromosome 4BL, Lr.locus-4B1 (70.9 cM) was mapped in the vicinity of Lr12, Lr31, and Lr49. The Lr.locus-4B1 was associated with response to race BBBQJ tested as seedling stage, which distinguishes them from the adult-plant resistance genes Lr12 and Lr49 (Dyck, 1991;Bansal et al., 2008). Likewise, Lr.locus-5B1 (2.7 cM) on 5BS was in close proximity to Lr52, which is closely linked to gwm443 and gwm234 and is also thought to confer leaf rust resistance in the Australian durum cultivar Wollaroi (Singh et al., 2010) and in accessions of the Watkins wheat collection (Dyck and Jedel, 1989). In addition, the locus Lr.locus-6A2 (88.2 cM) on chromosome 6A was located within the genomic region close to Lr64 (Kolmer et al., 2010;Marais et al., 2010, Kolmer, 2008. On chromosome 7AL, Lr.locus-7A3 (206.4 cM) appeared near the genomic region where the gene Lr20 (Watson and Luig 1966;Sears and Briggle 1969;Neu et al., 2002) was mapped.
Associated SNPs within the genomic region of Lr14a (7BL), which has been widely used, especially in the CIMMYT breeding programs, were not detected in the present study. This most probably is due to the rare frequency of Lr14a in the current germplasm collection as estimated based on the screening with the Lr14a diagnostic markers (gwm344 and gwm146) (data not shown). This agrees with the statement that GWAS has limited power to detect alleles that occur at low frequency in the germplasm (Myles et al., 2009;Brachi et al., 2011).
In summary, from this comparison, four of the currently identified loci are probably Lr31, Lr52, Lr64, and Lr20. However, allelism tests are necessary to determine Fig. 4. Continued. the relationship between these loci and the above-mentioned Lr genes.

Comparison to Previous Durum Wheat-Leaf Rust Association Mapping Study
The recently published tetraploid wheat consensus map (Maccaferri et al., 2015) containing different types of markers provided an opportunity to compare our GWAS results based on SNP markers with the only available AM of leaf rust response in durum wheat that was based on SSR markers (Maccaferri et al., 2010b) (Fig. 4).

Chromosome 1AS
The loci Lr.locus-1A1 (27.9 cM) and Lr.locus-1A2 (31.9 cM) were located close to wmc24 (28.1 cM). The latter was strongly associated with response to a durumspecific race from Italy and a bread wheat isolate from Poland. Two more loci, Lr.locus-1A3 (44.3 cM) and Lr.locus-1A4 (44.6 cM), were mapped close to wmc469 (46.1 cM). Marker wmc469 was previously reported to be associated with response to European and Mexican isolates tested at both adult plant and seedling stages.

Chromosome 2AL
Lr.locus-2A3 (107.6-107.7 cM) was located close to SSR marker gwm1045 (108.9 cM), which was shown to be associated with response to durum-wheat-type races from Italy both at adult plant and seedling stages.

Chromosome 3B
Lr.locus-3B2 (61.6 cM) was in the same genomic location as the SSR marker gwm779 (61.4 cM), which was associated with response at seedling stage to isolates collected from Poland, Italy, and the United Kingdom. Two other loci, Lr.locus-3B3 (100.2 cM) and Lr.locus-3B4 (100.7-101.6 cM), on the same chromosome, were close to barc164 (100.7 cM). The latter was associated with response at seedling stage to the durum-specific race BBBGJ collected from Italy.

Chromosome 4AL
Lr.locus-4A2 (71.2 cM) was closely mapped to barc155 (77.4 cM), which was previously associated mainly with response at seedling stage to BBBGJ from Italy.

Chromosome 4BL
Lr.locus-4B2 (70.9 cM) was in the vicinity of gwm1084 (78.5 cM). This SSR was associated with reaction at seedling stage to a triticale isolate (Poland) and two isolates from bread wheat that were collected from Poland and the United Kingdom.

Chromosome 5AL
Lr.locus-5A3 (90.6 cM) was in close proximity to gwm1236 (96.1 cM) and barc197 (96.3 cM). The marker gwm1236 was previously associated with leaf rust response at seedling stage to bread-wheat-type isolates collected from Poland and the United Kingdom, while barc197 was associated with response at adult stage to races in Italian field trials.

Chromosome 6A
Lr.locus-6A1 (5.9 cM) was in close proximity to SSR marker gwm459 (3.0 cM), which was associated with response at seedling stage to some bread-wheat-type isolates that originated from Poland and the United Kingdom. One more locus, Lr.locus-6A2 (88.2 cM), was proximal to wmc553 (90.3 cM) and gwm570 (90.5 cM). The wmc553 was associated with response at seedling stage to two bread-wheat-type races from the United Kingdom, while gwm570 was related to response at seedling stage to a bread-wheat-type isolate from the United Kingdom and a triticale isolate from Poland.

Chromosome 7AL
Lr.locus-7A3 (206.4 cM) was proximal to gwm344 (205.7 cM) and cfa2257 (204.2 cM). The SSR marker gwm344 was associated with response at seedling stage to some tested bread wheat isolates from the United Kingdom, a triticale isolate from Poland, and at adult-plant stage in field trials in Italy. The marker cfa2257 was related to response at seedling stage to bread wheat isolates from Poland and the United Kingdom and a triticale isolate from Poland.
Based on this comparison between the results of the current GWAS and the previous one performed by Maccaferri et al. (2010b), 17 associations were found in common between the two studies. Of these, 12 loci were associated with adult-plant response, while five loci were associated with seedling-stage response to leaf rust in the current study.

Comparison to Recent Bread Wheat-Leaf Rust Association Mapping Studies
To investigate possible cross-relationship between associated loci with leaf rust response found in durum and bread wheat, a comparison with GWAS on bread wheat germplasm was performed. Based on the genomic position of the SNPs in the tetraploid consensus map, eight identified loci in this study were previously reported in two AM analyses in bread wheat (Kertho et al., 2015;Gao et al., 2016).
For instance, the GWAS by Kertho et al. (2015) on winter wheat germplasm reported IWA3160 (50.1 cM) on 1A to be associated with response to the bread-wheat-type race TBDJ. This locus was proximal to the significant loci Lr.locus-1A3 (44.3 cM) and Lr.locus-1A4 (44.6 cM) for bread-wheat-type races in MN_StP_F. Similarly, IWA7429 (107.7 cM) on 2A and IWA5526 (135.1 cM) on 7A that were found to be in association with response to the breadwheat-type race TDBG in winter wheat were mapped close to Lr.locus-2A3 (107.6-107.7 cM) and Lr.locus-7A2 (133.7 cM), respectively. The marker IWA6244 (66.1 cM) on 3B that was significant for bread-wheat-type race MCDL was closely mapped to Lr.locus-3B2 (61.6 cM). Surprisingly, Lr.locus-2A3, Lr.locus-7A2, and Lr.locus-3B2 in the current work were associated with durum-type races in field experiments (MX_EB_F and MX_Ob_F).
In addition, a comparison with the GWAS on spring wheat germplasm by Gao et al. (2016) was conducted. This revealed that IWB74350/IWB73424 (2.4 cM) on 3B, IWB25253/IWB57347 (72.1-73 cM) on 4A, and IWB43173/IWB45939 (5.9 cM) on 6A, which were significant for North American bread-wheat-type races in spring wheat, were in close proximity to the Lr.locus-3B1 (4.2 cM), Lr.locus-4A2 (71.2 cM), and Lr.locus-6A1 (5.9 cM), respectively. The latter loci were also associated with response to bread-wheat-type races in field trials (MN_ Cr_F and MN_StP_F) in the present study. The results of this comparison suggest that some similarities between durum and bread wheat leaf rust resistance do exist.
Even though, some bread-wheat-type races (BBBDB, race mixture, and races in MN_StP_F and MN_Cr_F) used by Gao et al. (2016) were the same as those in the current work, there were very few similarities between loci reported. In addition, the race BBBQD (isolate collected from California) was also used by Gao et al. (2016); however, no common associations were observed with the present study. This suggests that leaf rust resistance in durum and bread wheat are under the control of different genetic loci.

Comparison with Quantitative Trait Loci Meta-Analysis for Leaf Rust Resistance in Durum and Bread Wheat
The identified SNPs in this work were also compared with the recent QTL meta-analysis study that was performed using 20 biparental mapping populations and 33 different parental lines of durum and bread wheat (Soriano and Royo, 2015). The comparison showed that nine identified loci in this study were previously reported in the QTL meta-analysis. For instance, two loci on 1A, Lr.locus-1A1 (27.9 cM) and Lr.locus-1A2 (31.9 cM), were in the genomic region of wmc95 (22.3 cM), which is a marker linked to resistant QTL in the bread wheat cultivars Apache and Sujata (Azzimonti et al., 2014;Lan et al., 2015). Another locus on 2A, Lr.locus-2A3 (107.6 cM), was positioned close to gwm339 (99.0 cM), which is one of the flaking marker of a QTL identified in the durum cultivar Creso (Marone et al., 2009). Likewise, Lr.locus-2B1 (94.2 cM) and Lr.locus-2B2 (95.2 cM) on 2B were closely mapped to gwm55 (94.0 cM) that is linked to QTL in bread wheat cultivars W-7984, Kariega, Avocet, and Carberry (Faris et al., 1999;Prins et al., 2011;Singh et al., 2014). The comparison also revealed that Lr.locus-3B2 (61.6 cM) on 3B was located close to gwm566 (59.0 cM), which is in the genomic region of QTL in the bread wheat cultivars TA4152-6 and Francolin#1 (Chu et al., 2009;Lan et al., 2014). Three loci, Lr.locus-2B4, Lr.locus-4B1, Lr.locus-5A1, that were associated with response to leaf rust at seedling stage in this study, were in the genomic regions of QTL in Apache and Cresso on 2B; Avocet, Creso, and Forno on 4B; and Oberkulmer on 5A; respectively (Messmer et al., 2000;Schnurbusch et al., 2004;William et al., 2006;Marone et al., 2009;Azzimonti et al., 2014). However, these previously identified QTL were associated with adult plant resistance.
In summary, the comparison between the 37 loci identified in the current study, characterized Lr genes, loci identified in previous GWAS, and QTL meta-analysis in biparental mapping populations in both durum and bread wheat revealed that 22 of the generated loci in the present work were previously reported. Consequently, may not be known to be previously associated with leaf rust response in durum and bread wheat. Five of those  are of special importance, as they are associated with leaf rust resistance to two or more Pt races. The identification of these loci through GWAS is a significant step in characterization of genes that may be used to broaden the genetic basis of leaf rust resistance in durum wheat germplasm.

Conclusions
This study identified durum wheat accessions with both race-specific and broad-spectrum resistance. The majority of durum accessions were resistant to races collected from bread wheat, while most were susceptible to races collected from durum wheat. Thirteen accessions showed resistance to races tested at seedling stage and to a race mixture in field experiments in both the United States and Mexico. Of these, eight accessions (PI 209274, PI 192051, PI 244061, PI 223155, PI 534304, PI 193920, PI 342647, and PI 195693) were also resistant at seedling stage to additional nine isolates, collected from Mexico, Argentina, France, Ethiopia, and Italy, and at adult stage in field trials in Ethiopia and Morocco. These broadspectrum resistant accessions could be a good leaf rust resistance source to introgress into locally adapted germplasm in breeding programs globally.
The GWAS revealed 88 SNPs representing 37 loci associated with leaf rust response across all durum wheat chromosomes except 1B and 7B. The comparison of their genomic regions with the known Lr genes, previous AM studies, QTL mapping in biparental populations in durum and bread wheat revealed 14 previously uncharacterized loci, of which, five were associated with leaf rust response to two or more Pt races tested at seedling stage or race mixture in field trials. The marker IWA754 (137.1 cM) on Lr.locus-2A4 was associated with response to all Pt races tested at seedling stage with the exception of BBBQJ. To validate the loci revealed by the GWAS study, biparental populations have been developed from selected accessions showing wide-spectrum leaf rust resistance. Furthermore, this study will facilitate the development of tightly linked markers for marker assisted selection in breeding programs.

Supplemental Information Available
Supplemental Figure S1. Heat map displaying the relationship matrix among durum accessions.
The red diagonal represents perfect relationship of each accession with itself. The symmetric off-diagonal elements represent relationship measures (Identity by state) for pairs of accessions. The blocks of warmer colors on the diagonal show clusters of closely related accessions.
Supplemental Table S1. Response of durum wheat accessions to Puccinia triticina races.
Supplemental Table S2. Marker proprieties of 3,659 SNPs used in the association analysis.
Supplemental Table S3. Chromosome locations and linkage disequilibrium between significant SNP-leaf rust associations.
Supplemental Table S4. List of resistant accessions and alleles of associated SNPs with resistance to each trait.
Supplemental Table S5. List of resistant accessions across experiments with additional tests.