Cross-Superfamily Molecular Convergence and Lineage-Specific Divergence in Plant-Parasitic Nematodes: Comparative Genomics Within Tylenchina

1. Introduction

Plant-parasitic nematodes (PPNs) are among the most destructive pathogens in global agriculture, causing direct economic losses of tens of billions of US dollars annually [1,2]. In the evolutionary history of the phylum Nematoda, the acquisition of the plant-parasitic habit was not a single event but occurred independently in multiple lineages. Within the suborder Tylenchina (order Rhabditida), several major superfamilies have evolved the ability to parasitize plants. Among these, the superfamily Aphelenchoidea encompasses numerous severely damaging plant pathogens; for example, Bursaphelenchus xylophilus, the causal agent of pine wilt disease, can kill host trees within weeks to months and poses a devastating threat to global pine forest ecosystems [3,4], while the congeneric Aphelenchoides besseyi is the causal agent of rice white-tip disease [5]. In addition, Ditylenchus destructor in the family Anguinidae (superfamily Sphaerularioidea) severely damages the yield and quality of tuber crops such as sweet potato and potato [6], whereas the superfamily Tylenchoidea comprises the agronomically devastating root-knot and cyst nematodes. Although these nematodes belong to different superfamilies and are only distantly related phylogenetically, each has independently evolved the capacity to penetrate host defenses, infect living plant tissues, and acquire nutrients, representing multiple parallel pathways in the transition of nematodes to plant parasitism [7,8].

Over the past two decades, research on the molecular mechanisms of plant–nematode interactions has overwhelmingly focused on species of the superfamily Tylenchoidea, particularly root-knot and cyst nematodes. These studies have established the classic molecular framework of nematode parasitism and pathogenesis: for instance, effector proteins secreted by the esophageal glands interfere with host immunity by targeting key immune nodes such as host calreticulin and chitinase [9]; cell wall-degrading enzymes acquired through horizontal gene transfer (e.g., cellulases and pectinases) facilitate penetration and migration within host tissues [8]; and significantly expanded detoxification gene families (e.g., cytochrome P450 and glutathione S-transferase) effectively counteract plant chemical defenses [10]. This research bias stems mainly from two factors. First, as the most widely distributed sedentary endoparasitic nematodes, root-knot and cyst nematodes of the Tylenchoidea cause the most severe economic losses. Second, they induce complex feeding structures such as giant cells or syncytia in the host and can be functionally validated efficiently using model plants such as Arabidopsis. By contrast, research on plant-parasitic nematodes of the Aphelenchoidea and of the family Anguinidae (Sphaerularioidea) has lagged markedly. This is partly attributable to differences in parasitic style—Anguinidae and Aphelenchoidea nematodes are mostly migratory endoparasites or ectoparasites that generally do not induce complex feeding structures—and, more importantly, to the long-standing lack of mature gene functional validation systems, which has greatly constrained the depth and scale of research on their parasitic mechanisms.

Approaching the problem from the perspective of comparative genomics and adaptive evolution offers an ideal entry point for deciphering the pathogenic mechanisms of non-model plant-parasitic nematodes. Phylogenetic analyses indicate that within the suborder Tylenchina, plant parasitism arose independently not only among the Aphelenchoidea, Tylenchoidea, and Anguinidae (Sphaerularioidea), but also multiple times within the Aphelenchoidea itself, where parasitic lineages are deeply nested among fungivorous and free-living relatives [11,12]. Such a phylogenetic topology indicates that the common ancestor of these nematodes was not a plant parasite; consequently, any parasitism-related molecular similarity among them cannot be simply explained as the inheritance of an ancestral trait, but more likely reflects innovative modifications and molecular convergence during adaptive evolution [13]. This makes the Aphelenchoidea and Anguinidae natural model systems for studying the evolution of molecular convergence. If distantly related or independently evolved plant-parasitic nematodes exhibit a high degree of consistency in terms of parallel gene family expansion and amino acid substitutions at homologous sites, this would provide incontrovertible genomic evidence for adaptive convergence and directly reveal the general rules by which a plant-parasitic lifestyle shapes the genome. Furthermore, broadening the perspective to the entire suborder Tylenchina and introducing the independently evolved parasitic lineage of the Tylenchoidea for cross-superfamily comparison can address a fundamental question in evolutionary biology: when different lineages within the Tylenchina unlocked the ability to parasitize plants, did they share a common set of molecular “rewiring rules”, or did they achieve the same outcome through entirely different evolutionary paths?

On this basis, we integrated representative species of the Aphelenchoidea, the Anguinidae (Sphaerularioidea), and the Tylenchoidea to conduct a systematic analysis of molecular convergent evolution. Our sampling design covered not only the multiple independent transitions to plant parasitism within the Aphelenchoidea but also included closely related non-plant-parasitic outgroups, establishing a multi-layered evolutionary comparison framework. This study aimed to address three core scientific questions: (i) What shared molecular innovations and convergent features did the Aphelenchoidea and Anguinidae evolve during their transitions to plant parasitism? (ii) Do these molecular convergences reflect a universal parasitic molecular response shared with the Tylenchoidea, or do they represent evolutionary paths unique to the Aphelenchoidea/Anguinidae? (iii) How do these convergent evolutionary features jointly shape the molecular mechanisms by which nematodes exploit their hosts, at both the genomic and protein-structural levels? By precisely dissecting the shared and unique genomic features across different parasitic transitions, this study aims to disentangle the core parasitic requirements determined by evolutionary constraint from lineage-specific adaptive mechanisms, to deepen our understanding of the parallel and convergent evolution of plant parasitism in nematodes, and to provide a solid theoretical foundation for developing broad-spectrum, durable molecular targets for nematode disease control.

2. Materials and Methods

2.1. Species Selection, Genome Data Acquisition, and Quality Control

To balance the precision of micro-scale molecular evolutionary rate comparisons with the statistical power of macro-scale genomic difference analyses, we designed a two-tier dataset framework comprising a “core species set” (Core Set, N = 7) and an “expanded species set” (Expanded Set, N = 16) (Table S1). The core set comprised three comparative pairs of plant-parasitic species and their closely related non-parasitic species: the first pair was Bursaphelenchus xylophilus BXYJ5 (plant-parasitic) and Bursaphelenchus okinawaensis (non-parasitic); the second pair was Aphelenchoides besseyi (plant-parasitic) and Aphelenchoides bicaudatus (non-parasitic); and the third pair was a cross-superfamily comparison of independently evolved plant-parasitic species: Ditylenchus destructor (plant-parasitic) from the superfamily Sphaerularioidea and Aphelenchus avenae (non-plant-parasitic) from the superfamily Aphelenchoidea. In addition, the free-living nematode Panagrellus redivivus (order Rhabditida, suborder Panagrolaimomorpha) was selected as an outgroup for ancestral state calibration. The core set was designed to yield a high-quality, quantitative set of 1:1 single-copy orthologous genes suitable for precise multiple sequence alignment, ancestral state reconstruction, convergent site identification, and molecular positive selection testing, while minimizing alignment gaps and trimming noise caused by the inclusion of too many distantly related species.

The expanded set additionally included Globodera pallida, Heterodera schachtii, and Meloidogyne graminicola (all important plant-parasitic nematodes) from the Tylenchoidea; Strongyloides ratti, Parastrongyloides trichosuri, and Steinernema carpocapsae (animal-parasitic and insect-parasitic taxa, respectively) from the superfamily Panagrolaimoidea; Bradynema listronoti (insect-parasitic) from the Sphaerularioidea; Acrobeloides nanus (free-living) from the superfamily Cephaloboidea; and Aphelenchoides fujianensis (non-plant-parasitic nematode) from the Aphelenchoidea. The expanded set was used specifically for non-parametric statistical tests of gene family copy number and KEGG pathway member counts, and—by increasing the sample size to provide adequate statistical power—for cross-lineage validation of the candidate convergent sites identified in the initial screen.

All genome sequences and their annotation files (including protein-coding sequences, CDS, and gene structural information in GFF3 format) were downloaded from the WormBase ParaSite database (https://parasite.wormbase.org/; version WBPS18) [14] and NCBI database. All selected core genomes were published chromosome-level or high-quality scaffold-level assemblies to ensure the reliability of subsequent comparative genomics and convergent evolution analyses.

2.2. Orthologous Gene Clustering and Phylogenetic Tree Construction

OrthoFinder (v2.5.4) was used with default parameters to cluster whole-genome protein sequences into orthologous gene sets to obtain orthologous gene sets among species [15]. All single-copy orthologous genes were extracted and aligned separately using MAFFT (v7.490) [16], and conserved regions were trimmed using trimAl (v1.4) [17]. To clarify the phylogenetic relationships among the seven nematode species, the amino acid sequences of all single-copy orthologous genes were concatenated into a supermatrix. ModelFinder was used to estimate the optimal amino acid substitution model based on the Bayesian Information Criterion (BIC) [18]. Maximum likelihood (ML) phylogenetic analysis was performed using RAxML-NG to construct the ML tree, with 1,000 bootstrap replicates to assess branch support [19].

2.3. Assessment of Conservation Between plant-Parasitic and Non-Plant-Parasitic Nematodes

To assess sequence conservation and selection pressure divergence between plant-parasitic and non-plant-parasitic nematodes, two species pairs were analyzed independently: B. xylophilus and B. okinawaensis, and A. besseyi and A. bicaudatus. Based on the whole-protein sequences of the two species, candidate 1:1 orthologous gene pairs were identified using the bidirectional best BLASTp approach (E-value < 10⁻⁵), with the additional requirements of sequence coverage ≥ 50% and identity ≥ 30%; only gene pairs meeting these criteria were retained. The corresponding CDS sequences were extracted, and protein alignments were mapped to codon alignments using Pal2NAL. The nonsynonymous substitution rate (Ka) and synonymous substitution rate (Ks) were then calculated for each gene pair using KaKs_Calculator (YN model). By comparing the differences in the Ka/Ks distributions between the two pairs [20], the degree of divergence in overall selection pressure between plant-parasitic and non-parasitic nematodes was assessed.

2.4. Screening of Amino Acid Convergent Sites

Amino acid convergent sites were screened in the single-copy orthologous gene set described above. The CSS (Convergence at Specific Sites) method published by Xu et al. was used to screen convergent sites [21], with the following criteria: when the amino acid state of the three non-plant-parasitic nematodes was fully consistent with that of the outgroup (P. redivivus), the ancestral state at these sites was considered identical to the outgroup state. Given the known ancestral state, if the amino acid states of the three plant-parasitic nematodes (B. xylophilus, A. besseyi, and D. destructor) were completely identical to one another and all differed from the ancestral state, the site was designated a potential convergent site. After obtaining genes containing potential convergent sites, the method of Zou and Zhang (2015) was used to calculate the expected number of convergent sites under random substitution for each target gene [13]; genes in which the number of convergent sites was significantly higher than expected (P < 0.05) were retained. Confirmation by alignment was then performed in the 16-species expanded validation set. When a convergent site within a gene maintained ≥ 80% conservation among the non-plant-parasitic nematodes of the expanded validation set and no convergent amino acid state identical to that of the plant-parasitic nematodes appeared, the site was recognized as a bona fide convergent site. If the convergent amino acid substitutions occurred only in the Aphelenchoidea and Anguinidae, the corresponding gene was designated an Aphelenchoidea/Anguinidae-specific convergent gene; if convergent amino acid substitutions in the same direction also occurred in the Tylenchoidea, the gene was designated a suborder-wide shared convergent gene.

2.5. Screening of Convergent Sites in Amino Acid Physicochemical Properties

To capture convergence signals at a higher level of physicochemical properties, four grouping schemes (GS1–GS4) were analyzed separately. Among these, GS1, GS2, and GS3 are coarse-grained classifications that partition the 20 amino acids into 6 or 7 categories based on polarity and volume, whereas GS4 is a fine-grained classification that partitions the amino acids into 12 categories based on side-chain structural similarity [22]. The grouping information is as follows:

· GS1 (6 groups): [PAGS], [CV], [T], [QNED], [LIFWMY], [HRK]

· GS2 (7 groups): [AGV], [ILFP], [YMTS], [HNQW], [RK], [DE], [C]

· GS3 (6 groups): [C], [AGPST], [NDQE], [RHK], [ILMV], [FWY]

· GS4 (12 groups): [G], [P], [F], [AILV], [ST], [DE], [C], [NQ], [RH], [K], [WY], [M]

For the amino acid alignment of each gene, if the amino acid states of the outgroup species and all non-plant-parasitic nematodes at a given site belonged to the same classification category, that site was defined as a conserved-property site. Among the conserved sites identified, each site was examined for the presence of a property convergence signal (i.e., a parallel shift by the plant-parasitic nematodes toward the same amino acid property category). After obtaining genes containing potential property convergent sites, genes in which the number of convergent sites was significantly higher than expected (P < 0.05) were retained. Finally, confirmation by alignment was performed in the 16-species expanded validation set; when a property convergent site within a gene maintained ≥ 80% category conservation among the non-plant-parasitic nematodes of the expanded validation set and no amino acid state of the same category as that of the plant-parasitic nematodes appeared, the site was recognized as a bona fide property convergent site. The classification rules were consistent with those for the amino acid sites described above.

2.6. Screening of Genes Under Shared Positive Selection in Plant-Parasitic Nematodes

To screen for genes showing positive selection signals in all three plant-parasitic nematode lineages, the codeml program in the PAML v4.9 package was used to perform branch-site model (Branch-site Model A, model = 2, NSsites = 2) analysis. Each plant-parasitic nematode (B. xylophilus, A. besseyi, and D. destructor) was independently designated as the foreground branch; in each run, the alternative model (fix_omega = 0, allowing the foreground branch ω2 > 1 to be freely estimated) was compared with the null model (fix_omega = 1, fixing the foreground branch ω2 = 1). The significance of the positive selection signal was assessed using the likelihood ratio test (LRT), with P < 0.05 as the threshold [23,24]. For genes passing the LRT, sites with a Bayes Empirical Bayes (BEB) posterior probability ≥ 0.95 were extracted as significant positively selected sites. Genes with significant positive selection signals detected in all three foreground-branch analyses were intersected and defined as the candidate gene set under shared positive selection in plant-parasitic nematodes.

2.7. Progressive Convergence Dynamics Analysis of Gene Family and KEGG Pathway Member Counts

To precisely identify the unique evolutionary features of the Aphelenchoidea and Anguinidae plant-parasitic nematodes and to evaluate their molecular convergence with the Tylenchoidea in plant-parasitic adaptation, we constructed a dual parallel fitting framework based on 16 representative nematode genomes. The analysis comprised two independent comparison groups: first, a core test set (Core Set) constructed from the Aphelenchoidea and Anguinidae plant-parasitic nematodes (N = 3, including B. xylophilus, A. besseyi, and D. destructor); and second, a reference set (Tylenchida Set) constructed from the Tylenchoidea plant-parasitic nematodes (N = 3, including G. pallida, H. schachtii, and M. graminicola). Both groups were independently compared against 10 non-plant-parasitic nematode outgroups (Non-PPN). Given that the lack of a well-established fossil record for nematodes makes accurate calibration of node divergence times difficult, traditional gene family dynamics software that relies on time-calibrated trees (such as CAFE) is limited in this system. We therefore introduced a statistical model that is insensitive to time calibration and can effectively handle data overdispersion (Quasipoisson GLM) to precisely capture convergent changes in pathway member counts.

After pre-filtering to remove low-variance terms in each comparison group (standard deviation SD < 0.5 or range < 2), a quasipoisson generalized linear model (Quasipoisson GLM) was used, given the discrete count nature of the data, to perform parameter fitting and Wald tests on gene family and KEGG pathway sizes. Taxa simultaneously satisfying an absolute log2 fold change |log2FC| ≥ 1 (corresponding to a ≥ 2-fold change in mean size between groups) and an uncorrected P < 0.05 were defined as significantly expanded or contracted.

Finally, by applying cross-logic discrimination to the two parallel fitting results, gene families and pathways showing significant differences were classified into three evolutionary paradigms: “Aphelenchoidea/Anguinidae-specific evolution”, showing significant changes only in the core test set; “suborder-wide shared universal response”, showing concordant significant changes in both the core set and the reference set; and “Tylenchoidea-specific paradigm”, showing significant changes only in the reference set.

3. Results

3.1. Phylogenetic Analysis and Independent Origins of the Parasitic Habit

The maximum likelihood (ML) phylogenetic tree constructed from the 1:1 single-copy orthologous genes of the seven core species showed a fully resolved topology, with all internal nodes receiving 100% bootstrap support (Figure 1a). The phylogeny clearly divided the core species into three main clades: clade I comprised the sister pair Aphelenchoides besseyi (plant-parasitic) and Aphelenchoides bicaudatus (non-plant-parasitic); clade II comprised the sister pair Bursaphelenchus xylophilus (plant-parasitic) and Bursaphelenchus okinawaensis (non-plant-parasitic); and clade III grouped and included Ditylenchus destructor (plant-parasitic), Aphelenchus avenae (non-plant-parasitic), and the free-living nematode Panagrellus redivivus, in which D. destructor and A. avenae formed a sister pair that subsequently clustered with the outgroup P. redivivus. In addition, the expanded-set phylogenetic tree comprising 16 species (Figure 1b) also achieved very high node support, and its topology was highly consistent with and robust relative to the seven-species tree in terms of core branching relationships, further confirming the reliability of the lineage divergence relationships described above.

This topology indicates that the three core plant-parasitic nematodes examined (A. besseyi, B. xylophilus, and D. destructor) do not form a monophyletic group; instead, they are dispersed among different evolutionary branches, and each plant-parasitic species forms a well-supported sister pair with a closely related non-plant-parasitic species (A. besseyi–A. bicaudatus; B. xylophilus–B. okinawaensis; D. destructor–A. avenae). This distribution pattern clearly indicates that the plant-parasitic habit has undergone at least three independent evolutionary origins among the taxa analyzed in this study. The well-supported topology provides a clear and reliable phylogenetic background for subsequent comparative genomic analyses and convergent site identification.

3.2. Characteristics of Orthologous Gene Clustering

The orthogroup clustering statistics showed that the proportion of clustered genes among the seven species ranged from 78.9% to 92.9%. The four species of the genera Aphelenchoides and Bursaphelenchus—A. besseyi (91.1%), B. xylophilus (92.1%), A. bicaudatus (92.9%), and B. okinawaensis (91.6%)—had the highest and highly similar clustering proportions. The non-plant-parasitic nematode A. avenae had a clustering proportion of 91.0%, comparable to the above species, but its proportion of species-specific genes was as high as 47.1%, markedly distinguishing it from the others. By contrast, the plant-parasitic nematode D. destructor (83.4%) and the free-living outgroup P. redivivus (78.9%) had noticeably lower clustering proportions, with correspondingly higher proportions of unassigned genes (16.6% and 21.1%, respectively), suggesting that the genomes of these two contain more lineage-specific or rapidly evolving gene sequences (Table 1). Within the Aphelenchoides branch (A. besseyi and A. bicaudatus), the clustering proportions of the two were similar (91.1% vs 92.9%), but the proportions of species-specific genes differed noticeably (8.4% vs 4.1%); within the Bursaphelenchus branch (B. xylophilus and B. okinawaensis), the proportion of specific genes in the plant-parasitic species B. xylophilus (8.9%) was higher than that in its non-parasitic sister species B. okinawaensis (5.2%). Taken together, global orthologous clustering did not reveal a uniform pattern of large-scale lineage-specific new gene expansion in plant-parasitic nematodes, and different lineages exhibited relatively independent evolutionary trajectories in terms of gene family conservation and specificity.

3.3. Global Divergence of Selection Pressure Between Plant-Parasitic and Non-Parasitic Lineages

After bidirectional best BLASTp alignment combined with strict filtering (coverage ≥ 50%, identity ≥ 30%), the combinations B. xylophilus–B. okinawaensis (Bxyl vs Boki) and A. besseyi–A. bicaudatus (Ab vs Abic) yielded 13,651 and 9,584 high-confidence 1:1 orthologous genes, respectively (Table 1). Because the divergence times between the two species pairs are relatively recent, a considerable proportion of orthologous gene pairs had a synonymous substitution rate Ks = 0, rendering the Ka/Ks ratio unestimable (65.7% and 64.3% in the two combinations, respectively). Therefore, subsequent selection pressure analyses were performed independently on the estimable effective subsets (n = 4,681 and n = 3,424, respectively).

The Ka/Ks distributions of both combinations showed a significant right-skewed pattern: for the Bx vs Boki combination, the mean Ka/Ks (0.122 ± 0.145) was markedly higher than the median (0.082), with 57.5% (2,690/4,681) of genes having Ka/Ks < 0.1; for the Abes vs Abic combination, the median was 0.156 (mean 0.180 ± 0.109), with 25.8% (884/3,424) of genes having Ka/Ks < 0.1 (Table 2). These results indicate that the genomes of both nematode pairs are overall dominated by strong purifying selection, and that the establishment of a plant-parasitic lifestyle did not fundamentally alter the nature of selection pressure at the genome-wide level.

At the same time, there was a systematic difference in the strength of selective constraint between the two combinations. The median Ka/Ks of the Abes vs Abic combination was approximately 1.9 times that of the Bxyl vs Boki combination (0.156 vs 0.082), and the proportion of genes with Ka/Ks ≥ 0.1 increased from 41.8% to 73.0%. Notably, very few genes in either combination satisfied Ka/Ks ≥ 1.0 (only 3 and 0, respectively, accounting for 0.06% and 0% of the effective subsets), suggesting that continuous amino acid substitutions driven by classic positive selection are extremely rare at the genome-wide level, and that the moderate elevation of ratios in most genes more closely resembles a relaxation of the precision of purifying selection.

3.4. Convergent Evolution at Amino Acid Sites and in Physicochemical Properties

Using the CSS method, 120 potential amino acid convergent sites (distributed among 112 orthologous genes) were initially identified among the three independently evolved plant-parasitic nematodes (B. xylophilus, A. besseyi, and D. destructor). After strict filtering using the 16-species expanded validation set, 10 genes with high-confidence convergent sites were ultimately identified (Figure 2). These comprised 6 suborder-wide shared convergent genes (alkaline ceramidase, asparagine synthetase, DOMON domain-containing protein, G protein-coupled receptor kinase, protein disulfide-isomerase, and Tr-type G domain-containing protein) and 4 convergent genes occurring specifically in the migratory endoparasitic lineage (ATP synthase subunit d, MARVEL domain-containing protein, meckelin, and Suf domain-containing protein). In the amino acid property convergence analysis (GS1–GS4 grouping), 2 property convergent genes occurring specifically in migratory endoparasitic nematodes were ultimately identified, namely zinc/iron permease and lipid-binding serum glycoprotein C-terminal domain protein (Figure 2). The neutral random substitution model test showed that the numbers of convergent sites in the 10 high-confidence amino acid convergent genes and the 2 physicochemical property convergent genes described above were all significantly higher than the random neutral expectation estimated from the phylogenetic tree topology and substitution matrix (Zou and Zhang model, P < 0.05), indicating that the convergent evolution at these sites was driven by strong natural selection rather than random mutational noise.

The positive selection test on foreground branches (Branch-site Model A) revealed significant differences in the number of positively selected genes detected among the three plant-parasitic nematode lineages: only 14 positively selected genes were identified in B. xylophilus, whereas 147 and 105 were identified in D. destructor and A. besseyi, respectively (Table S2–S4). Cross-species intersection analysis and KEGG enrichment analysis showed that there was no shared gene under positive selection common to all three plant-parasitic nematodes, and that none of the species-specific positively selected gene sets showed significantly enriched metabolic pathways. This negative result indicates that positive selection is highly dispersed in a lineage-specific manner, and that the molecular mechanisms driving the transition of nematodes to plant parasitism rely more on the precise convergence of specific key sites and the adjustment of macro-scale pathway sizes than on full-length positive selection across large numbers of genes genome-wide.

3.5. Convergent Evolution and Paradigm Divergence of KEGG Metabolic Pathways

Based on parameter fitting and dual parallel discrimination of the sizes of 568 KEGG pathways across 16 nematode genomes using the Quasipoisson GLM model, the pathways were classified into six evolutionary paradigms (Table S5 and Table S6). The results showed that the Aphelenchoidea and Anguinidae (Core Set) were relatively conserved in pathway size evolution, exhibiting overall local fine-tuning, whereas the Tylenchoidea (Tylenchida Set) underwent extensive and profound genome-wide pathway size reduction.

In terms of specific pathway features, only the unsaturated fatty acid biosynthesis pathway (ko01040) showed simultaneous significant contraction in both the Core Set and the Tylenchida Set (both satisfying |log2FC| ≥ 1, P < 0.05), indicating a high degree of evolutionary convergence in the reliance of different plant-parasitic nematode lineages on host-derived lipid resources and the streamlining of their own fatty acid synthesis metabolism. The Aphelenchoidea and Anguinidae-specific evolution was extremely conserved, with only 2 specifically contracted pathways (including the two-component system ko02020) and no significantly expanded pathways. By contrast, the Tylenchoidea exhibited large-scale functional reduction, with as many as 204 KEGG pathways showing specific contraction (accounting for 35.92%), covering numerous basic metabolic and protein degradation processes and reflecting the deep streamlining of genome-wide functions during its highly obligate parasitic specialization, while only 1 specifically expanded pathway (ko00524) was found (Table 3).

4. Discussion

4.1. Biological Causes and Evolutionary Logic Underlying the Weak Molecular Convergence Signal in Plant-Parasitic Nematodes

The molecular convergence signal detected in plant-parasitic nematodes in this study was generally weak: high-confidence convergence events were limited whether assessed by convergent substitutions at the amino acid level, convergence of selection pressure, or concordant changes in orthologous gene family sizes. This phenomenon does not stem from insufficient analytical sensitivity, but rather reflects, at a deeper level, the special evolutionary landscape of the transition of nematodes to plant parasitism.

On the one hand, the plant-parasitic habit did not originate from a single ancient monophyletic divergence within the phylum Nematoda, but rather evolved independently multiple times from non-plant-parasitic groups in different lineages [25]. The non-plant-parasitic control nematodes and the parasitic nematodes selected in this study are extremely closely related in evolutionary terms; for example, B. xylophilus and B. okinawaensis, and A. besseyi and A. bicaudatus, are sister species, having diverged relatively recently and sharing the vast majority of their ancestral genetic background. This high degree of background genomic similarity means that the very “opportunity” for parasitic species to undergo independent convergent substitutions at the sequence level is itself substantially compressed—the overwhelming majority of conserved genes retain their original purifying selection status, and the amino acid convergence signal that can be reliably captured by statistical models is therefore severely diluted [26].

On the other hand, from the perspective of life-history strategy, sedentary endoparasitic and migratory endoparasitic nematodes exhibit significant molecular divergence in host localization, tissue invasion strategies, the formation and maintenance mechanisms of feeding sites, and the composition of their effector repertoires. The systematic comparison by Qing et al. of the cell wall-degrading enzyme (CWDE) gene families across 8 nematode species provides direct evidence for this [27]. Notably, even D. destructor and A. besseyi, which are likewise migratory plant-parasitic nematodes, have undergone large-scale reorganization of their chromosome-level genome architecture relative to B. xylophilus, whereas B. xylophilus has maintained an ancestral architecture similar to that of free-living nematodes [6]. This implies that the three independently evolved parasitic lineages each explored parasitic strategies against distinctly different genomic structural backgrounds—if even chromosome-level genome architecture failed to converge, the scarcity of convergence signals at the protein sequence level is all the more logical.

This divergence pattern further suggests that much adaptive evolution tends to occur at higher-order functional levels such as gene expression regulatory networks, epigenetics, and metabolic pathway rewiring, rather than in the linear sequences of individual proteins [28]. Thus, the convergent evolution of plant-parasitic nematodes is more strongly reflected in convergence at the level of system function, and the scarcity of convergence at the sequence level is precisely consistent with their diverse mechanistic implementations.

4.2. A “General-Purpose Parasitic Molecular Arsenal” for Breaching Plant Defense Barriers and the Limits of Convergent Evolution

Although the Aphelenchoidea and Tylenchoidea represent completely independent origins of plant parasitism, the two groups appear to share a broadly predictable set of molecular features and evolutionary constraints for overcoming plant physical barriers, countering host innate immunity, and establishing basic parasitic relationships. At the macro-scale metabolic genomic level, analysis based on the Quasipoisson GLM model revealed that suborder-wide shared contracted pathways are extremely rare: among 568 KEGG pathways, only the unsaturated fatty acid biosynthesis pathway (ko01040) showed significant contraction in both the Aphelenchoidea and Tylenchoidea. Because the de novo synthesis of unsaturated fatty acids requires large amounts of reducing power (NADPH) and molecular oxygen, the concordant early loss of this energy-intensive pathway in parasitic nematodes raises the hypothesis that the transition to a plant-parasitic lifestyle may be accompanied by metabolic streamlining of endogenous lipid synthesis and a greater reliance on exogenous, host-derived lipid components [29]; this proposed lipid-acquisition scenario remains to be directly demonstrated.

To successfully infect their hosts, the two groups of nematodes underwent a certain number of site-level convergences at the micro-scale protein sequence level. Loss of GRK function in C. elegans causes chemosensory abnormalities, a mechanism involving inhibition of signal transduction by an inhibitory pathway [30]. Chemosensation is a core capability for host localization in plant-parasitic nematodes, and GRK convergence may reflect a convergent regulatory demand on the dynamics of sensory adaptation. The function of PDI in plant-parasitic nematodes has been experimentally validated. GpPDI1 was identified as an esophageal gland-secreted effector protein of G. pallida capable of triggering host defense-related cell death [31]. HsPDI is similarly localized to the esophageal gland in H. schachtii, and its silencing leads to significant ultrastructural changes in the syncytium, indicating that PDI is functionally required for feeding site establishment [32]. These findings support a secretory role for PDI during parasitism; we hypothesize that its convergent sites may influence both intracellular effector protein folding and extracellular host redox interference, although this remains speculative. Likewise, convergent variation in asparagine synthetase aligns with the physiological feature of plants using asparagine and glutamine as the principal nitrogen transport media [33], raising the possibility, as a candidate hypothesis, that such variation may relate to the ability of parasitic nematodes to assimilate and utilize specific amino acids in plant vascular sap. In addition, alkaline ceramidase participates in sphingolipid metabolism, disease resistance, and stress responses in plants [34]. Ceramide metabolism has been shown to affect the formation and stability of lipid raft microdomains. Convergent substitutions in nematode alkaline ceramidase might, as one hypothesis, relate to maintaining membrane microdomain stability against the host’s lipid-based toxic defense. Taken together, these site-level convergences are consistent with the hypothesis that, in crossing the major evolutionary transition to parasitism, nematodes may face broadly similar physiological barriers and thus converge along partly predictable evolutionary paths.

4.3. The “Dual-Track” Evolutionary Paradigm of Sedentary Endoparasitic and Migratory Parasitic Nematodes

After breaching the general parasitic barriers, parasitic nematodes, driven by differences in selection pressure, followed markedly divergent evolutionary paradigms along two dimensions—macro-scale genome size and micro-scale sequence adaptation—consistent with distinct adaptations of sedentary endoparasites and migratory parasites to their respective ecological niches and life-history strategies. In this study, as many as 204 KEGG pathways in the Tylenchoidea showed specific significant contraction, covering core pathways such as glycosaminoglycan degradation, RNA polymerase, and the proteasome, while only one specifically expanded pathway was found. After invading plant roots, these nematodes induce drastic reprogramming of host cells through the secretion of effector proteins, forming highly specialized feeding sites such as giant cells or syncytia [35]. These specialized cells serve as a stable nutrient sink continuously supplying free amino acids, nucleotides, and sugars, greatly reducing the selection pressure on nematodes to acquire nutrients autonomously. In this relatively stable host internal environment, the genomes of Tylenchoidea nematodes appear to have undergone large-scale reduction [36], which we suggest may have driven their genome-wide functions toward greater specialization and dependence on the host.

By contrast, only very few pathways underwent specific contraction in the Aphelenchoidea, whose metabolic pathway sizes exhibited extremely high conservation during evolution. This macro-scale functional retention may be closely related to their facultative and complex multi-host life history: Aphelenchoidea nematodes are thought to alternate between fungivory, infection of the aerial xylem or mesophyll of plants, and dispersal via vector insects, which may require retention of a broadly autonomous metabolic capacity to maintain plasticity [25]. While maintaining macro-scale metabolic flexibility, migratory endoparasitic nematodes (such as the Aphelenchoidea) appear to have remodeled stress-resistance and host-lipid-related functions at the micro-scale protein sequence level. We propose, as a candidate hypothesis, that under the reactive oxygen species (ROS) stress and toxic secondary metabolites encountered in plant tissues, the convergently substituted ZIP zinc/iron permease and Suf domain-containing protein may act together in mitochondrial stress mitigation: the former putatively transports transition metal ions [37], potentially supplying substrates for the Suf-mediated reassembly of iron-sulfur (Fe-S) clusters, which could help repair ROS-induced metalloprotein damage [38]. Together with the possible role of ATP synthase subunit d in maintaining cristae morphology and electron transport chain complex stability [39], these genes suggest a putative mitochondrial functional-protection module; however, direct functional evidence is currently lacking. At the level of lipid metabolism and immune regulation, we hypothesize that the specific convergence of the lipid-binding serum glycoprotein C-terminal domain protein may confer the potential to bind plant oxylipins and jasmonic acid (JA) precursors, which could in turn modulate host immune defense and/or host-lipid acquisition during life-history switching between fungal and plant parasitism [40]. Similarly, the putative strengthening of lipid raft microdomain transport by the MARVEL domain-containing protein [41] and the remodeling of the sensory neuron ciliary transition zone by the meckelin protein [42] together suggest a possible mechanism for chemotactic perception in complex host media; these remain candidate hypotheses pending functional validation.

In summary, this study suggests a dual-track pattern in the evolution of plant-parasitic nematodes: the Tylenchoidea appear to have followed a path of sedentary endoparasitism with extensive genome streamlining and dependence on host-specialized cells; whereas the Aphelenchoidea, while retaining broader metabolic versatility, may have adapted to high-stress aerial environments and multiple hosts through micro-scale modification of key mitochondrial, lipid, and sensory genes. These interpretations are offered as hypotheses to be tested in future functional studies.

4.4. Limitations and Future Perspectives

Although this study systematically elucidated the cross-superfamily molecular convergence and dual-track evolutionary paradigm of plant-parasitic nematodes, specific methodological and data bottlenecks remain. At the dataset level, owing to the current lack of high-quality closely related fungivorous or free-living sister-pair genomes for the Tylenchoidea, the amino acid site convergence (CSS) screening could not be carried out at an equivalent depth of paired comparison in this group, which to some extent limited the statistical detection power for Tylenchoidea-specific sequence convergent sites. At the statistical model level, in pathway evolution analysis based on the Quasipoisson GLM, the uneven distribution of non-parasitic outgroups on the phylogenetic tree introduced some interference into the test power for contraction and expansion of small-scale gene families. In addition, the current inferences rely mainly on static alignments of comparative genomics coding sequences and still lack support from transcriptomic and protein expression dynamic data across the spatiotemporal scale of infection.

To address these limitations, future research will combine single-cell transcriptomics with gene editing technologies for in-depth validation. On the one hand, following the single-nucleus transcriptome (snRNA-seq) and spatial transcriptome technical frameworks successfully established in recent years for the root-knot nematode Meloidogyne incognita and the pine wood nematode Bursaphelenchus xylophilus [43], we will precisely resolve the cell specificity and spatiotemporal expression dynamics during infection of core convergent genes such as PDI, GRK, and DOMON in nematode esophageal gland cells and head sensory neurons. On the other hand, using RNAi gene silencing and a CRISPR-Cas9-based site-directed mutagenesis system [25], we will perform in vivo functional knockout of key genes such as PDI, GRK, and ZIP zinc/iron permease as well as reversion validation of parallel mutation sites, to quantitatively assess their physiological functions in blocking host-seeking chemotaxis, effector protein folding, and antioxidant stress.

In terms of practical translation, the suborder-wide shared convergent sites identified in this study (such as the parallel variation of PDI and GRK) and the concordant contraction of the unsaturated fatty acid biosynthesis pathway may point to potential vulnerabilities associated with the transition to plant parasitism. Such shared targets, if functionally validated, might represent conserved nodes that are less easily bypassed through single-gene point mutations. We therefore suggest, as a long-term prospect rather than a demonstrated outcome, that broad-spectrum RNAi-based nematicides or small-molecule inhibitors targeting these sites and metabolic features could complement existing control strategies and help mitigate resistance issues in the control of cross-superfamily plant nematode diseases.

5. Conclusions

Through comparative genomics analysis of representative nematodes of the suborder Tylenchina, this study provides genomic evidence for cross-superfamily molecular convergence and lineage-specific divergence in plant-parasitic nematodes. The study shows that although the two groups of nematodes represent independent origins of parasitism, they appear to share a broadly predictable set of molecular features associated with crossing host defense barriers. The concordant contraction of the unsaturated fatty acid biosynthesis pathway suggests, as a hypothesis, a convergence toward greater dependence on exogenous lipids; at the functional-protein level, the site convergence observed in host recognition, effector protein folding, plant nitrogen assimilation, and interfacial stress perception may constitute part of the molecular framework underlying parasitic adaptation, and warrants functional confirmation. After crossing these basic barriers, the parasitic lineages appear to have diverged into a dual-track pattern shaped by life-history and niche differentiation. The Tylenchoidea, represented by sedentary endoparasites, showed evidence of large-scale genome streamlining and reduced autonomous metabolic functions, consistent with adaptation to the stable internal environment of specialized root feeding sites; whereas the migratory parasitic Aphelenchoidea retained broader metabolic versatility while exhibiting micro-scale remodeling of stress-resistance networks, which may contribute to adaptation to high-stress aerial environments and multiple hosts. This contrast between macro-scale reduction and micro-scale specialization may help inform understanding of the origin of parasitism and adaptive radiation in invertebrates, and offers candidate genomic targets that merit further validation for the development of broadly effective control strategies against nematodes.