Admixture and recombination among Toxoplasma gondii lineages explain global genome diversity

Abstract

Toxoplasma gondii is a highly successful protozoan parasite that infects all warm-blooded animals and causes severe disease in immunocompromised and immune-naïve humans. It has an unusual global population structure: In North America and Europe, isolated strains fall predominantly into four largely clonal lineages, but in South America there is great genetic diversity and the North American clonal lineages are rarely found. Genetic variation between Toxoplasma strains determines differences in virulence, modulation of host-signaling pathways, growth, dissemination, and disease severity in mice and likely in humans. Most studies on Toxoplasma genetic variation have focused on either a few loci in many strains or low-resolution genome analysis of three clonal lineages. We use whole-genome sequencing to identify a large number of SNPs between 10 Toxoplasma strains from Europe and North and South America. These were used to identify haplotype blocks (genomic regions) shared between strains and construct a Toxoplasma haplotype map. Additional SNP analysis of RNA-sequencing data of 26 Toxoplasma strains, representing global diversity, allowed us to construct a comprehensive genealogy for Toxoplasma gondii that incorporates sexual recombination. These data show that most current isolates are recent recombinants and cannot be easily grouped into a limited number of haplogroups. A complex picture emerges in which some genomic regions have not been recently exchanged between any strains, and others recently spread from one strain to many others.

Fig. 1.

Schematic of the four Toxoplasma crosses investigated in this study. Shown are four crosses (A–D) between strains (ME49, MAS, BOF, and P89) from four ancestral HGs (HGII, HG4, HG6, and HG9) that might have led to the creation of strains (TgCatBr5, GT1, VEG, and GUY-KOE) from HG8 (cross A), HGI (cross B), HGIII (cross C), and HG5 (cross D), respectively (9, 13). The fraction of the genome of TgCatBr5, GT1, VEG, and GUY-KOE that was almost identical to the putative parental strains is indicated below those strains.

Fig. 2.

Genomewide strain comparisons of SNPs identify regions with recent common ancestry. Nonexon SNPs between two distinct strains were identified for 50-kb bins along T. gondii chromosome VIII and the relative SNP rate was plotted (5-kb step size). The chromosome position in Mb (x axis) and SNPs/50 kb (y axis) are shown. Numbered arrows indicate regions with distinct pairwise SNP rates.

Fig. 3.

Phylogenetic analysis of haplotype blocks. (A) Two phylogenetic trees (1, 2) showing the relatedness for a group of 18 hypothetical strains (A–R) for distinct genomic regions. For most of the genome (tree 1), strains F and J are highly divergent, but for a particular haploblock (tree 2), they have very low rates of SNP diversity. To distinguish the directionality of transfer for that haploblock, a phylogenetic tree showing relatedness among strains is constructed. If the exchanged locus was derived from a J-like ancestor, we would expect the red tree, where both strains have the position of strain J in tree 2. Alternatively, if the locus was from an F-like ancestor, we would expect the green tree, where both strains are found in the F position. For each of the 157 haplotype blocks identified (Dataset S2B), phylogenetic trees were constructed based on the pairwise distance between strains (calculated from SNP frequencies) (Fig. S2). The y axis indicates relative SNP rate. (B) Pairwise distance tree for haploblock chrIb_7 and haploblock chrXI_1. Because of the low similarity among strains [except among strains from types I (RH, GT1), II (DEG, PRU, ME49), III (CEP, VEG), and HG12 (RAY, WTD3) lineages], we conclude that haploblock chrIb_7 is derived from different strains that existed before the recent recombinations here investigated. For haploblock chrXI_1, the position of the type III strains (CEP and VEG) and of TgCatBr44 has shifted, and for this haploblock, the SNP rate between types II and III strains and between TgCatBr44 and is very low. Colors match haploblocks as indicated in Dataset S2B. ME49.1 represents the ME49 genome data and ME49 represents the RNA-seq data. (C–E) Pairwise distance trees were constructed as described for the indicated regions. (F) Phylogenetic tree of the apicoplast sequence constructed using neighbor joining. Numbers at branches are bootstrap values (1,000 replicates). Branch lengths indicate number of nucleotide differences.