Figures
Abstract
Background
Worldwide the mosquito Aedes aegypti (L.) is the principal urban vector of dengue viruses. Currently 2.5 billion people are at risk for infection and reduction of Ae. aegypti populations is the most effective means to reduce the risk of transmission. Pyrethroids are used extensively for adult mosquito control, especially during dengue outbreaks. Pyrethroids promote activation and prolong the activation of the voltage gated sodium channel protein (VGSC) by interacting with two distinct pyrethroid receptor sites [1], formed by the interfaces of the transmembrane helix subunit 6 (S6) of domains II and III. Mutations of S6 in domains II and III synergize so that double mutants have higher pyrethroid resistance than mutants in either domain alone. Computer models predict an allosteric interaction between mutations in the two domains. In Ae. aegypti, a Ile1,016 mutation in the S6 of domain II was discovered in 2006 and found to be associated with pyrethroid resistance in field populations in Mexico. In 2010 a second mutation, Cys1,534 in the S6 of domain III was discovered and also found to be associated with pyrethroid resistance and correlated with the frequency of Ile1,016.
Methodology/Principal Findings
A linkage disequilibrium analysis was performed on Ile1,016 and Cys1,534 in Ae. aegypti collected in Mexico from 2000–2012 to test for statistical associations between S6 in domains II and III in natural populations. We estimated the frequency of the four dilocus haplotypes in 1,016 and 1,534: Val1,016/Phe1,534 (susceptible), Val1,016/Cys1,534, Ile1,016/Phe1,534, and Ile1,016/Cys1,534 (resistant). The susceptible Val1,016/Phe1,534 haplotype went from near fixation to extinction and the resistant Ile1,016/Cys1,534 haplotype increased in all collections from a frequency close to zero to frequencies ranging from 0.5–0.9. The Val1,016/Cys1,534 haplotype increased in all collections until 2008 after which it began to decline as Ile1,016/Cys1,534 increased. However, the Ile1,016/Phe1,534 haplotype was rarely detected; it reached a frequency of only 0.09 in one collection and subsequently declined.
Conclusion/Significance
Pyrethroid resistance in the vgsc gene requires the sequential evolution of two mutations. The Ile1,016/Phe1,534 haplotype appears to have low fitness suggesting that Ile1,016 was unlikely to have evolved independently. Instead the Cys1,534 mutation evolved first but conferred only a low level of resistance. Ile1,016 in S6 of domain II then arose from the Val1,016/Cys1,534 haplotype and was rapidly selected because double mutants confer higher pyrethroid resistance. This pattern suggests that knowledge of the frequencies of mutations in both S6 in domains II and III are important to predict the potential of a population to evolve kdr. Susceptible populations with high Val1,016/Cys1,534 frequencies are at high risk for kdr evolution, whereas susceptible populations without either mutation are less likely to evolve high levels of kdr, at least over a 10 year period.
Author Summary
Constant use of pyrethroid insecticides has driven mosquito populations to develop resistance. In Aedes aegypti, the primary mosquito vector of dengue, yellow Fever, and chikungunya viruses, pyrethroid resistance is primarily associated with mutations in the voltage-gated sodium channel protein. One mutation occurs in codon 1,016 and involves a replacement of valine with isoleucine (Ile1, 016), and a second located in subunit 6 of domain III in codon 1,534, replaces phenylalanine with cysteine (Cys1,534). In Mexico, we found that Cys1,534 was present in the same mosquito collections that were previously analyzed for Ile1,016. In this study, we performed a linkage disequilibrium analysis on both Ile1,016 and Cys1,534 in Mexican collections from 2000–2012. Our analysis suggests that pyrethroid resistance requires the sequential evolution of the two mutations and that Cys1,534 must occur first and appears to enable the Ile1,016 mutation to survive.
Citation: Vera-Maloof FZ, Saavedra-Rodriguez K, Elizondo-Quiroga AE, Lozano-Fuentes S, Black IV WC (2015) Coevolution of the Ile1,016 and Cys1,534 Mutations in the Voltage Gated Sodium Channel Gene of Aedes aegypti in Mexico. PLoS Negl Trop Dis 9(12): e0004263. https://doi.org/10.1371/journal.pntd.0004263
Editor: Audrey Lenhart, Centers for Disease Control and Prevention, UNITED STATES
Received: June 29, 2015; Accepted: November 4, 2015; Published: December 11, 2015
Copyright: © 2015 Vera-Maloof et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited
Data Availability: All relevant data are within the paper.
Funding: This work was supported by Fogarty Training Grant 2D43TW001130-08 “Training in Dengue Prevention and Control.” The study was also supported by the National Institutes of Health/National Institute of Allergy and Infectious Diseases (NIH/NIAID; International Collaborations in Infectious Disease Research Program U01-AI-088647). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Worldwide Aedes aegypti (L.) mosquitoes are the principal urban vectors of dengue, chikungunya, and yellow fever viruses. Approximately 2.5 billion people (40% of the human population) currently live with the risk of dengue transmission. In Mexico, Ae. aegypti is the primary vector of the four dengue virus serotypes (DENV1-4), the causative agents of dengue fever (DF), dengue hemorrhagic fever (DHF) and dengue shock syndrome (DSS). Mexico is severely affected by DF, DSS, and DHF because all four dengue serotypes co-occur in most states of Mexico. A recent review of dengue disease in Mexico [2] reported an increase in incidences from 1.72 per 100,000 in 2000 to 14.12 per 100,000 in 2011.
Currently the most effective means to reduce dengue transmission by Ae. aegypti is through reduction of larval and adult populations. In Mexico larval reduction is accomplished chiefly through the application of the organophosphate temephos to peridomestic larval breeding sites and through physical source reduction or alteration of potential water-holding containers. Following recommendations of the official Mexican policy for vector control, (NOM-032-SSA2-2002), pyrethroids were almost exclusively used to control adults in and around homes from 1999 to 2010.
Pyrethroid insecticides prolong the opening of the voltage gated sodium channel protein (VGSC) in insect nerves to produce instant paralysis and ‘‘knock-down.” The α-subunit of VGSC has four repeat domains, labeled I-IV, each of which contains six transmembrane helix segments, S1-S6. Pyrethroids preferentially bind to the open state of vgsc by interacting with two distinct receptor sites formed by the interfaces of the transmembrane helix S6 of domains II and III, respectively [1]. The original computer modeling studies [3] suggest that simultaneous binding of pyrethroids to S6 in both domains II and III is necessary to efficiently lock sodium channels in the open state. These models also predict that mutations in the S6 of domain III allosterically alter S6 in domain II via a small shift of IIS6 thus establishing a molecular basis for the coevolution of S6 mutations in domains II and III in conditioning pyrethroid resistance.
In 2006 we described a mutation, Ile1,016, in the S6 of domain II in Ae. aegypti that is associated with very high knock-down resistance (kdr) to the pyrethroid insecticide permethrin in mosquitoes homozygous for this mutation. We examined collections of Ae. aegypti from Mexico during 1996–2009 [4] and found that the overall Ile1,016 frequency increased from 0.1% in 1996–2000, to 2%–5% in 2003–2006, to 38.3%–88.3% in 2007–2009 depending upon collection location. In 2010 another vgsc mutation was described in the S6 of domain III in Ae. aegypti that was also strongly correlated with kdr and involved a cysteine replacement (Cys1,534Phe) [5–7]. A general trend in these studies was that Cys1,534 frequencies were generally higher and increased more rapidly than Ile1,016 frequencies in natural populations.
Based upon these observations and on the dual binding model [3], we analyzed fresly collected DNA from Ae. aegypti for Ile1,016 and Cys1,534 while DNA previously analyzed for Ile1,016 [4] were tested for the presence of Cys1,534. The purpose of this study was to test the hypothesis that mutations in the S6 of domains II and III coevolve in a dependent manner through various allosteric interactions as suggested by computer models [3, 8]. An analysis of linkage disequilibrium was performed on the two alleles in 1,016 (Val 1,016 (susceptible), Ile 1,016(resistant)) and on the two alleles in 1,534 (Phe 1,534 (susceptible), Cys1,534 (resistant)) to assess whether alleles at 1,534 and 1,016 evolve independently or in a correlated fashion through epistasis.
Materials and Methods
Mosquitoes
Larval mosquitoes were collected from the locations mapped in Fig 1 and listed in Table 1. At each collection site, we collected immatures from at least 30 different containers in each of three different areas located at least 100 m apart. This included water storage containers and discarded trash containers such as plastic pails, tires, and cans. Larvae were returned to the laboratory where they were reared to adults and then identified to species. The Viva Caucel collection was west of the city of Merida in Yucatán State (20.9979639°, 089.7174611°). The Vergel collection was from eastern Merida (Fig 1) (20.9575694°, -89.5886889°). Both were collected in 2011 by Universidad Autónoma de Yucatán. DNA was isolated from individual adult mosquitoes by the salt extraction method [9] and suspended in 150 mL of TE buffer (10 mM Tris-HCl, 1 mM EDTA, pH 8.0). The SNP identification, allele-specific polymerase chain reaction (PCR), melting curve conditions, and genotype readings followed published procedures [6, 10–12].
VV = Val1,016 homozygotes, VI = Val1,016/Ile1,016 heterozygotes, II = Ile1,016 homozygotes, FF = Phe1,534 homozygotes, FC = Phe1,534/Cys1,534 heterozygotes, CC = Cys1,534 homozygotes for Ae. aegypti in Mexico from 1996 to 2012.
Testing for associations between vgsc genotypes and kdr phenotypes
The F3 generation of the Viva Caucel and Vergel strains were exposed to 25 μg permethrin (Chem Service, West Chester, PA) coated 250 mL Wheaton bottles. In each bottle approximately fifty 3–4 days old mosquitoes were exposed for one hour. Active mosquitoes were transferred to cardboard cups and frozen at -80°C and formed the ‘alive’ group. Knocked down mosquitoes were transferred to a second cardboard cup and placed into an incubator at 28°C and 70% humidity. After four hours, newly recovered mosquitoes were aspirated, frozen and labeled as ‘recovered’. The mosquitoes that remained inactive were scored as ‘dead’.
Linkage disequilibrium analysis
There are four potential 1,016/1,534 dilocus haplotypes: Val1,016/Phe1,534 (VF), Val1,016/Cys1,534 (VC), Ile1,016/Phe1,534 (IF), Ile1,016/Cys1,534 (IC). The number of times (Tij) that an allele at locus i = 1,016 appears with an allele at locus j = 1,534 was determined by the program LINKDIS [13]. The program then calculated composite disequilibrium frequencies [14] because the phase of alleles at 1,016 and 1,534 are unknown in double heterozygotes. An unbiased estimate of the composite disequilibrium coefficient Δij [14, 15] was calculated as:
Where N is the sample size and pi and pj are the frequencies of alleles at locus i = 1,016 and locus j = 1,534 respectively. Bayesian 95% Highest Density Intervals (HDI) around pi and pj were calculated in WinBUGS[16].
A χ2 test was performed to determine if significant disequilibrium exists among all alleles at 1,016 and 1,534. The statistic was calculated and summed over all two-allele-interactions [15]:
The linkage disequilibrium correlation coefficient Rij [15] is distributed from -1 (both mutations trans) to 0 (1,534 and 1,016 mutations occur independently), to 1 (both mutations cis) and therefore provides a standardized measure of disequilibrium:
Where the Ci term corrects for departures from Hardy-Weinberg expectations:
where Hobs (i) is the observed frequency of i homozygotes. Departures from Hardy-Weinberg expectations were also expressed as Wright’s inbreeding coefficient (FIS) and calculated as 1- (Hexp/2p (1- p)) where Hexp is the observed frequency of heterozygotes. A χ2 test of the hypothesis FIS = 0 with one degree of freedom is:
Results
Testing for associations between vgsc genotypes and kdr phenotypes
The locations of all sampling sites are shown in Fig 1 and the latitude and longitude coordinates are provided in Table 1. The sample sizes and numbers of the nine dilocus genotypes (Three 1,534 genotypes x Three 1,016 genotypes) are listed in Table 1. From a total of 615 treated mosquitoes in Viva Caucel, 17.6% (n = 108) were scored as alive, 15.6% (n = 96) as recovered and 66.8% (n = 411) as dead (Table 2). Genotypes at 1,016 and 1,534 were identified in 95 randomly chosen individuals from each group. From a total of 337 treated Vergel mosquitoes, 48.1% (n = 162) were scored as alive, 20.5% (n = 68) as recovered and 31.5% (n = 106) as dead. We randomly chose 95, 68 and 95 Vergel individuals from each group, respectively to obtain the genotypes at 1,016 and 1,534 (Table 2).
In Viva Caucel, the frequency of the Ile1,106 allele was 0.746 and the frequency of the Cys1,534 allele was 0.926 (Table 3), while in Vergel Ile1,016 was at a slightly higher frequency of 0.80 while the Cys1,534 allele was close to fixation at 0.988. The Ile1,106 and Cys1,534 alleles were in positive disequilibrium in Viva Caucel, but were only marginally significant in Vergel.
Genotypes at the 1,016 and 1,534 loci were not independent, in agreement with the linkage disequilibrium analysis in Table 4. Table 5 is a three-way contingency analysis of genotypes at loci 1,016, 1,534 and numbers alive or dead individuals in Viva Caucel. Numbers of alive were not independent of genotypes at the 1,016 locus; specifically, numbers of alive were significantly greater in Ile1,016 homozygous mosquitoes than in heterozygotes or in Val1,016 homozygotes. Numbers of alive were also not independent of genotypes at the 1,534 locus; specifically, numbers of alive were significantly greater in Cys1,534 homozygous mosquitoes than in heterozygotes or in Phe1,534 homozygotes. In general, numbers of alive in the Viva Caucel strain were not independent of genotypes at either locus.
However, a problem with this analysis is that genotypes at the two loci are not independent. In this and previous studies [10, 11], Ile1,016 homozygous mosquitoes have the greatest survival, while few, if any heterozygotes or Val1,016 homozygotes survive. To evaluate Cys1,534 genotypes independently of Ile1,016 homozygous mosquitoes, we only compared the three Cys1,534 genotypes among Ile1,016 heterozygotes and Val1,016 homozygotes. A significantly larger proportion of Cys1,534 homozygotes survived.
Table 5 also shows the contingency analyses of Vergel mosquitoes. Genotypes at the 1,016 and 1,534 loci were not independent, while they were marginally significant in the linkage disequilibrium analysis in Table 4. Numbers of alive were not independent of genotypes at the 1,016 locus again because numbers of alive were significantly greater in Ile1,016 homozygous mosquitoes than in heterozygotes or in Val1,016 homozygotes. Numbers of alive were however independent of genotypes at the 1,534 locus; specifically because Cys1,534 was almost fixed in the Vergel strain.
Table 5 also shows the three-way contingency analysis between genotypes at loci 1,016 and 1,534 and the numbers recovered or dead in Viva Caucel. As in Table 4, genotypes at the 1,016 and 1,534 loci were not independent. The numbers of recovered mosquitoes were not independent of genotypes at the 1,016 locus; specifically-numbers recovered were significantly greater in Ile1,016 homozygous mosquitoes than in heterozygotes or in Val1,016 homozygotes. Numbers of recovered were also not independent of genotypes at the 1,534 locus; specifically, numbers of alive were significantly greater in Cys1,534 homozygous mosquitoes than in heterozygotes or in Phe1,534 homozygotes.
In general, numbers of recovered in the Viva Caucel strain were heavily dependent on genotypes at both loci. An interesting difference between the two loci is that 32% (28/88) of Ile1,016 heterozygotes recovered while only 3.6% (1/28) of Cys1,534 heterozygotes recovered. This difference was significant (χ2 = 7.59, df = 1, p-value = 0.006).
Table 5 also shows the same analysis of recovery but in Vergel mosquitoes. Genotypes at the 1,016 and 1,534 loci were not independent, while they were marginally significant in the linkage disequilibrium analysis in Table 4. Numbers of recovered were not independent of genotypes at the 1,016 locus, again because numbers of recovered were significantly greater in Ile1,016 homozygous mosquitoes than in heterozygotes or in Val1,016 homozygotes. However, numbers of recovered were independent of genotypes at the 1,534 locus; specifically because Cys1,534 was approaching fixation in the Vergel strain.
Spatial and temporal analysis of genotype frequencies
Table 6 contains the frequencies of Ile1,016 and Cys1,534 and their Bayesian 95% HDI. FIS was significantly greater than zero (heterozygote deficiency) in two of the 36 collections where Ile1,016 and Val1,016 alleles were segregating. In contrast, a significant heterozygote deficiency occurred in eight of the 53 collections where Cys1,534 and Phe1,534 were segregating and an heterozygote excess occurred in two collections.
FIS and the associated probabilities from the χ2 test to test whether FIS = 0.
The frequencies of the Ile1,016 and Cys1,534 alleles from 1999 to 2012 are plotted in Fig 2. The Cys1,534 allele appears sooner and increases more rapidly than Ile1,016. Only the states of Veracruz and Chiapas had sufficient samples over the years to compare the spatial distributions of Ile1,016 and Cys1,534 (Fig 3). It is very clear that Ile1,016 and Cys1,534 were increasing in frequency much earlier in Veracruz state in eastern Mexico than in Chiapas state in southwestern Mexico. It is also clear that in both states Cys1,534 was increasing in frequency much earlier than in Ile1,016. Starting in 2002, the frequency of Cys1,534 was greater than or equal to that of Ile1,016. In a yearly comparison of Ae. aegypti collection sites, 80 out of 87 sites (Table 6) had a frequency of Cys 1,534 being greater than the frequency of Ile1,016. In 6 of the 7 cases where the frequency of Ile1,016 exceeded that of Cys1,534, the difference was only from 1–2% and values were not different (overlapping 95% HDI). Only in Martınez de la Torre in 2002 was there a credible difference of 9%.
Linkage disequilibrium analysis
Linkage disequilibrium analysis can only be performed in datasets where alleles are segregating at both loci. There were 34 datasets that met this criteria of the 87 collections listed in Table 1. Table 7 lists the state, city and year of the 34 datasets along with linkage disequilibrium correlation coefficient Rij and its associated χ2 values and the probability of a greater χ2. Ile1,016 and Cys1,534 were in disequilibrium in the majority (21/34 = 62%) of datasets. For the most part, alleles in 1,534 and 1,016 were evolving in a correlated, dependent fashion. However, this analysis does not provide specific information about the four haplotypes.
The frequencies of the four potential dilocus haplotypes are plotted by year in Fig 4. The frequency of the susceptible Val1,016/Phe1,534 (VF) haplotype remained high from 1999–2003 (Fig 4A). No collections were made again until 2008, by which time frequencies had dropped to 0–0.6. Four years later, VF was approaching extinction in all collections. Fig 4B plots the frequency of the Val1,016/Cys1,534 (VC) haplotype. From 1999–2003, VC frequencies remained low (0–0.10). By 2008, frequencies had increased to 0.1–0.75. Four years later, VC was declining in frequency in two collections and was increasing in four collections. A very different trajectory occurred for Ile1,016/Phe1,534 (IF) (Fig 4C). From 1999–2002, the IF frequency remained low and only reached as high as 0.1 in two collections. By 2008 frequencies were approaching extinction and four years later similar trends were seen, even though VC and IC frequencies had increased dramatically. Fig 4D is a plot of the frequency of the resistant Ile1,016/Cys1,534 (IC) haplotype. From 1999–2002, the IC frequency was low and only reached 0.1 in one collection. By 2008 frequencies had increased dramatically in all collections and continued to increase in all collections up to 2012 when frequencies ranged from 0.5–0.9.
A) Frequency of the susceptible Val1,016/ Phe1,534 (VF) haplotype, B) Frequency of the Val1,016/Cys1,534 (VC) haplotype, C) Frequency of the Ile1,016/Phe1,534 haplotype and D) Frequency of the resistant Ile1,016/ Cys1,534 (IC) haplotype.
Discussion
The frequency of Cys1,534 has increased dramatically in the last decade in several states in Mexico including Nuevo Leon in the north, Veracruz on the central Atlantic Coast, and Chiapas, Quintana Roo and Yucatan in the south. The linkage disequilibrium analysis on the Ile1,016 and Cys1,534 alleles in Ae. aegypti collected in Mexico from 2000–2012 (Table 7) strongly supports statistical associations between 1,534 and 1,016 mutations in natural populations. Furthermore, the dynamics of haplotype frequencies during that time suggest pyrethroid resistance in the vgsc gene requires the sequential evolution of 1,534 and 1,016 mutations. Fig 4C suggests that the Ile1,016/Phe1,534 haplotype has a low fitness, even when pyrethroids are being released. For this reason Ile1,016 is unlikely to have evolved independently. Instead it is much more likely that the Cys1,534 mutation evolved first but conferred only a low level of resistance. This conjecture is strongly supported by the fact that in 80 of 87 collections (92%), the frequency of Cys1,534 was greater than the frequency of Ile1,016.
The findings of this study are different in many respects from those in a study of a Tyr1,575 substitution in Anopheles gambiae that occurs just beyond the S6 of domain III, within the linker between domains III and IV [17]. This linker contains a sequence of three amino acids (IFM) that close the sodium channel pore following activation, block the influx of sodium into the cell and restore the membrane resting potential. In contrast, Cys1,534 in Ae. aegypti occurs in the S6 of domain III. This is close to a Met1,524Ile substitution that has been associated with knockdown resistance in Drosophila melanogaster [18] and a Phe1,538Ile mutation that reduces sensitivity to deltamethrin in arthropods and mammals [19, 20].
Mutations in S6 of domain II, such as Phe1,014, Ser1,014 in An. gambie and Ile1,016 and Gly1,016 in Ae. aegypti are not directly in the binding pocket, but affect the resistance phenotype by preventing binding of insecticides and changing the conformation of the VGSC [3, 21]. In contrast, a binding site located in a hydrophobic cavity delimited by the IIS4-S5 linker and the IIS5/IIIS6 helices has recently been proposed [22] that it is accessible to the lipid bilayer and lipid-soluble insecticides. The methyl-cyclopropane (or equivalent structure) of pyrethroids and the trichloromethyl group of DDT appear to be critical features for the action of both pyrethroids and DDT. Both insecticides fit into a slot in a small pocket in the main hydrophobic cavity, flanked by Val1,529 and Phe1,530 on IIIS6. The binding site is formed upon opening of the sodium channel and is consistent with observations that pyrethroids bind preferentially to open channels. This binding pocket includes several known mutations in the S6 of domain III that reduce sensitivity to pyrethroids. Two nearby residues (Gly1,535 and Phe1,538) have been previously implicated in resistance in other insect species (23).
Study in which An. gambiae mosquitos were collected from a range of approximately 2000 km throughout West/Central Africa and had Tyr1,575 occurring at frequencies up to 30% in both M and S forms. Even though this mutation is seen over a large range of the continent, only a single Tyr1,575 haplotype occurred with a Phe1,014 haplotype background (possibly analogous in function to Ile,1016), which infers strong positive selection acting on a recent mutant [17]. In contrast to the present study, Phe1,014 is almost fixed in West Africa and the Tyr1,575 allele is increasing in frequency in M form but not in S form. Thus in contrast to the apparent evolution of Ile1,016 on a Cys1,534 background as reported here in An. gambiae, Tyr1,575 appears to have evolved on a Phe1,014 background. There are many potential reasons for this difference including the possibility that mutations within the S6 of domain III may produce a different resistance mechanism and have a different impact on fitness than mutations in the linker between domains III and IV. It is also possible that the specific changes of amino acids at these sites are unique and may confer different resistance phenotypes. In either case it seems likely that one of the mutations compensates for deleterious fitness effects of the other mutation and/or confers additional resistance to insecticides.
An interesting difference between the two mutations in the present study is that 32% of Ile1,016 heterozygotes recover from pyrethroid exposure but only 3.6% of Cys1,534 heterozygotes recover. Thus while Cys1,534 in synergy with Ile1,016 may confer greater survival following pyrethroid exposure, Ile1,016 may confer a greater ability to recover following knockdown in heterozygotes.
There was evidence of heterozygote deficiency in eight of the 53 collections and the average FIS among these eight collections was large and positive (0.580) while the average among all collections was 0.052. This suggests that the fitness of Phe1,534 and Cys1,534 homozygotes may be greater than the fitness of G/T heterozygotes (i.e. underdominance). While these parameters have been estimated at the 1,016 locus [23], no similar studies have involved the 1,534 locus and so the stability point beyond which the frequency of either allele would increase has not been determined. Since the Cys1,534 confers some degree of pyrethroid resistance (Tables 2–5), directional selection could increase the frequency of Cys1,534 beyond the underdominance stability point, at which stage the frequency of Cys1,534 would rapidly increase towards fixation.
Little is known of other mutations in the Ae. aegypti vgsc that may affect pyrethroid resistance. Codon 989 in the “super-kdr” region of domain II was assessed and no mutations were found [11]. Ile, Met and Val alleles occur at codon 1,011 [11] but these alleles were not associated with resistance in our initial survey of 1,318 mosquitoes from the 32 strains throughout Latin America [11]. The recombination dynamics of the Ae. aegypti vgsc are also poorly understood. Analysis of segregation between alleles at the 1,011 and 1,016 codons in F3 showed a high rate of recombination even though the two codons are only separated by a approximately 250 bp intron [11]. A maximum parsimony phylogeny of the intron spanning exons 20 and 21 in 88 mosquitoes with different genotypes in exons 1,011 and 1,016 indicated the presence of three clades with bootstrap support > 80%. These were arbitrarily labelled clades 1–3. The frequencies of Ile1,011, Met1,011, Val1,011, Val1,016, Ile1,016 and Gly1,016 (from Thailand only) were then compared among the three clades. The frequency of Ile1,011 was distributed independently among the three clades, as was Val1,011 and Met1,011. However, there was a very evident excess of Val1,016 alleles in clade 1 and an excess of Ile1,016 alleles in clade 2. Ile1,016 alleles occurred in disequilibrium with a large number of segregating sites in the intron and a large excess of Ile1,016 alleles were found to be associated with clade 2 in the phylogenetic analysis. This pattern is consistent with a hypothesis that a genetic sweep of the Ile1,016 allele and proximate intron sequences has occurred through DDT exposure and subsequently pyrethroid selection. Furthermore, the genetic sweep was recent enough that there has been insufficient time for recombination to disrupt the disequilibrium between the Ile1,016 allele and proximate intron sequences.
Recent work on the dual binding model may shed some light on the next steps in the evolution of pyrethroid resistance in the vgsc [8]. The Tyr1,575 mutation in An. gambiae was introduced alone into an Ae.aegypti sodium channel (AaNav1-1) [8] and then in combination with Phe1,014. Both substitutions were then functionally examined in Xenopus oocytes [8]. Tyr1,575 alone did not alter AaNav1-1 sensitivity to pyrethroids. However, the Tyr1575- Phe1014 double mutant was more resistant to pyrethroids than the Phe1014 mutant channel alone. Further mutational analysis showed that Tyr1,575 could also synergize the effect of Ser1,014 and Trp1,014, but not Gly1,014, or other pyrethroid-resistant mutations in subunit 6 of domains I or II. Computer modeling predicted that Tyr1,575 allosterically alters pyrethroid binding via a small shift of the subunit 6 of domain II. This establishes a molecular basis for the coexistence of Tyr1,575 with Phe1,014 in pyrethroid resistance, and suggests an allosteric interaction between IIS6 and Loop III/IV in the sodium channel.
The rapid increase in Cys1,534 (Fig 4B and 4D) cannot be the result of neutral forces such as genetic drift or founder’s effects. Parallel increases in Cys1,534 frequency occurred throughout Mexico. Even though the forces that caused an increase in the frequency of Cys1,534 are unclear, our results suggest that Ile1,016 in domain IIS6 arose from a Val1,016/Cys1,534 haplotype and was rapidly selected possibly because double mutants confer higher pyrethroid resistance. When combined with Phe1014, the Tyr1,575 mutation in An. gambiae increased resistance to permethrin and deltamethrin by 9.8- and 3.4-fold, respectively [8].
Fig 5 illustrates two models for the evolution of 1,534 and 1,016 mutations. Model 1 proposes that the 1,534 and 1,016 mutations occurred independently and became cis by crossing over. Model 2 instead proposes that 1,534 mutations occurred first because 1,016 mutations confer low fitness. Ile1,016 mutations then arose on a Val1,016/Cys1,534 background. These results suggest that knowledge of the frequencies of both 1,534 and 1,016 mutations are important to predict the potential of a population to evolve kdr. Obviously, the frequency of Ile1,016 by itself is a poor predictor (Fig 4C). Populations that are pyrethroid susceptible, but have high Val1,016/Cys1,534 frequencies, are at high risk for rapid kdr evolution. If our experience in tracking the frequencies of Ile 1,016 and Cys1,534 mutations over the past 15 years can be extended to other Ae. aegypti populations, then populations with intermediate to high frequencies of Cys1,534 might only be susceptible for 5–10 years. Conversely, pyrethroid susceptible populations without either mutation are unlikely to develop kdr quickly and might be susceptible for up to 10–15 years.
Model 1 proposes that the 1,532 and 1,016 mutations occurred independently and became cis through crossing over. Model 2 instead proposes that 1,532 mutations occur first because 1,016 mutations confer low fitness. Ile1,016 mutations then arise on a Val1,016/Cys1,534 background.
Author Contributions
Conceived and designed the experiments: FZVM KSR WCB. Performed the experiments: FZVM KSR AEEQ SLF. Analyzed the data: FZVM WCB. Contributed reagents/materials/analysis tools: WCB. Wrote the paper: FZVM KSR WCB.
References
- 1. Du YZ, Nomura Y, Satar G, Hu ZN, Nauen R, He SY, et al. Molecular evidence for dual pyrethroid-receptor sites on a mosquito sodium channel. P Natl Acad Sci USA. 2013;110(29):11785–90.
- 2. Dantes HG, Farfan-Ale JA, Sarti E. Epidemiological Trends of Dengue Disease in Mexico (2000–2011): A Systematic Literature Search and Analysis. Plos Neglect Trop D. 2014;8(11).
- 3. O'Reilly AO, Khambay BPS, Williamson MS, Field LM, Wallace BA, Davies TGE. Modelling insecticide-binding sites in the voltage-gated sodium channel. Biochem J. 2006;396:255–63. pmid:16475981
- 4. Garcia GP, Flores AE, Fernandez-Salas I, Saavedra-Rodriguez K, Reyes-Solis G, Lozano-Fuentes S, et al. Recent Rapid Rise of a Permethrin Knock Down Resistance Allele in Aedes aegypti in Mexico. Plos Neglect Trop D. 2009;3(10).
- 5. Harris AF, Rajatileka S, Ranson H. Pyrethroid Resistance in Aedes aegypti from Grand Cayman. Am J Trop Med Hyg. 2010;83(2):277–84. pmid:20682868
- 6. Yanola J, Somboon P, Walton C, Nachaiwieng W, Somwang P, Prapanthadara LA. High-throughput assays for detection of the F1534C mutation in the voltage-gated sodium channel gene in permethrin-resistant Aedes aegypti and the distribution of this mutation throughout Thailand. Tropical Medicine & International Health. 2011;16(4):501–9.
- 7. Alvarez LC, Ponce G, Saavedra-Rodriguez K, Lopez B, Flores AE. Frequency of V1016I and F1534C mutations in the voltage-gated sodium channel gene in Aedes aegypti in Venezuela. Pest Manag Sci. 2015;71(6):863–9. pmid:24935645
- 8. Wang LX, Nomura Y, Du YZ, Liu NN, Zhorov BS, Dong K. A Mutation in the Intracellular Loop III/IV of Mosquito Sodium Channel Synergizes the Effect of Mutations in Helix IIS6 on Pyrethroid Resistance. Mol Pharmacol. 2015;87(3):421–9. pmid:25523031
- 9.
Black WC, DuTeau NM. RAPD-PCR and SSCP analysis for insect population genetic studies. In: Crampton J BC, Louis C, editor. The Molecular Biology of Insect Disease Vectors: A Methods Manual. New York: Chapman and Hall 1997. p. 361–73.
- 10. Saavedra-Rodriguez K, Strode C, Suarez AF, Salas IF, Ranson H, Hemingway J, et al. Quantitative Trait Loci Mapping of Genome Regions Controlling Permethrin Resistance in the Mosquito Aedes aegypti. Genetics. 2008;180(2):1137–52. pmid:18723882
- 11. Saavedra-Rodriguez K, Urdaneta-Marquez L, Rajatileka S, Moulton M, Flores AE, Fernandez-Salas I, et al. A mutation in the voltage-gated sodium channel gene associated with pyrethroid resistance in Latin American Aedes aegypti. Insect Mol Biol. 2007;16(6):785–98. pmid:18093007
- 12. Urdaneta-Marquez L, Bosio C, Herrera F, Rubio-Palis Y, Salasek M, Black WC. Genetic relationships among Aedes aegypti collections in Venezuela as determined by mitochondrial DNA variation and nuclear single nucleotide polymorphisms. Am J Trop Med Hyg. 2008;78(3):479–91. pmid:18337347
- 13. Black WC, Krafsur ES. A Fortran Program for the Calculation and Analysis of 2-Locus Linkage Disequilibrium Coefficients. Theor Appl Genet. 1985;70(5):491–6. pmid:24253058
- 14. Cockerham CC, Weir BS. Digenic Descent Measures for Finite Populations. Genet Res. 1977;30(2):121–47.
- 15. Weir BS. Inferences About Linkage Disequilibrium. Biometrics. 1979;35(1):235–54. pmid:497335
- 16. Lunn DJ, Thomas A, Best N, Spiegelhalter D. WinBUGS—A Bayesian modelling framework: Concepts, structure, and extensibility. Stat Comput. 2000;10(4):325–37.
- 17. Jones CM, Liyanapathirana M, Agossa FR, Weetman D, Ranson H, Donnelly MJ, et al. Footprints of positive selection associated with a mutation (N1575Y) in the voltage-gated sodium channel of Anopheles gambiae. P Natl Acad Sci USA. 2012;109(17):6614–9.
- 18. Pittendrigh B, Reenan R, ffrenchConstant RH, Ganetzky B. Point mutations in the Drosophila sodium channel gene para associated with resistance to DDT and pyrethroid insecticides. Mol Gen Genet. 1997;256(6):602–10. pmid:9435785
- 19. He HQ, Chen AC, Davey RB, Ivie GW, George JE. Identification of a point mutation in the para-type sodium channel gene from a pyrethroid-resistant cattle tick. Biochem Bioph Res Co. 1999;261(3):558–61.
- 20. Wang SY, Barile M, Wang GK. A phenylalanine residue at segment D3-S6 in Nav1.4 voltage-gated Na+ channels is critical for pyrethroid action. Mol Pharmacol. 2001;60(3):620–8. pmid:11502895
- 21. Soderlund DM, Knipple DC. The molecular biology of knockdown resistance to pyrethroid insecticides. Insect Biochem Molec. 2003;33(6):563–77.
- 22. Davies TGE, Williamson MS. Interactions of pyrethroids with the voltage-gated sodium channel. Bayer CropScience Journal. 2009;62:159–77.
- 23. Barbosa S, Black WC, Hastings I. Challenges in Estimating Insecticide Selection Pressures from Mosquito Field Data. Plos Neglect Trop D. 2011;5(11).