Wheat genetic loci conferring resistance to stripe rust in the face of genetically diverse races of the fungus Puccinia striiformis f. sp. tritici

Analysis of a wheat multi-founder population identified 14 yellow rust resistance QTL. For three of the four most significant QTL, haplotype analysis indicated resistance alleles were rare in European wheat. Stripe rust, or yellow rust (YR), is a major fungal disease of wheat (Triticum aestivum) caused by Puccinia striiformis Westend f. sp. tritici (Pst). Since 2011, the historically clonal European Pst races have been superseded by the rapid incursion of genetically diverse lineages, reducing the resistance of varieties previously showing durable resistance. Identification of sources of genetic resistance to such races is a high priority for wheat breeding. Here we use a wheat eight-founder multi-parent population genotyped with a 90,000 feature single nucleotide polymorphism array to genetically map YR resistance to such new Pst races. Genetic analysis of five field trials at three UK sites identified 14 quantitative trait loci (QTL) conferring resistance. Of these, four highly significant loci were consistently identified across all test environments, located on chromosomes 1A (QYr.niab-1A.1), 2A (QYr.niab-2A.1), 2B (QYr.niab-2B.1) and 2D (QYr.niab-2D.1), together explaining ~ 50% of the phenotypic variation. Analysis of these four QTL in two-way and three-way combinations showed combinations conferred greater resistance than single QTL, and genetic markers were developed that distinguished resistant and susceptible alleles. Haplotype analysis in a collection of wheat varieties found that the haplotypes associated with YR resistance at three of these four major loci were rare (≤ 7%) in European wheat, highlighting their potential utility for future targeted improvement of disease resistance. Notably, the physical interval for QTL QYr.niab-2B.1 contained five nucleotide-binding leucine-rich repeat candidate genes with integrated BED domains, of which two corresponded to the cloned resistance genes Yr7 and Yr5/YrSp.


Introduction
Yellow rust (YR), caused by the biotrophic fungus Puccinia striiformis Westend f. sp. tritici (Pst), is a widespread pathogen of wheat (Triticum aestivum L.), and a substantial threat to global wheat production. Since the 2000s, a subset of genetically diverse and divergent Pst lineages have been responsible for recurrent YR epidemics in numerous wheatproducing regions (Ali et al. 2014. The rapid adaptation and subsequent spread of these lineages into previously hostile environments has given rise to more aggressive pathotypes that are generally better adapted to higher temperatures (Milus et al. 2009). In the United Kingdom (UK) and North-Western Europe, the historically clonal Pst populations have been largely displaced by a genetically diverse group of lineages (Hubbard et al. 2015a, b;Hovmøller et al. 2016). First detected simultaneously in 2011 across several European countries, the 'Warrior' and 'Kranich' Pst races Communicated by Beat Keller.
1 3 likely originated in the near-Himalayan region and rapidly spread throughout the continent as a group of genetically distinct lineages (Hovmøller et al. 2016). The ability of Pst to migrate over long distances, to locally adapt to new environments, and for new lineages to rapidly displace established populations, has had a notable impact on the resistant levels of wheat varieties in the UK (Hubbard et al. 2015a, b;Bueno-Sancho et al. 2017) and beyond (Wellings 2011;Hovmøller et al. 2016), highlighting the importance of the generation of genetic information and resources to support forward YR resistance breeding strategies.
The most efficient way to control the effects of wheat fungal diseases is via approaches that combine agricultural and agronomic practices, disease monitoring and genetic improvement of the wheat varieties grown . Resistance breeding focuses on two classes of rust resistance (R) genes. The first, termed 'seedling resistance' or 'all-stage resistance', confers qualitative resistance, typically to one or a low number of Pst isolates. The second, broadly termed 'adult plant resistance' (APR), typically provides quantitative resistance against multiple pathogen races (although APR specificity to Pst races does exist, e.g. Yr12 and Yr13. Johnson 1992;McIntosh et al. 1995), and are often effective against multiple biotrophic pathogens (e.g. Yr18/Lr34/Sr57/Pm38, conferring resistance to yellow rust, leaf rust, stem rust and powdery mildew; Krattinger et al. 2009). Over 300 genomic regions conferring YR resistance in wheat have been reported (Rosewarne et al. 2013;Wang and Chen 2017). Of these, approximately 80 are permanently designated yellow rust resistance (Yr) genes, recently summarised by Jamil et al. (2020). To date, 19 designated R genes controlling all-stage resistance to different wheat fungal pathogens have been cloned, and all but two (Sr60, a tandem kinase, Chen et al. 2020. Yr15/YrG303/YrH52, a kinase-pseudokinase protein, Klymiuk et al. 2018Klymiuk et al. , 2020 encode nucleotide-binding leucine-rich repeat (NLR) proteins: Lr1 (Cloutier et al. 2007), Lr10 (Feuillet el al. 2003), Lr21 (Huang et al. 2003), Lr22a (Thind et al. 2017), Sr13 (Zhang et al. 2017), Sr21 (Chen et al. 2018 (Marchal et al. 2018), Yr10 ; although see also Yuan et al. 2018 who question the result) and YrAS2388R (Zhang et al. 2019). Similarly, wheat R genes to other fungal pathogens are also encoded by NLR genes, such as the powdery mildew resistance genes Pm2 (Sánchez-Martin et al. 2016), Pm3b (Yahiaoui et al. 2004, Pm8 (Hurni et al. 2013), Pm21 (He et al. 2018) and Pm41 (Li et al. 2020). Host NLR proteins predominantly act by recognising the effector molecules that pathogens produce to inhibit host defence responses . To detect the large range of potential infecting pathogens, plant NLR gene families have radiated and diversified, commonly by localised gene duplication, unequal crossing-over and variation within the leucine-rich repeat (LRR) domains which bind pathogen effectors (Sarris et al. 2016). Some NLRs have evolved to include additional 'integrated' domains that may be involved in receptor activation or downstream signalling (Sarris et al. 2016), most commonly kinase and DNA-binding domains (Andersen et al. 2020;Steuernagel et al. 2020). Examples include the BED zinc finger domain in Yr7 and two alleles at a closely located paralogous gene termed Yr5 and YrSp (Marchal et al. 2018). NLRs can also act by monitoring the status of host receptor proteins, termed 'indirect recognition'. For example, the Arabidopsis NLRs RPM1 and RPS2 monitor RPM1-interacting protein 4 (RIN4) for cleavage or phosphorylation by several bacterial effectors (Mackey et al. 2002;Axtell et al. 2003;Andersson et al. 2006). Recently a wheat YR susceptibility gene encoding a branched-chain amino acid aminotransferase, TaBCAT1, has been identified via analysis of genes upregulated after Pst infection (Corredor-Moreno et al. 2021). In contrast to allstage resistance genes, each of the three map-based cloned YR adult plant resistance genes encode a different class of protein: Yr18/Lr34 an ABC transporter (Krattinger et al. 2009) involved in the translocation of abscisic acid (Krattinger et al. 2019), Yr36 a kinase-START domain protein (Fu et al. 2009) and Yr46/Lr67 a hexose transporter (Moore et al. 2015). Wheat breeding strategies historically focused on the use of qualitative yellow rust resistance genes, sometimes in isolation. Yr17 for example was a popular source of such resistance in North-Western Europe. The effectiveness of the genetic resistance conferred by Yr17 stopped soon after its deployment as a single resistance gene over a large wheat acreage (Bayles et al. 2000). Where possible, resistance breeding strategies now favour the more durable approach of combining race specific and non-race specific resistance genes to provide broad-spectrum resistance (Singh et al. 2005;Chen et al. 2014). Such use of genetic resistance will continue to be aided by the genetic and molecular characterisation of Yr genes and quantitative trait loci (QTL), understanding of how these genetic loci interact with each other in inbred and F 1 genetic backgrounds, availability of molecular markers to track favourable alleles and allelic combinations in breeding programmes, and the development and assessment of resistance gene cassettes containing multiple R genes (e.g. Luo et al. 2021).
YR resistance QTL and Yr genes have predominantly been identified in biparental (e.g. Rosewarne et al. 2013) and association mapping panels (e.g. Maccaferri et al. 2015;Zegeye et al. 2014;Kertho et al. 2015). The crossing of just two parents in a biparental population inevitably limits the number of resistance genes and alleles that can be investigated (Mackay 2001). While association mapping panels overcome this limitation by exploiting 1 3 historical recombination events in the germplasm collection used, population structure in such collections can lead to false marker-trait associations, and the effects of rare alleles will not be detected (Cockram & Mackay 2018). Multi-founder populations (reviewed by Cockram & Mackay 2018;Scott et al. 2020) now provide complementary resources for the genetic investigation of wheat disease resistance. Such designs include multiparent advanced generation inter-cross (MAGIC) populations, derived by inter-crossing all founders over multiple generations before the generation of inbred lines. As each MAGIC progeny line likely contains alleles from all founders dispersed throughout its genetic background, MAGIC populations allow the effects of multiple alleles to be assessed within a single unified population. We previously developed the eight-founder 'NIAB Elite MAGIC' population, estimated to capture > 80% of the genetic variation observed in UK wheat based on single nucleotide polymorphism (SNP) analysis (Mackay et al. 2014). Here we analyse this MAGIC population for resistance to YR at five trials conducted over two seasons and three sites in the United Kingdom (UK). Eight major adult plant resistance QTL were resolved, identified across growth seasons in which the new genetically diverse 'Warrior' group of Pst races were endemic in the UK. The gene space of these genetic loci was explored, candidate genes identified and genetic markers tagging resistant haplotypes validated for the four most significant QTL. Additionally, six minor effect QTL were identified. Collectively, this work highlights the role of combinations of adult plant YR resistance genes in the genetic control of the aggressive Pst pathotypes that have replaced the previous asexual populations of the pathogen.

Germplasm, trial design and field trials
The 'NIAB Elite MAGIC' wheat population (Mackay et al. 2014) consists of eight founders (the winter varieties Alchemy, Brompton, Claire, Hereward, Rialto, Robigus, Soissons and the facultative variety Xi19) inter-crossed over three generations, and the outputs of the crossing then selfed over multiple generations to produce > 1000 recombinant inbred lines (RILs). The population was grown in five field trials over two seasons (YR assessment seasons 2015 and 2016) at three sites in the UK: (1) NIAB-Cambridge trial ground (latitude 52.235010, longitude 0.097871), Osgodby (latitude 53.410161, longitude − 0.386770) and Rothwell (latitude 53.482597, longitude − 0.259779). All trials followed an incomplete randomised block design generated with the DEW experimental design software, formerly www. expde signs. co. uk but superseded by the R package 'blocksdesign' (Edmondson 2020). The numbers of MAGIC RILs and control lines assessed at each site and year combination are listed in Table 1 and Supplementary Table 1. Seed for each season's trials were sourced from nursery plots grown the preceding season. Accordingly, RIL 9 and RIL 10 seed were used for the 2015 and 2016 season trials, respectively. Standard agronomy practices were used for commercial wheat production at each location, but lacking application of chemical protection against YR. All trials were Autumn sown in the calendar year preceding YR assessment. Each line was sown in two 1-m rows, with six rows per plot. At the NIAB trial ground in Cambridgeshire, the central two rows were sown with the spreader wheat variety Vuka, known to be highly susceptible to all known UK races of YR, with RILs on either side. At the Osgodby and Rothwell sites in Lincolnshire, the central two rows were left empty. Instead, Vuka was present as a whole plot every three traverses. All trials apart from OSG16 were inoculated with a mixture of Pst races ('Solstice race' isolate 08/21 virulent on Yr 1, 2,3,4,6,9,17,25,32 and 'Warrior race' isolate 11/08 virulent on Yr 1, 2,3,4,6,7,9,17,25,32,Sp). Trial OSG16 was un-inoculated, and therefore exposed only to natural YR infection.

Phenotypic data
YR infection severity of the leaves of adult plants was assessed as the percentage of total leaf tissue with sporulating uredinia, estimated using the modified Cobb's scale (Peterson et al. 1948)

Trials analysis
A stepwise model selection approach was used to estimate the Best Linear Unbiased Estimators (BLUEs) for MAGIC RILs, integrating spatial and non-spatial mixed methods based on Restricted Maximum Likelihood (REML) implemented in Genstat, 18th edition (VSN International 2015). Three models were considered: Model 1-Blocking (genetic effects are estimated based only on the inter-and intra-block variation recovered from the model); Model 2-Spatial (only considers global and/or local field trends); Model 3-Spatial + blocking (combination of the above models). Initially, each model was optimised by including field trends running in either row or column direction or scoring order (the route used to score the trial), followed by between-model comparison to select the one that best fits the data. Further details of the three models used are given in Supplementary Text 1. The Akaike Information Criterion (AIC) was used as a measure for model selection (Akaike 1974). Each model was optimised using AIC as a measure of model fit improvement, and the model with the lowest AIC value selected. Natural log transformation of the disease severity data was performed in cases where residuals were observed not to be normally distributed (trials NIAB16 and OSG16), to improve the normality of the residuals. These were assessed visually with histograms illustrating the distribution of the residuals and Q-Q normality plots. The Shapiro-Wilk normality test was used to further support these observations. Phenotypic correlations were estimated among predicted means using Pearson's correlation coefficients and paired Wilcoxon signed-rank test using the Hmisc package (Harrell 2019). Broad sense heritability (h 2 ) was used as a measure of total phenotypic variation attributable to the genotypic effect. The VHERITABILITY function in Genstat (18th edition, VSN International 2015) was used to calculate h 2 for each trial and is based on the definition of heritability given by Cullis et al. (2006) and Piepho et al. (2007).

Genetic analysis
Genetic analysis was carried out using the 7369 SNPs mapped to unique positions on the MAGIC genetic map (Gardner et al. 2016). For all trials, two broad methods were used for genetic analysis, as previously described by Corsi et al. (2020).
(1) Single marker analysis (SMA): regression against allelic state at single markers using R/lme4 (Bates et al. 2015) in R (R Core Team 2020) using the following mixed model: where Y is the YR resistance value, μ is the adjusted YR score for MAGIC RIL x, G m is the fixed SNP marker effect, β is the population structure consisting of 'funnels' and 'plants within funnels' effects (from Mackay et al. 2014), and e as the residual error term. The model has one degree of freedom, since regression is carried out on binary allelic state. Multiple-test correction was carried out using R/ qvalue (Storey 2015), with a threshold of q < 0.05.
(2) Haplotype-based analysis, for which founder haplotype probabilities were calculated with the mpprob function in R/mpMap (Huang and George 2011) with a threshold of 0.5. Three types of haplotype-based analyses were conducted. Identity by descent (IBD): regression against haplotype probability estimates using R/qtl and the following mixed model: where Y, μ, β and e are as for SMA above, and Gp is the fixed term for founder probabilities. Here, the statistical model has up to seven degrees of freedom. A QTL significance threshold q < 0.05 was used. Interval mapping (IM): conducted in R/mpMap using the haplotype probability estimates. Composite interval mapping (CIM): conducted in R/ mpMap with 10 covariates using the haplotype probability estimates. Within mpMap, an automated forward selection process based on AIC values was used to select the best ten marker covariates for each MAGIC line. Significant QTLs were then selected in two stages. First, mpMap scans the 100 markers surrounding a particular marker location and selects QTLs based on a threshold of − log10(p) > 3. The number of significant QTL were then reduced by fitting a model with p < 0.05 and with phenotypic variance explained (PVE) > 0.5%. Significance threshold was estimated based on the simulation of the null distribution. Additionally, R/ mpMap outputs founder contributions for each significant QTL, computed using a regression approach at each marker location. To summarise the results of the four analyses, the p or q values from each mapping method for all adjusted and log transformed YR scores, and for all environments, were compiled into a single table and a 'consensus' peak marker for each QTL identified following the methods described in Supplementary Text 2. 'Major' and 'minor' QTL were defined as those explaining either ≥ 5% or < 5% of the phenotypic variance, respectively.

Bioinformatic analysis
The locations of the genetically mapped MAGIC SNPs on the physical map were determined using the SNP flanking DNA sequences ) as queries for BLASTn (Altschul et al. 1990) interrogation of the wheat reference genome assembly (RefSeq v1.0, IWGSC 2018). Where BLASTn hits of equal match were identified on more than one chromosome, genetic map position (Gardner et al. 2016) was used to assign BLASTn hits to chromosomes. Cloned wheat rust resistance R genes were used as queries for BLASTn searches of the wheat reference genome assembly (RefSeq v1.0) and associated gene annotation (RefSeq v1.1) and significant hits within QTL physical intervals listed. Additionally, analysis of the YR QTL physical intervals in the context of disease resistance gene density, candidate genes were first identified from the IWSGC RefSeq v1.0 assembly (IWSGC, 2018) high-and low-confidence gene functional annotations from v1.0 mapped on to gene model annotation RefSeq v1.1 (Alaux et al. 2018) where transcript.1 contained any of the following search terms: 'NB-LRR'; 'NBS-LRR'; 'NB-ARC'; 'TIR-NBS'; 'LRR family protein'; 'Leucine-rich repeat domain'; 'Plant disease resistance response'; 'Nucleotide-binding site leucinerich repeat'; 'disease resistance protein (TIR class)'. The percentile ranks of all gene counts and resistant gene counts for each QTL were calculated using windows of the same lengths (9.85-133.91 Mbp) sampling the genome (all 21 chromosomes plus the unassigned chromosome) every 100 bp. Plots were generated using the Circos visualisation tool (Krzywinski et al. 2009). QTL resistant gene enrichment p values were calculated using the binomial cumulative probability function where the chance of success is 0.0365 (genome-wide, resistant genes / total genes), the number of tests is the number of genes within the QTL and the number of successes is the number of resistant genes within the QTL.

Analysis of YR infection in the MAGIC population
To assess YR resistance at the adult plant stage (Fig. 1), infection was assessed in the MAGIC population grown across five trials in the UK, with phenotyping conducted at two time-points per site. Three of the trials were conducted in 2015 (sites NIAB15, OSG15 and ROTH15) and two in 2016 (NIAB16 and OSG16). In 2016, YR scores were skewed towards the resistant end of the scoring scale, and a normal distribution was not observed. Normality for the 2016 season trials was improved using the natural log transformation, and so was used in all subsequent analyses for that year. Normality tests on the 2015 season trials found the residuals to all be normally distributed, and so no further transformation was required. To account for variation within the field, three linear model approaches were applied to the 2015 raw data and the 2016 log transformed data, and the best models selected via Akaike Information Coefficient (AIC) (Supplementary  Table 3). YR infection scores in the MAGIC founders spanned from 0 to 100%, with similar trends in overall founder ranking observed between sites and years (Fig. 2a). Soissons (0-3% infection) and Robigus (70-100% infection) were the most resistant and susceptible parents, respectively, regardless of year, site and scoring time-point. While the remaining parents varied somewhat in ranking dependent on year and location, two main founder groups were observed: relatively resistant (percentage infection range: 3-33; mean = 10.4: Alchemy, Hereward and Xi19) or relatively susceptible For the MAGIC RILs, YR scores were similarly distributed at the different sites within each year, but differed considerably between years (Fig. 2b). Nevertheless, correlation of RIL YR scores between test environments was high, ranging from 0.80 to 1.00 (p < 0.001). Transgressive segregation in the RILs was observed in 2015 and 2016, both above and below that observed for the most susceptible and resistant founders, respectively (Fig. 2b). In 2015, as YR progressed through the season and susceptibility increased among the MAGIC population, a subset of 31 RILs remained highly resistant (%inf < 1) at all three sites. Intermediate YR scores (10-80%inf) were normally distributed, with this trend more evident at OSG15 and ROTH15, where Pst developed more gradually compared to at NIAB15. In 2016, percent YR infection was more skewed towards resistance, compared to the previous year. High levels of resistance (%inf < 1) were maintained by 32 RILs, of which eight were the same as in 2015. Overall, 121 and 171 MAGIC lines exhibited a resistant response (%inf < 10) throughout the 2015 and 2016 scoring seasons, respectively, at all three locations.
Pathotype data for 30 Pst isolates sampled from the trial sites and surrounding regional areas in the 2015 and 2016 seasons showed the predominant races to be Warrior 3/Old European and Warrior 4 in 2015 and Warrior 4 and Warrior 1 in 2016 (Supplementary Table 5). To support the assumption that these quantitative sources of resistance represented APR mechanisms, rather than all-stage resistance, the eight MAGIC founders were phenotyped for YR infection at the seedling stage to four Pst isolates collected in the geographical regions in which the trials were undertaken in 2015 and  Table 6). Given that Xi19 was not found to confer resistance alleles at any of the identified loci in the absence of resistance at the same locus from other founders, and that resistance sources were quantitative in nature, our working assumption is that the YR QTL we identified represent APR mechanisms.    Table 8).

Pairwise QTL interactions
Next, using the wheat reference genome assembly and gene model annotation, we investigated the genomic locations of the 14 YR QTL in respect to overall gene density, and to the density of genes annotation as being involved in disease resistance, including those belonging to the NBS-LRR gene family (Supplementary Fig. 1; Supplementary  Table 9). Seven of the 14 QTL were located in genomic regions that contain very high numbers of genes (> 95th percentile, genome-wide). Additionally, four QTL were located in regions with significantly higher number of 'resistance genes', based on protein families associated with resistance (> 95th percentile, genome-wide), although these were dominated by NBS-LRR genes typical of all-stage resistance. Furthermore, six QTL had a significantly higher number of 'resistance genes' than expected from the QTL gene count, compared to as expected genome-wide, calculated using the binomial cumulative probability function (p > 0.05). Details of potential candidate genes within the QTL physical intervals, based on gene model functional annotations and predicted protein domains, as well as sequence searches using the 22 wheat rust R genes cloned to date (Supplementary Table 10), are detailed in Supplementary Table 9. In addition to NLR genes characteristic of all-stage resistance genes, the colinear QYr.niab-2A.2 and QYr.niab-2D.1 were notable for the prominence of receptor-like kinase candidate genes, which encode proteins with protein-kinase and LRR domains. Similarly, the broadly colinear QYr.niab-6A.2 and QYr.niab-6B.1 QTL both had high numbers of F-box/LRR domain encoding genes in additional to multiple copies of canonical NLRs, with QYr.niab-6B.1 also containing a PK-NLR (TraesCS6B02G099900; although its 6A homoeologue TraesCS6A02G093400LC is a pseudogene and so not annotated as a PK-NLR). QYr.niab-2B.1 was found to span the wheat reference genome gene models reported to be homologous to the cloned wheat YR resistance genes Yr7 and Yr5/ YrSP, represented in the reference genome assembly of cv. Chinese Spring by gene models TraesCS2B02G488000 and TraesCS2B02G488600/TraesCS2B02G488700, respectively. Both Yr7 and Yr5/SP encode NLRs with an integrated BED domain, with three additional BED-NLR candidate genes located very close-by (TraesCS2B02G488400, TraesCS2B-02G734100LC and TraesCS2B02G48900). With the exception of Yr7 and Yr5/YrSP, none of the remaining 13 cloned wheat R genes, or genes with high sequence similarity to these R genes, were located in the QTL intervals identified here.

Development of KASP markers tagging resistant haplotypes
For the four most robust YR QTL investigated above, we used SNPs from the 90 k SNP array able to discriminate resistant and susceptible haplotypes or SNPs in the MAGIC founders and converted these to the KASP genotyping platform (Supplementary Table 2; Supplementary Fig. 2). The utility of these markers was further assessed in terms of ability to discriminate haplotypes within haploblocks identified at the QTL peaks of each of the four loci in a panel of 403 predominantly winter wheat varieties previously genotyped with the 90 k SNP array ( Fig. 6; Supplementary Table 11). For MAGIC QTL QYr.niab-1A.1, four haplotypes were defined in the variety panel within a large haploblock of 55 SNPs (174.98-187.84 cM) that spanned the majority of the MAGIC QYr.niab-1A.1 genetic interval (177.50-194.67 cM) (Fig. 6a). Haplotype 1A-hap1 was identified in just 10 varieties (2.5%), including Hereward which was identified as the only MAGIC founder to carry a resistant allele at the QTL. The most common haplotypes, 1A-hap2 (53%, including MAGIC founder Alchemy) and 1A-hap3 (40%, including MAGIC founders Brompton, Claire, Rialto, Robigus and Xi19), differed by just one of the 55 SNPs within the interval, and were both designated susceptible haplotypes, based on the allele effects at this locus identified in the MAGIC population. MAGIC founder Soissons, which carried a susceptible allele at the 1A locus, was one of the 15 varieties suspected to carry a genetic recombination within the haploblock, combining a 'susceptible' 1A-hap2/1A-hap3 region with a 'resistant' 1A-hap1 region. The co-dominant KASP markers developed for this locus were located on chromosome 1A at 178.013 cM (BobWhite_c44164_402) and 185.804 cM (tplb0021i12_383), and within this interval discriminated the haplotype associated with resistance (1A-hap1) from the two haplotypes associated with susceptibility (1A-hap2/1A-hap3) (Fig. 6a). The fourth haplotype was only identified in a single variety (Isidor, 1A-hap4). As this haplotype was not captured in the MAGIC founders, it was not possible to predict which class of YR allele it carried.
For QYr.niab-2A.1, co-dominant KASP markers were developed for two SNPs (BS00062679_51 and BS00022641_51) which together were able to discriminate the resistant Rialto/Xi19 allele from the more susceptible alleles from the remaining six founders. Analysis of linkage disequilibrium in the variety panel at the peak of the MAGIC QTL identified seven haplotypes, of which four were relatively common (frequency > 10%). These common haplotypes were found in one or more of MAGIC founders carrying susceptible alleles at the locus and were therefore assigned as likely representing susceptible alleles (Fig. 6b). Haplotype 2A-hap1 was present at low frequency in the variety panel (5%), including the two MAGIC founders carrying the resistance allele at the locus: Rialto and Xi19. The two KASP markers (BS00062679_51 and BS00022641_51) were able to discriminate Rialto/Xi19 alleles from the remaining six susceptible MAGIC founders. While analysis in the variety panel found these SNPs to be located in the adjacent haploblock, this was still within the MAGIC QTL interval, and together the two SNPs differentiated the resistant class with all but 12 accessions in the variety panel (Fig. 6).

Pst infection at the trial sites
The spread of genetically diverse exotic Pst races to European environments and their displacement of the previously clonal Pst races has resulted in sudden shifts in wheat resistance ratings. The first of these Pst races to be detected in Europe is termed the 'Warrior' race (www. wheat rust. org) and was characterised by several notable traits, including relatively large reductions in resistance in varieties that previously carried effective long-term adult plant resistance ) and high production of sexual stage spores (teliospores) (Rodriguez-Algaba et al. 2014) indicative of evolution from a sexual population. In 2015, a second genetically diverse Pst pathotype called was detected for the first time in the UK, termed the 'Kranich' race (UKCPVS 2016). This race is broadly related to the 'Warrior' group had been previously detected in continental European countries, and both races are thought to have originated from sexually recombining populations in the near-Himalayan region in Asia (Hovmøller et al. 2016). Four of our trials were inoculated with a mixture of 'Solstice' and 'Warrior' Pst races, while OSG16 was naturally infected. However, natural infection was especially high in both the 2015 and 2016 seasons, with infection pressure being very high even at the time of trial inoculation. Testing of Pst isolates from the two locations used in our field trials shows these locations were dominated by natural infection by the Warrior group of races in 2015 and 2016, as was the case within the wider Cambridgeshire and Lincolnshire regions (Supplementary Table 5), and as mirrored across most of the UK (UKCPVS 2012(UKCPVS , 2016. Thus, natural infection was most likely predominant in all five trials for the following reasons: (1) the high natural YR infection levels in 2015 and 2016 at these trial sites, as well as across the UK. (2) At our Lincolnshire trial sites OSG and ROTH, which were included in these surveys, the Warrior group of Pst races dominated. (3) Our own isolate pathotyping from the NIAB15 and ROTH15 trials. (4) Our observation of natural YR infection prior to trial inoculation, and mirroring regional reports of very high YR disease pressure starting from autumn 2015 when mild conditions allowed rapid early spread (UKCPVS 2012). Therefore, we conclude that the resistance QTL identified in the MAGIC population most likely conferred adult plant resistance to genetically diverse Pst races that characterise the recent rapid shifts in the population type of this pathogen that from 2011 replaced the previous clonal Pst forms.

MAGIC yellow rust resistance genes
The YR resistance loci identified at the adult plant stage controlled a large percentage of the phenotypic variance, predominantly accounted for by four loci on chromosomes 1A, 2A, 2B and 2D. Our finding that in almost all cases, combining two or more of these QTL resulted in increased resistance provide experimental evidence that such stacking should be effective in providing strong genetic control for YR. While good resistance was generally provided by stacking two QTL, in practice the use of more loci would provide increased security against future partial or full breakdown of any single resistance locus. While there was no evidence for any of the 'major' QTL breaking down over the two seasons investigated, subsequent reports indicate that the resistance conferred by alleles from the MAGIC founders Alchemy and Claire at QYr.niab-2D.1 has been overcome (Simon Berry, personal communication), indicating it may represent an all-stage resistance locus. Gradual degradation of the effective adult plant YR resistance historically conferred by Claire since its release in 1999 (Powell et al. 2013) can first be traced to the period following the incursion of the 'Warrior' group of Pst races, when in 2011-2012 Claire went from having the highest resistance score of 9 down to an intermediate score of 6 (UKCPVS 2016). Thus, while Claire's YR resistance had partially broken down, it nevertheless contained sources of resistance unaffected by at least some of the new Pst races present in our 2015 and 2016 seasons, but which was further eroded by the subsequent breakdown of QYr.niab-2D.1. Interestingly, four adult plant YR resistance QTL have been identified on chromosomes 2B, 2D and 7B in field-grown trials of a Claire × Lemhi (C × L) population grown in trials between 2003 and 2007 (Powell et al. 2013). The trials were conducted before the introduction of the Warrior type races in 2011, and Claire alleles conferred resistance at both loci. Of these, two QTL co-located with the adult plant resistance loci identified in the MAGIC population. CxL QTL QYr.niab-2B overlapped with MAGIC QTL QYr.niab-2B.1 (based on the CxL physical interval defined by markers wPt-0950 and wPt-9190: 685.047-750.121 Mbp). However, the MAGIC founder Soissons confers the allele with the highest resistance at the QYr. niab-2B.1 locus, indicating that the underlying loci at these chromosome 2B QTL are most likely different, or that alleles with higher resistance than that conferred by Claire are present. CxL QTL QYr.niab-2D.2 co-located with our MAGIC QTL QYr.niab-2D.1 (based on CxL markers EST18a and wmc817a: 637.651-645.997 Mbp), for which resistance alleles in our MAGIC population were conferred by Claire and Alchemy. Notably, Alchemy has Claire in its pedigree [Alchemy = Clare × (Consort × Woodstock)] (Fradgley et al. 2019). Indeed, the haplotypes of Claire and Alchemy, based on the 90 k SNP data, are identical across the QTL confidence interval (from SNP RFL_Contig1128_620 to BS00010685_51, 188.01-198.86 cM), indicating the chromosome 2D resistance alleles carried by MAGIC founders Claire and Alchemy are identical by descent. The absence of the chromosome 2D and 7B CxL QTL in our MAGIC population may reflect changed virulence profiles of the genetically diverse Pst races prevalent in our 2015/2016 season trials compared to the CxL trials that preceded the incursion of exotic Pst races. Alternatively, it is possible that none of the MAGIC founders carry the susceptible allele that originated from Lemhi, which is a US variety. Interestingly, a recent study using the 'NIAB Elite MAGIC' population found a robust QTL conferring resistance to the necrotrophic fungal pathogen Parastagonospora nodorum (the causal agent of Septoria nodorum blotch, SNB) located in the genomic region as the major YR resistance QTL QYr.niab-2A.1 . However, comparison of the predicted allelic effects at these QTL finds the resistant alleles at the YR QTL carried by the founders Rialto and Xi19 to be associated with susceptibility for SNB. Analysis of multiple traits in the same population allows such correlations to be identified and further investigated. Indeed, as noted by Scott et al. (2020), the comparatively high genetic diversity and genetic recombination captured by multi-founder populations makes them well suited for the genetic analysis of multiple traits within a single experimental population, so maximising the chances of identifying potential trade-offs between traits.

Future exploitation of the major YR resistance loci
We identified sources of YR resistance originating from different subsets of the population founders, highlighting the benefit of conducting analyses in multi-founder populations. For the four most significant MAGIC QTL, analysis of the haplotypes present in the varietal panel of 403 accessions indicated between two-to-four 'common' haplotypes (defined here as occurring at a frequency of > 10%) were present per locus. The relatively low haplotype diversity at these loci in north-western wheat is in agreement with that identified across the genome in a recent analysis of 16 wheat varieties selected to maximise genetic diversity within this geographic region (Scott et al. 2021). Given such findings, it was perhaps unexpected that for three of our four large effect YR resistance loci, the resistant haplotype was rare in the variety panel (1A-hap1 = 2%, 2A-hap1 = 5%, 2B-hap1 = 7%), further highlighting their possible usefulness as immediate targets for genetic improvement of YR resistance. There could be several reasons for the low frequency of these resistance haplotypes. For example, in the case of the resistance allele conferred by MAGIC founder Hereward at the chromosome 1A QTL, while Hereward was well known and widely grown in the UK for its good grain quality (Mackay et al. 2014), it is the parent of very few subsequent varieties in the wheat pedigree (Fradgley et al. 2019), so reducing the likelihood of further use of the resistance allele in subsequent varieties. Resistance alleles may also be genetically linked to chromosomal regions conferring reduced performance for other agronomically important traits, and so be selected against. Indeed, our chromosome 1A resistance locus was also notable in that the entire QTL interval of ~ 11 cM identified in MAGIC was nested within a wider haploblock of ~ 13 cM in the varietal panel. Given the importance of YR resistance for breeder selection, the rarity of the 1A-hap1 haplotype from Hereward and the relatively large genetic interval of the haploblock in the varietal panel (spanning a chromosomal region that the MAGIC population shows undergoes frequent genetic recombination), it is possible a selective sweep for other agronomically important traits has occurred across this region. Indeed, analysis of changes in allelic diversity across the genome previously identified a SNP adjacent to our peak SNP for the chromosome 1A YR resistance QTL as being under strong selection (located on chromosome 1A at ~ 557 Mbp, based on SNP BS00032825_51; Fradgley et al. 2019). Furthermore, genome-wide association scans in a subset of the varietal panel we use here identified a nearby genetic locus controlling grain yield (at ~ 544 Mbp; White et al. 2021). The peak of our chromosome 1A YR resistance QTL lies towards the distal end of the haploblock, at 568 Mbp. As the KASP markers developed here distinguish the resistant 1A-hap1 haplotype from the remaining susceptible haplotypes, these could be used to help combine a distal region containing the YR resistance haplotype 1A-hap1 with a proximal region consisting of the common 1A-hap2 or 1A-hap3 haplotypes, and investigate the consequences on agronomic traits and disease resistance. Thus, consideration of our chromosome 1A YR resistance locus illustrates how complementary sources of genetic, genomic, molecular and phenotypic data and resources might be integrated to effectively exploit sources of genetic resistance.

Analysis of YR QTL physical intervals
All QTL identified were outside of the pericentromeric regions that are characterised by low genetic recombination and encompass approximately half of the wheat genome (Fig. 5). This, combined with high heritability of the phenotype and the strong contrasting phenotypic effect of the majority of the eight 'major' resistance QTL means that it could be useful to exploit the residual heterozygosity in the MAGIC RILs to create near isogenic lines for specific QTL for subsequent fine-mapping studies. Ultimately, combining such knowledge with a better understanding of the genes conferring YR susceptibility, such as the recently identified branched-chain amino acid aminotransferase gene TaBCAT1 that modulates amino acid metabolism (Corredor-Moreno et al. 2021), as well as gene editing (reviewed by Kumar et al. 2019), may allow future informed design and use of resistance alleles to enhance durable resistance to YR and other wheat pathogens. To date, just three YR adult plant resistance genes have been map-based cloned: Yr18/Lr34 (Krattinger et al. 2009), Yr36 (Fu et al. 2009) and Yr46/Lr67 (Moore et al. 2015), and no homologues of these genes were found within the physical intervals of any of the YR resistance QTL we identified here. However, the allstage resistance gene Yr5 is somewhat atypical in that it is known to confer resistance to a broad range of Pst isolates worldwide (Marchal et al. 2018). Yr5 is located in a region of chromosome 2B known to contain numerous yellow rust resistance genes (Feng et al. 2015;Luo et al. 2008), including Yr7 and YrSP. Here we find the physical interval of our MAGIC YR resistance locus QYr.niab-2B.1 to overlap with the Yr7 and Yr5/YrSp loci. While Yr7 and YrSP no longer provide adequate resistance in the field, Yr5 continues to remain effective against a range of isolates worldwide (Marchal et al. 2018). Yr5 (represented in the wheat reference genome of cv. Chinese Spring by gene model TraesCS2B01G488700) and Yr7 (for which no equivalent Chinese Spring gene model is thought to be present; Marchal et al. 2020) are paralogous BED-NLR genes with ~ 78% DNA identity across their coding regions (Marchal et al. 2018). YrSP is a truncated allele of Yr5 resulting in the loss of most of the LRR coding region (Marchal et al. 2018). Resistant Yr5 and YrSP alleles are not reported to be present in European wheat, with diagnostic molecular tests further confirming their absence in all eight MAGIC founders (Marchal et al. 2018). Similarly, molecular analysis has shown Yr7 to be absent from the MAGIC founders (Marchal et al. 2018). In our MAGIC population, resistance at QTL QYr.niab-2B.1 came from the Soissons allele. Collectively, this indicates that while Yr7 and Yr5/YrSP do not themselves confer resistance at QYr.niab-2B.1, it is possible that different mutations at their underlying genes, or allelic variation at nearby paralogous BED-NLRs in the vicinity, could underly resistance at this locus. Indeed, recent analysis of wheat genome assembly sequences has found copy number variation for BED-NLRs within the vicinity of the Yr5 and Yr7 locus (Marchal et al. 2020). Ultimately, further investigation of these, and other, candidates at the QTL is required to narrow down the underlying causative gene and polymorphism. This process will be greatly aided by the finalisation of genome sequence assemblies for all eight MAGIC founders, two of which are currently available (Walkowiak et al. 2020) and the remaining six close to completion, as well as investigation of wheat gene expression resources (e.g. Borrill et al. 2016), including samples collected across different stages of wheat Pst infection (Adams et al. 2021).