Local adaptation and vector-mediated population structure in Plasmodium vivax malaria.

Plasmodium vivax in southern Mexico exhibits different infectivities to 2 local mosquito vectors, Anopheles pseudopunctipennis and Anopheles albimanus. Previous work has tied these differences in mosquito infectivity to variation in the central repeat motif of the malaria parasite's circumsporozoite (csp) gene, but subsequent studies have questioned this view. Here we present evidence that P. vivax in southern Mexico comprised 3 genetic populations whose distributions largely mirror those of the 2 mosquito vectors. Additionally, laboratory colony feeding experiments indicate that parasite populations are most compatible with sympatric mosquito species. Our results suggest that reciprocal selection between malaria parasites and mosquito vectors has led to local adaptation of the parasite. Adaptation to local vectors may play an important role in generating population structure in Plasmodium. A better understanding of coevolutionary dynamics between sympatric mosquitoes and parasites will facilitate the identification of molecular mechanisms relevant to disease transmission in nature and provide crucial information for malaria control.

Introduction

Worldwide, 30–40 Anopheles mosquito species are capable of transmitting malaria (Kiszewski et al. 2004), but the susceptibility of a mosquito to malaria infection and its ability to successfully transmit the disease is highly variable among species (Billingsley and Sinden 1997; Alavi et al. 2003) and even among genotypes within a species (Lambrechts et al. 2005). Coadaptation of vectors and parasites has been proposed to explain these differences (Cohuet et al. 2006). Local adaptation, whereby parasites exhibit greater infectivity on local as opposed to foreign hosts, is predicted when the parasite has an evolutionary advantage over the coevolving host, such as a shorter generation time, larger population size, or higher recombination rate (Kaltz and Shykoff 1998; Thompson 1999; Kawecki and Ebert 2004; Brandt et al. 2007). In order to be transmitted to humans, the malaria parasite must complete a series of developmental stages within the mosquito vector, where it encounters a number of population bottlenecks (reviewed in Abraham and Jacobs-Lorena 2004). Parasite adaptations that optimize transmission success on local hosts might contribute to the mosaic of malaria parasite–vector compatibilities observed in nature.

In southern Mexico, Plasmodium vivax parasites exhibit different compatibilities with respect to 2 local mosquito vectors, Anopheles pseudopunctipennis and Anopheles albimanus. Previous studies (Gonzalez-Ceron et al. 1999; Rodriguez et al. 2000) have suggested an association between mosquito compatibility and the central repeat motif of the malaria circumsporozoite (csp) gene that encodes a major surface protein of malaria sporozoites. The csp variants differ from each other in the amino acid composition of the central repetitive region (Rosenberg et al. 1989; Kain et al. 1992). Two variants, termed VK210 and VK247, are circulating in southern Mexico where their distributions roughly mirror those of the vectors: VK210 is found with An. albimanus on the coast and VK247 co-occurs with An. pseudopunctipennis in the foothills (Rodriguez et al. 2000). The vectors have largely nonoverlapping geographic distributions (Manguin et al. 1995; Rubio-Palis and Zimmerman 1997), although both mosquitoes and csp variants are found in the town of Tapachula, situated between the 2 ecoregions (elevation: 120–220 m).

Despite the apparent association between csp variants and mosquito compatibility (Gonzalez-Ceron et al. 1999; Rodriguez et al. 2000), the circum-sporozoite protein (CSP) is unlikely to control the observed compatibility phenotypes. Subsequent studies indicate that destruction of the incompatible parasite occurs prior to the expression of the CSP: VK210 parasites were unable to successfully traverse the An. pseudopunctipennis peritrophic membrane (Gonzalez-Ceron et al. 2007), and VK247 parasites in An. albimanus were destroyed during midgut invasion and early oocyst development (Gonzalez-Ceron et al. 2001). Further, in Columbia the reverse pattern was observed: An. albimanus fed on blood infected with VK247 parasites had a higher number of salivary gland sporozoites than those fed on VK210-infected blood (Gonzalez et al. 2001).

It is well known that underlying genetic population structure can result in spurious associations between candidate loci and a phenotype (Falush et al. 2003). If P. vivax in southern Mexico is divided into genetically isolated populations, then the association between csp variants and mosquito compatibility may simply reflect this. In this study, we surveyed P. vivax using a panel of recently developed microsatellites (MSs) to assess the presence of population structure. We found that P. vivax in our sample area comprised 3 genetic populations and that each population optimally infected the sympatric mosquito species. This study highlights the potentially misleading effect of population structure on association analysis and the importance of employing a coevolutionary framework to study epidemiologically relevant traits such as infectivity.

Materials and Methods

DNA Samples, Genome-wide Amplification, and Microsatellite Typing

Plasmodium vivax DNA samples were extracted from preexisting anonymous blood spots on Whatman No. 2 filter paper collected from Chiapas State, Mexico, between 1997 and 2005 (Gonzalez-Ceron et al. 1999). A total of 230 DNA samples were analyzed for the current study. MS loci were selected from approximately 160 genome-wide MSs identified by the Carlton laboratory. We used 16 Mexican isolates to subscreen 50 of these loci. Twenty-six loci were rejected either because they were monomorphic or because results were difficult to interpret due to the presence of stutter peaks and/or imperfect repeats. Primers and a description of the remaining 24 MS loci are given in supplementary table 1 (Supplementary Material online).

The small amount of parasite DNA contained in a single blood spot was not sufficient for the multimarker typing performed here. Therefore, we used whole-genome amplification with the Qiagen Repli-g Kit (Qiagen, Valencia CA) to increase genomic DNA yield. Whole-genome–amplified DNA was washed with ethanol to remove extraneous peaks observed in preliminary analyses. Typically, 2 μl of stock DNA extracted from filter paper resulted in a 50-μl sample with a DNA concentration of 200–300 ng/μl. Loci were amplified in 10-μl reactions (1.5 mM MgCl2, 2 μM of each primer, and 0.5 U Taq polymerase) following standard protocols (Schlötterer 1998). Polymerase chain reaction (PCR) parameters were as follows: 2 min initial denaturation at 94 °C, then 42 cycles of 20 s at 94 °C, 10 s at 52 °C, 10 s at 47 °C, and 30 s at 60 °C, followed by 5 min at 60 °C. Length variation of PCR products was measured on a 16 capillary ABI3100 DNA Genetic Analyzer (Applied Biosystems, Foster City, CA). Fragment size polymorphism was analyzed using GENEMAPPER software, also from Applied Biosystems.

Estimate of Genetic Diversity

The predominant allele present in each PCR product was used to estimate allele frequencies at each locus. MS variability was measured by expected heterozygosity (He) (Nei 1978) and variance in repeat number (Vr) (Moran 1975), using the program MICROSATELLITE ANALYZER v4.05 (Dieringer and Schlotterer 2003). He and Vr were corrected for small sample sizes by multiplying n by n/(n − 1), where n is the number of chromosomes analyzed. To assess the level of multiple infections in our sample, we counted the proportion of isolates with a single allele at all 24 loci. We also measured the mean number of alleles per locus per infection.

Determination of Population Structure

We used the Bayesian model–based clustering method implemented in STRUCTURE v 2.1 (Pritchard et al. 2000) to test for the presence of population structure in our sample. The number of genetic partitions (K) was inferred by calculating the probability of the data given a certain prior value of K, [P(D|K)], over a large number of Markov chain Monte Carlo (MCMC) iterations. Genetic partitions were characterized by different allele frequencies, and individuals were probabilistically assigned to 1 or more partition according to their allele distribution. MCMC searches consisted of 500,000 “burn-in” steps followed by 1,000,000 iterations. To check for convergence of P(D|K), we performed 10 replicate runs at each K value from 1–10 under the admixture model with correlated allele frequencies. To aid in determining the most appropriate number of partitions, we used an ad hoc quantity that evaluates the second-order rate of change of the likelihood function with respect to K (Evanno et al. 2005).

We estimated population pairwise FST values for the partitions defined in the STRUCTURE analysis using θ (Weir and Cockerham 1984) in Arlequin 3.11 (Excoffier et al. 2005). The significance of pairwise FST estimates was calculated using permutation tests (N = 10,000) and Bonferroni corrections for multiple comparisons.

Principle Coordinate Analysis (PCA)

Relationships among genetic partitions were visualized using principal coordinate analysis (PCA; Gower 1966) performed in Genalex 6 (Peakall and Smouse 2006). PCA generates a set of orthogonal axes for which each successive dimension maximizes the remaining variance in the data. Input for PCA analysis consisted of individual pairwise genetic distance matrices of 1) the proportion of shared alleles (Dps,) and 2) Nei's standard genetic distance corrected for sample size (D), calculated using MsAnalyzer (Dieringer and Schlotterer 2003).

Multivariate Spatial Genetic Autocorrelation

We used multivariate spatial genetic autocorrelation (Smouse and Peakall) to test the hypothesis that P. vivax genetic structure has a spatial component that reflects the distributions of the 2 vectors, An. pseudopunctipennis and An. albimanus. We measured correlations between P. vivax pairwise geographic and pairwise squared genetic distance matrices using the program Genalex (Peakall and Smouse 2006). The autocorrelation coefficient r provides a measure of pairwise genetic similarity among isolates within a particular distance size class. We used progressively larger size classes starting with 5 km and increasing up to the maximum distance within the study area (70 km). Ten thousand random permutations of genotypes among geographic locations were used to generate a null distribution under the hypothesis of no genetic structure (Smouse and Peakall 1999; Peakell et al. 2003). In addition, 10,000 bootstrap replicates of pairwise comparisons were used to generate 95% confidence intervals (CIs) for the autocorrelation coefficient r (Peakell et al. 2003). We performed the analysis twice, both with and without samples from Tapachula.

Mosquito Infectivity Experiments

Data from 183 laboratory colony mosquito feeding experiments conducted previously at the Centro de Investigacion de Paludismo in Tapachula between 1997 and 2005 (Gonzalez-Ceron and Rodriguez 1991; Ramsey et al. 1994; Gonzalez-Ceron et al. 1999, 2007) were reanalyzed for associations between parasite population and mosquito infectivity. Only feeding experiments with no missing data were included. Both the prevalence (proportion of infected mosquitoes) and the infection intensity (mean number of oocysts) were measured.

For the analysis of infection prevalence, we used log-linear models on the number of infected mosquitoes with an offset for total number of mosquitoes. We expressed results using the relative rates (RRs) of infecting An. pseudopunctipennis compared with An. albimanus and the ratio of those relative rates (RRR). For example, to examine the effects of csp variant, we used the RRR for VK247 versus VK210:

where P(AP,VK247) is the probability of infecting An. pseudopunctipennis when the csp repeat variant was VK247, and the other notation is analogous. Values of RRR >1 in this example indicate that the RR of infecting An. pseudopunctipennis compared with An. albimanus is larger when the csp repeat variant is VK247. We explored several of these models to describe this RRR with and without additional effects for parasite genetic population. Further, we used the log-linear models to estimate several different RRRs to examine the RR of infecting An. pseudopunctipennis compared with An. albimanus for different pairs of genetic populations. To control for correlations within parasite isolate, we used the generalized estimating equation (GEE) with the working independence correlation structure. We performed Wald-type tests with t distributions and K–p degrees of freedom (dfs), where K is the number of unique parasite isolates and p is the number of parameters (Lipsitz et al. 1994).

For the analysis of mean oocyst number, we also used a GEE with a working identity matrix in order to account for the correlation within parasite isolate. As with infection prevalence, we used the t distribution on the Wald test with K–p dfs. We compared means 1) among genetic population within each mosquito and 2) between mosquitoes for each genetic population. To test for significant interaction between mosquito species and parasite genotype, we again used the GEE method with an F distribution and K–p dfs. GEE analyses were carried out in R 2.5.1 using the gee package (Carey 2007; R Development Core Team 2007).

Results

Parasite Samples and Genetic Diversity

Parasite samples from a total of 98 localities within a 100-km radius in Chiapas, Mexico, were included in this study (supplementary table 2, Supplementary Material online). MSs loci were located on 9 of the 14 chromosomes (supplementary table 2, Supplementary Material online). Multiple alleles were scored at a given locus if minor peaks were more than one-third the height of the predominant peak (Imwong et al. 2007). Using this criterion, the majority of our samples appear to be single infections. The mean number of alleles per locus per infection was 1.01 (standard deviation ± 0.01). Out of 234 typed samples, 197 (84%) carried a single allele at each locus. Among the 37 putatively mixed infections, 33 (89%) showed multiple alleles at only 1–2 loci out of 24, suggesting that they represent mutations within infection rather than multiple infections (Imwong et al. 2007). The 4 isolates with >2 loci showing multiple alleles were assumed to be mixed infections and were removed from further analyses. He for the 24 MS loci ranged from 0.12 to 0.87. Three markers with low He (0.12–0.21) were retained on the assumption that even markers with low mutation rates would contribute information for the estimation of demographic processes.

Highly Differentiated Plasmodium vivax Populations

Our data showed that P. vivax in southern Mexico exhibits population structure and restricted gene flow in an area of less than 100 km2. STRUCTURE analysis indicated that P. vivax comprised 3 genetic subpopulations. The estimated logarithm of the probability of the data [P(D|K)] plateaued at K = 5–6 (fig. 1). However, the K estimator derived from the rate of change of K (ΔK) (Evanno et al. 2005) indicated that the best estimate of K was 3. One population, hereafter termed c1, includes the majority of coastal samples (88%), as well as outgroup isolates 045V from Brazil and Sal1 from El Salvador (fig. 1). The other 2 populations, hereafter termed f1 and f2, predominate in the foothills and included a number of individuals with mixed ancestry (fig. 1).

F

Population structure of Plasmodium vivax in southern Mexico. (A) Plot of the log probability of the data [ln P(D)] and the second-order rate of change constant (ΔK) for K = 1–10. A peak for ΔK at K = 3 suggests that 3 populations best fit the data. (B) Proportional membership of isolates in the 3 clusters identified by the STRUCTURE program. Each isolate is represented by a thin vertical line. Colors identify genetic clusters: yellow, c1; blue, f1; red, f2. Isolates are organized by geographical origin.

Significant genetic differentiation was detected among all 3 populations (FST = 0.279–0.364). He, median Vr, number of alleles, and number of private alleles were all highest in c1 (table 1), even though this sample size was smaller than the other 2 (58 vs. 83 and 89). Populations f1 and f2 showed similar levels of genetic diversity for all measures. These results suggest that the coastal parasite population may be either older or larger than the 2 foothill populations. The distributions of both mosquitoes extend beyond our study area (Manguin et al. 1995; Molina-Cruz et al. 2004), and it is likely that the parasite populations are more extensive than what we have sampled here.

Table 1

Variability Measures for 24 MS Loci in Plasmodium vivax from Southern Mexico

Population	n	He	Vr	Mean No. of Alleles	Mean No. of Private Alleles
c1	59.2	0.71	15.7	7.9	3.6
f1	82.7	0.32	6.7	4.3	0.3
f2	88.4	0.43	9.8	4.4	0.4
Total	230.4	0.61	14.4	9.1	—

NOTE.—n, mean effective number of chromosomes; He, mean heterozygosity among 24 MS markers; Vr, median variance in repeat number.

We observed temporal fluctuations in the proportional representation of the P. vivax populations in our sample (fig. 2). The proportion of c1 isolates was significantly lower during the period of 2002–2005 as compared with 1997–2001 (z = −2.62, P = 0.008). The opposite pattern was seen for f2, which comprised a significantly higher proportion of the sampled isolates during 2002–2005 as compared with 1997–2001 (z = −3.14, P = 0.0009). The proportion of f1 isolates remained relatively unchanged throughout the sampling period.

F

Temporal variation in the incidence of Plasmodium vivax populations c1, f1, and f2 in parasite samples from Chiapas, Mexico, during 1997–2005.

Restricted Gene Flow Between Plasmodium vivax Populations

To better characterize the degree of genetic relatedness among P. vivax populations, we performed principal coordinates analysis (PCA). Results based on Dps and D did not differ greatly from each other; therefore, we present results for Dps only (fig. 3). The 3 genetic populations showed little overlap. Axis 1 primarily separated c1 from f1 and f2, whereas axis 2 further separated f2 from f1. We observed a cline rather than discreet clusters for f1 and f2 isolates. The dispersion of samples along a line suggests mixed ancestry (Patterson et al. 2006), and this is confirmed by the central location of admixed f1/f2 individuals along the cline (fig. 3).

F

PCAs of the Plasmodium vivax populations. Plot of the first 2 PCA axes based on the proportion of shared alleles (Dps). The percentage of variation explained by each axis is indicated. Individuals with mixed f1/f2 ancestry (membership >0.20 in each) are shown as gray.

We detected significant positive spatial structure extending up to 40 km (fig. 4). Results did not differ greatly for analyses including Tapachula isolates; therefore, we present only the results from the non-Tapachula samples here. In the presence of positive structure, r, the autocorrelation coefficient, should decrease with increasing distance class size (Peakell et al. 2003) and the distance class for which r is no longer significant can be taken as the upper limit of positive structure (Double et al. 2005). There was no evidence of significant structure within either the coastal or the foothill ecoregions when analyzed separately (data not shown) despite the presence of 2 parasite populations in the foothills. The mean distance between coastal and foothill sites is 34.3 km (SD ± 10.2), and the observed positive autocorrelation extending to 40 km is likely to reflect restricted gene flow between parasites collected from the coastal and foothill ecoregions.

F

Spatial genetic structure of Plasmodium vivax in southern Mexico. (A) Study area in southern Chiapas State, Mexico. Circle size is proportional to frequency in the sample, with the exception of Tapachula (n = 135). Individuals with mixed ancestry are shown as a single color corresponding to the population with the highest membership proportion. Yellow, c1; blue, f1; red, f2. (B) Spatial genetic autocorrelation (r) across multiple distance classes. Samples from the town of Tapachula were not included. Ninety-five percent CIs for estimates of r based on 1,000 bootstrap replicates are shown as error bars, and 95% CIs based on 1,000 random permutations for the null hypothesis of no spatial autocorrelation are shown by the bold dashes surrounding r = 0.

Enhanced Infectivity of Sympatric Mosquito Species

There was a clear phenotypic separation between c1 and f1/f2 populations with respect to mean oocyst number in the 2 mosquito species, with each parasite population being more compatible with the sympatric mosquito (fig. 5). The main effects of mosquito species and parasite population were both significant (P = 8.32 × 10−04 and P = 0.0016, respectively). The mean number of oocysts was significantly different for c1 compared with f1/f2 within each vector. No differences between f1 and f2 were observed in either vector. Most importantly in terms of coevolutionary dynamics, the interaction between mosquito species and parasite genetic population had a highly significant effect (P = 1.191 × 10−08; fig. 5).

F

Comparison of the infection intensities of 3 parasite populations in 2 mosquito species. Anopheles pseudopunctipennis (AP) and Anopheles albimanus (AA) were challenged with blood infected with Plasmodium vivax isolates belonging to 1 of 3 populations: c1, f1, and f2. Error bars represent 95% CIs. Data are from 183 laboratory colony feeding experiments conducted between 1997 and 2005. A strong interaction between parasite genetic background and mosquito species can be seen.

We also examined the RR of infecting An. pseudopunctipennis compared with An. albimanus for the CSP VK247 variant over the analogous RR for the VK210 variant, termed the ratio of relative rate (RRR) of VK247 versus VK210. We first modeled this RRR without accounting for parasite population, giving an RRR of 8.82 (95% CI = 5.33–14.62), indicating elevated infectivity in An. pseudopunctipennis for VK247 parasites. A separate analysis for each parasite population (c1, f1, f2) shows that this effect is confined to c1 (RRR = 106.74); however, the CI is unreliable as there is only one c1/VK247 isolate in our sample (table 2). Removing this isolate does not greatly alter the overall RR (8.75, 95% CI = 5.31–14.43). Further examination of the RR of infecting An. pseudopunctipennis compared with An. albimanus for each parasite population/CSP variant combination revealed that the overall RRR is being driven by the low RR for c1/VK210 (table 2). Within populations f1 and f2, CSP variant had no effect on infection prevalence.

Table 2

The RR of Infecting Anopheles pseudopunctipennis Compared with Anopheles Albimanus (e.g., an RR of 5 Means that the Rate of Infection is 5 Times Larger for An. pseudopunctipennis than for An. albimanus)

Genetic Population/csp Variant	Isolates (n)	Mosquitoes (n)		95% CI
			RR	Lower	Upper
c1/VK247	1	35	12.500	NaN	NaN
c1/VK210	36	1081	0.117	0.056	0.244
f1/VK247	66	2006	5.025	4.050	6.235
f1/VK210	7	200	4.500	1.295	15.636
f2/VK247	54	1546	6.269	4.627	8.494
f2/VK210	19	616	4.495	2.472	8.174
c1/ALL	37	1116	0.170	0.075	0.384
f1/ALL	73	2206	4.872	3.921	6.055
f2/ALL	73	2162	5.679	4.294	7.512

NOTE.—NaN, not a number.

When CSP variant is excluded from the model, the RRR for f1 versus c1 is 28.64 (95% CI = 12.33–66.55) and for f2 versus c1 is 33.39 (95% CI = 14.12–78.99). In contrast, the RRR for f1 versus f2 does not differ from 1 (RRR = 0.858, 95% CI = 0.602–1.222). These data support the hypothesis that c1 is better adapted to infect An. albimanus and f1 and f2 are better adapted to infect An. pseudopunctipennis.

Discussion

Adaptation of a parasite to its local host is predicted when the parasite is under stronger selection pressure and/or has a greater evolutionary potential (e.g., larger population sizes, shorter generation times, and higher mutation rates) than the host (Gandon et al. 1996; Gandon and Michalakis 2002; Brandt et al. 2007). We have identified distinct P. vivax genetic populations in southern Mexico and have demonstrated that they differ markedly from each other in their ability to successfully infect local mosquito vectors. Further, we showed that sympatric mosquito–parasite combinations are more compatible than allopatric ones, constituting evidence for local adaptation (Woolhouse et al. 2002). Our data support the view that asymmetries in the coevolutionary arms race between mosquito and parasite favor local adaptation of the parasite. Selection pressure is predicted to be greater for the parasite than for the mosquito because successful host exploitation is essential for parasite survival, whereas the cost of infection incurred by the mosquito (Lyimo and Koella 1992; Hogg and Hurd 1997) must be balanced against the cost of maintaining resistance mechanisms (Hurd et al. 2006). Additionally, in low transmission regions such as southern Mexico, where infection incidence among wild-caught mosquitoes is <0.01%, most hosts may never encounter a parasite.

Although 3 parasite populations were identified, 2—f1 and f2—are not distinguishable from each other with regard to mosquito infectivity, geographic distribution, or genetic diversity. One possible explanation for the identification of 2 populations in the foothills is that one may constitute a more recent immigration to the region and that not enough time has elapsed for complete admixture to have occurred between them. The distribution of f1/f2 haplotypes in the PCA plot as 2 widely separated foci connected by a smaller number of haplotypes is suggestive of this scenario (fig. 3), as is the absence of f2 in the 1997 collection (fig. 2), although we cannot completely rule out the possibility that f1/f2 constitute a single parasite population.

Both ecological and evolutionary factors associated with the mosquitoes may act to promoted divergence and restrict gene flow among parasite lineages adapted to different vectors. Anopheles albimanus and An. pseudopunctipennis differ in their ecological requirements, resulting in minimal overlap in both space and time. Anopheles albimanus is coastal and is found mainly at elevations <100 m. Anopheles pseudopunctipennis prefers higher elevations and is often the only vector above 600 m (Rubio-Palis and Zimmerman 1997). Although their distributions overlap in and around Tapachula where all 3 parasite populations were collected, we found little evidence for recombination between c1 and f1/f2 isolates (fig. 1). Anopheles albimanus densities peak during the rainy season (May–November) coinciding with a peak in coastal malaria infections. In contrast, An. pseudopunctipennis densities and the number of malaria infections in the foothills peak during the dry season (December–April) and are associated with pools forming in dry river beds (Savage et al. 1990) that are eliminated by heavy rains during the rainy season (Fernandez-Salas et al. 1994).

Anopheles albimanus and An. pseudopunctipennis are not closely related but rather belong to nonsister subgenera Nyssorhynchus and Anopheles, respectively (Krzywinski et al. 2001). They likely diverged >90 to 106 MYA, the estimated divergence time for sister subgenera Cellia and Anopheles (Krzywinski et al. 2006). Thus, these 2 mosquitoes present potentially divergent environments to developing malaria parasites. As such, adaptations for successful transmission in one species may not be advantageous in the other. The strong interaction between parasite genetic background and mosquito species that we observed in oocyst number supports this view.

A limitation of this study is its dependence on laboratory mosquito colonies that are likely to have greatly reduced genetic diversity compared with natural mosquito populations, due to selection and genetic drift (Norris et al. 2001). On the other hand, we have shown that the infectivity phenotype is consistent across a range of parasite genotypes collected during a 7-year period. Further infection experiments with recently established mosquito colonies are needed to directly address the contribution of the mosquito genotype.

Although the effects of artificial selection (i.e., drug pressure) on Plasmodium genetic diversity have been extensively studied (Wootton et al. 2002; Roper et al. 2004; Nash et al. 2005; Laufer et al. 2006; Nair et al. 2007), this is to our knowledge the first study to examine the effects of natural selection. Here we have identified strong local adaptation of malaria parasites to their vectors on a small geographic scale, resulting in distinct parasite populations. The previously observed association between csp alleles and mosquito species in southern Mexico can be attributed to the confounding effects of population structure. The presence of only a single c1/VK247 isolate in our sample was surprising given that CSP repeat variants are segregating freely in f1 and f2 and that c1 is more genetically diverse than f1 or f2. Nonetheless, the overall geographical distribution of CSP variants strongly suggests that CSP does not play a role in adaptation to these local vectors. Vector-mediated population structure may also explain the differences in susceptibility to csp variants observed for Anopheles aquasalis and Anopheles darlingi in Brazil (Silva et al. 2006).

Local adaptation and vector–parasite compatibilities may play an important role in generating population structure in Plasmodium, especially in Central and South America, where extensive genetic subdivision on small geographic scales has been observed for both Plasmodium falciparum and P. vivax. (Anderson et al. 2000; Machado et al. 2004; Ferreira et al. 2007). Local adaptation to vectors offers a possible explanation for the paradoxical finding of high genetic diversity coupled with strong linkage disequilibrium in P. vivax in Brazil (Ferreira et al. 2007), as well as the absence of a relationship between geographic and genetic distance observed for P. falciparum, also in Brazil (Machado et al. 2004). Spatial heterogeneity in the prevalence and distribution of lineages has also been reported for avian malaria on an even smaller geographic scale than this study (Wood et al. 2007). We have shown that genetically distinct P. vivax populations differ from each other in terms of vector compatibility, but the mechanism and genes responsible for these differences remain unknown. We have, however, identified the relevant components for a genetic cross (c1 vs. f1/f2) and subsequent genome-wide scan for loci under selection, a necessary next step to identify parasite genes contributing to vector–parasite compatibility.

Supplementary Materials

Supplementary tables 1 and 2 are available at Molecular Biology and Evolution online (http://www.mbe.oxfordjournals.org/).