Akebia trifoliata (Thunb.) Koidz, belonging to the family Lardizabalaceae (Liu et al., 2007), presents great potential for use as a fruit and oil crop and a medicinal and ornamental plant (Li L. et al., 2010; Chen et al., 2017). After extensive artificial cultivation, wild plants will inevitably suffer from various diseases, such as kiwifruit bacterial canker (Scortichini et al., 2012), sweet potato root rot disease (Ma et al., 2020) and apple fruit anthracnose (Zhou and Zhou, 2017); this is especially true for perennials, such as black cottonwood, which often suffer attacks from various pathogens or herbivores before reaching the reproductive stage because of their relatively long life cycle (Tuskan et al., 2006). Likewise, A. trifoliata is a typical perennial that is attacked by various pathogens, such as Alternaria tenuissima, causing leaf spot (Ye et al., 2013; Zhang et al., 2015), Colletotrichum gloeosporioides, causing leaf anthracnose (Cheng et al., 2020; Pan et al., 2020), Nigrospora sphaerica, causing fruit shrinkage (Hong et al., 2021), and Oidium subgenus Pseudoidium, causing powdery mildew (Garibaldi et al., 2004), which has resulted in large yield losses and significantly reduced the economic value of the species.

During the arms race between hosts and their pathogens, great numbers of variations arise in pathogens due to their short growth cycle and simple genomic sequences, requiring plants (especially perennials) to evolve a set of sophisticated and effective defense systems to combat them (Tuskan et al., 2006). Beyond the first line of defense against pathogens imposed by the thickened cell wall, plants usually build the second defense system by producing resistance proteins, mainly encoded by NBS genes, to resist invading pathogens after the failure of the cell wall defense system (Andersen et al., 2018). By directly or indirectly recognizing pathogen-secreted effectors, these proteins confer resistance to various pathogens, including fungi, bacteria and viruses, by initiating a series of defense responses, such as hypersensitive responses, activating signaling pathways and consequently inhibiting the plant infection process (Andersen et al., 2018).

The NBS gene family, encoding proteins containing nucleotide binding sites (NBSs) and C-terminal leucine-rich repeats (LRRs), is the largest family of plant resistance (R) genes, accounting for over 60% of detected and cloned R genes in all plant species (Kourelis and van der Hoorn, 2018). The NBS domain can bind ATP/GTP, resulting in phosphorylation to transmit disease resistance signals downstream, which plays a key role in combating almost all pathogens (Bent, 1996). According to the N-terminal structure of NBS genes, scientists have commonly divided these genes into three subfamilies: (TIR)-NBS-LRR (TNL) genes, containing TIR domains with homology to Toll/interleukin-1 receptors, and (CC)-NBS-LRR (CNL) and (RPW8)-NBS-LRR (RNL), characterized by CC and RPW8 domains, respectively (Xiao et al., 2001; Collier and Moffett, 2009; Wang et al., 2013). The RNL clade is composed of two different lineages: the Nicotiana benthamiana N-required gene 1 (NRG1) and Arabidopsis activated disease resistance gene 1 (ADR1) lineages. To clearly indicate the absence of the TIR domain compared with TNLs, we refer to CNLs and RNLs as the non-TNL (nTNL) subclass. Generally, TNL and CNL proteins are mainly responsible for recognizing specific pathogens, while RNL proteins may play an auxiliary role in downstream defense signal transduction.

With the completion of the sequencing of many plant genomes in succession, genome-wide analysis has increasingly become an important tool for studying the genetic diversity and evolution of NBS resistance genes (Meyers et al., 2003; Xue et al., 2019; Zhang et al., 2020). The number of NBS genes in a plant genome ranges from dozens to more than 2,000 (Shao et al., 2019; Andersen et al., 2020), leading to the discussion of whether the NBS gene number is related to genome size and, if so, whether the relationship is positive or not (Liu et al., 2014; Zhang et al., 2016; Goyal et al., 2020). In addition, the composition of the three subclasses is usually not the same or is even highly lopsided. For instance, Dioscorea rotundata possesses 166 CNLs, only one RNL and no TNLs (Zhang et al., 2020), and Brassica napus contains 461 TNLs, 180 CNLs and no RNLs (Alamery et al., 2018), suggesting a major discrepancy between the different subclasses. The differences in the composition of NBS genes are considered to be responsible for the diversity and specificity of resistance to various pathogens (Jupe et al., 2012). NBS genes generally arise via various types of duplications, mainly including tandem and dispersed duplications, which have a critical influence on the different arrangements of NBS genes, such as the occurrence of singletons or clustered loci on chromosomes (Meyers et al., 2003; Leister, 2004).

The exploration of the characteristics of NBS genes usually involves phylogenetic and structural analyses in diverse plant species, which has greatly accelerated the identification and utilization of functional R genes. Although there has been no research progress on A. trifoliata NBS resistance genes to date, the available genome and transcriptome data of A. trifoliata not only enrich the genomic resources available for this species but also create favorable conditions for the systematic study of NBS genes at the genomic level. In this study, we comprehensively illustrated the NBS gene profile of A. trifoliata using multiple analysis tools. A collection of R gene resources valuable in production were screened, providing meaningful data for further identifying functional R genes and accelerating the future breeding of resistant cultivars of A. trifoliata.

NBS genes are the largest disease R gene family that exists in all plants, and broad genome-wide analyses of this family have been carried out in recent years (Jones et al., 2016). Despite the existence of recent publications on the A. trifoliata genome (Huang et al., 2021), the corresponding genomic data files are still unavailable online. Therefore, a high quality A. trafoliata genome assembled and submitted by our group to NCBI was employed to identify and characterize NBS candidates (see detailed data in section “Materials and Methods”), including analyses of their numbers, structures, distribution, duplications and gene expression profiles.

Previous studies have shown that the number of NBS genes varies greatly among different species and usually reaches hundreds in a single plant species, to confer protection against diverse and rapidly evolving pathogens. There are over 2,000 NBS genes in the extremely large wheat genome, but these genes are extremely scarce in C. papaya, Thellungiella salsuginea and some orchids (Shao et al., 2016; Xue et al., 2019; Andersen et al., 2020). Here, we identified 73 NBS genes, accounting for only 0.3% of the 24,138 annotated genes of A. trifoliata (Table 2), which suggested that the number of NBS genes in its genome was small. In previous studies, some authors have suggested that a positive correlation exists between NBS gene number and genome size within five samples, where the number increases with increasing genome size (Goyal et al., 2020), while the others have inferred that there is no significant correlation between them also in five samples (Zhang et al., 2016). In this study, total 13 representative species including four basal eudicots, four core eudicots and five monocots (listed in Table 1), were employed to address this issue mainly because of the evolutionary relationship with A. trifoliata, available knowledge about NBS genes, and economic importance. The results showed a significant correlation coefficient (0.96) between NBS gene number and genome size among the 13 species; however, the value decreased to −0.14 after we removed wheat as the out-group (Table 1). Therefore, the correlation between NBS gene number and genome size could depend on both the size and components of the genome sample employed.

Despite the small number of NBS genes in A. trifoliata, the types of these genes are diverse, including 9 subclasses within 3 classes (Table 2), and these classes are distinct in composition (a 3:1 ratio of nTNL:TNL). Previous studies have shown that the nTNL:TNL ratio varies widely, ranging from 1:3 in B. oleracea (Liu et al., 2019) to an extremely high ratio of 18:1 in Capsicum annuum (Qian et al., 2017) among dicots, while monocots only contain the nTNL subclass because of the loss of the TNL subclass after the divergence of dicots and monocots (Qian et al., 2017). These large disparities in composition not only indicate that distinct NBS subgroups have experienced different evolutionary events, such as different types of genomic duplication and magnitudes of natural selection, but also provide key clues for establishing the relationships between speciation and the growth environment.

We found that there were markedly more CNL genes (68.5%) in A. trifoliata than TNL (26.0%) and RNL (5.5%) genes, which showed conserved evolution without large-scale expansion (Table 1). Moreover, the four RNLs formed an extremely conserved branch that did not deviate from the CNLs in Figure 5, and similarly low proportions of RNLs have been reported in some orchids (Xue et al., 2019) and D. rotundata (Zhang et al., 2020), in distinct contrast to the large numbers of CNLs and TNLs. A reasonable explanation for this phenomenon is that RNLs only assist other functional NBS genes as downstream genes with conserved functions, as reported for “Helper NBS-LRRs” (Bonardi et al., 2011; Collier et al., 2011). In fact, the number of NBS subfamilies retained in the genome usually depends on the genomic duplication type and the selection pressure exerted during evolutionary history (Mun et al., 2009). The analysis of the Ka/Ks ratio to explore the evolutionary difference between CNL and TNL showed that both CNLs and TNLs have experienced purifying selection, and the significantly lower Ka/Ks value of CNLs than TNLs at the p = 0.01 level (Supplementary Table 5) demonstrated that greater selective constraints and purifying selection in the CNL group and relaxation of purifying selection in the TNL group might have resulted in ancient diversification of NBS genes (Mun et al., 2009).

As the primary source of raw genetic material, WGD is widely involved in plant evolutionary processes and is especially common in angiosperms (Akoz and Nordborg, 2019), but only two NBS genes from the ancestral WGD event survived in A. trifoliata (Supplementary Table 1). This was probably because the synteny of NBS genes resulting from WGD was destroyed in the process of long-term evolution or because the number of NBS genes might have been reduced after WGD via rapid species-specific gene loss (Die et al., 2018a; Song et al., 2019). Recently, some studies have suggested that tandem and dispersed duplication events could be a major driving force of NBS gene expansion (Baumgarten et al., 2003; Qian et al., 2017); however, the main duplications in each subclass of NBS are still unclear. In the present study, we found that in the A. trifoliata genome, CNLs were mainly expanded by the dispersal of individual or small groups of duplicated genes to unlinked loci, while TNLs were mainly expanded by tandem duplication, which produced closely related genes in the same family that clustered together (Supplementary Table 1; Leister, 2004). Additionally, the chromosomal distribution of NBS genes in the A. trifoliata genome (Figure 2) agreed well with the differences in gene expansion, which could also explain the greater number of CNLs relative to TNLs.

The other interesting characteristic involved the relationship between exon number and gene length in the A. trifoliata genome, where a significant positive relationship existed between exon number and gene length among the 73 NBS genes (Table 2). We further found that TNLs showed larger values than CNLs for both exon number and gene length in all other six dicots on our list (Supplementary Table 7), including Raphanus sativus (Ma et al., 2021), A. thaliana, B. oleraces, B. rapa (Mun et al., 2009), B. napus (Alamery et al., 2018), Helianthus annuus (Neupane et al., 2018), and Manihot esculenta (Lozano et al., 2015). All results from the above species indicated that independent gains or losses of exons occurred in distinct NBS subclasses, and the greater exon number and longer gene lengths of TNLs could result from a more ancient evolutionary history than CNLs, which evolved more recently (Meyers et al., 2003).

Structurally, the NBS domains of the characterized R genes shared eight common conserved motifs (Meyers et al., 2003). Through MEME analysis, we found that 45.2% of the identified genes contained all eight motifs, and P-loop (66), GLPL (68) and kinase-2 (70) motifs were the most common ones (Supplementary Table 3). Both the P-loop and kinase-2 motifs function in ATP/GTP binding (Traut, 1994; Meyers et al., 2003), and the existence of these typically conserved NBS motifs in R proteins is essential for protein function. Although the specific function of GLPL is still poorly understood, its importance in disease resistance has been demonstrated by site mutation analysis in the context of flax rust (Dodds et al., 2001). Thus, the higher the frequency distribution of a motif is, the more conserved it is. In addition, the order of these eight motifs was shown to be relatively inflexible, and only two (AtNBS19 and AtNBS63) of the 73 identified NBS genes in A. trifoliata actually showed an order change (Figure 1 and Supplementary Table 2). We inferred that the conserved order of these motifs could be more important than the sequences themselves.

We further observed that the last residue of the kinase-2 motif exhibited the main difference between TNL and nTNL: 88.9% of nTNLs contained W as the last residue of the kinase-2 motif, while it was replaced with D in 68.4% of TNLs (Supplementary Table 3), which was consistent with previous works (Meyers et al., 2003; Die et al., 2018b). Initial reports suggested that the conserved aspartate of kinase-2 is indispensable for the concentration-dependent binging of Mg²⁺ to ATP, which is essential for its function in plant defense signaling (Tameling et al., 2002). Thus, we inferred that the TNL proteins with more conserved aspartate residues may be more functional ATP-binding proteins with ATPase activity.

Homologous proteins with sequence similarity usually retain similar functions in evolution, and consequently we aimed to select some key and well known NBS resistance genes from other species, which was conducive to the identification of close functional homologs in A. trifoliata (Lozano et al., 2015). Then, the TNL (Figure 4) and CNL (Figure 5) trees were separately reconstructed by adding reference TNLs and CNLs/RNLs. Many TNL candidates clustered with functional genes for resistance to bacterial or fungal root rot diseases in TNL-B, and five CNL candidates clustered with the Dm3 gene for resistance to downy mildew in clade CNL-2b under 100% bootstrap support, while only CNL candidates were included in clade CNL-4, indicating that these genes might provide resistance to unknown A. trifoliata pathogens or play a role in non-host resistance responses. The widespread A. trifoliata NBS genes in the phylogenetic trees suggested that the different NBS subclasses could provide broad resistance to various pathogen species, while a few NBS genes would confer resistance to diverse races of the same pathogen. Furthermore, CNLs showed a higher level of topological complexity than TNLs in Figure 3, which might be helpful for broadening the resistance spectrum to achieve a broadly suitable distribution and could explain the long evolutionary history of this basal eudicot (Sun et al., 2016; Qian et al., 2017).

Comparative transcriptomic analysis can provide important insight into gene structures, functions and application prospects (Lowe et al., 2017). In the present study, we found that 90.4% of NBS genes were expressed in various fruit tissues at different developmental stages, and most of them were expressed at low levels (Figure 6A), which indicated that they have broad functions under field conditions. The relatively large number of NBS genes (46 of 66) that were mainly expressed in the rind at higher expression levels, especially in later stages, also indicated that most of these genes functionally conferred broad horizontal resistance to various pathogens in the field because only the rind tissue of A. trifoliata is exposed to pathogens. In addition, rind is usually considered the first barrier to pathogen infection, and hence the higher NBS gene expression compared with both flesh and seed in the later stages is reasonable. Only AtNBS51 (RNL) was identified as a gene with intermediate or high expression, and 27 NBS genes were identified as genes with low expression (Supplementary Table 6), which indicates that no pathogen can escape the defense exerted by a few NBS genes under a given condition. In addition, the high proportion (75%) of stably expressed RNLs despite the small proportion (5.5%) of RNLs among all 73 genes further demonstrated that with the help of RNLs, many CNLs and TNLs of A. trifoliata can confer broad resistance to various pathogen species and/or different races of the same species. Together, all evidence suggests that the small number of NBSs with structural diversification in A. trifoliata can provide defense against attacks by all types of pathogens that may be encountered. From the perspective of plant pathology, this provides a reasonable explanation for the wide distribution of wild A. trifoliata in various geographical and ecological regions and its continued healthy growth over many years.

In conclusion, only 73 non-redundant NBS genes were identified in the A. trifoliata genome in this work, but they showed high diversity, including three classes (CNL, TNL and RNL) with nine subclasses. Evolutionarily, although almost all NBS genes were retained by purifying selection, the duplication style was different: CNLs were mainly expanded by dispersed duplication, while TNLs were mainly expanded by tandem duplication. Structurally, the order of the eight conserved motifs in the NBS domain was more conserved than their amino acid sequences. Functionally, differential expression analysis showed that many of NBS genes could play an actual role in fighting pathogens by conferring broad-spectrum resistance, which may have contributed to the genetic basis favoring adaptation to various environments during the evolutionary of A. trifoliata. Therefore, this study paves the way for improving disease resistance during the domestication of A. trifoliata and further enriches the available information about NBS genes, especially in basal eudicots.

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

XY and PL conceived and designed the project. SZ obtained the data. XY analyzed the data and wrote the manuscript. CC, WC, HuY, HaY, JG, FT, TR, JS, MZ, and PF participated in the data analysis and discussion. PL revised the manuscript. All authors contributed to the discussion of the results, reviewed the manuscript and approved the final article.

WC and HaY were employed by the company Sichuan Akebiatrifoliata Biotechnology Co., Ltd. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.