Gene gain and loss during evolution of obligate parasitism in the white rust pathogen of Arabidopsis thaliana

Abstract

Biotrophic eukaryotic plant pathogens require a living host for their growth and form an intimate haustorial interface with parasitized cells. Evolution to biotrophy occurred independently in fungal rusts and powdery mildews, and in oomycete white rusts and downy mildews. Biotroph evolution and molecular mechanisms of biotrophy are poorly understood. It has been proposed, but not shown, that obligate biotrophy results from (i) reduced selection for maintenance of biosynthetic pathways and (ii) gain of mechanisms to evade host recognition or suppress host defence. Here we use Illumina sequencing to define the genome, transcriptome, and gene models for the obligate biotroph oomycete and Arabidopsis parasite, Albugo laibachii. A. laibachii is a member of the Chromalveolata, which incorporates Heterokonts (containing the oomycetes), Apicomplexa (which includes human parasites like Plasmodium falciparum and Toxoplasma gondii), and four other taxa. From comparisons with other oomycete plant pathogens and other chromalveolates, we reveal independent loss of molybdenum-cofactor-requiring enzymes in downy mildews, white rusts, and the malaria parasite P. falciparum. Biotrophy also requires "effectors" to suppress host defence; we reveal RXLR and Crinkler effectors shared with other oomycetes, and also discover and verify a novel class of effectors, the "CHXCs", by showing effector delivery and effector functionality. Our findings suggest that evolution to progressively more intimate association between host and parasite results in reduced selection for retention of certain biosynthetic pathways, and particularly reduced selection for retention of molybdopterin-requiring biosynthetic pathways. These mechanisms are not only relevant to plant pathogenic oomycetes but also to human pathogens within the Chromalveolata.

Figure 1. Genomic sequencing data and coverage of A. laibachii Nc14 and Em1 assemblies.

(A) Reads generated for A. laibachii Nc14 using Illumina genome analyzer version 1 (GA1) or version 2 (GA2). (B) Distribution of genomic coverage. Grey fields indicate the total amount of sequence represented by the 100-bp windows with corresponding coverage. (C) Reads generated for A. laibachii Em1. (D) Distribution of genomic coverage showing Em1 reads aligned to the Nc14 genome using MAQ aligner. Nc14 and Em1 show a major peak at 226× and 43× coverage, respectively. A second peak is detected at 112× or 22×, showing half the coverage of the main peak, indicating highly heterozygous regions that were not merged in the assembly or hemizygous regions.

Figure 2. Distribution of contig length, nucleotide coverage, SNP frequency, and cDNA coverage in the A. laibachii assembly.

(A) Genomic distribution of contig length (N length) versus contig number (N number). N lengths were calculated by ordering all sequences according to their length and then adding the length from longest to shortest until the summed length exceeded 10% (N10), 20% (N20), etc., up to 100% (N100) of the assembled contigs (32.7 Mbp). Blotting the N length versus the N number (number of contigs in each N category) indicates that 90% of the assembled genome show high continuity, while the last 10% are highly fragmented. (B) Average coverage for each category for Nc14 (red) and Em1 (green). In all, 90% of the genome shows low variation, consistent with 210–240× coverage for Nc14 and 40–50× for Em1. The last 10% show highly elevated coverage, indicating unresolved highly repetitive regions present in Nc14 and Em1. (C) Distribution of heterozygosity in each N category shows elevated levels in the set of short contigs. Heterozygous positions were accepted only if coverage was >180× and <350× for Nc14 (red) or >27× and <80× for Em1 (green). SNPs between Nc14 and Em1 were calculated ignoring heterozygous positions (lilac). (D) Alignment of Nc14 cDNA and summing up all regions showing >2× coverage indicate that the more continuous part of the genome contains more transcribed regions than the highly repetitive regions of the genome (in the histogram, N length and N number are cumulative while read depth, SNP frequency, and cDNA coverage are presented as binned data).

Figure 3. Repeats identified in the A. laibachii Nc14 contigs.

Initial run of RepeatScout produced a library of 1,252 consensus repetitive sequences that include transposable elements, recently duplicated paralogous genes, and other dispersed duplicated regions. (A) Inset: The distribution of lengths of the identified repeats versus their frequency in the genome is shown. The majority of repeats fell into the category of short and rare in the assembly. The primary plot in (A) shows that the majority of the longest and most frequent repeats in the genome are transposon elements (shown in green and Table 1), while Albugo-specific repetitive sequences are mostly short (shown in red). (B) Summary of the proportion of the repetitive sequences (percent) in the A. laibachii Nc14 genome.

Figure 4. A. laibachii has a compact genome with expression clusters.

(A) Synteny between A. laibachii, Py. ultimum, H. arabidopsidis, and P. infestans. The region shown is an example of the dense clustering of genes in the pentafunctional AROM polypeptide and a P-type ATPase. The AROM polypeptide comprises five enzymes of the shikimate pathway in one enzyme. With increasing genome size the distance between both genes increases and re-organisations occur (red, synteny without inversion; blue, inverted regions). (B) Plotting the distance between transcriptional islands based on the 5′ to 3′ orientation of the forward strand reveals that transcriptional regions are clustered close together. The maximum peak reflects the average intron size. Regions with no 3′ but with 5′ distance and vice versa reflect overlapping 3′ and 5′ non-coding regions of genes. Analysing the distance distribution between transcriptional units reveals a median distance between genes of 45 bp, showing that within transcribed regions, nearly all the DNA sequence corresponds to genes. (C) Plotting the 5′–3′ distance for all genes from ATG to stop to the next gene confirms the gene clustering. Only 10.8% of all genes have a distance to the next gene or the end of the contig greater than 3 kb. Summing the distance between these genes contributes to only 10.9 Mbp of the genome because of the close packaging, while summing the distance of the few genes that are not in clusters contributes to 8.4 Mbp of the genome.

Figure 5. Molecular divergence between A. laibachii and other species based on pairwise comparisons.

(A) Molecular divergence based on all pairwise comparisons of the one-to-one orthologues. In the figure, the cumulative frequencies of amino acid identity across each set of potential orthologous pairs is presented, indicating that although H. arabidopsidis and A. laibachii are both biotrophs, H. arabidopsidis is less diverged from P. infestans than it is from A. laibachii (e.g., in the H. arabidopsidis–A. laibachii comparison, ∼22% of all orthologues show an amino acid identity of <50%, while only ∼14% in a Py. ultimum–A. laibachii comparison show an amino acid identity of <50%). A. laibachii shows the highest amino acid identity to Py. ultimum. (B) Molecular divergence between A. laibachii and other species based on the subset of core eukaryotic genes to show stability of the test. Results are consistent with the one-to-one orthologue analyses although differences between A. laibachii, P. infestans, H. arabidopsidis, and Py. ultimum are less obvious, indicating the lack of selection pressure on the core eukaryotic genes . For comparative reasons, a tree using ITS2 sequences is added. The represented tree is a maximum likelihood tree produced with PhyML.

Figure 6. Validation and identification of potential delivery motifs.

To identify potential effector delivery motifs we analysed RXLR (A), RXLQ (B), and the new CHXC (C) effector candidates for enrichment in the secretome. Motif shuffling was used to identify background levels. Calculating the cumulative hypergeometric probability to analyse the enrichment of secreted proteins (red) over non-secreted proteins (blue) for each of the permutated motifs reveals a significant enrichment of CHXCs in the secretome (p[X≥x] = 2.19×10⁻²²). None of the RXLR or RXLQ motifs or permutations shows significant enrichment. There is also enrichment for CXHC proteins in the secretome. Except for one CHHC protein, there is no overlap between the two motif classes. The logo blot (D) clearly indicates that RXLR-containing proteins are conserved only within the selected amino acids, while for CHXCs it is not only the motif but also sequences C-terminal to it are conserved, including conserved glycine, leucines, and a tyrosine.

Figure 7. P. capsici test for delivery motif.

To identify known and new classes of transfer motifs in A. laibachii, the P. capsici–N. benthamiana translocation assay was used. This test system is based on Avr3a-mediated avirulence in plants carrying R3a . (A) virulence assay to show that transgenic P. capsici is not impaired in growth; (B) Hypersensitive response (HR) assay on R3a-carrying plants for delivery assay. The Avr3a RxLR translocation domain is replaced by the N-terminus of different A. laibachii effector candidates. For the assay, RXLR1, CRN3, and CHXC9 carrying the newly identified CHXC motif were used. Our results validate that known motifs like the CRN motif are functional while the selected RXLR shows low delivery efficiency. CHXC9 shows the same efficiency as Avr3a does, and dependency of the CHXC motif could be identified (statistical analyses using the Tukey test; means with the same letter are not significantly different; error bars denote standard error of the mean). wt, wild type.

Figure 8. Candidate A. laibachii Nc14 effectors contribute to Pst DC3000 virulence.

(A) Arabidopsis plants (4- to 5-wk-old) were spray inoculated with 5×10⁸ CFU Pst DC3000 lux harbouring candidate effectors cloned in pEDV6. Bacterial growth was measured as an increase in luciferase photon emission per gram fresh weight per second (photon/g[fw]/sec). The histogram represents the log median of photon emission of three independent experiments, each with four technical replicates. Error bars denote standard error of the mean. Two-way ANOVA: #, p<0.001; **, p<0.01; *, p<0.05 from AvrRps4(AAAA). (B) Plants 4- to 5-wk-old were infected with 5×10⁸ CFU of Pst DC3000 ΔAvrPto/ΔAvrPtoB harbouring candidate effector cloned in pEDV6. Bacterial populations were sampled 4 d post-inoculation. The histogram represents the median colony count of two independent experiments, each with more than four technical replicates. Error bars denote standard error of the mean.

Figure 9. Result of neighbour-joining analyses using N-termini of all predicted CHXCs or CXHCs from the genomes of P. infestans, Py. ultimum, H. arabidopsidis, T. pseudonana, Pl. falciparum, E. siliculosus, C. merolae, Ch. reinhardtii, V. carteri, S. parasitica, as well as Ar. thaliana.

The outer ring summarises clades with N-termini predominantly carrying CHXC or CXHC motif or mixed clades (CXHC/CHXC) into classes. A. laibachii CHXCs are mainly clustered in the CHXC class (green), containing besides A. laibachii distantly related CHXCs from S. parasitica, V. carteri, Ch. reinhardtii, and Ar. thaliana. CHXCs are distant from endoplasmic reticulum proteins like disulphide isomerases that predominantly carry the CXHC motif and are grouped within the CXHC class (red). Between the CHXC class and the CXHC class, mixed clades contain protease and defensin homologues (orange) or Ar. thaliana cystein-rich proteins (violet). (Names in green indicate A. laibachii CHXCs and in yellow, A. laibachii CXHCs. Blue indicates CHXCs from other species; magenta indicates CXHCs from other species; 16 amino acids before and 45 amino acids after the CHXC or CXHC motif in the N-terminus were used. The tree is midpoint rooted. All bootstrap counts refer to 1,000 replications.). Ath, Ar. thaliana; Cla, Ch. reinhardtii; Cme, C. merolae; Ect, E. siliculosus; Hpa, H. arabidopsidis; Pfa, Pl. falciparum; Pin, P. infestans; Pul, Py. ultimum; Spr, S. parasitica (Spr); Tps, T. pseudonana; Vca, V. carteri.

Figure 10. Gain and loss of genes and pathways for selected Chromalveolata in comparison to A. laibachii.

It was hypothesized that the last common ancestor of Chromalveolata was a brown-alga-like organism with genes from green and red algae integrated into the nuclear genome after primary and secondary endosymbiosis ,. While some heterokonts kept their secondary endosymbiont and, in the case of diatoms, acquired a silicated bipartite cell wall , others lost their secondary endosymbiont. We postulate that after the loss of the endosymbiont, convergent evolution led to effector proteins like PEXEL , and RXLR precursors. PEXEL effectors might have enabled Pl. falciparum to achieve more complex interactions with its host and establish intercellular growth. In addition to the RXLR effector proteins, oomycetes acquired or evolved another class of effectors, the CRNs and a secreted invertase that allows use of sucrose from host plants . Oomycetes that are biotrophs or hemibiotrophs today lost their thiamine biosynthesis pathway and, in the case of A. laibachii, evolved a new “CHXC” effector class. After taking up the biotroph lifestyle, the genomes of Pl. falciparum, H. arabidopsidis, and A. laibachii started a gene reduction that is exemplified by looking at enzymes that require molybdenum cofactors and the molybdopterin biosynthesis pathway. Hemibiotroph P. infestans instead shows a strong genome expansion . In this context, H. arabidopsidis showed a genome expansion and acquired biotrophy late, based on the loss of only one molybdenum-dependent enzyme.