Multiple populations of artemisinin-resistant Plasmodium falciparum in Cambodia

Abstract

We describe an analysis of genome variation in 825 P. falciparum samples from Asia and Africa that identifies an unusual pattern of parasite population structure at the epicenter of artemisinin resistance in western Cambodia. Within this relatively small geographic area, we have discovered several distinct but apparently sympatric parasite subpopulations with extremely high levels of genetic differentiation. Of particular interest are three subpopulations, all associated with clinical resistance to artemisinin, which have skewed allele frequency spectra and high levels of haplotype homozygosity, indicative of founder effects and recent population expansion. We provide a catalog of SNPs that show high levels of differentiation in the artemisinin-resistant subpopulations, including codon variants in transporter proteins and DNA mismatch repair proteins. These data provide a population-level genetic framework for investigating the biological origins of artemisinin resistance and for defining molecular markers to assist in its elimination.

Figure 1. Population Structure in the samples set analysed

(a) Neighbour-joining (NJ) tree of all 825 samples, based on a pairwise distance matrix. The West African (WAF) samples are very clearly separated from South-East Asian (SEA) one. Within SEA, we observe that Thailand samples form a distinct cluster from Cambodian samples. Even within Cambodia, we see clear separation between a NE Cambodia/Vietnam group and a western Cambodian cluster. Sample in western Cambodia form a common clade, in spite of discernible population structure. (b) Three-dimensional plot of Principal components analysis (PCA) using all 825 samples, showing the first three components (PC1, PC2 and PC3). PC1 separates between continental clusters in Africa and Asia, while PC2 and PC3 largely account for the variability in western Cambodian samples (blue circles). (c, d) PCA plots using only SEA samples, showing (c) PC1 vs. PC2, and (d) PC2 vs. PC3. Once African samples are removed, all first three components are driven by variability in western Cambodia. Accordingly, we identified three “outlier” clusters in Cambodia (labelled KH2, KH3 and KH4), whose initial composition was determined by dividing the plot into quadrants. A “core” group (labelled KH1) is composed of north-eastern and western Cambodian samples, which cluster close to other SEA countries. Between the core and outlier group, several samples (labelled KHA) have intermediate PC values, suggesting admixture.

Figure 1. Population Structure in the samples set analysed

(a) Neighbour-joining (NJ) tree of all 825 samples, based on a pairwise distance matrix. The West African (WAF) samples are very clearly separated from South-East Asian (SEA) one. Within SEA, we observe that Thailand samples form a distinct cluster from Cambodian samples. Even within Cambodia, we see clear separation between a NE Cambodia/Vietnam group and a western Cambodian cluster. Sample in western Cambodia form a common clade, in spite of discernible population structure. (b) Three-dimensional plot of Principal components analysis (PCA) using all 825 samples, showing the first three components (PC1, PC2 and PC3). PC1 separates between continental clusters in Africa and Asia, while PC2 and PC3 largely account for the variability in western Cambodian samples (blue circles). (c, d) PCA plots using only SEA samples, showing (c) PC1 vs. PC2, and (d) PC2 vs. PC3. Once African samples are removed, all first three components are driven by variability in western Cambodia. Accordingly, we identified three “outlier” clusters in Cambodia (labelled KH2, KH3 and KH4), whose initial composition was determined by dividing the plot into quadrants. A “core” group (labelled KH1) is composed of north-eastern and western Cambodian samples, which cluster close to other SEA countries. Between the core and outlier group, several samples (labelled KHA) have intermediate PC values, suggesting admixture.

Figure 1. Population Structure in the samples set analysed

(a) Neighbour-joining (NJ) tree of all 825 samples, based on a pairwise distance matrix. The West African (WAF) samples are very clearly separated from South-East Asian (SEA) one. Within SEA, we observe that Thailand samples form a distinct cluster from Cambodian samples. Even within Cambodia, we see clear separation between a NE Cambodia/Vietnam group and a western Cambodian cluster. Sample in western Cambodia form a common clade, in spite of discernible population structure. (b) Three-dimensional plot of Principal components analysis (PCA) using all 825 samples, showing the first three components (PC1, PC2 and PC3). PC1 separates between continental clusters in Africa and Asia, while PC2 and PC3 largely account for the variability in western Cambodian samples (blue circles). (c, d) PCA plots using only SEA samples, showing (c) PC1 vs. PC2, and (d) PC2 vs. PC3. Once African samples are removed, all first three components are driven by variability in western Cambodia. Accordingly, we identified three “outlier” clusters in Cambodia (labelled KH2, KH3 and KH4), whose initial composition was determined by dividing the plot into quadrants. A “core” group (labelled KH1) is composed of north-eastern and western Cambodian samples, which cluster close to other SEA countries. Between the core and outlier group, several samples (labelled KHA) have intermediate PC values, suggesting admixture.

Figure 1. Population Structure in the samples set analysed

(a) Neighbour-joining (NJ) tree of all 825 samples, based on a pairwise distance matrix. The West African (WAF) samples are very clearly separated from South-East Asian (SEA) one. Within SEA, we observe that Thailand samples form a distinct cluster from Cambodian samples. Even within Cambodia, we see clear separation between a NE Cambodia/Vietnam group and a western Cambodian cluster. Sample in western Cambodia form a common clade, in spite of discernible population structure. (b) Three-dimensional plot of Principal components analysis (PCA) using all 825 samples, showing the first three components (PC1, PC2 and PC3). PC1 separates between continental clusters in Africa and Asia, while PC2 and PC3 largely account for the variability in western Cambodian samples (blue circles). (c, d) PCA plots using only SEA samples, showing (c) PC1 vs. PC2, and (d) PC2 vs. PC3. Once African samples are removed, all first three components are driven by variability in western Cambodia. Accordingly, we identified three “outlier” clusters in Cambodia (labelled KH2, KH3 and KH4), whose initial composition was determined by dividing the plot into quadrants. A “core” group (labelled KH1) is composed of north-eastern and western Cambodian samples, which cluster close to other SEA countries. Between the core and outlier group, several samples (labelled KHA) have intermediate PC values, suggesting admixture.

Figure 1. Population Structure in the samples set analysed

(a) Neighbour-joining (NJ) tree of all 825 samples, based on a pairwise distance matrix. The West African (WAF) samples are very clearly separated from South-East Asian (SEA) one. Within SEA, we observe that Thailand samples form a distinct cluster from Cambodian samples. Even within Cambodia, we see clear separation between a NE Cambodia/Vietnam group and a western Cambodian cluster. Sample in western Cambodia form a common clade, in spite of discernible population structure. (b) Three-dimensional plot of Principal components analysis (PCA) using all 825 samples, showing the first three components (PC1, PC2 and PC3). PC1 separates between continental clusters in Africa and Asia, while PC2 and PC3 largely account for the variability in western Cambodian samples (blue circles). (c, d) PCA plots using only SEA samples, showing (c) PC1 vs. PC2, and (d) PC2 vs. PC3. Once African samples are removed, all first three components are driven by variability in western Cambodia. Accordingly, we identified three “outlier” clusters in Cambodia (labelled KH2, KH3 and KH4), whose initial composition was determined by dividing the plot into quadrants. A “core” group (labelled KH1) is composed of north-eastern and western Cambodian samples, which cluster close to other SEA countries. Between the core and outlier group, several samples (labelled KHA) have intermediate PC values, suggesting admixture.

Figure 1. Population Structure in the samples set analysed

(a) Neighbour-joining (NJ) tree of all 825 samples, based on a pairwise distance matrix. The West African (WAF) samples are very clearly separated from South-East Asian (SEA) one. Within SEA, we observe that Thailand samples form a distinct cluster from Cambodian samples. Even within Cambodia, we see clear separation between a NE Cambodia/Vietnam group and a western Cambodian cluster. Sample in western Cambodia form a common clade, in spite of discernible population structure. (b) Three-dimensional plot of Principal components analysis (PCA) using all 825 samples, showing the first three components (PC1, PC2 and PC3). PC1 separates between continental clusters in Africa and Asia, while PC2 and PC3 largely account for the variability in western Cambodian samples (blue circles). (c, d) PCA plots using only SEA samples, showing (c) PC1 vs. PC2, and (d) PC2 vs. PC3. Once African samples are removed, all first three components are driven by variability in western Cambodia. Accordingly, we identified three “outlier” clusters in Cambodia (labelled KH2, KH3 and KH4), whose initial composition was determined by dividing the plot into quadrants. A “core” group (labelled KH1) is composed of north-eastern and western Cambodian samples, which cluster close to other SEA countries. Between the core and outlier group, several samples (labelled KHA) have intermediate PC values, suggesting admixture.

Figure 2. Chromosome Painting

A high-level view of chromosome painting across all 14 chromosomes is shown here. Each line represents a sample; samples are organized by subpopulation, in the order KH1, KHA, KH2, KH3 and KH4, starting from the bottom. The colours of segments indicate the cluster found to be the most probable donor of each segment (purple=KH1, blue =KH2, red=KH3, green =KH4). Callouts (a-c) show details of some important features of the painting plots. (a) KH1 samples contain many short segments of DNA sequence labelled as KH2/3/4, consistent with KH1 representing the ancestral population from which KH2/3/4 are derived. (b) KH2/3/4 samples comprise very long segments labelled as belonging to other populations, suggesting recent admixture. (c) KHA samples contain both short and long chunks of DNA labelled as KH2/3/4, consistent with the hypothesis that it is an admixed population, with contribution from these ART-R subpopulations.

Figure 2. Chromosome Painting

A high-level view of chromosome painting across all 14 chromosomes is shown here. Each line represents a sample; samples are organized by subpopulation, in the order KH1, KHA, KH2, KH3 and KH4, starting from the bottom. The colours of segments indicate the cluster found to be the most probable donor of each segment (purple=KH1, blue =KH2, red=KH3, green =KH4). Callouts (a-c) show details of some important features of the painting plots. (a) KH1 samples contain many short segments of DNA sequence labelled as KH2/3/4, consistent with KH1 representing the ancestral population from which KH2/3/4 are derived. (b) KH2/3/4 samples comprise very long segments labelled as belonging to other populations, suggesting recent admixture. (c) KHA samples contain both short and long chunks of DNA labelled as KH2/3/4, consistent with the hypothesis that it is an admixed population, with contribution from these ART-R subpopulations.

Figure 3. Genetic and Phenotypic differentiation between Cambodian subpopulations

(a) Ancestry analysis of the 293 Cambodian samples, based on chromosome painting analysis of all Cambodian samples (see Methods). In this figure, each vertical bar represents a single samples, and the bar is coloured according to the proportion of the genome that was determined to have originated in each of the four KH1-KH4 clusters, as defined by PCA analysis (KH1=purple, KH2=blue, KH3=red, KH4=green). Samples were ordered according to these proportions, and any sample with ≥80% of the genome originating from a single cluster was assigned the label of that cluster. Any sample that was not labelled as belonging to one of the resulting four subpopulations was assigned to the KHA group. The final partitioning of samples into subpopulations is shown in the figure. Most of the samples classified as KHA appear to be composed of a mixture of KH1 and one or more of the remaining clusters, consistent with the hypothesis of a recombinant population. (b) Boxplot showing the distribution of artemisinin response phenotypes (parasite clearance half-life) in the Cambodian parasite subpopulations thus identified (n=212). Each box represents the interquartile range of values, split at the median; whiskers extend to the furthest point that is within 1.5 times the length of the box. (b) An analogous plot derived from an analysis of Western Cambodian parasites only (i.e. excluding the Ratanakiri site; n=166). A table of parasite clearance half-life with p-values is given in Supplementary Table 1.

Figure 3. Genetic and Phenotypic differentiation between Cambodian subpopulations

(a) Ancestry analysis of the 293 Cambodian samples, based on chromosome painting analysis of all Cambodian samples (see Methods). In this figure, each vertical bar represents a single samples, and the bar is coloured according to the proportion of the genome that was determined to have originated in each of the four KH1-KH4 clusters, as defined by PCA analysis (KH1=purple, KH2=blue, KH3=red, KH4=green). Samples were ordered according to these proportions, and any sample with ≥80% of the genome originating from a single cluster was assigned the label of that cluster. Any sample that was not labelled as belonging to one of the resulting four subpopulations was assigned to the KHA group. The final partitioning of samples into subpopulations is shown in the figure. Most of the samples classified as KHA appear to be composed of a mixture of KH1 and one or more of the remaining clusters, consistent with the hypothesis of a recombinant population. (b) Boxplot showing the distribution of artemisinin response phenotypes (parasite clearance half-life) in the Cambodian parasite subpopulations thus identified (n=212). Each box represents the interquartile range of values, split at the median; whiskers extend to the furthest point that is within 1.5 times the length of the box. (b) An analogous plot derived from an analysis of Western Cambodian parasites only (i.e. excluding the Ratanakiri site; n=166). A table of parasite clearance half-life with p-values is given in Supplementary Table 1.

Figure 3. Genetic and Phenotypic differentiation between Cambodian subpopulations

(a) Ancestry analysis of the 293 Cambodian samples, based on chromosome painting analysis of all Cambodian samples (see Methods). In this figure, each vertical bar represents a single samples, and the bar is coloured according to the proportion of the genome that was determined to have originated in each of the four KH1-KH4 clusters, as defined by PCA analysis (KH1=purple, KH2=blue, KH3=red, KH4=green). Samples were ordered according to these proportions, and any sample with ≥80% of the genome originating from a single cluster was assigned the label of that cluster. Any sample that was not labelled as belonging to one of the resulting four subpopulations was assigned to the KHA group. The final partitioning of samples into subpopulations is shown in the figure. Most of the samples classified as KHA appear to be composed of a mixture of KH1 and one or more of the remaining clusters, consistent with the hypothesis of a recombinant population. (b) Boxplot showing the distribution of artemisinin response phenotypes (parasite clearance half-life) in the Cambodian parasite subpopulations thus identified (n=212). Each box represents the interquartile range of values, split at the median; whiskers extend to the furthest point that is within 1.5 times the length of the box. (b) An analogous plot derived from an analysis of Western Cambodian parasites only (i.e. excluding the Ratanakiri site; n=166). A table of parasite clearance half-life with p-values is given in Supplementary Table 1.

Figure 4. Evidence for multiple founder effects in Cambodian subpopulations

(a) Minor allele frequencies (MAF) spectra for 7 populations. For each population, SNPs with MAF > 0 were assigned to one of 5 equally-sized MAF bins, and we show counts in these bins. A sub-sampling procedure, averaged over 1000 iterations, was used to minimize confounding by differences in sample size. MAF frequency spectra of KH1 and KHA are consistent with the rest of SEA, with a high proportion of low-MAF SNPs. In contrast, the ART-R subpopulations (KH2, KH3 and KH4) have considerably fewer SNPs at low frequency, and more homogeneous frequency spectra. (b) Genome-wide decay of Linkage Disequilibrium (LD), expressed using a classic LD measure (r²). SNPs with MAF >0.15 were organized into equal-size allele frequency bins. Each bin was analyzed independently and corrected for effects of sample size and population structure, before bin results are aggregated (see Supplementary Note). African populations exhibit rapid LD decay, with no detectable LD within 1 kbp. In Cambodian ART-R subpopulations, LD decays over much greater distances than in other groups (including KH1 and KHA), suggesting highly reduced haplotype diversity. (c) Genome-wide distribution of haplotype homozygosity in different populations. Each box represents the interquartile range of values, split at the median; whiskers extend to the furthest point that is within 1.5 times the length of the box. A sliding window of width 201 SNP was used to construct local haplotypes for each sample. While “core” populations of South-East Asia (SEA and KH1) have a level of haplotype diversity consistent with that seen in Africa (WAF), the ART-R subpopulations (KH2, KH3 and KH4) exhibit much higher homozygosity, indicating reduced overall haplotype diversity.

Figure 4. Evidence for multiple founder effects in Cambodian subpopulations

(a) Minor allele frequencies (MAF) spectra for 7 populations. For each population, SNPs with MAF > 0 were assigned to one of 5 equally-sized MAF bins, and we show counts in these bins. A sub-sampling procedure, averaged over 1000 iterations, was used to minimize confounding by differences in sample size. MAF frequency spectra of KH1 and KHA are consistent with the rest of SEA, with a high proportion of low-MAF SNPs. In contrast, the ART-R subpopulations (KH2, KH3 and KH4) have considerably fewer SNPs at low frequency, and more homogeneous frequency spectra. (b) Genome-wide decay of Linkage Disequilibrium (LD), expressed using a classic LD measure (r²). SNPs with MAF >0.15 were organized into equal-size allele frequency bins. Each bin was analyzed independently and corrected for effects of sample size and population structure, before bin results are aggregated (see Supplementary Note). African populations exhibit rapid LD decay, with no detectable LD within 1 kbp. In Cambodian ART-R subpopulations, LD decays over much greater distances than in other groups (including KH1 and KHA), suggesting highly reduced haplotype diversity. (c) Genome-wide distribution of haplotype homozygosity in different populations. Each box represents the interquartile range of values, split at the median; whiskers extend to the furthest point that is within 1.5 times the length of the box. A sliding window of width 201 SNP was used to construct local haplotypes for each sample. While “core” populations of South-East Asia (SEA and KH1) have a level of haplotype diversity consistent with that seen in Africa (WAF), the ART-R subpopulations (KH2, KH3 and KH4) exhibit much higher homozygosity, indicating reduced overall haplotype diversity.

Figure 4. Evidence for multiple founder effects in Cambodian subpopulations

(a) Minor allele frequencies (MAF) spectra for 7 populations. For each population, SNPs with MAF > 0 were assigned to one of 5 equally-sized MAF bins, and we show counts in these bins. A sub-sampling procedure, averaged over 1000 iterations, was used to minimize confounding by differences in sample size. MAF frequency spectra of KH1 and KHA are consistent with the rest of SEA, with a high proportion of low-MAF SNPs. In contrast, the ART-R subpopulations (KH2, KH3 and KH4) have considerably fewer SNPs at low frequency, and more homogeneous frequency spectra. (b) Genome-wide decay of Linkage Disequilibrium (LD), expressed using a classic LD measure (r²). SNPs with MAF >0.15 were organized into equal-size allele frequency bins. Each bin was analyzed independently and corrected for effects of sample size and population structure, before bin results are aggregated (see Supplementary Note). African populations exhibit rapid LD decay, with no detectable LD within 1 kbp. In Cambodian ART-R subpopulations, LD decays over much greater distances than in other groups (including KH1 and KHA), suggesting highly reduced haplotype diversity. (c) Genome-wide distribution of haplotype homozygosity in different populations. Each box represents the interquartile range of values, split at the median; whiskers extend to the furthest point that is within 1.5 times the length of the box. A sliding window of width 201 SNP was used to construct local haplotypes for each sample. While “core” populations of South-East Asia (SEA and KH1) have a level of haplotype diversity consistent with that seen in Africa (WAF), the ART-R subpopulations (KH2, KH3 and KH4) exhibit much higher homozygosity, indicating reduced overall haplotype diversity.