This study used a panel of 296 bread wheat genotypes from released varieties and advanced breeding lines developed primarily over the last twenty years from 25 research centres across India. DNA extraction and DArT sequencing of the genotypes were done per the protocol demonstrated by Li et al. (2016). The wheat accessions were sequenced with the DArTseq^® technology at the Genetic Analysis Service for Agriculture (SAGA) at CIMMYT, Mexico. Markers with a minor allele frequency of less than 10% (2804 markers) or more than 30% missing data points (96 markers) were excluded from further analysis. A total of 9668 SNPs were finally used for the GWAS analysis, and their physical positions on the reference whole genome sequence (IWGSC: Chinese Spring RefSeq v1.0, International Wheat Genome Sequencing Consortium et al., 2018) were acquired from the database https://wheat-urgi.versailles.inra.fr/). The SNP markers were given names based on their chromosome location followed by clone ID, e.g., 5BL:3955588.

The panel was evaluated for 3 years (2019, 2021, 2022) for SNB and 2 years for TS (2019, 2021) in greenhouse at the seedling stage. Mexican P. tritici-repentis (Ptr) isolate MexPtr1 and P. nodorum isolate MexSn4 were used for resistance screening against TS and SNB, respectively. Both isolates are ToxA producers based on inoculation experiments using differential genotypes, infiltration experiments, and PCR with the ToxA-specific marker (data not shown). The isolates were grown on V8-PDA media (Lamari and Bernier, 1989), and conidiospore concentrations for inoculation were adjusted to 4 × 10³ spores mL⁻¹ (MexPtr1) and 1 × 10⁷ spores mL−1 (MexSn4) (Singh et al., 2007; Singh et al., 2016). The response of TS and SNB was tested on seedlings in a greenhouse at 22°C day and 18°C night temperatures with a 16-h photoperiod. Experiments were set up in a randomized complete block design with two replicates, with four plants grown in plastic containers as experimental units to derive mean values for further analysis. Erik and Glenlea were used as resistant and susceptible controls, respectively. Inoculations were performed when the second leaf was fully expanded at around 14 days after planting. The inoculum was applied to seedlings with a hand sprayer until runoff (about 0.5 mL inoculum per plant). The trays were transferred to a humid chamber (RH 100%, 20°C) once the leaves were dry to promote infection, and the plants were returned to the greenhouse bench after 24 hours. Both diseases were rated on a linear scale of 1–5 at seven days after inoculation (dpi) (Feng et al., 2004; Hu et al., 2019).

All 9668 SNP markers were used to calculate a kinship matrix, clusters among individual genotypes, and a heat map using the traditional Van Raden (2008) equation using TASSEL v5 (http://www.maizegenetics.net, accessed on 25 Oct 2022). Linkage disequilibrium (LD) was analysed with R² among SNP markers plotted against the physical distances in mega base pairs (Mb) across the 21 wheat chromosomes. The genotypic data was numerically transformed for population structure analysis utilizing XLSTAT (v. 2022.1). Structure 2.3.4 software was used to obtain the population structure (Pritchard et al., 2000). The admixture model was adjusted with a 100,000 burn-in period followed by 500,000 marker chain Monte Carlo (MCMC) iterations. The subpopulation test range was maintained at K1 - K5, with five iterations. The actual subpopulations were assessed using the DK approach (Earl and vonHoldt, 2012), which was verified using the Structure Harvester program (Web v0.6.94, Earl and vonHoldt, 2012) as per the method described by Evanno et al. (2005). The output summary calculated the standard deviation and average logarithm of the probability of the observed likelihood [LnP(D)]. The log-likelihood of the data was computed for each class (K = 1 to 5) to determine LnP(D) for each MCMC step. A neighbour-joining tree was created using TASSEL 5.0 (Bradbury et al., 2007) and visualized using the iTOL website (Letunic and Bork, 2021).

For GWAS, three models (MLM, MLMM and FarmCPU) implemented in the R package GAPIT3 (Wang and Zhang, 2021) were used. The first was the mixed linear model (MLM=K+Q), which was based on the kinship matrix (K) and the principal component (PC) (Yu et al., 2006). The second method was the Multiple Loci Mixed Linear Model (MLMM), which employs forward-backward stepwise linear mixed-model regression to include associated markers as covariates. In contrast to MLM and MLMM, the third model FarmCPU (fixed and random, circulating probability unification); used all linked markers in a fixed-effect model and optimised the linked markers in a separate random-effect model, reducing false positives and false negatives and enabling quick computation. Population structure, principal components, and kinship were used as covariates in MLM and MLMM, and their QQ plots were compared. Our phenotypic data fit the FarmCPU model better, and QTL showed higher significance; therefore, only data analyzed by the FarmCPU model was chosen for additional examination. The Bonferroni correction (α=0.1) was used for an exploratory significance threshold to uncover putative QTL (Chan et al., 2010). QTL was recognized as robust when linked markers reached the strict -log10(p) criteria of 3.0 in at least being detected across-year by either threshold adopted in this investigation and consistent across the models. The Quantile-Quantile (QQ) plots were also examined to determine the point at which the observed p-values diverge from those predicted by the null hypothesis. SNP marker annotations were obtained using the databases http://www.cerealsdb.uk.net and https://triticeaetoolbox.org. The physical positions of markers on the reference whole genome sequence (IWGSC: Chinese spring RefSeq v1.0, International Wheat Genome Sequencing Consortium et al., 2018) were acquired from the database https://wheat-urgi.versailles.inra.fr/ and https://plants.ensembl.org/Triticum_aestivum/Info/Index. If two significant markers shared a 10-Mbp interval or had substantial LD (R²>0.8) with one another, they were regarded as belonging to the same QTL region.

Kompetitive Allele-Specific PCR (KASP) marker Ktsn1 specific for Tsn1 (Dreisigacker et al., 2016) was used for genotyping the panel. Genomic DNA was isolated from 21-day-old seedling leaves, DNA concentration was measured with a Nanodrop 8000 spectrophotometer (Thermo Scientific), and the final concentration of 10 ng/ml was obtained by diluting the DNA with sterile PCR-grade water. KASP genotyping was undertaken using the PACE master mix per the manufacturer’s guidelines (3CR biosciences). In brief, the assay mix for 25 samples was prepared as 3µl VIC primer, 7.5 µl common primer, and 14.5 µl of PCR grade H2O. Finally, for reaction assembly, individual reactions contained 5 µl 2X PACE, 0.138 µl assay mix, and 5 µl genomic DNA. The samples were run on a CFX96 Real-time PCR system (Bio-Rad). The PCR steps included 94°C for 15 min, 10 cycles of touch down program (94°C-0.20 min, 65 to 57°C-1 min), followed by 30 cycles of (94°C-0.20 min, 57°C-1 min). The results were analyzed in Bio-Rad CFX manager 3.1, and output was obtained as.csv files. The results obtained were compared against the linked SNP calls from the SNP array.

For haplotype analysis, stable MTAs on chromosomes 1B, 2A, 5B, and 6B across experiments were selected. Eight markers that fitted the requirement that they were over the -log10(p) 3.0 threshold across the environments were chosen for the haplotype analysis: 5BS:1102120, 2AS:7487614, 2AS:1094287, 5BL:5324846, 5BL:3955588, 6B:1085698, 1B:1129298, and 5BL: Ktsn1. For TS, three markers, 5BL:5324846, 5BL:3955588, and 5BL : Ktsn1 were selected for the analysis. The corrected disease severities between haplotypes were compared using the Wilcoxon test implemented in R package ggpubr (Kassambara, 2020). The phenotypic data available for the seedling resistance to SNB in 2019, 2021, 2022, and TS in 2019 and 2021 was used for the analysis.

Data for three years (2019, 2021, and 2022) were used for this analysis. The corrected mean disease index from BLUP was used for haplotype analysis. Marker trait associations for SNB and TS were chosen from various models (MLM, MLMM, FarmCPU) to investigate the effect of stacking resistance alleles ( Table S1 ). Resistant alleles were determined by mean comparison of corrected disease severity between alleles based on the Wilcoxon test, using the R package ‘ggpubr’ (Kassambara, 2020). Wheat lines were classified according to the number of resistant alleles present. The t-test (p < 0.05) was implemented in the R package ‘multcomView’ (Graves et al., 2015) to differentiate the significant groups.

To assess the gene effect of Tsn1 and Snn3, two makers, 5BL:5324846 and 5BL:3955588, linked to Tsn1 and data obtained from KASP assay Ktsn1 was compared with SNP 5BS:1102120 linked with Snn3. The interaction was discovered by comparing the mean between alleles using the t-test in the R package ‘ggpubr’ (Kassambara, 2020).

The corrected disease index, genotypic and phenotypic variance, and heritability estimates were obtained from META-R software (Alvarado et al., 2020). Additionally, mean comparisons with the Wilcoxon test, t-test, and data visualization were executed using packages `ggpubr`, ‘dplyr’, and ‘multcomView’ using R software (Version 4.2.1, R Core Team, 2022).

The population displayed significant phenotypic variation for SNB and TS resistance with a skewed distribution towards the lower disease in all experiments ( Figure S1 ). A high proportion of genotypes exhibited high resistance to both diseases (See Table S1 for the disease scores) and the top 20 entries resistant to both diseases are presented in Table 1 , including PBW 277, PBW 771, and HD 3319.

Analysis of variance revealed significant effects on genotype and genotype-by-environment for both diseases ( Table 2 ). SNB showed higher overall heritability (0.87) than TS (0.75), with heritability in individual experiments ranging from 0.90 to 0.96.

Among the 9668 SNP markers selected for GWAS in the 296 bread wheat genotypes, 21.78% were from the A genome, 28.63% from the B genome, 12.03% from the D genome and 21.33% from unknown chromosomes. Population structure analysis based on the K means cluster approach divided the population into 4 major clusters that were depicted as a neighbour-joining tree ( Figure 1 ). The kinship analysis also distinguished the population into 4 major groups and presented it as a kinship matrix-based heat map ( Figure 2 ). The average extent of LD considered physical distance taken for the decay of R2 to reach a critical value of 0.10 across the genome, was approximately 10 Mb ( Figures 2 , S2, S3 ).

The exploratory -log10(p) threshold for the panel ranged from 2.92 to 11.53 for SNB and 2.96 to 15.86 for TS ( Tables S1 , S2 ). In total, 49 marker-trait associations (MTAs) were detected on various chromosomes for SNB ( Table S1 ; Figures 3 , S4 ). Six QTL on chromosomes 1BL, 2AS, 5BL, and 6BL were particularly significant for SNB over all three years and were considered robust QTLs. The QTL Q.CIM.snb.1BL on chromosome 1BL was identified with FarmCPU, and its peak marker 1BL:1129298 (p = 2.86E-04 to 8.97E-06) was located at 450.50Mbp on 1BL ( Table 3 ). Two QTLs on 2AS, Q.CIM.snb.2AS1 and Q.CIM.snb.2AS.2, were detected in all models with their peak markers being 2AS:7487614 (59.43Mbp) and 2AS:1094287 (88.18Mbp), respectively. Another robust QTL, Q.CIM.snb.5BS, was identified at 5BS:1102120 (p=2.31E-05 to 5.99E-12) and was associated with Snn3-B1. Two Tsn1-associated markers 5BL:5324846 and 5BL:3955588 were detected across years and models on chromosome 5B at 566.04 and 588.4 Mbp, respectively ( Table 3 ). The QTL Q.CIM.snb.6BL (p= 4.16E-04 to 1.61E-09) was detected using farmCPU, on the long arm of the chromosome 6B with peak marker 6BL:1085698 (669.13 Mbp).

For TS, 79 MTA (-log10(p) >3) were detected on various chromosomes ( Table S2 ; Figure 3 ). Two Tsn1-associated markers 5BL:5324846 and 5BL:3955588 were significantly detected across years. The SNP 5BL:5324846 (p=3.51E-05 to 5.41E-06) was detected across the three models (MLM, MLMM, FarmCPU). While the SNP 5BL:3955588 (p=1.33E-05 to 3.63E-13) was detected with a negative effect across the years (2019, 2021) from models MLM and MLMM only ( Table 3 ).

Out of the 49 MTA identified for SNB resistance, 15 were selected to analyse their association with predicted genes and the corresponding biological functions. The criteria for shortlisting MTAs was that they were detected in at least two models or two experiments (See Table 4 ). Two Tsn1-associated markers, 5BL:5324846 and 5BL:3955588, were found to be associated with the transcripts TraesCS5B02G409000, TraesCS5B02G387500, and TraesCS5B02G387000. These transcripts are known to possess biological function pathogen stress response (PADRE) domain, Serine/threonine-protein kinase and Ubiquitin-like proteins respectively. Another marker, 2AS:7487614, has been found to be associated with two genes, namely TraesCS2A02G107100 and TraesCS2A02G107200. These genes are known to have biological roles as a potassium transporter and a plant peroxidase, respectively. As for other MTAs, several noteworthy biological functions have been identified, including cysteine-type endopeptidase inhibitor activity (6B:1085698), NAD(P)-binding domain (4A:1082366), Cytochrome c oxidase subunit 5c, and Serine/threonine-protein kinase, active site (both for 3B:2251334). These functions are believed to play an active role in disease resistance ( Table 4 ).

The KASP assay revealed that 83 genotypes (25.85%) carry the susceptible allele Tsn1, whose SNB indices ranged mostly between 2.4 to 4.2 ( Table S3 ). The allelic combinations at Tsn1 and Snn3 were tested using the results obtained from the Ktsn1 marker and two Tsn1 linked markers, 5B:5324846 and 5B:3955588, along with Snn3-B1 linked marker 5BS:1102120. The comparison revealed that the allelic combination of recessive alleles tsn1/snn3 imparts best resistance to SNB, followed by the tsn1/Snn3 combination, whereas the allele combinations Tsn1/snn3 and Tsn1/Snn3 showed similar levels of resistance, indicating the higher phenotypic effects of tsn1( Figure 4 ).

There were significant differences in the corrected SNB index between resistant and susceptible haplotypes. When compared to the resistance haplotypes “AA/TT” for tsn1 (mean SNB index 1.69), the susceptibility haplotype “GG” for Tsn1 had a higher disease index (3.4) ( Figure 5 ). Furthermore, the haplotypes were tested for the four significant SNB QTLs, Q.CIM.snb.1BL, Q.CIM.snb.2AS.1, Q.CIM.snb.2AS.2, and Q.CIM.snb.6BL ( Figure S5 ). Similarly, there was a significant difference between susceptible and resistant haplotypes for TS index (3.2 for the susceptibility haplotypes and 1.57 for the resistance haplotypes; Figure 5 ).

The stacking of resistance alleles for SNB and TS exhibited a clear trend that the number of resistance alleles is positively correlated with disease resistance ( Figure 6 ). For example, cultivars DBW 242 and DBW 246, each possessing 35 resistant alleles, exhibited very high SNB resistance, whereas HD 2932 and Kharchiya 65 having the least number of resistance alleles were highly susceptible. A similar trend was observed for TS resistance.