Evidence of strain structure in Plasmodium falciparum var gene repertoires in children from Gabon, West Africa

Abstract

Existing theory on competition for hosts between pathogen strains has proposed that immune selection can lead to the maintenance of strain structure consisting of discrete, weakly overlapping antigenic repertoires. This prediction of strain theory has conceptual overlap with fundamental ideas in ecology on niche partitioning and limiting similarity between coexisting species in an ecosystem, which oppose the hypothesis of neutral coexistence. For Plasmodium falciparum, strain theory has been specifically proposed in relation to the major surface antigen of the blood stage, known as PfEMP1 and encoded by the multicopy multigene family known as the var genes. Deep sampling of the DBLα domain of var genes in the local population of Bakoumba, West Africa, was completed to define whether patterns of repertoire overlap support a role of immune selection under the opposing force of high outcrossing, a characteristic of areas of intense malaria transmission. Using a 454 high-throughput sequencing protocol, we report extremely high diversity of the DBLα domain and a large parasite population with DBLα repertoires structured into nonrandom patterns of overlap. Such population structure, significant for the high diversity of var genes that compose it at a local level, supports the existence of "strains" characterized by distinct var gene repertoires. Nonneutral, frequency-dependent competition would be at play and could underlie these patterns. With a computational experiment that simulates an intervention similar to mass drug administration, we argue that the observed repertoire structure matters for the antigenic var diversity of the parasite population remaining after intervention.

Keywords: Gabon; Plasmodium falciparum; parasite diversity; strain structure; var genes.

Fig. 1.

Frequency distribution and diversity of DBLα types in the population. (A) For the Bakoumba isolates, clustering of the 13,058 DBLα sequences in the dataset at 96% nucleotide identity resulted in 6,404 DBLα types. As depicted in the frequency distribution of DBLα types for Bakoumba, 61.2% DBLα types (3,917) were encountered in only 1 of the 200 isolates, whereas 0.2% DBLα types (13) were found in ≥30 of the 200 isolates with the maximum DBLα type frequency being 76. Identical or near-identical matches to the 10 most “common” DBLα types were found in the public sequence databases (Table S3). (B) The Bakoumba type accumulation curve for the sampling reveals a large number of DBLα types in the population. Despite in-depth sampling of the local population, the up-sloping curve demonstrates that more DBLα types are to be found with further sampling of the population.

Fig. S1.

Effect of sample size on Chao1 richness estimate.

Fig. S2.

Effect of nucleotide identity threshold for complete linkage clustering on rarefaction curve. High diversity of DBLα types sequenced is evident throughout a range of identity thresholds. In particular, the shape of the rarefaction curve and the total number of DBLα types identified does not change markedly within the nucleotide identity range of 0.90–0.98. Below 0.88, the slope of the curve drops more quickly, and then stabilizes again around 0.82.

Fig. 2.

Diversity profiles. Hill diversity numbers of the observed data show high levels of unevenness in the frequency of DBLα types (black). The x axis is the index of the Hill number (see text), and the y axis is the associated level of diversity for each of the measures. Examples of more common forms of diversity profiles are illustrated in gray. The calculation of the Hill numbers follows the procedure of Chao et al. (34) (Materials and Methods).

Fig. 3.

Pairwise comparisons of repertoire overlap. The pairwise type-sharing (PTS) scores (representing the proportion of DBLα types shared between two isolates) was calculated for all possible pairwise comparisons in the studies conducted in Bakoumba. A score of 0 represents no shared DBLα types between two isolates; a score of 1 represents complete identity of all DBLα types between two isolates. For Bakoumba, there were 19,900 possible pairwise comparisons of the 200 isolates in the study. (A) PTS scores for Bakoumba ranged from 0 to 0.984, with median and mean scores of 0.017 and 0.025, respectively. Overall, there was little to no overlap between most pairs of isolates, with 83.5% of scores being less than or equal to 0.05. (B) Bakoumba PTS scores depicted on a heat map reveal this low overlap among most pairs of isolates (light shading). A minority of isolates, however, appear more closely related to each other based on relatively higher levels of DBLα repertoire overlap (>0.50 PTS, darker shading). These clusters are apparent in the dendrogram in Fig. S4.

Fig. S4.

Dendrogram of Bakoumba isolates based on DBLα repertoire overlap.

Fig. S3.

Frequency distribution of sampled repertoire sizes. The repertoire sizes (number of DBLα types per isolate) for the isolates in this study ranged from 1 to 258, with a median repertoire size of 56.

Fig. 4.

Distribution of median and maximum PTS isolate similarity scores. Histograms showing the distribution of the median and maximum scores are found in a random ensemble of 1,000 samples. (A) The mean of the median scores in the sample is 0.023, with SD (std) = 0.0004, whereas the observed median in the data is 0.0169 (indicated by the red arrow). This is significantly lower than expected by chance (p_v = 0.002, and z score = −16.70). Similarly, (B) the mean of the maximum scores found in the sample is 0.199, with std = 0.088, whereas the observed maximum in the data is 0.983 (indicated by the red arrow). This is significantly higher than expected by chance (p_v = 0.017 and z score = 9.286). See Materials and Methods for information on the randomization process.

Fig. S5.

Pianka isolate similarity scores. (A) Distribution of Pianka scores for all isolate couples in the observed data. Similar to the distribution of PTS scores (as shown in Fig. 3A), the majority of the isolates have very low overlap, whereas a minority of couples appear more closely related to each other. Based on the same ensemble of 1,000 random samples used in Fig. 4, (B) the expected mean of the median scores is 0.026 with SD (std) = 0.0005, whereas the observed median (indicated by the red arrow) is 0.019; this is significantly lower than expected by chance (p_v = 0.002 and z score = –16.117); and (C) the expected mean of the maximum scores is 0.24 with std = 0.09, whereas the observed maximum (indicated by the red arrow) is 0.98; this is significantly higher than expected by chance (p_v = 0.018 and z score = 8.275). See Materials and Methods for information on the randomization process.

Fig. 5.

An illustration of the co-occurrence matrix of DBLα types in the data. The DBLα types in the matrix are sorted according to their observed frequency (i.e., the number of isolates they were observed in). Sum “A” represents the total number of DBLα couples, each of which was observed exactly once, but both in the same isolate. These couples were observed a significantly higher number of times than expected by chance using our random ensemble of 1,000 matrices (Materials and Methods). The z score equals 12.99 [the number of SDs the observed data departs from the mean; z score = (observed − mean)/std].

Fig. 6.

Distribution of Hill measures for different removal percentiles. Histogram distributions of marginal means of the random ensemble in gray for Hill measures: (A) H₀ and H₁, and (B) H₂ and H₃ (Materials and Methods), for removal of 40%, 60%, and 80% of the isolates. The corresponding values observed in the data are marked in red (indicated by the red arrows).