Anti-folate combination therapies and their effect on the development of drug resistance in Plasmodium vivax.

Can we predict the rise and spread of resistance to multi-drug therapy in a more predictable manner? We raise this question after analyzing over 500 Plasmodium vivax isolates collected from different, geographically isolated regions of China for sequence variation in and around the dhfr and dhps genes. We find: that resistance lineages have arisen at least once in each region; that there appears to have been little movement of parasite populations between these areas; and that highly resistant parasites contain dhfr and dhps alleles that are in linkage disequilibrium. We show a direct relationship between this linkage disequilibrium and a parasite's fitness in the absence of drug pressure. Such fitness would increase the spread of drug resistant phenotypes and is thus a selectable trait. These conclusions raise questions about the appropriate use of some other drug combinations to prevent and treat infection.

Pyrimethamine has been used alone and in combination with sulphadoxine (SP) as a broad spectrum treatment for microbial infections which includes Plasmodium. Despite the fact that these drugs are rarely recommended for the treatment of P. vivax, selective pressure on the P. vivax dhfr and dhps genes is commonly observed in endemic areas12. This is of particular concern because, SP treatment remains the safest option for intermittent preventive treatment (IPT) of malaria among pregnant women and infants34.

Pyrimethamine is a competitive inhibitor of dihydrofolate reductase (dhfr) and sulphadoxine of dihydroteroate synthetase (dhps). Both dhfr and dhps are enzymes in the folate pathway, and inhibition of either affects the synthesis of dTMP and methionine5. Mutations associated with resistance have been identified in the dhfr and dhps genes6789. Epidemiological studies, underpinned by in vitro analyses, suggest that drug pressure leads to the stepwise acquisition of mutations along with increasing levels of resistance51011. A moderate level of resistance is observed in 58R+117N double mutants of pvdhfr. An approximately 40-fold increase in resistance appears yeast carrying corresponding mutations12. This level is further augmented by mutations at codons 57 and 61, resulting in a 500-fold increase in resistance to pyrimethamine12. SP resistance is associated with mutations in codons 382, 383 and 553 of pvdhps and either double or triple mutant alleles of pvdhfr5. Clinical investigations most often indicate that treatment failure is associated with these changes213.

The use of drugs in combination to treat a pathogen serves to lengthen the effective life span of each individual drug14. It has been suggested that process of evolution of resistance to SP has its own modality because these drugs target two genes in the same pathway15. We ask whether the unique constraint on the evolution of SP resistance also alters the selective process itself.

Results

The geographical distribution of dhfr and dhps mutation haplotypes

Our study involves parasites from three different regions of China. Pyrimethamine and sulphadoxine were used in Yunnan Province and Hainan Island for a period of over 20 years starting in the early 1960s1617. In contrast, Central China used pyrimethamine during the same period1819.

We present here sequence comparisons for the pvdhfr gene from 561 isolates and for the pvdhps from 534 isolates. Variation at positions 13, 57, 58, 61, 117 and 173 of the pvdhfr and 382, 383, 512, 553 and 585 of pvdhps genes is shown in Table 1a and 1b. We refer to each haplotype based upon variation in these 6 positions for dhfr (i.e., IFSTSI is the pvdhfr wild type) and these 5 positions for dhps (i.e., SAKAV is the pvdhps wild type). We found additional mutations at position 99 of the dhfr gene and 399 of the dhps gene but these have not been shown to be involved in anti-folate resistance1220.

Table 1

Distribution of mutations in pvdhfr among field isolates from different area

Sampling areas	pvdhfr haplotypes
	13	57	58	61	117	173	n	(%)
Yunnan	I	L	R	M	T	I	80	36.36%
	I	I	R	M	T	I	52	23.64%
	I	L	S	M	T	I	12	5.45%
	I	I	S	M	T	I	7	3.18%
	I	L	R	T	T	I	2	0.91%
	I	L	R	T	T	F	1	0.45%
	I	L	S	T	T	I	1	0.45%
	I	F	S	T	T	I	1	0.45%
	I	F	R	T	N	I	18	8.18%
	I	F	S	T	N	I	11	5.00%
	I	F	R	T	N	L	2	0.91%
	L	F	R	T	N	I	1	0.45%
	I	L	R	T	N	I	1	0.45%
	I	F	R	M	N	I	1	0.45%
	I	F	S	T	S	I	19	8.64%
	I	L	R	T	S	I	4	1.82%
	I	F	R	T	S	I	4	1.82%
	I	L	S	T	S	F	3	1.36%
Hainan	I	F	R	T	N	I	58	100%
Central China	I	F	R	M	T	I	2	0.71%
	I	L	R	M	T	I	2	0.71%
	I	F	S	T	N	I	137	48.41%
	I	F	S	T	S	I	142	50.18%

Codon mutations are indicated in bold type;

n = the number of isolates;

% = the percent of total isolates;

Regional summaries

Yunnan: In this province SP was widely used and we know of no regions where pyrimethamine was used alone. Non-synonymous mutations resulted in amino acid changes in the dhfr protein in the following positions: 13L, 57L, 58R, 61M, 117N/T and 173L/F to the extent of 0.45%, 74.09%, 75.45%, 69.09%, 85.91% and 3.18%, respectively. Synonymous mutations were observed in codon 15, 38, 69 and 134. A small number of the isolates had the ancestral wild type allele (8.64%). The majority of isolates showed the presence of quadruple mutations at residues 57, 58, 61 and 117 (60.00%). Alleles with triple mutations and double mutations accounted for 11.82% of the total.

In pvdhps, mutation in codon 382, 383 and 553 occurred most frequently (20.67%, 82.21% and 57.21%). One non-synonymous mutation, 512T/M, and 3 synonymous mutations (codons 372, 495, 561) were also found. The most prevalent haplotype is double mutant allele SGKGV (47.60%), followed by wild type (17.79%). Single mutant alleles SGKAV, double mutant alleles AGKAV and triple mutant alleles AGKGV were identified in isolates from Yunnan as well (12.50%, 12.50% and 8.17%).

Hainan: SP was used on this island and we were unable to find reports of pyrimethamine being used by itself. Among Hainan isolates mutations at 58R and 117N in pvdhfr were found in 100% of the isolates. All samples displayed the identical nucleotide sequence and had no indels in the tandem repeat region. A majority of isolates had the wild type allele SAKAV (70.97%) in pvdhps. Alleles with 383G mutation were identified (29.03%).

Central China area: Pyrimethamine alone was used in this region. Sulphadoxine was not introduced. The ancestral wild type and the single mutant haplotype at codon 117 were prevalent (50.18% and 48.41%, respectively) for the pvdhfr gene among the 283 isolates. Mutations in codon 57, 58 and 61 were also observed (0.71%, 1.41% and 1.41%). The triple and quadruple mutant alleles were found in 2 isolates respectively. In pvdhps, all 264 examined isolates were wild type.

We found striking regionally centered differences in the dhfr and dhps genes of P. vivax. A total of 17 resistant haplotypes were identified for pvdhfr in Yunnan whereas only 1 and 3 resistant halplotypes were found in Hainan and Central China, respectively. In addition, more than 60% of the isolates in Yunnan had the highly resistant allelic types with quadruple mutations for pvdhfr gene while only 2 of 283 isolates had the allele in Central China area. For the pvdhps gene, all the isolates from Central China were the wild type whereas the majority of the isolates from Yunnan (82.21%) had alleles with various mutations.

Evaluation of linkage disequilibrium between paired alleles of pvdhfr and pvdhps

To investigate whether mutations of pvdhfr and pvdhps, which located different chromosome, are linked, we examined linkage disequilibrium (LD) for each of the 21 pairs of SNPs in genes of isolates from Yunnan. As shown in Table 2, trans-chromosome LD was observed between 117 and 383 (r2 = 0.2479, D' = 0.5996, p<0.001), 61 and 383 (r2 = 0.2078, D' = 0.6564, p<0.001), 57 and 383 (r2 = 0.1618, D' = 0.5160, p<0.001). Intragenic LD was also measured between 57 and 61, 57 and 117, 58 and 61, 58 and 117, 61 and 117, 383 and 553.

Table 2

Linkage disequilibrium values D' and r2

	dhfr				dhps
	57	58	61	117	382	383	553
57	*	0.4084	0.9711	0.6672	0.269	0.516	0.2658
58	0.1492	*	0.4961	0.6220	0.1045	0.2931	0.2767
61	0.7482	0.1746	*	1	0.3824	0.6564	0.3502
117	0.1865	0.1813	0.3325	*	0.8134	0.5996	0.5004
382	0.0069	0.0009	0.0176	0.0264	*	1	0.2992
383	0.1618	0.0584	0.2078	0.2479	0.0579	*	1
553	0.0328	0.0318	0.0718	0.0487	0.0361	0.2822	*

|D′|, above diagonal;

r2, below diagonal;

Microsatellite diversity flanking dhfr alleles

We determined microsatellite polymorphism on chromosome 5 in a region spanning approximately ±100 kb around the pvdhfr gene for 369 parasite isolates. We used He as the measure of variation. Six microsatellites were polymorphic with 57 different alleles identified. Isolates bearing sensitive alleles were highly variable with a mean He of 0.759, similar to that seen for neutral markers. We observed a clear diminution in He value at the microsatellites nearest the dhfr gene of the isolates with the resistant alleles as shown in Supplementary Table S4 and Fig. S1. He values at markers situated more than 37 kb in the downstream direction and more than 38 kb in the upstream direction were 0.535 ~ 0.840.

Multiple origins of resistant dhfr alleles

To better understand the origin and spread of pyrimethamine resistance in China, the microsatellite data set was analyzed using the Network. Where multiple haplotypes are present as a result of independent mutations, the haplotypes appear far from one another and to be separated by a large number of mutational steps on the diagram. In cases where multiple haplotypes are present as a result of mutations accumulating in one genetic background, the haplotypes will cluster and be separated by a limited number of mutational steps. The 369 P. vivax isolates were classified into 338 haplotypes based on the 6-loci, microsatellite profile around dhfr gene. Compared with the scatter pattern of dhfr wild type, the microsatellite haplotypes of the 117N alleles in Central China Area and 58R/117N alleles in Hainan were clustered, respectively. The microsatellite haplotypes in 57I(L)/58R/61M/117T alleles in Yunnan were more dispersed but still clustered (Fig. 2a). Hence, we were able to determine that pyrimethamine resistance had arisen independently in the Yunnan, Hainan and Central China areas.

Figure 1

(a) Geographic distribution of P.vivax isolates collected in China. (b) Diagram of pvpppk-dhps, pvdhfr and the flanking microsatellites, and the strategy of amplifying.

Figure 2

(a) Network diagram depicting three independent origins for different mutant dhfr alleles in China. The network diagram showing the origins for mutant pvdhfr alleles in China based on the 6-loci microsatellite profile around dhfr. The size of the circle is proportional to the number of isolates harboring particular haplotype. The proportion of dhfr alleles on that haplotype background are indicated in pie charts. The lines represent evolutional steps. The gray dots represent hypothetical median vectors. (b) Network diagram assessing the evolution of the resistance-conferring mutations in P.vivax dhps-pppk genes in Yunnan. The network diagram showing the evolution of the mutations in pvdhps-pppk in Yunnan based on the mutation in dhps-pppk gene and intron region. The size of the circle is proportional to the number of isolates showing particular allele. The lines represent evolutional steps connecting haplotypes.

Network diagram of P. vivax dhps and intron region

We examined the pppk gene and intron region based on the mutation analysis of the dhps coding region to assess the evolution of these resistance-conferring mutations. Besides the mutations in pvdhps coding region, we found other mutations in intron regions (intron1 G141T/A262T, intron2 G2327A/G2389T/C2410T) and codon M205I in pppk gene (Supplementary Table S5). In addition, indels existed in the repetitive domain at the pvpppk-dhps C-terminal region.

All the mutations (except the indels) were taken into account in the network diagram describing variation in Yunnan (Figure 2b). We suggest that mutation A383G appeared first, and was then followed by either of the double mutant alleles S382A/A383G or A383G/A553G. Either double mutant can further evolve to S382A/A383G/A553G. The other triple mutant alleles A383G/K512M/A553G, A383G/K512T/A553G derived from the double mutant A383G/A553.

Discussion

Many pathogens are dependent upon their own metabolic pathway to produce folate while humans depend entirely on salvage pathways to meet their metabolic needs. This makes the folic acid synthetic pathway an attractive target for treating humans for various infections. While pathogens that are resistant to single drugs often arise quickly, the combining of drugs for the treatment of infections lengthens the useful life span of each individual drug. Sulphadoxine and pyrimethamine, which target separate enzymes in this pathway, are often used in combination to treat a variety of infections. We ask whether mutants that have evolved in an environment that targets multiple enzymes in the same pathway have been selected solely on the basis of individual drug-enzyme interaction. Our results indicate that this is not the case. Extended drug pressure acts to favor resistant parasites on the basis of the association of particular combinations or pairs of enzymes. This raises certain concerns about the way multi-drug therapy (MDT) is applied.

If not for geographical isolation or adaptive disadvantages, the earliest resistance-conferring alleles would have swept from a single origin into other endemic regions in China. This did not happen. As shown in Table 1a and 1b, each region has its own distinctive pattern of dhfr and dhps variation. We analyzed over 500 samples from Yunnan, Hainan Island and Central China area by monitoring the pvdhfr and pvdhps genes and their surrounding microsatellite sequences. We found that each region had a distinctive pattern of dhfr variation indicative of past anti-malarial drug use (Table 1a).

We therefore monitored microsatellite polymorphism in the regions surrounding the dhfr genes for evidence of recent origins2122. Three distinct groups appeared distant from one another on the Network diagram as described in the Results section (Fig. 2a). This result supports the idea that resistance has multiple origins although we are aware that constraints on allele size at microsatellite loci can be a confounding factor23. We proceeded to dissect the population distribution of each area.

Our analysis of samples from Central China area provided an interesting data set for comparisons. Pyrimethamine was added to table salt as a broad-spectrum antibiotic for over 20 years, continually exposing parasites to a higher dose of the drug here than in the other provinces. Sulphadoxine was never introduced into the area.

A survey of recently collected isolates in the area indicated that today's population is made up of nearly equal numbers of wild-type parasites and those that carry a single mutation in their dhfr gene (Table 1a). We also found just two isolates (0.71% of the isolates) with quadruple dhfr mutant (ILRMTI). Microsatellite analysis has shown that these have arisen in Central China. They represent a small percentage of the population relative to wild type alleles suggesting that they carried a higher cost to fitness. Below, we contrast this with their abundance in the isolates of Yunnan Province.

The absence of mutation in the dhps gene from Central China is consistent with the fact that sulphadoxine was not introduced into the area. This absence also indicates that movement of resistant strains into Central China has not occurred to any significant extent.

The use of SP treatment of P. falciparum on Hainan Island was limited to several villages in the southern and southwestern regions of the island from between 1968 and 1972. It was also used as a partner drug with artemisinins or primiquine for the treatment of chloroquine-resistant P. falciparum. Although the drug's use on the island was limited, we found that the P. vivax isolates from Hainan Island all carried the pvdhfr double mutant, S58R/S117N, indicating that all samples had a history of SP exposure. Only 29.03% of the isolates however, carried a mutation associated with sulphadoxine resistance. Although parasites with resistance to both drugs had a selective advantage during the period of drug usage, we presume that some reversion to the wild type dhps occurred subsequently.

In Yunnan, SP was used frequently as a prophylactic treatment for malaria between the mid 1960s and the early 1990s. The prophylactic use of SP meant that the parasite's exposure to pyrimethamine was very similar to that found in Central China, except for the presence of sulphadoxine. This makes the comparison of data from the two regions very interesting. Our analysis of 220 isolates showed that quadruple mutants at the dhfr represented over 60% of the population while wild type alleles represented less than 9%. There were 18 different dhfr alleles and 7 different dhps alleles. We noted that 36.36% of the population had the same quadruple mutation (ILRMTI) that had arisen in Central China but there represented less than 1% of isolates. The absence of sulphadoxine pressure in Central China, during the period of drug pressure, is the most obvious difference that could be associated with changes in the dhfr locus. How this could affect the dhfr locus remained unclear. We looked for other factors associated with dhfr resistance that might explain the difference.

The presence of extensive polymorphism in Yunnan Province allowed us to address some factors that are involved in the evolution of drug resistance with quantifiable data. Here we are concerned with the population diversity after the drug pressure has been removed. Factors affecting such diversity are reviewed by Anderson et al and include: no cost mutations, compensatory evolution, and genetic linkage or co-selection between the resistance markers and other selected markers242526.

For example, linkage disequilibrium is a measure of co-segregation of genes. Genes from different chromosomes, such as dhfr and dhps, would be expected to segregate independently unless there is selective pressure to favor the offspring with particular pairs or combinations of genes. In the case of SP pressure there would a selective force to favor the association of the resistant alleles during their evolution but, after the drug pressure had been removed, the most obvious reason for the association would be gone. Survival of mutant forms of dhfr and dhps in the population would then relate to their relative cost in fitness with respect to each other and wild type parasites.

We showed that the highly resistant dhfr and dhps alleles have remained in a state of disequilibrium. The fact that they are in disequilibrium indicates that a selective force favors their association. In Yunnan, the association is apparent at only some of the positions related to drug resistance. We have shown that mutations in positions 57, 61, 117 of pvdhfr (chromosome 5) and position 383 of pvdhps (chromosome 14) of the isolates from Yunnan are significantly associated with linkage disequilibrium.

We investigated other available data sets in an attempt to find some common features that would lead us to a better understanding of this phenomenon1. In Thailand where a similar degree of polymorphism also exists among highly resistant triple or quadruple mutations as Supplementary Table S6, the pvdhfr and pvdhps alleles also exist in a state of linkage disequilibrium. Analysis of the isolates from Thailand indicates that position 57 and 61 of dhfr and 553 of dhps are strongly linked. We note that the association of genes exists but is centered on different codons here than in Yunnan. Further these are neither the codons most frequently associated with resistance nor those that independently result in resistance. The data from Thailand also indicates that the situation in Yunnan is not the result of a rare concurrence of multiple factors (“the perfect storm”). The single operative appears to be the protocol of drug delivery. The outcome is the development of highly drug-resistant pathogens that successfully compete with wild type parasites in the absence of drug pressure.

The above indicates that cooperative interaction between the genes is a factor in the evolution of the highly resistant forms of resistance. The association between disequilibrium and the longevity of particular alleles in the absence of drug pressure (e.g. the difference in representation of wild type dhfr and the ILRMTI mutant in Central China and Yunnan) is apparent. The rapid spread of quadruple mutation throughout Asia and Africa is consistent with this interpretation. Our data indicate that sulphadoxine pressure affects the course of the development of pyrimethamine resistance.

This report has implications for the use of multi-drug therapy and in particular for the control of malaria: (1) Although SP therapies have been replaced as the first-line malaria treatment, the intermittent preventive treatment (IPT) is still recommended for infants and pregnant women. The continued use of SP as a broad-spectrum antibiotic serves to increase the frequency of highly resistant phenotypes of malaria and the risk of fatal infection for pregnant women. (2) The use of SP can affect the evolution of the malaria parasites even when not prescribed for their treatment. (3) Other anti-folate combination, such as co-trimoxazole, that inhibits successive steps in the folate synthesis pathway should be noted for the selection to mildly resistant dhfr alleles in malaria. In addition, co-trimoxazole is widely used as a broad-scale treatment for pathogen infections and as a prophylactic for people living with HIV2728. Mutations of the Pneumocystis dhfr and dhps genes are observed and associated with prophylaxis failures2930. Highly resistant mutations in both genes can be expected when the drug pressure is increased, as observed with malaria. (4) Although our observations relate to SP treatment, we must question whether the implications extend beyond SP treatment. Selection by other antimalarial combination therapy, such as artemisinin-based combination therapy (ACT), needs to be monitored similarly. As is the case described here there may be a drug level or condition at which selection is so intense that it alters the course of the evolution of resistance. A different form of association between the genes may be expected to occur, which accelerates the spread of highly drug resistant phenotypes. We therefore suggest that MDT, despite its great promise, has potentially dangerous downsides.

Methods

Study sites and sampling information

P. vivax isolates were collected from Yunnan, Hainan Island and Central China area (Figure. 1a). Yunnan province borders Burma, Laos, and Vietnam. Details on sites of study are presented in the Supplementary information.

Isolates were collected from symptomatic patients that were positive for P. vivax by microscope. Samples were collected after informed consent and the clinical protocol was approved by the Internal Review Board of Second Military Medical University. We obtained 277 samples from different areas along the Nu River in Yunnan (June 2006 to January 2009), 322 samples from Central China area (July 2007 to August 2008) and only 35 samples from Hainan island (July 2007 to September). We also included isolates collected previously (August 2003 to October) for this study to compensate sample insufficiency. The filter papers were thoroughly dried and stored in a sealed plastic bag with desiccant at −20°C.

Genotyping of Pvdhfr and Pvdhps Genes

DNA was extracted from blood spots on filter paper using QIAmp Mini kits (Qiagen, Valencia, CA, USA). The chromosomal localization of pvdhps and pvdhfr gene were shown in Fig. 1b. The pvdhfr and pvdhps genes were amplified and sequenced. Primer sequences and details of the PCR conditions are given in Supplementary Table S1. The nested products were sequenced on an ABI 3730XL Genetic Analyzer (Applied Biosune, Shanghai City, CN). Amino-acid sequences were compared with wild-type sequences (GenBank accession no. X98123 for pvdhfr and Genbank accession no. AY186730 for pvdhps), using MEGA software. All samples were subjected to neutral microsatellite (15 loci) analysis and multiply infected samples were excluded. In addition, isolates contained mixed alleles at one or more point mutation in codons based on multiple peaks in the sequencing were excluded in this study.

Microsatellite selection and typing

Six microsatellite markers were selected in the regions near the pvdhfr gene (±100 kb) (Figure 1b). The total length of each microsatellite fragment was greater than 100 bp with the reiterating sequence being as short as 2 nucleotides and as long as 8 nucleotides. The description of microsatellite loci, primers and thermal cycling parameters are given in Supplementary Table S2 and S3. Length variation of PCR products was measured on QIAxcel (QIAGEN, Switzerland). Fragment size polymorphism was analyzed using the Bio Calculator software.

Statistical analysis

Linkage disequilibrium was assessed for resistant forms of pvdhfr and pvdhps genes in Yunnan isolates. Pairwise LD was measured with Pearson's correlation (r2) and Lewontin's coefficient (D'). The calculation was performed using Arlequin 3.01 with Markov chain steps of 10,000 and 1000 dememorisation steps. The significance of associations was estimated using Chi-square test for a significance level of 0.001.

The genetic diversity for each microsatellite locus was measured by calculating the expected heterozygosity (He) and the number of alleles per microsatellite locus (A) using GenALEx 6. He was calculated for each locus as follows: , where n is the number of isolates genotyped for that locus and pi is the frequency of the ith allele. The sampling variance for He was calculated as . For He analysis samples were grouped by alleles. The allele types with over 40 samples from a single region are plotted in Supplementary Fig. S1.

To examine the probable number of origins of pyrimethamine resistant pvdhfr alleles in China, a median-joining network was constructed using 6 loci haplotypes in the software Network 4.6.0.0 (http://www.fluxus-engineering.com/sharenet.htm). This can be also used in analysising the relationships between putative sulphadoxine drug resistance–conferring mutations. Networks were chosen for analysis because of their emphasis on multi-furcating relationships and low divergence, in the better analysis the biological events of repetitive elements or genes subjected to selective pressures by constructing the phylogeny of regions with reticulate evolution.

Author Contributions

S.D., D.Z. and W.P. conceived and designed this study. S.D., R.Y. and D.Z. performed the experiments. S.D., R.Y., T.M. and W.P. analyzed the data. X.S. and H.Z collected the samples, S.D., T.M. and W.P. drafted the manuscript. All authors contributed to the interpretation of the study.