Comparative genomic analysis reveals the difference of NLR immune receptors between anthracnose-resistant and susceptible sorghum cultivars

Abstract

Anthracnose, caused by Colletotrichum sublineola, is a prevalent disease that noticeably affects global sorghum production. Nucleotide-binding leucine-rich repeat receptors (NLRs) are crucial for disease resistance. Here, we report the differences in the number, expression profile, and gene structure of NLRs between the anthracnose-resistant and susceptible sorghum cultivars. Through a systematic anthracnose disease assay on 365 sorghum accessions, we identified the American improved cultivar BTx623 as the resistant and the Chinese improved glutinous cultivar Guojiaohong1 (GJH1) as the susceptible cultivar. Then we sequenced the genome of GJH1 and identified 239 NLRs, substantially fewer than the 302 in BTx623. Although the collinear NLRs are highly conserved between GJH1 and BTx623, more than half of the non-collinear NLRs showed notable mutations or structural variations. During C. sublineola infection, BTx623 exhibited a higher number of highly expressed and inducible NLR genes than GJH1 did. Moreover, we identified some candidate anthracnose resistance genes that are potentially valuable for disease-resistant breeding. Therefore, our data provide genetic resources for developing disease-resistant glutinous sorghum.

Background

Sorghum [Sorghum bicolor (L.) Moench] is the fifth most widely produced cereal worldwide (Khoddami et al. 2023), with an annual yield exceeding 68 million tons across approximately 44 million hectares (Zheng et al. 2011; Abreha et al. 2021). Despite its agricultural significance, sorghum productivity is severely constrained by anthracnose, a fungal disease caused by Colletotrichum sublineola. This pathogen can cause yield losses up to 36% in susceptible genotypes, and as high as 86% in certain inbred lines under optimal conditions (Cota et al. 2017; Fu et al. 2020). The impact of anthracnose is particularly severe in regions such as southwest China, where the high temperatures and humidity exacerbate the susceptibility of local glutinous sorghum. This increased vulnerability results in reduced yields, diminished quality, and significant challenges for local Baijiu-brewing industry.

The defense mechanisms of sorghum against pathogens are categorized into two main types: pathogen-associated molecular pattern-triggered immunity (PTI) and effector-triggered immunity (ETI) (Sun et al. 2020; Abreha et al. 2021; Ngou et al. 2022). PTI is a basal immune system that responds to common pathogenic patterns, such as bacterial flagellin and lipopolysaccharides. In contrast, ETI represents a more specialized defense mechanism mostly involving interactions between plant nucleotide-binding leucine-rich repeat (NLR) proteins and pathogen avirulence (Avr) effectors (Narusaka et al. 2009; Li et al. 2016; Saur et al. 2019; Van de Weyer et al. 2019). Although ETI often triggers a hypersensitive response to restrict pathogen spread, the co-evolutionary dynamics between hosts and pathogens frequently result in the breakdown of this defense. Consequently, there is a critical and ongoing need to identify novel NLR genes to enhance and sustain plant resistance (Tamborski and Krasileva 2020; Barragan and Weigel 2021).

The deployment of NLR genes in crop breeding has notably enhanced their disease resistance (Liu et al. 2020; Sun et al. 2020; Martin et al. 2022). For example, the rice NLR gene Pi9 confers broad-spectrum resistance against Magnaporthe oryzae, the causative agent of rice blast disease. Its introduction into elite cultivars has greatly improved rice yield in disease-prone areas. Similarly, the soybean NLR gene Rpp1 provides resistance to Asian soybean rust (Phakopsora pachyrhizi), a global threat to soybean production (Wei et al. 2023). In potatoes, the NLR gene Rpi-blb1, originally identified in a wild Solanum species, has been pivotal in breeding varieties resistant to late blight caused by Phytophthora infestans (Xie et al. 2015). Recent advances in gene-editing technologies have further expanded the potential of NLR proteins. For example, targeted editing of the rice NLR gene Pik-1 has enhanced resistance to M. oryzae without compromising yield or quality (Białas et al. 2021). These successes underscore the critical role of the NLR genes in addressing escalating disease pressures and ensuring food security. Notably, recent studies have identified four NLR genes that confer resistance to sorghum anthracnose (Lee et al. 2022; Mewa et al. 2022; Habte et al. 2023), highlighting their importance in managing sorghum disease. However, the differences in NLR genes between resistant and susceptible sorghum cultivars remain to be elucidated.

In this study, we conducted anthracnose disease assay on 365 sorghum accessions and identified BTx623 and Guojiaohong1 (GJH1) as the representative resistant and susceptible genotypes, respectively. BTx623, a modern non-glutinous variety developed in the United States, offers the most comprehensive genome assembly and annotation for sorghum. In contrast, GJH1 is a renowned glutinous sorghum cultivar widely used for Baijiu-brewing in southwest China. Through a systematic comparative analysis, we identified candidate NLR proteins with potential for breeding programs aimed at enhancing anthracnose resistance in glutinous sorghum.

Results

Identification of sorghum accessions resistant and susceptible to anthracnose

To identify anthracnose-resistant germplasm, we conducted spray and punch inoculation disease assays on 365 sorghum accessions, classified into six subpopulations (Pop1-Pop6) (Additional file 1: Table S1). Our results revealed that improved American grain sorghum accessions (Pop1) exhibited the highest resistance, followed by wild accessions (Pop5) and sweet accessions (Pop6). In contrast, glutinous sorghum accessions (Pop3) were the most susceptible to anthracnose, and Southwestern landrace sorghum (Pop4) and Northern grain accessions (Pop2) showing moderate susceptibility (Fig. 1a and Additional file 1: Table S1). Spray inoculation grading results further revealed that a higher proportion of Pop1 accessions were highly resistant, whereas Pop3 accessions were predominantly susceptible (Fig. 1b and Additional file 1: Table S1). Therefore, we selected BTx623 and GJH1 as the representative resistant and susceptible accessions, respectively, due to the availability of comprehensive genomic information for BTx623 and GJH1 as a widely used glutinous sorghum cultivar for Baijiu-brewing in Southwest China.

Fig. 1

Anthracnose resistance assessment of 365 sorghum accessions. a Mean incidence disease area of punch inoculation among six sorghum subpopulations. Different letters indicate significant differences at P < 0.05 according to one-way ANOVA followed by Tukey’s HSD test for multiple comparison. b Incidence grade of spray inoculation among six sorghum subpopulations. R, resistance; S, susceptibility. Pop1, American improved grain varieties; Pop2, Improved grain varieties of northern China; Pop3, glutinous sorghum of southwest China; Pop4, landrace sorghum of southwest China; Pop5, wild sorghum accessions; Pop6, sweet sorghum accessions. c BTx623 exhibited more resistance to these six isolates than GJH1. Top panel: GJH1, and Bottom panel: BTx623. The six C. sublineola strains (ZG-FS-3, ZG-DA-3, ZG-DA-20, HB-LC-7, CQ-YC-6, and YB-NX-1) were isolated across China. d Quantification of sorghum leaf lesions after 7 dpi following infection with six strains. Disease lesion areas were determined from 6 cm leaf sections using Image J and analyzed according to one-way ANOVA followed by Tukey’s HSD test for multiple comparison (P < 0.05). Values are means ± s.d. calculated from four biological replicates. Different letters above bars indicate significant differences. **P < 0.01; ***P < 0.001

To further confirm the disease phenotypes of BTx623 and GJH1, we conducted punch inoculation assays with six C. sublineola strains collected from diverse sorghum-producing regions across China. Disease lesions were significantly smaller in BTx623 than in GJH1 for all tested strains (Fig. 1c, d). These results underscore the markedly higher anthracnose resistance of BTx623 compared to GJH1.

Characterization of sorghum NLR genes

To compare NLR genes between BTx623 and GJH1, we sequenced, assembled, and annotated the GJH1 genome. These data have been deposited in the Genome Warehouse at the National Genomics Data Center under accession number PRICA026490. To determine the number of NLRs in GJH1 and BTx623 genome, we performed HMMSEARCH using the NB-ARC Pfam model PF00931, followed by manual curation. We identified 239 NLR genes in GJH1, whereas we identified 302 NLRs in the genome of the anthracnose-resistant cultivar BTx623 (Additional file 1: Table S2), a number higher than a previous annotation (Yang and Wang 2021). This finding indicates a substantial genetic variation in the NLR gene family between the two cultivars.

Then, we analyzed the physicochemical properties of the 302 NLR proteins in BTx623 using ProtParam, including protein length, molecular weight (MW), grand average of hydropathicity (GRAVY), and isoelectric point (pI) (Additional file 1: Table S2). The NLR proteins were generally large, averaging 1022.46 amino acids, with the longest (Sobic.005G192400) comprising 2730 amino acids. The GRAVY values ranged from 0.721 (Sobic.002G004500) to 0.082 (Sobic.003G329300). The pI values ranged from 4.91 to 9.39, with the lowest and highest values observed in Sobic.005G192400 and Sobic.003G279700, respectively. The MW ranged from 15.8 (Sobic.009G093301) to 271.55 kDa (Sobic.010G205600).

Based on the typical NLR domains found in Poaceae, i.e., CC, NBS, and LRR, we classified the sorghum NLRs into four categories: CC-NBS (CN), CC-NBS-LRR (CNL), NBS-LRR (NL), and solely NBS (N). Most NLRs belonged to the CNL group, totaling 187, followed by the CN and NL group, with 62 and 35 proteins, respectively. Eighteen proteins consisted solely of NBS domains (Fig. 2a). Among the 302 NLRs, 20 contained at least one atypical NLR domain or integrated domain (ID), representing 13 distinct Pfam domains, such as Pkinase_Tyr, WD40, FNIP, and WRKY (Fig. 2b). Both WD40 and FNIP were found in tandem duplicated NLRs. These findings underscore the rich diversity of NLR genes in sorghum.

Fig. 2

Classification and chromosome distribution of NLR genes in S. bicolor. a Number of NLR genes in different structural classes. CNL, CC-NBS-LRR; CN, CC-NBS; N, NBS; NL, NBS-LRR. b Number of NLR genes based on the inclusion of integration domains. c The 302 NLR genes identified in BTx623 were mapped to the ten chromosomes of S. bicolor. Chromosome numbers are indicated at the top of each bar, and the sizes of chromosomes are represented by the vertical scale in megabases (Mb). d Number of clustered, paired, and dispersed NLR genes in the genome of S. bicolor. Clustered NLR genes are those located within 200 kb of each other in the genome. Paired NLR genes are a subset arranged in a head-to-head orientation, positioned oppositely and transcribed divergently

Next, we generated detailed chromosomal maps of the NLR genes using the high-quality genome annotation of BTx623 because such an information is not available in the literature. The NLR genes were unevenly distributed across the 10 chromosomes (Fig. 2c). Chromosome 5 contained the most NLR genes (98 NLRs, 32.45%), followed by chromosomes 2 and 8, with 46 and 44 NLR genes, respectively. In contrast, chromosome 4 had the fewest, with only 7 NLR genes (Fig. 2c). NLR genes were predominantly located on the chromosome arms and telomeric regions, with fewer in the pericentromeric and centromeric regions.

To further analyze NLR clustering, we defined an NLR cluster as a group of genes located within 200 kb of each other on the genome. Using this definition, we identified 213 NLR genes within such clusters (Fig. 2d). We also observed paired NLR, which refer to two adjacent NLR genes that often work together to detect pathogen effectors and activate immune responses (Saur et al. 2019; Narusaka et al. 2009). In BTx623 genome, we identified 11 pairs of such NLR genes, with five pairs located on chromosome 5 and four pairs on chromosome 2.

To clarify the evolutionary relationships of the NLRs, we used IQ-TREE (v. 2.2.0) to construct a maximum likelihood (ML) evolutionary tree of the 302 sorghum NLR proteins. The phylogenetic tree divided the NLR proteins into six distinct groups (Fig. 3). Groups I and II had fewer members, with 24 and 25 NLR proteins, respectively. In group I, apart from Sobic.008G020900 and Sobic.008G112550, all NLR proteins contained the three typical structural domains: CC, NBS, and LRR. This group showed a close phylogenetic relationship with the rice disease resistance protein RPM1 (Additional file 1: Table S3). The anthracnose-resistant NLR ARG1 was located in group II, within the same clade as the rice NLR protein Pit (Additional file 1: Table S3). Groups III (57 members) and VI (56 members) contained proteins with a tyrosine kinase (Pkinase_Tyr) domain, showing similarity to rice NLR proteins Pi56 and Pi5. Group IV included a high proportion of NLR proteins with incomplete domain architectures, highly homologous to rice NLR proteins RGA2, Pi37, and Pish (Additional file 1: Table S3), and included the sorghum NLR protein Sobic.008G166400 (ARG4) and Sobic.008G177900 (ARG5) (Fig. 3). Group V had the largest number of NLR proteins (77 members), followed by group VI with 63 members, suggesting that gene duplication evens likely contributed to their expansion. Group V also included the NLR protein ARG2, which confers anthracnose resistance in sorghum (Mewa et al. 2022). Additionally, groups V and VI contained NLR genes homologous to rice blast resistance proteins, such as Pigm, SCR8, Pik, and RGA5 (Fig. 3 and Additional file 1: Table S3).

Fig. 3

Phylogenetic relationship of the 302 NLR proteins. Each clade in the phylogenetic tree is represented by a different color. Red arrowheads indicate the anthracnose-resistant NLR genes that have been cloned in sorghum: Sobic.007G085400 (ARG1), Sobic.005G047700 (ARG2), Sobic.008G166400 (ARG4), and Sobic.008G177900 (ARG5)

Previous studies have demonstrated that NLRs contain 20 conserved motifs, four of which are specific to TIR-type NLRs (TNLs) (Jupe et al. 2012). Using MEME, we identified 16 conserved motifs among the 302 NLR proteins, each occurring more than 20 times and ranging from 11 to 40 amino acids in length (E-values < 2.3e−230) (Additional file 2: Figure S1). Motif 1 appeared 1000 times, followed by motif 12, which appeared 636 times, and both are conserved in the LRR domain. Motifs 6 and 16, which are conserved in the CC domain, appeared 251 and 154 times, respectively. This number also closely matched the 249 predicted NLR proteins containing CC domain. The NBS domain was the most conserved across the NLR proteins. We identified 11 conserved motifs within the NBS domain, including p-loop (motif 2), RNBS-A (motif 14), Walker-B (motif 7), RNBS-B (motif 9), RNBS-C (motif 5), GLPL (motif 8), RNBS-D (motif 4), and MHDV (motif 13). Additionally, three unreported motifs (motifs 3, 10, and 11) were relatively conserved, appearing 54, 260, and 242 times, respectively (Additional file 2: Figure S1). By integrating evolutionary tree with gene structure, we observed that each subgroup had similar motif compositions and exon–intron organization. This provided additional support for the taxonomic classification and phylogenetic relationships of the NLRs (Additional file 2: Figure S2).

Collinearity and genetic variation of the NLR genes between BTx623 and GJH1

To analyze the collinearity between the genomes of BTx623 and GJH1, we identified 27,177 genes that display strict collinear relationships, constituting 79.63% of the total genes (Additional file 1: Table S4). Notably, large segment inversions were identified near the centromeres of chromosomes 5 and 7 (Additional file 2: Figure S3). A recent study attributed the inversions to assembly errors in BTx623 (Wei et al. 2024). This finding further emphasizes the high-quality genome assembly of GJH1. Then, we identified 195 NLR genes with strict collinearity between BTx623 and GJH1, with 70.85% of them showing identical coding regions (Fig. 4a and Additional file 1: Table S4). In contrast, the 103 non-collinear NLR genes in BTx623 displayed a high proportion (57.94%) of substantial mutations, including large effect mutations, structural variations, and presence-absence (Table 1). Altogether, 113 NLR genes in BTx623 exhibited major genetic variation compared to those in GJH1 (Table 1). The substantial variation in NLR genes between the two cultivars highlights their potential for identifying NLR genes conferring anthracnose resistance.

Fig. 4

Duplication events of NLR Genes in S. bicolor. a Collinearity analyses of NLR genes between BTx623 and GJH1. Red and blue represent BTx623 and GJH1 chromosomes, respectively. Gray lines in the background show the collinear relationship of the whole gene pairs, while the purple line highlights the collinear NLR gene pairs between BTx623 and GJH1. b Collinearity blocks containing whole genome duplication gene pairs, with red lines indicating NLR gene pairs. c Ks values of different types of duplicated NLR gene pairs. d Ka/Ks values of different types of duplicated NLR gene pairs

Table 1 Variations within NLR family genes between BTx623 and GJH1 genomes

Duplication of the NLR genes in sorghum

To explore the duplication events of NLR genes in sorghum, we used the MCScanX tool to categorize their duplication types. We successfully identified two pairs of NLR genes from whole-genome duplication (WGD) events: Sobic.005G219300 paired with Sobic.008G131900, and Sobic.003G329300 paired with Sobic.009G180900 (Fig. 4b). We also identified 90 NLR genes in tandem duplications and 82 in proximal duplications, consistent with the cluster formation patterns typically observed in NLR genes (Additional file 2: Figure S4). This pattern suggests that asymmetrical crossing over during meiosis is the primary driver of NLR gene expansion. Other NLR genes were predominantly found in dispersed duplications, likely driven by transposon duplication and natural selection. Only four singleton NLR genes were identified (Additional file 2: Figure S4).

Next, we calculated the synonymous substitution rates (Ks) of the duplicated gene pairs to estimate their divergence time. NLR genes duplicated through WGD exhibited high Ks values (Fig. 4c), aligning with their origin from an early genome duplication event in the Poaceae. In contrast, most NLR gene pairs arising from tandem and proximal duplications showed low Ks values, indicating the recent duplication events (Fig. 4c). To understand selection pressures on NLR genes, we calculated the Ka/Ks ratios of the duplicated NLR gene pairs. Only two tandem-duplicated NLR gene pairs had Ka/Ks ratios greater than 1, indicating prominent positive selection during domestication (Fig. 4d). Among them, AGR4 was identified as a key gene for resistance to anthracnose. The Ka/Ks ratios of other NLR gene pairs were less than 1, suggesting that most of the NLR genes evolved under purifying selection (Fig. 4d). These results suggest that purifying selection plays pivotal roles in driving adaptation and resistance during sorghum domestication.

Expression of the NLR genes in response to C. sublineola

To explore the differences in the expression of NLR genes in GJH1 and BTx623 in response to C. sublineola, we analyzed transcriptomic data collected at 0-, 12-, 24-, and 48-h post-inoculation (hpi) (Fig. 5a and Additional file 1: Table S5). In BTx623, 105 NLR genes were highly expressed, including 53 that were expressed at very high levels, compared to 95 and 42 in GJH1, respectively (Fig. 5b). This indicates a stronger basal expression of NLR genes in BTx623 than in GJH1. At 24 and 48 hpi, the number of induced NLR genes in BTx623 reached 100 and 97, respectively, surpassing the 91 and 22 in GJH1. Among them, 12 NLR genes were significantly induced at both 24 and 48 hpi of the two materials (Fig. 5c). These findings imply differential response between BTx623 and GJH1, and a sustained and robust defense responses in BTx623, underpinning its resistance to C. sublineola.

Fig. 5

Transcriptional profile analysis during S. bicolor infection by C. sublineola. a Comparison of sorghum anthracnose symptoms at 5 DPI using spray infection on sorghum cultivars BTx623 and GJH1. The former is disease-resistant, while the latter is susceptible. b Number of NLR genes with different expression levels in the two cultivars. Very high expression: FPKM > 10; High expression: 1 < FPKM < 10; Low expression: 0.1 < FPKM < 1; No expression: FPKM < 0.1. c Venn diagram comparing differentially responsive NLR genes to C. sublineola infection in BTx623 and GJH1 at 24 dpi and 48 dpi. d Temporal co-expression analysis revealing 14 modules of sorghum gene expression during C. sublineola infection. e Bioinformatic analysis of the transcriptional regulatory networks of NLR responses to early C. sublineola infection. Red squares represent NLRs, and green circles represent transcription factors (TFs). Green arrowheads indicate the key TFs

To further understand the defense response of BTx623 to C. sublineola, we conducted a weighted gene co-expression network analysis (WGCNA) to identify key NLR genes associated with the infection. The analysis grouped all expressed sorghum genes into 14 modules. Among these, the royalblue and black modules were strongly associated with the C. sublinola infection at 24 and 48 hpi, respectively. The royalblue module included one NLR gene, whereas the black module contained 42 NLR genes (Fig. 5d and Additional file 1: Table S6). Using these genes and their corresponding transcription factors (TFs), we constructed a core co-expression regulatory network that is responsive to C. sublineola infection. Key TFs such as Sobic.001G120900 and Sobic.009G157500, along with crucial NLR genes like Sobic.005G047700 (ARG2) and Sobic.005G218400, were identified as central components of this network (Fig. 5e).

In addition to genetic variation and differential responses to C. sublineola, notable differences exist in NLR gene expression between BTx623 and GJH1 (Fig. 6a). A detailed comparison revealed that 46 NLR genes in BTx623 exhibited significantly higher expression than their counterparts in GJH1, including 25 NLR genes that were either unexpressed or expressed at very low levels in GJH1. Examples include Sobic.005G183500, Sobic.006G050700, Sobic.005G167500, and Sobic.005G167400 (Fig. 6a). This distinct expression profile suggests the presence of unique NLR genes in BTx623 that may enhance its disease resistance. Among the 46 highly expressed genes in BTx623, 30 responded to C. sublineola infection, further implicating their role in anthracnose resistance. Conversely, 18 NLR genes exhibited higher expression in GJH1 than in BTx623 (Fig. 6a), implying their differential roles in response to C. sublineola.

Fig. 6

Differential expression profiles and RT-qPCR analysis of key NLR genes in BTx623 and GJH1. a Validation of RNA-seq data through qRT-PCR analyses of 12 selected NLR genes in response to C. sublineola. Error bars represent mean ± SE (n = 3). Asterisks denote significant differences: *P ≤ 0.05; **P ≤ 0.01; ***P ≤ 0.001. b Differential expression of NLR genes at different infection stages between BTx623 and GJH1. Red color represents high expression of the NLR genes, and white color represents low-level expression at the corresponding infection times

To validate the 12 NLR genes expression profiles that were notably induced at 24 and 48 hpi in both GJH1 and BTx623, we conducted qRT-PCR (Fig. 6b). The qRT-PCR results closely mirrored the transcriptomic data, confirming their responses to C. sublineola infection. The NLR genes such as Sobic.008G177900, Sobic.005G053900, and Sobic.010G205700 showed the highest expression. Specifically, the expression of Sobic.008G177900, Sobic.008G133200, and Sobic.005G053900 were significantly higher in BTx623 than in GJH1, whereas Sobic.009G138000 and Sobic.010G205700 were more highly expressed in GJH1 than in BTx623 (Fig. 6b). Additionally, Sobic.007G085400 (ARG1) and Sobic.008G177900 (ARG5), two previously identified NLR genes involved in anthracnose resistance, displayed similar expression patterns with previous reports (Lee et al. 2022; Habte et al. 2023). These results further support their important role in sorghum anthracnose resistance.

Candidate NLRs for anthracnose-resistance breeding in Chinese glutinous sorghum

To predict the candidate NLRs that potentially contribute to anthracnose resistance, we selected those with the genetic variation and differential responses to C. sublineola between BTx623 and GJH1. We identified a cohort of NLR candidates that may play potential roles of anthracnose resistance in Chinese glutinous sorghum (Fig. 7). Among them, 14 NLR genes exhibited genetic variation, were induced by C. sublineola infection, and showed higher expression in BTx623 than in GJH1. Notably, five of these NLR genes overlap with previously reported genes or QTLs associated with anthracnose resistance, including ARG5, Sobic.009G013300, Sobic.009G013100, Sobic.005G075600, and Sobic.005G076100. Thirteen NLR genes displayed genetic variation and differential expression between BTx623 and GJH1, whereas 16 NLR genes were specifically induced by C. sublineola infection and differentially expressed between the two cultivars. Additionally, 12 NLR genes with genetic variation were also induced by C. sublineola infection, including AGR2. The 52 NLR genes that responded to C. sublineola infection may also contribute to anthracnose resistance, including ARG1 and AGR4 (Fig. 7). These candidate NLR genes offer opportunities for further functional characterization.

Fig. 7

The candidate NLR gene for anthracnose resistance. A comprehensive analysis identified a total of 14 NLR genes with genetic variation, response to C. sublineola, and differential expression between BTx623 and GJH1, 13 NLR genes with genetic variation and differential expression between BTx623 and GJH1, 16 NLR genes responsive to C. sublineola and differentially expressed between BTx623 and GJH1, and 12 NLR genes with genetic variation and response to C. sublineola. Red indicates NLR genes that have been cloned or overlap previously identified QTLs associated with anthracnose resistance, blue represents the anthracnose resistant NLR gene that has been cloned in sorghum

Discussion

NLR genes have been extensively studied in various species, including Arabidopsis, rice, and wheat (Wang et al. 2019; Van de Weyer et al. 2019; Dinh et al. 2020; Li et al. 2021). These genes play a critical role in breeding disease-resistant crops by providing effective defense against various pathogens. Here, we compared and analyzed the differences in NLR genes between BTx623 and GJH1. BTx623 has been reported to be susceptible to several C. sublineola isolates in the U.S. (Lee et al. 2022; Habte et al. 2023). However, our results indicate that BTx623 exhibits resistance to six virulent isolates collected from China. This observation aligns with a previous report that BTx623 showed the resistance to C. sublineola isolates from China (Xu et al. 2020). We speculate that this discrepancy is presumably due to the differences between U.S. and Chinese strains. In this study, we presented the genome-wide repertoire of NLR differences between BTx623 and GJH1, which could aid in further identification of the resistance genes for disease-resistance breeding of glutinous sorghum in Southwestern China.

The question of whether resistant cultivars possess more NLRs than susceptible ones remains debated (Jupe et al. 2012; Thatcher et al. 2023). In this study, we identified 302 NLRs in BTx623 and 239 in GJH1, suggesting that a higher number of NLRs may contribute to disease resistance in sorghum. This observation aligns with a recent report indicating that a higher proportion of NLR genes in plant genomes often mobilizes substantial resources to sustain defense mechanisms (Giolai and Laine 2024). In addition to the number of NLR genes, we also identified substantial genetic differences between the two varieties by collinearity analysis (Table 1). The non-collinear NLR genes exhibit high rates of significant mutations, structural variations, and presence-absence variations (Table 1). Chinese sorghum, known as 'Kaoliang', has undergone distinct evolutionary trajectories and artificial domestication compared to U.S. grain sorghum (Cisse and Ejeta 2003; Zhang et al. 2023). Therefore, the modern varieties from the U.S. can be used to improve the resistance of glutinous sorghum in China.

The NLR genes, which are often clustered in genome (Fig. 2d), enable the development of new resistance specificities through dynamic processes like gene duplication and unequal crossing over (Huang et al. 2022). Understanding the factors that influence NLR genes diversity within sorghum is essential for developing effective strategies to enhance disease resistance (Arun et al. 2017). The higher number of NLR genes in sorghum compared to maize suggests a significant expansion during domestication, primarily driven by historical tandem and proximal duplication events (Fig. 4 and Additional file 2: Figure S4). This explains why most novel NLR genes often emerged following the whole-genome duplication, particularly during the domestication process (Tamborski and Krasileva 2020). During domestication, outcrossing between sorghum varieties of diverse origins and races contributed to the duplication events. Similar patterns of NLR gene duplication have been observed in other crops, such as rice and wheat (Arun et al. 2017; Van de Weyer et al. 2019; Wang et al. 2019; Kim et al. 2021; Wang et al. 2022). For example, the Pit1 and Pit2 genes in rice originated from an ancestral Pit gene through tandem duplication (Li et al. 2024).

In sorghum, purifying selection plays a critical role in eliminating disadvantageous mutations, thereby maintaining the stability and functionality of essential genes. This is reflected in our study (Van de Weyer et al. 2019), where the calculated Ka/Ks ratio of most duplicated NLR gene pairs was less than 1 (Fig. 4c), indicating that negative selection has constrained functional divergence after gene duplication. Despite the extensive genetic diversity within the NLR gene family, amino acid sequences have remained stable through sub-functionalization, a process that preserves essential functions while allowing evolutionary refinement (Cusack and Wolfe 2007). To fully leverage these findings, future research should prioritize the functional validation of candidate NLR proteins using advanced gene-editing technologies like CRISPR/Cas9, elucidating their roles in disease resistance pathways and enabling their application in breeding programs.

Comparative genomics and bioinformatics tools accelerate the characterization of NLR genes in sorghum by leveraging functional predictions, thereby reducing the time and resources needed for experimental validation (Mosquera et al. 2009; Han et al. 2023; Yan et al. 2023). During C. sublineola infection, nearly half of the NLR genes were noteworthily modulated in both GJH1 and BTx623. BTx623 exhibited a more robust response, evidenced by a higher number of induced NLR genes and sustained expression over a longer period, compared to GJH1, particularly at 48 hpi (Fig. 5). This differential expression profile suggests that BTx623 may possess a more effective NLR-mediated defense against C. sublineola than GJH1. Furthermore, analyzing the NLR gene expression profiles can aid in identifying key genes and elucidating their roles in defense against C. sublineola (Denison et al. 2011; Brown and Hudson 2015; Lee and Zhang 2015; Ramos et al. 2021; Xu et al. 2022). Through co-functional network analysis, we identified key NLR genes and transcription factors, including AGR2 and the Pik-2 homolog Sobic.005G218400 (Fig. 6b). These candidate genes represent promising targets for functional validation.

Conclusions

This study provides a panoramic view of NLRs in sorghum, classifying them into six distinct groups based on phylogenetic relationships. The tandem and proximal duplication of NLR genes is the primary driver of their expansion and appears to be largely influenced by out-crossing. By comparing genetic variation, response to C. sublineola, and expression profile between BTx623 and GJH1, we identified candidate NLRs that likely contribute to anthracnose resistance in glutinous sorghum. In conclusion, these candidate NLRs, which exhibit substantial genetic differences between resistant and susceptible cultivars, present promising targets for functional validation in the future.

Methods

Sorghum accessions and anthracnose disease assays

A total of 365 sorghum accessions were classified into six subpopulations based on their genetic relationship, including BTx623 and GJH1, and were grown under controlled greenhouse conditions with temperatures of 28/25°C (light/dark), 12/12 h (light/dark) photoperiod and 70% humidity. A mycelial block was placed in 50 mL of liquid PDA medium and incubated at 28°C with shaking at 180 rpm for 4 days to culture C. sublineola. Plants were inoculated following a standardized protocol using C. sublineola spores at a concentration of 1 × 10⁵ spores/mL, applied at the five-leaf stage (Zhong et al. 2018). Mock-inoculated controls were included for comparison.

Anthracnose resistance was assessed using a grading scale modified from previous reports on sorghum (Xu et al. 2020; Lee et al. 2022). Plants with no visible symptoms or lesions were classified as R + , those with a small number of chlorotic flecks and less than 5% of the leaf area covered by lesions were classified as R. Plants with reddening or hypersensitive spots and no acervuli present, covering between 5 and 20 percent of the leaf area, were classified as R-. Necrotic lesions containing acervuli and covering between 20 and 50 percent of the leaf area were classified as S-, while coalescing necrotic lesions with acervuli covering 50%–80% of the leaf area were classified as S. Finally, plants that were nearly dead with 80%–100% of the leaf area covered by coalescing necrotic lesions were classified as S + .

Six C. sublineola strains (ZG-FS-3, ZG-DA-3, ZG-DA-20, HB-LC-7, CQ-YC-6, and YB-NX-1), isolated from diverse sorghum-producing regions across China, were used to evaluate the anthracnose resistance of BTx623 and GJH1.

Genome-wide analyses of the NLR genes in sorghum

The NLR proteins were identified using Hidden Markov Model (HMM) searches that targeted the conserved NB-ARC domain (Pfam: PF00931). Proteins containing the NB-ARC domain were further analyzed to confirm their classification by querying the Pfam (http://pfam.xfam.org/), CDD (https://www.ncbi.nlm.nih.gov/cdd), and SMART (http://smart.embl.de/) databases (Letunic et al. 2021; Mistry et al. 2021; Wang et al. 2023). The BTx623 genome data, including the sequences of genomic DNA, predicted coding regions, and protein sequences, were retrieved from the Sorghum bicolor v.3.1.1 genome in Phytozome 13 (Paterson et al. 2009; McCormick et al. 2018). The NLR genes identified were then classified based on their domain architecture and chromosomal location. The physicochemical properties of the NLR proteins that were identified, such as protein length, MW, GRAVY, and pI, were analyzed using the ProtParam tool (https://web.expasy.org/protparam/) (Hassan et al. 2022).

The motif discovery tool MEME v5.5.0 (http://meme-suite.org/tools/meme) was used to identify the conserved motifs within the NLR proteins (Buske et al. 2010). The structural features of the genes, including exon–intron organization and conserved motifs, were visualized using TBtools (https://github.com/CJ-Chen/TBtools) (Chen et al. 2020). The MapGene2chromosome v2 (MG2C) software tool (http://www.mg2c.iask.in/mg2c_v2.0/) was used to map the chromosomal distribution of the NLR genes based on their genomic coordinates (Chao et al. 2015).

Phylogenetic analysis of the NLR genes

The phylogenetic relationships of the NLR proteins were elucidated using multiple sequence alignments of full-length NLR sequences with the software MUSCLK (Edgar 2004). Phylogenetic trees were constructed using maximum likelihood methods in IQ-TREE2 (v.2.0.6), with robustness assessed by bootstrap analysis. The trees were visualized using the Interactive Tree Of Life (iTOL) platform (https://itol.embl.de/itol.cgi) (Letunic and Bork 2021).

Analysis of collinearity and gene duplication

The collinearity and gene duplication within the sorghum genome were analyzed using MCScanX (Wang et al. 2024). The types of duplicated genes were categorized based on sequence similarity and physical proximity. The ratio of nonsynonymous to synonymous substitutions (Ka/Ks) was calculated using the KaKs_Calculator to infer evolutionary pressures (Zhang et al. 2022). The divergence times were estimated using the formula, T = Ks/(2 × 6.1 × 10⁻⁹) × 10⁻⁶ million years ago (MYA), which aligned with previous genome duplications in the Poaceae (Kim et al. 2013; Lee et al. 2020).

RNA extraction, RNA-seq, and qRT-PCR analysis

To elucidate the transcriptional response of sorghum NLR genes following infection with C. sublineola, an RNA-seq data analysis was conducted. The total RNA was extracted from leaf tissues collected at various time points of post-inoculation using the TRIzol reagent (BioTeke, Beijing, China). RNA-seq libraries were prepared following established protocols, and the mRNA transcriptome was sequenced on an Illumina HiSeq 2500 platform at Beijing Biomarker Technologies Co., LTD (Beijing, China). HISAT2 alignment software was utilized to quality filter the raw reads and align to the sorghum reference genome. The differentially expressed genes (DEGs) were analyzed using DESeq2 to identify the responsive genes.

The transcriptional data were validated by real-time quantitative reverse transcription PCR (qRT-PCR) that was performed on 12 selected NLR genes that responded to infection with C. sublineola. The primers were designed with Primer5 (Additional file 1: Table S7). The cDNA was synthesized from 2 μg of total RNA using NovoScript Plus All-in-one 1st Strand cDNA Synthesis SuperMix (E047-01A, Suzhou, China). The sorghum SbEIF4a gene served as an internal control to normalize the data for expression, which ensured the accuracy of quantitative analyses (Du et al. 2022).

Statistical analysis

A one-way analysis of variance (ANOVA) was performed using SPSS 20.0 (IBM, Inc., Armonk, NY, USA), followed by Tukey's HSD test for multiple comparisons. Significant differences were defined at P < 0.05.