Genome-wide identi cation and characterization of NBS-encoding genes in Raphanus sativus L. and their roles related to Fusarium oxysporum resistance

Yinbo Ma Chungnam National University College of Agriculture and Life Sciences Sushil Satish Chhapekar Seoul National University College of Agriculture and Life Sciences Lu Lu Chungnam National University College of Agriculture and Life Sciences Sangheon Oh Chungnam National University College of Agriculture and Life Sciences Sonam Singh Chungnam National University College of Agriculture and Life Sciences Chang Soo Kim Chungnam National University College of Agriculture and Life Sciences Seungho Kim Neo Seed Co. Gyung Ja Choi Korea Research Institute of Chemical Technology Yong Pyo Lim Chungnam National University College of Agriculture and Life Sciences Su Ryun Choi (  srchoi@cnu.ac.kr ) Chungnam National University College of Agriculture and Life Sciences https://orcid.org/0000-00021401-9875

pests [22]. Speci cally, Fusarium wilt caused by the soil-borne fungal pathogen Fusarium oxysporum, can severely damage radish and numerous other types of vegetables [23]. Because of its broad host range, F. oxysporum can survive at relatively high soil temperatures (> 24 °C) and remain viable even in the absence of any host plant, making it di cult to control [24,25]. In addition to the TNL-type resistance genes FOC1 and FocBo1, which are responsible for the resistance of Brassica oleracea to F. oxysporum, many NBS-LRR genes related to Fusarium wilt resistance have been recently identi ed in diverse plant species [26][27][28].
Other genes belonging to the NBS-LRR family have been detected in various plants, including Arabidopsis thaliana [17,29], Brassica rapa [30], chickpea [31], and Gossypium species [32]. Unfortunately, to the best of our knowledge, the potential roles of radish NBS-LRR genes related to disease resistance have not been investigated. The available radish genome sequence is a useful resource for the whole-genome identi cation of transcription factor families [33,34]. However, the effects of F. oxysporum infections on radish NBS-LRR genes and their families remain unexplored. Therefore, combining bioinformatics and gene expression analyses to systematically study the evolution, expression, and potential functions of NBS-LRR genes may help to improve our understanding of the regulatory networks involved in radish plant growth and in response to F. oxysporum.
In this study, a genome-wide analysis of the radish genome identi ed 225 NBS-encoding genes (99 full NBS-LRR and 126 partial NBS genes) divided into two subclasses (CNL and TNL). These genes were further characterized regarding chromosomal locations, structures, and duplications. Additionally, with a focus on the NBS-LRR genes, we examined the encoded conserved domains as well as phylogenetic relationships and synteny with genes from A. thaliana, B. rapa, and B. oleracea. Furthermore, the pathogen-induced NBS-LRR gene expression pro les indicated that some R genes are differentially expressed in genetically different germplasm of radish (resistant and susceptible). Our results provide crucial insights into the evolution of this gene family in the radish genome. Moreover, the extensive R gene data presented herein may be useful for accelerating future breeding efforts aimed at improving the disease resistance of radish and other Brassicaceae crops.

Results
Identi cation and classi cation of NBS-encoding genes in R. sativus To comprehensively identify potential NBS-encoding genes in radish, the hidden Markov model (HMM) pro le NB-ARC (Pfam: PF00931) from the Pfam database was used to screen the protein sequences encoded in the radish genome [36]. A total of 488 gene candidates with an NBS-LRR domain were identi ed. These candidate NBS-encoding genes were manually screened and functionally annotated according to the closest A. thaliana homolog. Finally, 225 non-redundant NBS-encoding R gene candidates were identi ed in the Rs1.0 genome ( Table 1, Additional le 1: Table S1). The NBS-containing candidate proteins were classi ed into the TNL or CNL subfamilies based on their HMM pro les and the NCBI Conserved Domain Database, with relationships visualized in a phylogenetic tree (Fig. 1, Table 1).
The phylogenetic tree clearly distinguished the TNL and CNL genes in two separate clades (Fig. 1). Gene names were assigned based on their domain type and chromosomal position. The TNL subfamily included 80 genes with full-length domains (TIR, NBS, and LRR) as well as 54 TN (TIR and NBS, but no LRR), 25 NL TNL (NBS and LRR, but no TIR), and 15 N TNL (no TIR or LRR) genes. The remaining 51 genes belonging to the CNL subfamily included 19 with full-length domains (CC, NBS, and LRR) as well as 10 CN (CC and NBS, but no LRR), 9 RN (RPW8 and NBS, but no LRR), 2 NL CNL (NBS and LRR, but no CC), and 11 N CNL (no CC or LRR) genes. Genes encoding RPW8-NBS-LRR proteins were not detected in the radish genome. We also analyzed the well-characterized NBS-encoding genes from the genetically closely related plant species A. thaliana, B. rapa, and B. oleracea (Table 1) using the same method. A total of 164 (A. thaliana), 212 (B. rapa), and 244 (B. oleracea) NBS-encoding genes were identi ed. Genomic distribution among radish chromosomes Of the 225 NBS-encoding genes, 202 were mapped onto nine radish chromosomes, whereas the other 23 genes were located on several scaffolds which were not mapped on the chromosome of R. sativus (Fig. 2, Additional le 1: Table S1). We determined the distribution of CNL, TNL, and partial NBS-encoding genes on different chromosomes (Additional le 2: Fig. S1). Moreover, the TNL genes were almost uniformly distributed on chromosomes R01 to R09, whereas the CNL genes were not detected on chromosomes R03 and R06. Furthermore, chromosome R09 had the most NBS-encoding genes (41) There is no known pattern in the chromosomal distribution of NBS-encoding genes, with most of them detected in clusters. This distribution may facilitate sequence exchanges via recombination mispairing. We identi ed NBS-encoding gene clusters based on previously established criteria that an NBS-encoding cluster should have two or more genes separated by fewer than 200 kb, with no more than eight non-NBSencoding genes in between [35]. On radish chromosomes, 146 (72%) NBS genes were mapped in 48 clusters, whereas the remaining 55 genes were detected as singletons ( Fig. 2 and Additional le S2). Our analysis revealed that chromosome R09 has the most NBS genes (41; 20.30% of the mapped genes) distributed in eight clusters, in addition to nine singletons. Cluster sizes varied across the genome (2-11 genes). Cluster 44 was the largest, with 11 genes belonging to the TNL subfamily.
To clarify the evolutionary relationships among genes, we analyzed the distribution of NBS-LRR genes among the crucifer blocks in the radish genome. Of the 24 identi ed blocks (Fig. 2), the U block on chromosomes R02, R04, and R08 was the largest (48 NBS-encoding genes), followed by the R block (24 genes), D block (23 genes), E block (23 genes), and F block (17 genes). The U block may be one of the most important in terms of NBS-encoding genes.
Gene characteristics and structure The lengths of the genomic and coding sequences of the 225 NBS-encoding genes as well as the length, molecular weight (MW), and isoelectric point (pI) of the corresponding proteins were comprehensively analyzed (Additional le 1: Table S1).

Cis-element analysis
Cis-elements in promoters are usually involved in gene regulation. Thus, to further investigate the potential regulatory networks of the NBS-encoding genes, the cis-elements in the sequences 2 kb upstream of the start codon of 225 NBS-encoding genes were analyzed with PlantCARE (Additional le 1: Table S3). A total of 70 cis-elements were detected, including 10 hormone-responsive elements, 33 lightresponsive elements, 4 elements related to abiotic stress, 8 elements associated with tissue-speci c expression, and 15 other elements. The cis-element DRE, which is a common cis-acting element in promoter and enhancer regions, was detected in the promoter region of 218 NBS-encoding genes.
Additionally, an AT-rich fragment, which is a core promoter element located approximately 30 bp upstream of the transcription start site, was detected in 216 NBS-encoding genes. This sequence was considered to be an essential element in the promoter of the NBS-encoding genes. The promoter region of 159 NBS-encoding genes included the CGTCA-motif, which is a cis-acting regulatory element associated with methyl jasmonate-responsiveness. The abscisic acid-responsive element in uencing abscisic acid responsiveness was detected in the promoter region of 165 genes, suggesting that NBS-encoding genes are involved in plant responses to pathogen infections. Moreover, a TCA-element, which is involved in salicylic acid responses, was detected in 97 gene promoters, whereas the TGA-element related to auxin responses was identi ed in 82 gene promoters. Furthermore, the P-box and GARE-motif associated with gibberellin-responsiveness were present in 55 and 40 NBS-encoding gene promoters, respectively. These results may be relevant for developing a method for identifying candidate genes related to disease resistance.

Conserved motifs and phylogenetic relationships among TNLs and CNLs
To investigate the TNL and CNL gene and protein structures, we built a phylogenetic tree based on the full-length amino acid sequences encoded by the R. sativus TNL and CNL genes. The phylogenetic analysis indicated that the radish NBS-LRRs can be divided into two large groups, namely CNL and TNL ( Fig. 4 and Additional le 1: Table S4). The RsTNL group comprised four subgroups. Of the 80 RsTIR proteins, 49, 10, 5, and 16 belonged to subgroups RsTNL-1, RsTNL-2, RsTNL-3, and RsTNL-4, respectively. These results were identical to those of phylogenetic analyses of R. sativus and A. thaliana (Additional le Fig. 2).
To further elucidate the potential functions and diversi cation of the TNL and CNL genes in R. sativus, 20 encoded conserved motifs were identi ed and numbered 1-20 based on the MEME program (Additional le 1: Tables S4 and S5). The TIR domain was detected in all RsTNLs. Additionally, the RsTNL-1 subgroup members had most of the motifs. The RsTNL-2 and RsTNL-3 proteins lacked motifs 20 and 11, respectively. The proteins in subgroup RsTNL-4 were missing motifs 11, 12, and 16 ( Fig. 3). The RsCNL group members mostly had only six motifs, including the CC domain. Common motif compositions were revealed within subgroups. However, regarding the motif types and numbers, there was considerable diversity among subgroups. This suggests the proteins within subgroups are functionally similar.

Tandem duplication and synteny analyses of NBS-encoding genes
Whole-genome and tandem duplications are critical events for enhancing genome complexity and evolutionary novelty. In the R. sativus genome, 34 of 225 NBS-encoding genes (15.11%) were associated with tandem duplications and were distributed in 15 tandem arrays of 2-5 genes (Additional le 1: Table  S6). Our data also revealed variability in the number of duplicated genes per tandem duplication event and an uneven distribution of these duplications on ve of nine chromosomes. Genes encoding domains (e.g., CNL, TNL, and TN) were present in the tandem arrays. The 14 tandem duplication events detected for chromosome R09 involved ve RsTN genes. In contrast, single tandem duplication events occurred on chromosomes R02, R06, R07, R08, and R09, each involving two genes. Chromosome R08 had the most tandem arrays (six tandem groups containing 13 genes), re ecting a hot spot for the distribution of NBSencoding genes. An analysis of our data according to BLASTP and MCScanX methods identi ed 20 segmental duplication events involving 32 NBS-encoding genes ( Fig. 4 and Additional le 1: Table S7).
The Ka/Ks values (Additional le 1: Table S8) of the pairs of segmentally duplicated genes were less than 1, indicating these genes evolved under negative selection. These results suggest that both tandem and segmental duplication events were a major driving force for the evolutionary expansion of the NBSencoding genes in the radish genome.

Comparative synteny analyses of orthologous pairs of NBS-encoding genes
To further investigate the phylogenetic relationships among the radish NBS-encoding genes, we constructed three synteny maps comparing radish with A. thaliana, B. rapa, and B. oleracea. A total of 209 pairs of NBS-encoding genes (Additional le 1:  Table S9). There were about 97 homologous gene pairs between R. sativus and B. oleracea, among which 47, 19, and 4 NBS-encoding genes in radish were retained as one, two, and three copies, respectively, in B. oleracea, with the remaining genes lacking syntenic relationships. The partial NBS-encoding genes RsTN11, RsTN31, RsNL07, and RsN19 were syntenic with complete TNL genes (Bo2g010720, Bo8g104700, Bo9g061200, and Bo9g029350) in the B. rapa genome. Additionally, two RsCNL genes (RsCNL03 and RsCNL05) as well as four RsTNL genes (RsTNL09, RsTNL14, RsTNL62, and RsTNL72) corresponded to partial NBS genes (Bo6g020950, Bo1g048080, Bo3g154220, Bo6g089230, Bo2g126980, and Bo3g006960) in the B. rapa genome. We detected a total of 70 radish NBS-encoding genes syntenic with 67 NBS-encoding genes in the B. oleracea genome (Additional le 1: Table S9). An analysis of the synteny among the radish, A. thaliana, B. rapa, and B. oleracea genomes revealed the complete NBS-LRR genes were more syntenic than the partial genes.
Finally, a comparative analysis of the orthologous pairs of NBS-LRR genes among four species (R. sativus, A. thaliana, B. rapa, and B. oleracea) revealed 22 R. sativus NBS-encoding genes with corresponding copies in the A. thaliana, B. rapa, and B. oleracea genomes ( Fig. 6 and Additional le 1: Table S10). There was a one-to-one relationship among the NBS-encoding genes between the R. sativus and A. thaliana genomes. Additionally, half of the radish genes had single-copy syntenic genes in the B. rapa and B. oleracea genomes. Furthermore, the B. rapa and B. oleracea genomes had two and three copies of RsTN38, respectively, as well as three copies of RsTN01. Accordingly, these genes may have been important for the evolution of the NBS-LRR gene family.
The Ka/Ks ratio indicates the selective pressure on genes during evolution. We examined the Ka/Ks ratio for the orthologous gene pairs to determine the evolutionary selection patterns of NBS-LRR genes among R. sativus, A. thaliana, B. rapa, and B. oleracea. The Ka/Ks ratio re ects the number of non-synonymous substitutions per non-synonymous site (Ka) and the number of synonymous substitutions per synonymous site (Ks). The ratios for the orthologous gene pairs were estimated for each branch of the phylogenetic tree using the KaKs Calculator. The segmentally and tandemly duplicated NBS-LRR gene pairs as well as all orthologous NBS-LRR gene pairs had a Ka/Ks ratio < 1 (Additional le 1: Table S11), suggesting that the radish NBS-LRR gene family might have experienced strong purifying selection pressure during evolution.
Radish NBS-LRR gene expression pro les in response to F. oxysporum f. sp. raphanin 59 (FOR59) To identify NBS-LRR genes responsive to a FOR59 infection and determine their spatiotemporal expression patterns, we analyzed the transcriptome data for all NBS-LRR genes. The transcriptome data were generated for whole plant seedlings of the 'YR4' and 'YR18' genotypes infected with FOR59; the whole seeding of 'YR4' and 'YR18' were collected on days 0, 1, 3, and 6 after the inoculation of F. oxysporum for the subsequent sequencing on the Illumina platform. Furthermore, we extracted the NBS-LRR genes from the generated RNA-seq data. Transcriptome data were obtained for 171 of 225 NBSencoding genes based on the conserved domains and gene IDs. The remaining unidenti ed genes may be unexpressed in response to a FOR59 infection. Among these 171 NBS-encoding genes, 29 were not differentially expressed between 'YR4' and 'YR18', and were minimally expressed at all time intervals. Thus, these genes were excluded from the gene expression analysis. The fragments per kilobase of exon per million fragments mapped (FPKM) values obtained for the 142 remaining NBS-LRR genes varied from 0 to 278.884. A total of 80 NBS-encoding genes had an FPKM value less than 1 FPKM (low expression), whereas there were 84 genes with more than 1 FPKM value, and less than 10 (high expression). Extremely high expression represented in 13 genes with more than 10, and less than 50; 1 gene with more than 50, less than 100; and 1 gene with more than 100, less than 300. The differentially expressed genes following the FOR59 inoculation were analyzed with the R package edgeR [36], which revealed 36 genes that were more highly expressed in the resistant 'YR4' line than in the susceptible 'YR18' line at various time-points (Fig. 7). Seven genes (RsCN10, RsN21, RsN26, RsNL19, RsTN04, RsTNL15, and RsTNL38) were upregulated only in 'YR18', whereas six genes (RsCN05, RsCNL17.3, RsN22.1, RsTN31, RsTN32, and RsTN43) were up-regulated in 'YR4'. Most importantly, the expression of ve genes (RsN25.2, RsNL07, RsTN18, RsTNL09, and RsTN17) gradually increased over time in the 'YR4' plants, but remained relatively stable in the 'YR18' plants, implying they may be crucial for the resistance to Fusarium wilt. However, the RsTNL51 expression level was higher in 'YR18' than in 'YR4', and gradually increased during the infection period (Fig. 7g). These ndings suggest a total of 75 NBS-LRR genes contribute to the resistance of radish to Fusarium wilt, with six genes (RsN25.2, RsNL07, RsTN18, RsTNL09, RsTN17, and RsTNL51) potentially crucial for the resistance.

Responses to FOR59 infections
The expression patterns of NBS-LRR-encoding genes in 'YR4' and 'YR18' plants were con rmed by a quantitative real-time (qRT)-PCR assay. Speci cally, we validated the expression of only the CNL and TNL genes. Of the 19 RsCNL and 80 RsTNL genes in the radish genome, approximately 40 genes encoded proteins with amino acid differences between 'YR4' and 'YR18'. The expression patterns of these genes in 'YR4' and 'YR18' at 0, 1, 3, 6, 9, and 12 days after inoculation (DAI) were determined, with the 18S rRNA gene used as an internal control. The RsTNL03 (Rs093020) and RsTNL09 (Rs042580) genes on chromosome R02 were signi cantly more highly expressed in the resistant 'YR4' plants than in the susceptible 'YR18' plants ( Fig. 8), indicating these two genes are related to the Fusarium wilt resistance of radish. Interestingly, both genes were similarly expressed in 'YR4' and 'YR18' on day 0, but the infection rapidly up-regulated the expression levels in 'YR4' plants. The expression levels continued to increase until 9 DAI, further suggesting their potential role in Fusarium wilt resistance. Homologs of RsTNL09 (Rs042580) were identi ed in A. thaliana (AT4G36150), B. rapa (Bra010552 and Bra011666), and B. oleracea (Bo3g154220 and Bo7g117810).

Discussion
The NBS-LRR genes form a large family of stress resistance genes that are ubiquitous in all plant species. Genome-wide analyses of NBS-LRR gene families have been conducted for numerous species with sequenced genomes [17,20,37]. In the current study, we identi ed candidate NBS-LRR genes and studied their distribution, structure, clustering, duplication, synteny, and conservation. We identi ed 225 non-redundant NBS-encoding R gene candidates in the radish genome (  [17]. This distinct pattern of TNL gene abundance suggests that the TIR domain is more functionally active than the CC domain in Brassicaceae species [38]. The NBS-LRR genes were unevenly distributed across the radish genome, with chromosomes R09 and R03 having the most and fewest genes, respectively. Additionally, 23 of the genes were located on diverse scaffolds (Fig. 2). Similar to other species, the NBS-LRR genes in radish mainly exist in clusters because of their rapid evolution [39,40]. About 72% of the radish NBS-LRR genes were detected in 48 clusters, which is higher than the corresponding percentages in A. thaliana (61.7%), B. rapa (59.4%), and B. oleracea (60.3%) [17]. Moreover, 11 genes were included in cluster 44 of the radish genome. Gene families may expand because of polyploidizations, and tandem and local duplications are the most commonly evaluated mechanisms underlying gene family expansions [41]. We identi ed 66 NBS-encoding genes (29.33%) in radish that underwent tandem and segmental duplications, which is lower than the duplication rates for A. thaliana (55.7%), B. rapa (47.1%), and B. oleracea (43.3%) [17]. The TN genes were predominantly present in tandem arrays, and only three pairs of CN genes were detected in tandem arrays 2 (RsCn02 and RsCN03), 4 (RsCN07 and RsCN08), and 8 (RsCNL11 and RsCNL12). Interestingly, RsCNL14 and RsTNL59, which belong to two distinct subgroups (CNL and TNL), had undergone segmental duplications. These two genes were located in cluster 40 of the radish genome. We speculate that the NBS-LRR genes (especially TNL subfamily genes) in the radish genome had undergone inter-and intraspeci c replications. Because radish and A. thaliana are Brassicaceae species, we investigated the crucifer blocks in the radish genome containing the identi ed NBS-LRR genes. A total of 45 NBS-LRR genes were located in the U block distributed on different chromosomes (R02, R04, and R08). The genes were also relatively abundant in the R (24 genes), D (23 genes), E (23 genes), and F (17 genes) blocks (Fig. 2). These observations suggest that the R. sativus genome structure arose following the rearrangement and divergence from a common ancestor with A. thaliana [42,43]. All tandemly duplicated genes were found in the same clusters (i.e., highly similar), whereas most of the segmentally duplicated genes did not form clusters and were located on different chromosomes. The relatively few tandem and segmental duplications of NBS-encoding genes may help to explain why radish has evolved more slowly than other Brassicaceae species [44].
Gene duplication is considered as a major force for evolution [45]. The frequent sequence variation that can occur in the duplication process, can lead to domain's structural variation, directly involved in protein functional role [46,47]. In addition to this expansion, variation leads to neofunctionalization in family members. The neofunctionalization has been reported in KCS gene family and MADS box gene of Arabidopsis [48,49]. According to the genome comparison in Brassica family members, we identi ed that the structural variation happened in many NBS-encoding genes which in uenced the domain structure (Additional le 1: Table S12). For example, RsNL01 has the homolog gene in A. thaliana, B. rapa, and B. oleracea, but the number of exons was reduced, additionally LRR domain was added and TIR domain was lost in the R. sativus. Taken together of this results, we possibly suggest that these variation in the NBS-encoding genes might result in neofunctionalization or loss of function during the radish evolution.
The conserved structural domains of the radish TNL and CNL proteins were examined in this study. The NBS-LRR genes encoding the TIR domain were more common than the genes encoding the CC domain in the analyzed species (R. sativus, A. thaliana, B. oleracea, and B. rapa). However, the functions of several partial NBS-encoding genes lacking one or two domains (TIR, CC, and LRR) remain unknown, but their presence in various plant species imply they are important [17]. In the radish genome, the average number of exons was higher for TNL genes than for CNL genes, with half of the CNL genes containing only one exon (Fig. 4). This difference between gene types may be due to the conservative nature of CNL gene replications, which involve many regulatory components. This result is also consistent with those reported for A. thaliana, B. rapa, and B. oleracea [17,20,37]. Although the TNL and CNL genes are related to signaling and resistance speci city during pathogen recognition, the genes in these two subfamilies vary regarding sequences and associated signaling pathways and they cluster separately according to phylogenetic analyses [50]. Previous studies proved that the TNL group forms four phylogenetic clades in B. rapa [18] and B. oleracea [37]. Our motif analysis with the MEME program uncovered diverse motif compositions in the different RsTNL subfamilies. For example, the RsTNL-4 subgroup members had lost motifs 11, 12, and 16, whereas these motifs are present in the RsTNL-1 subgroup members. These differences may contribute to the functional divergence among the TNL proteins in radish.
Radish is an agronomically important root vegetable crop. Its genome, which contains triplicated segments, has intermediate characteristics between the Brassica A/C and B genomes, suggesting radish originated from a Brassica species [34]. According to a comparison with the Brassica A (Br), B (Bn), and C (Bo) genomes, the radish genome has been positioned between the Brassica A/C and B genomes [34].
Thus, an analysis of the synteny among the NBS-LRR genes of radish, A. thaliana, B. rapa, and B. oleracea may lead to novel insights into the evolutionary characteristics of RsNBS-encoding genes as well as the phylogenetic relationships with the genes in the other three species. In the current study, 39, 73, and 97 orthologous pairs were respectively identi ed between radish and A. thaliana, B. rapa, and B. oleracea (Fig. 5), implying that R. sativus is more closely related to B. oleracea than to A. thaliana or B. rapa.
Moreover, we identi ed genes in the A. thaliana, B. rapa, and B. oleracea genomes that correspond to 22 R. sativus NBS-encoding genes (Fig. 6). Additionally, the orthologous gene pairs may have existed before the ancestral divergence, with important roles related to the evolution of the NBS-LRR gene family [51].
Furthermore, we compared the Ka/Ks values of the orthologous gene pairs between R. sativus and A. thaliana, B. rapa, and B. oleracea lineages. For the CNL and TNL gene types, there were no signi cant differences in the orthologous gene pairs among the three species (Additional le 1: Table S11). We speculate that the NBS-LRR genes in R. sativus, A. thaliana, B. rapa, and B. oleracea may have been exposed to diverse selection pressures to develop genetic resistance to the same pathogen. Additionally, some pathogens may be speci c to certain Brassica species.
In a previous study regarding radish resistance to Fusarium wilt, we isolated 68 NBS-RGAs and 46 SRLK-RGAs from two Fusarium wilt-resistant radish inbred lines, and the genetic diversity was analyzed with RGA-speci c primers [52]. We also identi ed a major QTL (Fwr1) containing the ORF4 gene, encoding a serine/arginine-rich protein kinase that may mediate Fusarium wilt resistance [53]. To further explore the mechanism underlying Fusarium wilt resistance, we pro led the expression of NBS-LRR genes by screening transcriptome data for the resistant ('YR4') and susceptible ('YR18') lines at different timepoints during an infection (Fig. 7). A total of 75 NBS-LRR genes were identi ed as likely involved in the Fusarium wilt resistance of the 'YR4' line, including RsNL07, RsTNL09 (chromosome R02), RsTN17 and RsTN18 (chromosome R04), and RsN25 (scaffold), which exhibited gradually up-regulated expression in infected 'YR4' plants. Accordingly, the genes belonging to the same cluster were similarly expressed. The identi ed genes and expression pro les based on transcriptome data were veri ed by a qRT-PCR analysis of the 'YR4' and 'YR18' lines at several time-points after plants were inoculated (Fig. 8). Of the analyzed genes, RsTNL03 and RsTNL09 were more highly expressed in 'YR4' (resistant) than in 'YR18' (susceptible) at different time-points. The RsTNL09 (Rs042580) gene is homologous to AT4G36150 in A. thaliana, Bra010552 and Bra011666 in B. rapa, and Bo3g154220 and Bo7g117810 in B. oleracea. The Bra010552 gene encodes a TNL protein, and is one of the candidate genes for clubroot resistance in B. rapa [54]. The RsTNL03 gene, which is a homolog of AT4G16890 in A. thaliana, is likely involved in salicylic acid-dependent early plant defense responses [38][39][40]. Our RNA-seq and qRT-PCR data also indicated that RsTNL06 (Rs053740) was up regulated in 'YR18' than in 'YR4' after the FOR59 infection, taken together with the morphological characteristics (Additional le 1: Table S13) of 'YR18' (severe infected) suggests that RsTNL06 might be potential negative regulators of Fusarium wilt resistance in radish. Interestingly, RsTNL06 was identi ed as an ortholog of the B. oleracea gene Bo7g106630, and may contribute to plant responses to Fusarium species [37]. The potential positive and negative regulators of Fusarium wilt resistance that were identi ed following the comprehensive analyses of RNAseq and qRT-PCR data will need to be thoroughly functionally characterized to assess their utility for enhancing the Fusarium wilt resistance of radish.

Conclusions
The NBS-LRR genes are important for the resistance of radish plants to various pathogens. A comprehensive analysis of the radish genome revealed 225 NBS-encoding genes, including 99 NBS-LRR genes and 126 partial NBS genes. Additionally, of the 202 NBS-encoding genes mapped onto nine chromosomes, 72% were located in clusters. Another 23 genes were located on different scaffolds. An analysis of syntenic and phylogenetic relationships among the NBS-LRR genes from A. thaliana, B. rapa, and B. oleracea provided valuable insights into the evolutionary characteristics of radish NBS-LRR genes. Moreover, 99 complete NBS-LRR genes were grouped in two main families (TNL and CNL), with the TNL genes further divided into four subgroups with highly similar exon-intron structures and motif compositions. In this study, we identi ed 75 NBS-LRR genes that protect radish plants from Fusarium wilt based on the transcriptome data. We also determined that RsTNL03 (Rs093020) and RsTNL09 (Rs042580) positively associated with the resistance to F. oxysporum, whereas RsTNL06 (Rs053740) negatively associated with Fusarium wilt resistance in radish. The phylogenetic and gene expression data presented herein may help to clarify NBS-LRR gene functions.

Identi cation of NBS-encoding genes
To identify the NBS-encoding genes, we downloaded the whole-genome sequence from the radish database (http://radish-genome.org/). We applied hmmsearch of the HMMER (version 3.3) program [55] based on the HMM corresponding to the Pfam NBS (NB-ARC) domain (PF00931) to screen and identify the NBS-encoding genes in the radish genome. We also selected proteins with an E value < 1e −20 for sequence alignments with ClustalW [56]. We constructed the radish-speci c NBS HMM with the hmmbuild module of HMMER (version 3.3) to rescan the radish protein database, and proteins with an Evalue less than 0.01 were selected for further analyses [17,57].
The N-terminal of NBS-containing proteins usually include the TIR and CC domains, whereas the Cterminal contains the RPW8 and LRR domains. We used the CLC main workbench (version 7.9) [58] and the Pfam (version 32) program [59] to detect domains in the NBS-containing proteins. These results were con rmed with the NCBI Conserved Domains Tool [60] and MEME (Multiple Expression motifs for Motif Elicitation) [61]. The CC domain in the protein sequences was identi ed with Paircoil2 [62], with a P score cut-off of 0.025. The NBS-encoding genes in the genomes of A. thaliana (Araport11), B. rapa (version 1.5), and B. oleracea (version 1.1) were also identi ed using the same method. The A. thaliana genome was download from the TAIR database (www.arabidopsis.org), whereas the B. rapa and B. oleracea genomes were downloaded from the database on the BRAD website (http://brassicadb.org/brad/).

Multiple sequence alignment and phylogenetic analysis of NBS-encoding genes
These analyses con rmed the segregation between the two major NBS-encoding subpopulations (TNL and CNL) in the radish genome and clari ed the phylogenetic relationships among the genes on the major branches. Multiple full NBS-containing protein sequences were aligned with the MUSCLE [63] program. The MEGA-X [64] program was used for the phylogenetic analysis, which was completed with the maximum likelihood method based on the Whelan and Goldman model [65]. Finally, the maximum likelihood method involving the LG model with 500 bootstrap replicates was used to construct the phylogenetic tree. Two additional phylogenetic trees were constructed with the same method. One tree was used for analyzing the relationships between the CNL and TNL genes in radish, whereas the other trees were used to verify that the radish CNL and TNL genes are consistent with the corresponding A. thaliana genes.

Chromosomal locations, structures, and duplication events of NBS-encoding genes
The Mapchart (version 2.2) software was used to map all identi ed NBS-encoding genes onto R. sativus chromosomes based on their physical positions indicated in the R. sativus genome database. Several genes were organized in diverse NBS-LRR clusters, in which at least two NBS-encoding genes were localized in a 200-kb region and were separated by a maximum of eight non-NBS-LRR genes [35].
The tandem and segmental duplications of the NBS-encoding genes in the radish genome were analyzed with the default parameters of MCScanX [67] and TBtools [68]. Gene duplications were con rmed based on the following two criteria: (a) the length of the shorter aligned sequence covered > 70% of the longer sequence; (b) the similarity of the two aligned sequences was > 70% [69]. The non-synonymous (Ka) and synonymous (Ks) substitutions for each duplication event were calculated with the KaKs Calculator (version 2.0) [70].
The predicted TNL and CNL proteins were analyzed with MEME [71], with default iterative cycles and the maximum number of motifs set to 20. Additionally, TBtools was used to visualize the structures of the TNL and CNL genes according to the genomic and coding sequences.
Analysis of the orthologous gene pairs between R. sativus and three Brassicaceae species Orthologous genes are important for investigating the evolutionary associations among diverse species.
In this study, we identi ed the orthologous gene pairs between the R. sativus and A. thaliana, B. rapa, and B. oleracea genomes with the MCScanX program. Speci cally, the following parameters were used: e = 1e −20 , u = 1, and s = 5 [72]. The orthologous pairs of NBS-encoding genes were extracted and the maps for analyzing synteny were prepared with TBtools. The 24 genome building blocks from the ancestral karyotype were assigned to the radish genome as previously described [34].

Plant growth and Fusarium oxysporum inoculation
The 'YR4' and 'YR18' lines included in this study are highly resistant and susceptible to F. oxysporum, respectively. To do the sampling, pathogen inoculation was used the modi ed root-dipping methods following the method by Baik et al. [73]. F. oxysporum f. sp. raphanin 59 was supported from the Korea Research Institute of Chemical Technology (Daejeon, Korea). Ten-day-old seedling of germplasm were grown at control condition at 25 °C and further were subjected to water rinse to remove the soil from the root, and inoculated in spore suspension of F. oxysporum with concentration of conidia . 3×10 6 [74].
Previously generated transcriptome data for the FOR59-inoculated 'YR4' and 'YR18' lines at four timepoints (0, 1, 3, and 6 DAI) were analyzed (Accession number PRJNA643982) (Additional le 1: Table S14). The growth conditions of plant materials were mentioned earlier and two biological replications were used for RNA-Seq. The brief summary of RNA sequencing and library preparation is as follows. Firstly rRNA was removed from total RNA using RIBO COP rRNA depletion kit (LEXOGEN, Inc., Austria). Total RNA was used to prepare cDNA libraries using the NEBNext Ultra II Directional RNA-Seq Kit (NEW ENGLAND BioLabs, Inc., UK) as per manufacture's instruction. Subsequently, libraries were checked using the Agilent 2100 bioanalyzer (DNA High Sensitivity Kit) to evaluate the mean fragment size.
Quanti cation was performed using the library quanti cation kit using a StepOne Real-Time PCR System (Life Technologies, Inc., USA). High-throughput sequencing was performed as paired-end 100 sequencing using Hiseq X10 (Illumina, Inc., USA). Total RNA-Seq reads were mapped using TopHat software tool [75] in order to obtain bam le (alignment le). The alignment les were used for assembling transcripts, estimating their abundances and detecting differential expression of genes using cu inks. We used the FPKM (fragments per kilobase of exon per million fragments) as the method of determining the expression level of the gene regions. The FPKM data were normalized based on Quantile normalization method using EdgeR within R [36]. Additionally, heat maps were produced with TBtools. Availability of data and material

List Of Abbreviations
All the RNA-seq data of radish inbred lines (YR4 and YR18) has been deposited at NCBI (Accession number PRJNA643982) and provided in the article and in the form of tables as additional le.

Competing interests
The authors declare that they have no competing interests.