The study used 368 UCR minicore genotypes described by (Muñoz-Amatriaín et al., 2021). The UCR minicore contains diverse landraces and breeding materials from 50 countries covering Africa, Asia, North and South America, Europe, and Australia. The minicore population, previously genotyped with high-density SNPs, represents the existing genetic and phenotypic diversity of cultivated cowpea while maintaining a sample size that can be managed by most researchers and breeders for evaluating traits of interest (Muñoz-Amatriaín et al., 2021). For marker development, the present study used 60 diverse cowpea genotypes along with 46 F₁ progenies for SNPs technical verification and another independent set of 20 genotypes consisting of some known resistant and susceptible sources commonly used as standard checks in breeding programs for marker validation.

The study was conducted at three sites in Northern Nigeria where there is a widespread of Striga gesnerioides. The first site was a Striga hot spot at IITA Minjibir Research Farm, Kano, Latitude 12° 14′ 35.30″ N and Longitude 8° 66′ 62.10″ E. The second Striga hot spot site was at Malam Madori, Jigawa, located at Latitude 12° 33’ 36.32” N and Longitude 9° 59’ 9.56” E. Previously, the samples of Striga races collected from these two sites have been characterized as belonging to race SG3, which is a dominant S. gesnerioides race in Nigeria among other less frequent races labeled as SG1, SG5 and SG6 (Ohlson and Timko, 2020). The third site designated for validation screening is situated at the IITA station in Kano, located at latitude 11°58’50.0”N and longitude 8°33’26.8”E. This facility includes a nursery that has been deliberately infested with Striga gesnerioides race SG3 to facilitate controlled research on cowpea resistance to this parasitic weed.

The minicore genotypes were evaluated in the two locations in Northern Nigeria for two years in Minjibir and one year in Malam Madori. The genotypes were planted in a two-row plot of 4 meters long and at a spacing of 0.75m between the rows and 0.2m within the row. The trials were established as alpha lattice designs with two replications. Variation in Striga infestation was assessed using two criteria: First, the Striga score was evaluated on a scale of 0-3, where 0 = no Striga emergence; 1 = few Striga emerged; 2 = moderate Striga emergence; 3 = heavy Striga infestation The second method was the presence-absence rating, where ‘0’ indicates that Striga is absent in the plot and ‘1’ indicates Striga is present. In addition, a validation set of 20 cowpea genotypes was screened in an artificially infested nursery using a randomized complete block design experiment. For the past several years, this nursery has been dedicated to Striga screening following heavy inoculation of the soil with seeds of Striga collected from infested fields at Malam Madori, Jigawa State, Nigeria. Prior to soil inoculation, Striga seeds were preconditioned through surface-sterilization using 10% sodium hypochlorite for ~10 minutes followed by incubation in the dark at 29°C for 10 days. The soil was then inoculated with Striga seeds at a rate of 2 g/m², evenly distributing the seeds across the soil surface. The 20 cowpea genotypes were each planted in a 2-row plot of 4m and spaced at 0.75m by 0.2m between and within rows, respectively, and the experiment had two replications. The presence or absence of Striga was assessed for each plot, including the number of Striga plants that emerged, and the data was analyzed with Analysis of Variance (ANOVA) and means compared using LSD test.

The SNP data was obtained from (Muñoz-Amatriaín et al., 2017). The genotyping was done at the University of Southern California Molecular Genomics Core lab (Los Angeles, California, USA) following the procedure described by (Muñoz-Amatriaín et al., 2021). Briefly, genomic DNA was extracted from the young leaves of individual plants using DNeasy Plant Kit (Qiagen, Valencia, California, USA). The Cowpea iSelect Consortium Array including 51,128 SNPs (Muñoz-Amatriaín et al., 2017) was used for genotyping each DNA sample. SNPs were called in GenomeStudio (Illumina Inc., San Diego, California, USA) and manually curated to remove those with high levels (>20%) of missing data and/or heterozygous calls. The genotype data contained 47129 SNPs after removing the contigs with no chromosomal assignment. Further filtering was conducted, setting the minimum minor allele frequency at 0.05 and maximum heterozygous proportion at 0.1, resulting in 41510 SNPs for downstream analyses.

The phenotypic data were analyzed using a general linear model tailored for an alpha lattice design, implemented through the agricolae package in R (de Mendiburu, 2020). The model applied is represented as:

Where: Y _ i j k denotes the observed phenotypic value for the i ^th treatment in the k ^th block within the j ^th replication., μ is the general mean, τ _ i represents the fixed effect of the i ^th treatment, γ _ j signifies the effect of the j ^th replication, and ρ _ k ( j )   accounts for the effect of the k ^th incomplete block nested within the j ^th replication, ϵ _ i j k is the random error term associated with the observation. The mean values for each minicore accession were extracted from the model using the LSD.test () function and later utilized in GWAS analysis.

GWAS was conducted using the 41510 filtered SNPs in TASSEL v 5.2.20 (Bradbury et al., 2007) and rMVP package v 1.03 (Yin et al., 2021), utilizing the Mixed Linear Model (MLM) (Zhang et al., 2010) and the Fixed and random model Circulating Probability Unification (FarmCPU) (Liu et al., 2016) on the following traits: Striga score and Striga presence-absence rating. The FarmCPU model uses a multiple loci linear mixed model (MLMM) and incorporates multiple markers simultaneously as covariates in a stepwise MLM to partially remove the confounding between testing markers and kinship (Liu et al., 2016). A genomic PCA matrix (P) and kinship matrix (K) were used to capture the population structure and relatedness among individuals in the panel (Kang et al., 2008). In TASSEL, kinship, and principal components were computed and fitted together with the SNPs in the MLM model to account for relatedness and population structure, respectively. The MLM statistics were also utilized separately in CMplot package 3.1.3 (Yin, 2016) to generate customized GWAS Manhattan and QQ plots. Decisions on significant GWAS signals were based on both the conservative Bonferroni (Henry, 2015) and a less conservative false discovery rate (FDR) (Glickman et al., 2014; Ongom et al., 2022, 2024) corrections of multiple statistical tests to reduce the risk of a type I error.

Candidate genes linked to Striga resistance were pinpointed by aligning the positions of significant SNPs from the GWAS with the cowpea reference genome (version 1.1) available on Phytozome 13 Genome Browser (https://phytozome-next.jgi.doe.gov/info/Vunguiculata_v1_1, accessed on 26 January 2025). A length of 200 kb was added or removed from either end of the significant marker to locate potential regions for comparison based on the LD rate of the cowpea minicore population (Ongom et al., 2022). In selecting candidate genes, the following criteria were used: (i) genes of known function in cowpea related to the trait under study, (ii) genes with function-known orthologs in Arabidopsis related to the trait under study, and (iii) genes pinpointed by the peak SNPs. Putative candidate genes were subsequently researched in the literature for verification. We then used chromoMap v4.1.1 package (Anand and Rodriguez Lopez, 2022) to visualize the distribution of selected genes on target chromosomes where GWAS signals were discovered.

Technical verification of KASP assay: Seventeen SNP markers spanning Striga resistance gene regions were selected to develop KASP assays. The design sequences of the selected SNPs are presented in Table 1 . The validity of these KASP assays was verified using 60 cowpea genotypes and 46 heterozygous F₁ progenies. The samples for this SNP technical verification were collected in three replicates per genotype. The leaf samples and SNP design sequences were sent to Intertek Lab Sweden for assay development. Leaf sampling followed the standard procedure previously described (Ongom et al., 2021). Genotyping data were visualized using SNPviewer version 5.1.1.27582 (LGC Genomics, 2025) and SNP cluster plots were generated.

Testing the marker-trait association: To validate the performance of these markers for tracking Striga resistance, we deployed the candidate markers in a panel of 20 cowpea validation genotypes consisting of some known resistant and susceptible sources. These cowpea genotypes were independently phenotyped in a an artificially Striga-infested nursery to confirm their reaction status. The number of Striga plants that emerged in each plot was counted, and the data was used to categorize the 20 genotypes into either resistant (R) or susceptible (S) with zero Striga count = R and any count greater than zero = S. The 20 cowpea genotypes were then genotyped with the 17 candidate KASP markers. A single marker analysis was then conducted to verify the consistency of these candidate markers in tracking Striga resistance in cowpea. A chi-square test of independence and a t-test were performed in R for each marker with the null hypothesis that each marker genotype is independent of the observed variation in the phenotypic status of the cowpea genotypes. The results of these tests were visualized using plots generated by ggplot2 and ggstatsplot packages.

Phenotypic assessment

The phenotypic data from both Minjibir and Malam Madori sites depicted skewed distributions for Striga ratings, with a greater proportion of the minicore genotypes showing susceptibility to Striga. This is portrayed by the stacked bar chart in Figure 1A , where only 2.7 to 6.5% of the minicore showed resistance across the two testing sites. In contrast, a mean Striga score of 2.35 (scale 1-3) was registered for the combined data across the two locations ( Figure 1B ). The susceptible genotypes triggered extensive growth of Striga plants, causing severe yellowing and death of cowpea plants. At the same time, no emergence of Striga was observed in plots planted with resistant genotypes ( Figure 1C ).

GWAS discovered one major association signal for Striga resistance on chromosome Vu11 that was consistently significant in both locations and five other minor signals on chromosomes Vu02, Vu03, Vu07, Vu09, and Vu10 ( Figure 2A ). The major association signal on chromosome Vu11 and two minor regions on chromosomes Vu07 and Vu10 were detectable in both Minjibir and Malam Madori locations, while the regions on Vu03, Vu02, and Vu09 were significant only in one location ( Figure 2A ). The QQ plot shows that the observed p-values largely align with the expected distribution under the null hypothesis, indicating no evidence of inflation. However, the upward deflection at the tail suggests the presence of a small number of SNPs exhibiting stronger association signals than expected by chance.” ( Figure 2B ). It was also evident from the QQ plot that the GWAS signals were stronger in Malam Madori than in Minjibir and the combined location data set for Striga scores ( Figure 2B ).

A total of 309 significant SNPs with the Boneforrini threshold of -log10(p) > 5.92 were found to be associated with Striga resistance ( Supplementary Table 1 ). Out of these, seventeen (17) representative SNPs that had R² values above 7.5% were selected for KASP assay development ( Table 2 ). The strongest GWAS signal on chromosome Vu11 was represented by five significant SNPs, accounting for 15.8% to 19.8% of phenotypic variance ( Table 2 ). The other minor signals on chromosomes Vu02, Vu03, Vu07, Vu09, and Vu10 explained 6.8% to 14.1% of phenotypic variance ( Supplementary Table 1 ).

Searching candidate genes within 200kb of peak GWAS signals identified 178 genes on four target chromosomes, with 28 on Vu11, 38 on Vu10, 56 on Vu07, and 57 on Vu02 ( Supplementary Table 2 ). Sixty-four of the 178 genes were clustered closely to the peak SNPs on the target chromosomes ( Figure 3 ). Further examinations revealed 20 unique annotated proteins associated with the 64 genes ( Figure 3 ), given that most of these gene loci encoded similar functional proteins. Of the 64 candidate genes identified, 11 were positioned within approximately -167461 bp to 191285 bp of the lead SNP on Chromosome Vu11, five genes were within ~ -139927 bp to 83997 bp of the association peak on Vu10, 36 genes were within ~ -190061 bp to 202188 bp of the signal on Vu07 and 12 were within ~ -145716 bp to 125678 bp of the association locus on Vu02 ( Table 3 ). Several of these clusters contain genes with related biological functions. For instance, up to five genes on Chromosome Vu02 belonged to the Pentatricopeptide repeat (PPR) superfamily protein. On chromosome Vu07, up to 28 genes had functions related to Cysteine-rich receptor-kinase-like protein, while 10 genes on chromosome Vu11 had functions related to LEUCINE-RICH REPEAT-CONTAINING PROTEIN ( Table 3 ). The rest of the genes on each chromosome had unique functions, including, among others, C2H2-like zinc finger protein, myb transcription factor, MADS-box transcription factor family protein, ETHYLENE-RESPONSIVE TRANSCRIPTION FACTOR ERF003, and WRKY DNA -binding domain (WRKY) ( Table 3 , Figure 3 ).

KASP assays were successfully developed for the 17 representative SNPs having strong associations with Striga resistance in cowpea. The SNP calls were verified using a unique technical validation set of cowpea genotypes, including homozygous genotypes and highly heterozygous F₁. The raw genotype calls displaying the KASP assay performance of the 17 candidate SNPs in the cowpea technical validation set are provided in Supplementary Table 3 . The KASP assay technical verification results revealed 12 SNPs as having good quality assays, allowing easy scoring of the alleles ( Figure 4A ). Two of the SNPs were of medium quality, with some ambiguity in differentiating between homozygotes and heterozygotes ( Figure 4B ). Two other SNPs were rated as having inconclusive quality due to their sensitivity to DNA concentration and the resultant difficulty in scoring the marker genotypes ( Figure 4C ), while one SNP did not form scorable clusters and was regarded as having a bad quality ( Figure 4D ). A detailed report on the quality assessment of each SNP is presented in Table 4 . SNPs that were easy to score formed three genotype clusters representing the two Mendelian homozygotes and one heterozygote. Those with inconclusive results formed either one or two genotype clusters, while the bad-quality SNPs did not create any meaningful clusters and were considered failed SNPs ( Table 4 ). Overall, the performance of 14 out of the 17 SNPs was deemed acceptable and selected as candidates for further validation.

The validation exercise involved screening 20 cowpea genotypes in an artificially Striga-infested nursery and comparing the Striga phenotypic data with the marker genotypes on these cowpea genotypes. The phenotypic distribution of the 20 cowpea genotypes in response to Striga infestation is presented in Figures 5A, B . The resistant cowpea genotypes had zero Striga emergence, as portrayed by zero mean and median. In contrast, the susceptible cowpea genotypes showed dispersed distributions with the mean and median Striga count clearly above zero ( Figures 5A, B ). Consequently, the counts of Striga emergence as influenced by the 20 cowpea genotypes were significantly different, with each genotype recording varying mean Striga count ( Table 5 ). Five cowpea genotypes were entirely immune to Striga (no Striga emerged in these plots), while the other genotypes differed in the average number of Striga ranging from 3-10.5, indicating different response levels to Striga ( Table 5 ). As expected, an old variety named Achishiru registered the highest number of Striga, followed by the landraces TVu-867 and Vita7, among the susceptible genotypes compared to low Striga counts registered in some of the released varieties like IT07K-318-33, IT98K-1111-1, IT08K-150-12, IT89KD-288 and IT98K-589-2 ( Table 4 ), suggesting that these are Striga tolerant genotypes.

Using a chi-square test of independence, the study compared the association between each of the 14 candidate SNPs with the phenotypic status of the 20 cowpea genotypes ( Table 6 ). Five SNPs in close proximity to the major association signal on chromosome Vu11 displayed highly significant deviations (P ≤ 0.01) from the null expectation of no marker–phenotype correlation, supporting their strong association with genomic regions influencing Striga resistance in cowpea ( Table 6 ). Five other SNPs positioned close to each other on chromosome Vu02 showed moderate statistical significance (P ≤ 0.05). The rest of the SNPs on Chromosomes Vu07 and Vu10 revealed no significant chi-square values ( Table 6 ).

Further efforts to group the phenotypic reactions by marker alleles revealed a strong marker-phenotype association for the candidate markers on chromosome Vu11 and a weak relation for the other SNPs ( Figure 6 , Supplementary Figures 1A–D ). For instance, the GG allele of the marker snpVU00077 on Vu11 was responsible for 100% of the cowpea genotypes that had a Striga-resistant phenotype, while the alternative allele CC was carried by 94% of phenotypically susceptible genotypes, with only 6% phenotypic misclassification ( Figure 6 ). This observation was supported by a highly significant chi-square test and high Cramer’s correlation ( V ^ _Cramer =0.86). Similar results were observed for all five proximal SNPs on Vu11 ( Supplementary Figure 1D ), suggesting that these markers are strongly associated with genes underlying Striga resistance in cowpea. However, a relatively high percentage of misclassifications and V ^ _Cramer less than 0.5 were observed for the remaining candidate markers on chromosomes Vu02, Vu07, and Vu10 ( Figure 6 , Supplementary Figures 1A–C ), suggesting weak marker-trait association.

Similarly, a t-test was conducted to compare the differences in the mean Striga count between cowpea genotypes carrying the two alleles of each marker. Results for four markers representing the candidates on chromosomes Vu02, Vu07, Vu10, and Vu11 are presented in Figure 7 . Highly significant differences (p ≤ 0.001) between the mean of the allelic groups were found only for the markers on chromosome Vu11, while the other markers were not significant, strongly supporting the chi-square test results. The marker alleles conferring resistance to Striga had significantly lower mean Striga counts than the alternative alleles. The results of the t-test for all 14 markers specifying the alleles underlying resistance and susceptibility and how they differ have been presented in Supplementary Figures 2A–D .

The present study deciphered genomic regions controlling S. gesnerioides resistance in cowpea by utilizing a high-density single nucleotide polymorphism marker and diverse cowpea minicore genotypes. The phenotypic data from the two test sites revealed significant genetic variation in resistance to Striga gesnerioides with skewed distribution, indicating that major-effect genes likely govern resistance. The finding aligns with prior studies, which suggested that resistance to S. gesnerioides in cowpea is often controlled by specific resistance (R) genes exhibiting major effects. Notably, the cowpea genotypes B301 and IT82D-849 exhibited resistance controlled by a single dominant gene effective against Striga races SG1, SG2, and SG3 found in Nigeria (Singh and Emechebe, 1990; Atokple et al., 1993; Li et al., 2009). In a cross of cowpea genotypes HTR (from Niger) and Wango-1 (from Burkina Faso), it was reported that resistance to S. gesnerioides race SG1 in HTR was controlled by one or two dominant genes that are non-allelic to the genes in B301 and IT82D-849 (Timko et al., 2007; Li et al., 2009). In addition, monogenic but recessive inheritance of resistance to S. gesnerioides race SG3 from Niger was also reported (Touré et al., 1997). However, the interaction between cowpea and S. gesnerioides is complex due to the existence of multiple parasite races. Seven distinct races (SG1 to SG7) were identified, each capable of overcoming specific resistance genes in cowpea (Li et al., 2009; Ohlson and Timko, 2020). This variability necessitates discovering key loci involved, followed by stacking multiple resistance genes to achieve broad-spectrum and durable resistance.

A genome-wide association analysis revealed a major association signal on chromosome Vu11, alongside f significant minor signals on chromosomes Vu02, Vu07, and Vu10 linked to Striga resistance. The prominent Vu11 signal aligns closely with the a coiled-coil nucleotide-binding site leucine-rich repeat (CC-NBS-LRR) protein, encoded by the RSG3–301 gene, was implicated in race-specific resistance to S. gesnerioides in cowpea (Li and Timko, 2009). This proximity suggests that similar R-genes may underlie the Vu11 association signal.​ Additionally, a BLAST analysis using primers from earlier marker technologies (e.g., AFLP, SSR, and SCAR) returned best matches within the same Vu11 region ( Table 7 ), further supporting the colocalization of historical markers with our GWAS-identified signals. The previously reported genes spanning the same area on chromosome Vu11 included Rsg3 (conferring resistance to Striga race SG3 and SG5), Rsg2–1 (effective against Striga race SG1), and Rsg4–3 (effective against Striga race SG3) (Ouédraogo et al., 2001; Essem et al., 2019). Earlier linkage mapping studies identified resistance loci to Striga in cowpea via bi-parental populations. Ouédraogo et al. (2002) mapped a race-specific resistance gene (designated Rsg3 and Rsg994) to linkage group 1 (equivalent to chromosome Vu10 in the reference genome), conferring resistance to Striga race SG1. Boukar et al. (2004) mapped another resistance locus, Rsg1, to linkage group 1 (chromosome Vu08), conferring race-specific resistance to SG3 The studies above used traditional linkage mapping in F₂ bi-parental populations to identify the Striga resistance loci. In comparing these classical mapping results with our GWAS findings, we observed a strong correspondence between previously reported loci and our newly identified association signals (see Table 7 ). This concordance supports the involvement of these genomic regions, now validated by independent GWAS in mediating Striga resistance.

Through our candidate gene analysis, we identified 20 unique proteins that have been implicated in plant defense and immune response signaling, collectively annotated to 64 cowpea genes spanning the GWAS-identified association intervals. Notably, among these is a leucine-rich repeat (LRR)-containing protein mapped near the major association peak on chromosome Vu11.LRR domains are well-established as pivotal in pathogen recognition and activation of plant defense responses, as seen in both extracellular receptor kinases and intracellular NB-LRR immune receptors (Jones and Jones, 1997; McHale et al., 2006). An earlier study in cowpea pinned leucine-rich repeat (CC-NBS-LRR) protein to a gene-for-gene resistance mechanism in the interactions between Striga and cowpea, with a corresponding gene in cowpea named RSG3-301 (Li and Timko, 2009). The authors also unlocked the hypothesis of race specificity resistance by silencing the RSG3–301 gene in cultivar B301, which rendered it susceptible to race RG3 but remained resistant to races SG2 and SG5. Another gene PTHR22952:SF183 - TRANSCRIPTION FACTOR TGA5 on Vu11, has been found to work along with TGA2 and TGA6, playing a crucial role in plant systemic acquired resistance (SAR), a broad-spectrum defense mechanism induced after a local infection by avirulent pathogens (Zhang et al., 2003). On chromosome Vu10, the likely genes involved in plant defense were ETHYLENE-RESPONSIVE TRANSCRIPTION FACTOR ERF003, WRKY DNA-binding domain, and PF03763 - Remorin, C-terminal region (Remorin_C). Ethylene-responsive transcription factors (ERFs) are integral components of the APETALA2/ERF superfamily, playing pivotal roles in regulating plant responses to biotic and abiotic stresses (Müller and Munné-Bosch, 2015). ERFs regulate molecular response to pathogen attack by binding to specific cis-acting elements in the promoters of stress-responsive genes, thereby modulating their expression (Müller and Munné-Bosch, 2015). WRKY transcription factors, on the other hand, have been associated with responses to various pathogens in cowpea, particularly exhibiting differential expressions when the plants were challenged with Fusarium oxysporum, a pathogen responsible for Fusarium wilt (Hao et al., 2024). This finding underscores the significant role of WRKY genes in cowpea’s defense mechanisms against biotic stresses. Given the conserved nature of WRKY transcription factors across plant species, it is plausible that manipulating specific WRKY genes in cowpea could enhance resistance to Striga. However, targeted studies are required to identify which WRKY genes are involved and to elucidate their mechanisms of action in the context of Striga resistance. In the context of plant defense, remorins have been observed to interact with receptor-like kinases (RLKs) and pathogen effectors, suggesting a role in the early stages of immune signaling (Yu, 2020). Chromosome Vu07 harbored genes encoding MADS-box transcription factors, cysteine-rich receptor-like kinases, MYB-like proteins, and GATA-type zinc finger transcription factors, all of which play roles in plant development and stress responses (Abdullah-Zawawi et al., 2021; Zhang et al., 2023a, b). The genes on chromosome Vu02, included, among others, pentatricopeptide repeat (PPR) proteins, C2H2-like zinc finger proteins, glucuronosyltransferase PGSIP8-like, receptor-like kinases, and lipid transfer proteins that contribute to the complex network of plant defense mechanisms, each playing distinct roles in responding to biotic and abiotic stresses (Goff and Ramonell, 2007; Finkina et al., 2016; Zhang et al., 2020; Huang et al., 2022; Meng et al., 2024).

Seventeen markers tagging GWAS-identified association intervals for Striga resistance were converted into KASP assays and independently validated across 20 cowpea genotypes with known resistance profiles. Five of these markers, all situated near the major association signal on chromosome Vu11, were confirmed to be associated with resistance. KASP is a globally recognized technology for SNP genotyping, known for being user-friendly, relatively cheap, and can be automated, ensuring high throughput genotyping (Ongom et al., 2021; Dipta et al., 2024). The five Striga-associated SNP markers designated by Intertek KASP assay IDs; snpVU00075, snpVU00076, snpVU00077, snpVU00078, and snpVU00079 are all positioned proximally to each other within the major QTL region on Vu11. Previously reported marker systems for Striga resistance were mostly AFLPs, SCARs, and SSRs (Ouédraogo et al., 2001, 2002; Boukar et al., 2004; Essem et al., 2019). Although some of these old marker technologies were reportedly effective in marker-assisted selection (MAS) for Striga resistance (Larweh et al., 2017; Omoigui et al., 2017), they have limitations, including labor-intensive gel electrophoresis or capillary electrophoresis, lack of automation and cumbersomeness for high-throughput applications, high cost per data point, and low reproducibility and accuracy (Mut et al., 2008; Dipta et al., 2024). The present study underscores the potential of new KASP-based SNP markers for efficient MAS in cowpea breeding programs aimed at developing Striga-resistant varieties.