Hybridization and polyploidy enable genomic plasticity without sex in the most devastating plant-parasitic nematodes

Abstract

Root-knot nematodes (genus Meloidogyne) exhibit a diversity of reproductive modes ranging from obligatory sexual to fully asexual reproduction. Intriguingly, the most widespread and devastating species to global agriculture are those that reproduce asexually, without meiosis. To disentangle this surprising parasitic success despite the absence of sex and genetic exchanges, we have sequenced and assembled the genomes of three obligatory ameiotic and asexual Meloidogyne. We have compared them to those of relatives able to perform meiosis and sexual reproduction. We show that the genomes of ameiotic asexual Meloidogyne are large, polyploid and made of duplicated regions with a high within-species average nucleotide divergence of ~8%. Phylogenomic analysis of the genes present in these duplicated regions suggests that they originated from multiple hybridization events and are thus homoeologs. We found that up to 22% of homoeologous gene pairs were under positive selection and these genes covered a wide spectrum of predicted functional categories. To biologically assess functional divergence, we compared expression patterns of homoeologous gene pairs across developmental life stages using an RNAseq approach in the most economically important asexually-reproducing nematode. We showed that >60% of homoeologous gene pairs display diverged expression patterns. These results suggest a substantial functional impact of the genome structure. Contrasting with high within-species nuclear genome divergence, mitochondrial genome divergence between the three ameiotic asexuals was very low, signifying that these putative hybrids share a recent common maternal ancestor. Transposable elements (TE) cover a ~1.7 times higher proportion of the genomes of the ameiotic asexual Meloidogyne compared to the sexual relative and might also participate in their plasticity. The intriguing parasitic success of asexually-reproducing Meloidogyne species could be partly explained by their TE-rich composite genomes, resulting from allopolyploidization events, and promoting plasticity and functional divergence between gene copies in the absence of sex and meiosis.

Fig 1. Distribution of CDS mapping one to several loci in the Meloidogyne genomes.

Occurrences (y axis) of Meloidogyne CDS mapping at minimum 95% identity on minimum 2/3 of their length to one, two, three, four or up to 15 loci (x axis) in their respective genomes. In M. hapla (red), >89% of the CDS map to one single locus while >85% of the CDS map to multiple loci in the apomicts M. incognita (green), M. javanica (blue) and M. arenaria (violet).

Fig 2. A genome structure consistent with absence of meiosis.

Five pairs of duplicated homologous collinear regions co-occur on a same scaffold in M. incognita (A) and five in M. arenaria (B). All curves show connections between collinear gene pairs used by MCScanX to define duplicated regions (blue curves show tandem duplications and purple curves show palindromic duplications). Grey lines represent the gene density on the scaffolds.

Fig 3. Example of structural and evolutionary relationship between pairs of duplicated regions.

A. Circos [31] plot showing the collinear gene pairs (forming homologous regions) that were used for phylogenetic analyses (units = kb). All curves show the connections between the collinear gene pairs used by MCScanX to define segmental duplications. In each Circos plot, color codes are as follows. Collinear orthologs between M. hapla and any of the three asexuals species are in grey. Collinear ‘homoeologs’ within asexual species are in purple. Collinear orthologs between M. arenaria and M. javanica are in green. Collinear orthologs between M. arenaria and M. incognita are in yellow. Collinear orthologs between M. incognita and M. javanica are in red. The outer blue lines represent the gene density on the scaffolds. B. Maximum-likelihood phylogeny of concatenated alignments of collinear protein-coding genes used to form blocks with SH-like branch support. Topologies identical to the mitochondrial phylogeny were considered to represent the maternal contribution to the nuclear genome. The other topologies were considered as representative of paternal contributions.

Fig 4. Phylogenetic relationships between duplicated regions in the genomes of apomictic Meloidogyne.

The three possible topologies for bifurcating trees separating M. incognita, M. javanica, M. arenaria and their sexual relative M. hapla are represented as (1), (2) and (3) and their relative observed frequencies are indicated in the associated pie chart. The frequencies were calculated from the 40 phylogenetic trees containing at least one monophyletic clade with the 3 apomictic Meloidogyne that were constructed from the concatenated alignments covering a total of 2,202 protein-coding genes.

Fig 5. Combinations of topologies for duplicated regions in the three parthenogenetic Meloidogyne.

(A) Schematic representation of two duplicated regions each containing 3 collinear genes (a, b, c and a’, b’, c’). (B) ML trees combining twice the same topology among (1), (2) and (3), for the two duplicated regions (further detailed in Fig 4). These trees suggest the two duplicated regions have the same evolutionary history. (C) ML trees combining two different topologies for the duplicated regions. These trees suggest the two duplicated regions have different evolutionary histories. The relative frequencies of trees combining twice the same (red) or two different (blue) topologies are indicated in the big central pie chart. Relative frequencies within the red and blue categories are indicated by small pie charts next to the corresponding schematic tree. (D) The most frequently observed ML trees consist in a combination of topologies (1) and (2).

Fig 6. Allele-like versus homeologous relationship between genes in more than two copies.

(A) An example of homoeologous relationship where each M. arenaria gene copy (arrows) clusters with the copy of another species. (B) An example of allele-like relationship where two of the three M. incognita gene copies (arrows) are more similar to one another than they are to a copy from another species.

Fig 7. Mitochondrial phylogeny.

Consensus phylogeny obtained using ML and Bayesian analyses on 14 concatenated mitochondrial genes (12 protein-coding and 2 rRNA). Posterior probability (above) and bootstrap (below) support values are given at each corresponding branches. Phylogenies were rooted with Pratylenchus vulnus as an outgroup. The tree is represented without taking branch lengths into an account for better visibility of phylogenetic relations; the same tree with actual branch lengths is available in S6 Fig. Species with a (L) suffix indicate sequences coming from the genome sequencing effort undertaken as part of this paper. Species with a (D) suffix indicate sequences coming from other databases. Full species names and accession numbers for the sequences coming from external databases are given in S1 Table.

Fig 8. Distribution of the Ka / Ks ratio for pairs of collinear genes in asexual Meloidogyne.

Histogram displaying the distribution of the ratio of rates of non-synonymous (Ka) / rates of synonymous (Ks) nucleotide substitutions for pairs of duplicated collinear genes in M. incognita, M. javanica and M. arenaria. Pairs with a Ka / Ks > 1 indicate positive selection between genes.

Fig 9. Distribution of M. incognita genes and copies in the expression clusters.

For each of the 24 expression clusters, the information ‘n = a/b/c’ indicates. a: the total number of M. incognita genes having this expression pattern. b: the number of times the two genes from a same pair have this same expression pattern. c: the number of genes part of a pair having this pattern, while the other cognate gene has a different expression pattern (i.e. belongs to another expression cluster). Hence, c refers to gene pairs with diverged expression patterns. The expression ranks across the four developmental life stages (eggs, J2 infective juveniles, J3-J4 larval stages and adult females) are represented in colors as follows rank1: lowest expression (blue), rank2: second lowest expression (green), rank3: second highest expression (orange) and rank4: highest expression (red).

Fig 10. Pfam domain abundance in asexual Meloidogyne as a function of the abundance in M. hapla.

The x-axis represents abundance of Pfam domains in M. hapla and the y-axis represents abundance of the same domains in M. incognita (red diamonds) M. javanica (green squares) and M. arenaria (yellow triangles). Linear regressions of M. incognita, M. javanica and M. arenaria Pfam domain abundances as a function of the abundance in M. hapla are plotted alongside their respective equations and correlation coefficients. Pfam domains associated to transposable elements (TE) are represented as crossed red, green and yellow squares in M. incognita, M. javanica and M. arenaria, respectively.