# Forty years of erratic insecticide resistance evolution in the Mosquito Culex pipiens

Pierrick Labbé, Claire Berticat, Arnaud Berthomieu, Sandra Unal, Clothilde

Bernard, Mylène Weill, Thomas Lenormand

► To cite this version:

Pierrick Labbé, Claire Berticat, Arnaud Berthomieu, Sandra Unal, Clothilde Bernard, et al. Forty years of erratic insecticide resistance evolution in the Mosquito Culex pipiens. PLoS Genetics, 2007, 3 (11), pp.2190-2199. 10.1371/journal.pgen.0030205. halsde-00452664

HAL Id: halsde-00452664

https://hal.science/halsde-00452664

Submitted on 30 Nov 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Forty Years of Erratic Insecticide Resistance Evolution in the Mosquito *Culex pipiens*

Pierrick Labbé¹,²*, Claire Berticat¹, Arnaud Berthomieu¹, Sandra Unal¹, Clothilde Bernard¹, Mylène Weill¹, Thomas Lenormand²

1 Equipe Génétique de l'Adaptation, Institut des Sciences de l'Evolution, CNRS Université Montpellier 2, Montpellier, France, 2 Centre d'Ecologie Fonctionnelle et Evolutive, UMR 5175, Montpellier, France

One view of adaptation is that it proceeds by the slow and steady accumulation of beneficial mutations with small effects. It is difficult to test this model, since in most cases the genetic basis of adaptation can only be studied a posteriori with traits that have evolved for a long period of time through an unknown sequence of steps. In this paper, we show how *ace-1*, a gene involved in resistance to organophosphorous insecticide in the mosquito *Culex pipiens*, has evolved during 40 years of an insecticide control program. Initially, a major resistance allele with strong deleterious side effects spread through the population. Later, a duplication combining a susceptible and a resistance *ace-1* allele began to spread but did not replace the original resistance allele, as it is sublethal when homozygous. Last, a second duplication, (also sublethal when homozygous) began to spread because heterozygotes for the two duplications do not exhibit deleterious pleiotropic effects. Double overdominance now maintains these four alleles across treated and nontreated areas. Thus, *ace-1* evolution does not proceed via the steady accumulation of beneficial mutations. Instead, resistance evolution has been an erratic combination of mutation, positive selection, and the rearrangement of existing variation leading to complex genetic architecture.

Citation: Labbé P, Berticat C, Berthomieu A, Unal S, Bernard C, et al. (2007) Forty years of erratic insecticide resistance evolution in the mosquito *Culex pipiens*. PLOS Genet 3(11): e205. doi:10.1371/journal.pgen.0030205

## Introduction

Adaptation is often envisioned as a slow and regular improvement, a view embodied by Fisher's geometrical model of adaptation, whereby mutations fix if they bring the current phenotype closer to an optimum [1,2]. However, the trajectory towards the optimum can take a truly tortuous path [3], as adaptation uses almost any allele that brings it closer to the optimum, regardless of its negative side effects. "Evolution is a tinkerer," as emphasized by Jacob [4], may indeed be more than a pretty metaphor. That beneficial mutations often have deleterious pleiotropic effects is well established (i.e., they generate a fitness cost [1,5-8]). Pleiotropy causes an "evolutionary inertia" [9] whereby beneficial mutations often only ameliorate the side effects of the last beneficial mutation. This process of "amelioration" [10] can follow a scenario a la Fisher [11] with modifiers or compensatory mutations occurring at different loci, or with allele replacement at the same locus (Haldane [12]). Both cases have been reported (e.g., [10,13-16]). However, it is perhaps less often appreciated that pleiotropy can have more dramatic consequences. First, the pleiotropic effects of beneficial mutations may be more complex than simply deleterious. Second, they may trigger the evolution of the genetic architecture and gene number. This paper illustrates these two aspects with the tortuous path taken during the evolution of insecticide resistance in *Culex pipiens* mosquitoes in southern France.

In natural populations of the mosquito *Culex pipiens*, various resistance genes have been selected over the course of ~40 y of control using organophosphorous (OP) insecti-

cides (see [17] for a review). OP insecticides kill by inhibiting acetylcholinesterase (AChE1) in the central nervous system. The genetic basis of OP resistance involves two main loci, the super-locus *Ester* and the locus *ace-1*, both of which have major resistance alleles. Strong direct and indirect evidence indicate that the resistance alleles have a "fitness cost"; i.e., they are selected against in the nontreated areas (see reviews in [17,18]). At the *Ester* locus, the different resistance alleles are not identical and one has slowly replaced the other [15,16].

At the *ace-1* locus, the resistance allele present worldwide, *ace-1*<sup>R</sup>, displays a single amino acid mutation (G119S), changing the glycine at position 119 into a serine. This large effect point mutation confers high resistance towards OP insecticides due to lower affinity and has arisen independ-ently in several mosquito species [19-21]. However, this mutation also exhibits a strong deleterious side effect: G119S is less efficient at degrading acetylcholine (ACh) than the susceptible variant of AChE1. Overall, G119S causes a more than 60% reduction in enzymatic activity in the absence of insecticide [22], which alters the optimal functioning of

**Editor:** David L. Stern, Princeton University, United States of America

**Received** May 16, 2007; **Accepted** October 3, 2007; **Published** November 16, 2007

**Copyright:** © 2007 Labbé 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.

**Abbreviations:** AChE1, acetylcholinesterase; OP, organophosphorous

* To whom correspondence should be addressed. E-mail: Pierrick.Labbe@ed.ac.uk

<figure><img src="image_2.png" /></figure>
## Author Summary

Adaptation is not always a straightforward process, and often results from natural selection tinkering with available variation. We present in this study just such a tortuous natural selection pathway, which allows the mosquito *Culex pipiens* to resist organophosphorous insecticides. In the Montpellier area, following the use of insecticide to control mosquito populations, a high-resistance allele of the insecticide target enzyme appeared. But this allele also displayed strong deleterious side effects. Recently, several duplicated haplotypes began to spread in natural population that put in tandem a susceptible and a resistant allele. We show that the duplicated haplotypes actually display reduced side effects compared to the resistant allele when in the heterozygous state, but also new and strong costs in the homozygote. This pattern leads to an unexpected equilibrium between four different alleles across treated and nontreated areas. The story of resistance in *C. pipiens* is indeed far from a slow progression toward a “perfect” adaptation. Rather, selection for resistance to insecticide is a long process of trial and error leading to an uncommon genetic architecture.

cholinergic synapses of the central nervous system and probably causes the various developmental and behavioral problems that have been identified in individuals carrying *ace-1*<sup>R</sup> [23–25]. Field surveys also quantified the overall fitness cost of this allele in the wild. Fitness values are 1 and 0.89 for a susceptible and resistant homozygote, respectively, in the absence of insecticide (i.e., when only the cost applies [26,27]). However, in a treated habitat, only *ace-1*<sup>R</sup>-carrying individuals survive.

A few years after the appearance and spread of *ace-1*<sup>R</sup>, several duplications of the *ace-1* gene appeared in wild populations, each involving a resistant and a susceptible copy of the *ace-1* gene on the same chromosome [28–30]. We refer to these new duplicated haplotypes as *ace-1*<sup>D</sup>. The *ace-1* copies involved in these duplicated haplotypes are barely distinguishable from local single-copy alleles. In particular, within a population, the sequence of the local *ace-1*<sup>D</sup> resistant copy is always identical to the sequence of the single-copy *ace-1*<sup>R</sup> allele [30]. In the Montpellier area of Southern France, the appearance of a duplicated haplotype was indirectly inferred and traced back to at least 1993 [29]. More recently, sequence data of *ace-1* confirmed the presence of duplication, but surprisingly also revealed the presence of two distinct duplicated haplotypes in Montpellier area (*ace-1*<sup>D2</sup> and *ace-1*<sup>D3</sup>) [30]. Other independent duplications were also found in the Caribbean (Martinique, *ace-1*<sup>D1</sup>) and Philippines (Palawan, *ace-1*<sup>D4</sup>) [30].

Gene duplication is an important type of mutation because it fosters the evolution of new functions [12,31–33]. In the case of *ace-1*, a strong genetic constraint drives resistance evolution, as the degree of resistance and the ability to degrade ACh trade off; OP and ACh molecules compete for the same active site of AChE1. Duplication may be a way to disentangle the two functions, i.e., by improving synapse signalling and mosquito's fitness while maintaining resistance.

Thus, our working hypothesis is that the spread of *ace-1* duplications is driven by the increased AChE1 activity in individuals carrying *ace-1*<sup>D</sup> compared to those carrying *ace-1*<sup>R</sup>. In the last survey eight y ago [29], indirect estimation indicated a rapid initial increase of *ace-1*<sup>D</sup>, suggesting a

strong selective advantage of *ace-1*<sup>D</sup> over *ace-1*<sup>R</sup> (in the range of 3%–6% of *ace-1*<sup>R</sup> fitness [29]). Based on these previous findings, we expected the duplicated haplotype to quickly replace *ace-1*<sup>R</sup>. The aim of this study was to test this prediction and to understand more precisely the fitness relationship between the various genotypes (single resistance or susceptible alleles, duplicated haplotypes and their heterozygotes) in treated and nontreated areas. Surprisingly, we show that *ace-1*<sup>D</sup> did not replace *ace-1*<sup>R</sup>. Several experiments indicated that duplications, while favored when heterozygous, exhibit very strong deleterious effects when homozygous. Our work thus highlights a novel scenario of adaptive evolution and trade-offs that hinge on the genetic architecture underlying the expression of resistance variation.

## Results

In the following, and according to the nomenclature used in Labbé et al [30], the duplicated haplotypes *ace-1*<sup>D2</sup> and *ace-1*<sup>D3</sup> will be denoted D₂ and D₃, susceptible copies being denoted D₂(S) and D₃(S), and resistant copies, D₂(R) and D₃(R), respectively. D refers to either D₂ or D₃. The single copy resistant allele (*ace-1*<sup>R</sup>) will be denoted as R. Finally, single susceptible alleles (*ace-1*<sup>S</sup>) will be denoted S. The only S allele present in the strains used in laboratory experiments and originating from the susceptible reference strain SLAB [34] will be denoted S<sub>SLAB</sub>.

We first analyzed the frequency variation of R and D (= D₂ + D₃) across treated and nontreated areas around Montpellier from 1986 to 2002, using a purely descriptive model (see Methods). As in previous analyses [26,27,29], we found that their frequency showed clinal variation: R and D are more frequent in the treated than the nontreated area (Table 1). The straightforward interpretation of these clines is that R and, to a lesser extent, D have a selective advantage in the presence of insecticide over S, but that they have a fitness cost in its absence. Migration leads to the smooth decline of frequency in the nontreated area and prevents fixation of resistance in the treated area. We also confirm that D was rare until around 1993, when it started to rapidly increase in frequency. Surprisingly, D frequency stopped increasing around 1996 and reached a plateau at about 22% (this holds when only summers clines are considered, analysis not shown) (Figure 1; Table 2).

However, this analysis does not distinguish D₂ and D₃ duplications. To assess the frequency and the distribution of D₂ and D₃ along the cline, five populations were sampled in 2005, on the same transect as for cline analysis, two in the treated area (Maurin and Lattes), one at the limit between the two areas (Distill), and two in the nontreated area (Viols and Ganges, Figure 2). Molecular tests available allow detection of D₂ and D₃ only when S<sub>SLAB</sub> is the only S allele present. Resistant females from these samples were thus mated with SLAB-TC males (S/S) and their offspring analyzed (see Methods). Despite the large crossing (700 females from each sample), the number of egg rafts was very low (this is usually observed for the *C. p. pipiens* subspecies, P. Labbé and M. Weill, personal observation), and some females were able to lay eggs several times, so that no reliable frequencies could be estimated. Nevertheless, molecular analyses revealed that

<figure><img src="image_3.png" /></figure>
<table><caption>**Table 1.** Parameters of the Best-Fitted Clines</caption><thead><tr><th rowspan="2">Season</th><th rowspan="2">Year</th><th colspan="2">R</th><th colspan="2">D</th><th rowspan="2">%TD</th></tr><tr><th><i>h</i></th><th><i>a</i> × 10<sup>-4</sup></th><th><i>b</i> × 10<sup>-2</sup></th><th><i>h</i></th><th><i>a</i> × 10<sup>-2</sup></th></tr></thead><tbody><tr><td rowspan="10">Summer</td><td>1984</td><td>0.71</td><td>—</td><td>—</td><td><i>0.01</i></td><td>—</td><td>—</td></tr><tr><td>1986</td><td>0.52</td><td>101.7</td><td>0</td><td><i>0.06</i></td><td>0</td><td>0.00</td></tr><tr><td>1991</td><td>0.86</td><td>3.0</td><td>1.4</td><td><i>0.05</i></td><td>0</td><td>58.20</td></tr><tr><td>1995</td><td>0.58</td><td>2.1</td><td>2.5</td><td>0.35</td><td>0</td><td>92.16</td></tr><tr><td>1996</td><td>0.52</td><td>0</td><td>3.2</td><td>0.22</td><td>0</td><td>79.43</td></tr><tr><td>1999</td><td>0.87</td><td>0</td><td>5.9</td><td>0.13</td><td>0</td><td>84.90</td></tr><tr><td>2001</td><td>0.76</td><td>0</td><td>4.5</td><td>0.24</td><td>0</td><td>84.89</td></tr><tr><td>2002</td><td>0.79</td><td>3.8</td><td>2.5</td><td>0.21</td><td>0</td><td>3.0</td></tr><tr><td>2005</td><td>0.81</td><td>0</td><td>3.8</td><td>0.19</td><td>0</td><td>83.31</td></tr><tr><td>Winter</td><td>1995</td><td>0.39</td><td>1.6</td><td>0.3</td><td>0.20</td><td>9.1</td><td>81.73</td></tr><tr><td>1996</td><td>0.34</td><td>2.9</td><td>0</td><td>0.19</td><td>0</td><td>5.0</td></tr><tr><td>1999</td><td>0.69</td><td>2.7</td><td>1.1</td><td>0.11</td><td>0</td><td>1.0</td></tr><tr><td>2002</td><td>0.54</td><td>2.8</td><td>0.2</td><td>0.34</td><td>61.3</td><td>3.2</td></tr><tr><td rowspan="3">Spring</td><td>1993</td><td>0.48</td><td>0</td><td>4.2</td><td><i>0.04</i></td><td>0</td><td>0</td></tr><tr><td>1996</td><td>0.74</td><td>0.4</td><td>3.3</td><td>0.26</td><td>0</td><td>86.05</td></tr><tr><td>2000</td><td>0.59</td><td>11.8</td><td>0</td><td>0.20</td><td>1.6</td><td>68.30</td></tr></tbody></table>

The estimated value of parameters $h_{ij}$, $a_{ij}$, and $b_{ij}$ (see Materials and Methods, Equation 1) of the best model after likelihood ratio tests simplification procedure are given for each allele *i* and each sample *j*. Values not significantly different from 0 are in italics. R frequency decreases from the coast inland (i.e., slope parameters *a* and/or *b* different from zero in the best model) for each year and each season. D is detected significantly (i.e., *h* > 0) only since summer 1995, in agreement with Lenormand et al. [29]. Over-dispersion is generally low (<1.5) around fitted models (although higher for spring and summer clines). All fitted models explain more than 50% of the total observed deviance (%TD, as in [16]), except for spring 1993 (35%, probably because of the poor sampling). doi:10.1371/journal.pgen.0030205.t001

both duplications are now segregating in these populations (Table 3).

This plateau frequency pattern strongly suggests that equilibrium has been reached among the four alleles segregating in these populations, S, R, $D_2$, and $D_3$. The next

<table><caption>**Table 2.** Time Variation of Resistance Clines</caption><thead><tr><th>Model</th><th>Dev</th><th>df</th><th>F</th><th>p-Value</th><th>%TD</th></tr></thead><tbody><tr><td>Complete</td><td>18,342.17</td><td>63</td><td>—</td><td>—</td><td>74.92</td></tr><tr><td>α<sub>R</sub> = β<sub>R</sub></td><td>18,406.30</td><td>62</td><td>25.71</td><td>&lt; 0.001</td><td>72.26</td></tr><tr><td>α<sub>D</sub> = β<sub>D</sub></td><td>18,444.83</td><td>62</td><td>41.16</td><td>&lt; 0.001</td><td>70.67</td></tr><tr><td>α<sub>R</sub> = 0</td><td>18,381.46</td><td>62</td><td>15.75</td><td>&lt; 0.001</td><td>73.29</td></tr><tr><td>α<sub>D</sub> = 0</td><td>18,385.22</td><td>62</td><td>17.26</td><td>&lt; 0.001</td><td>73.13</td></tr><tr><td>β<sub>R</sub> = 0</td><td>18,604.97</td><td>62</td><td>105.36</td><td>&lt; 0.001</td><td>64.04</td></tr><tr><td>β<sub>D</sub> = 0</td><td>18,342.65</td><td>62</td><td>0.19</td><td>0.6615</td><td>74.90</td></tr></tbody></table>

The different model simplifications for alleles' frequencies evolution over the period 1986–2002 (see Materials and Methods, Equation 2) are presented in the Model column. The corresponding deviance (Dev) and degrees of freedom (df) are indicated. Each simplification is tested with likelihood ratio test (see Methods): the *F* value and the associated *p*-value are indicated. The percentage of the total deviance (%TD, as in [16]) explained by the corresponding model is also indicated. The model simplification indicates that (i) *R* changes in frequency at different rates ($\alpha_R$ different from $\beta_R$) before and after ~1995 ($t^* = 9.0$, support limits [8.6, 9.3]): it first decreases in frequency ($\alpha_R = -0.13$) and then increases ($\beta_R = 0.14$), (ii) *D* changes in frequency at different rates before and after ~1995 ($\alpha_D \neq \beta_D$): it first increases in frequency ($\alpha_D = 0.55$) and then remains constant ($\beta_D$ not different from zero). The fitted frequency variations are illustrated in Figure 1. doi:10.1371/journal.pgen.0030205.t002

step in the analysis was therefore to infer, from the clinal patterns observed in the last years (1999–2002), the different genotypic fitnesses (in treated and nontreated areas) that could correspond to this equilibrium. With four alleles, there are ten diploid genotypes and therefore nine relative fitnesses to estimate in each habitat. Since we do not have access to the respective frequencies of $D_2$ and $D_3$, some of these parameters cannot be estimated separately (see Methods). We report fitness estimates that assume that the two duplications have identical fitness effects. We developed a maximum-likelihood analysis by combining exact simulations with an optimization routine. We used previous estimation for migration: *C. pipiens*

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Time Variation of the Different Alleles' Maximum Frequency</figcaption></figure>

The evolution of the maximum frequency (coastal frequency) for each *ace-1* resistance allele estimated from the descriptive model is presented. The blue triangles and red squares represent, respectively, the maximum frequency of R and D (= $D_2 + D_3$) estimated for each cline independently (blue and red bars represent the confidence intervals, respectively), the blue and red lines representing the frequency estimated from the complete dataset conjointly (see text). Note that no tools are currently available to distinguish $D_2$ and $D_3$ in the field, only their total frequency being thus estimated (see Materials and Methods). doi:10.1371/journal.pgen.0030205.g001

<figure><img src="image_5.png" /></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Sample Site Locations in the Northwest Southeast Transect in the Montpellier Area</figcaption></figure>

Samples are indicated with black circles. The dashed line represents approximately the border between treated and untreated areas (modified from [16]). *C. pipiens* is present in the whole area. doi:10.1371/journal.pgen.0030205.g002

displays a migration of 6.6 km.generation⁻¹/² [26]. This analysis indicated that the stable cline we observe could be expected at equilibrium between migration and selection. To achieve this equilibrium, the different genotypes should display relative fitness, as reported in Table 4. Since the equilibrium situation seems to have been reached around 1995–1996, it therefore suggests that both duplications spread before that date. Importantly, this equilibrium requires a situation of double

overdominance: heterozygotes involving duplication, i.e., (D₂/S), (D₃/S), (D₂/D₃), (D₃/R), and (D₂/R), must have a higher fitness than R, D₂, and D₃ homozygotes in absence of insecticide. It also requires that R and S homozygotes are the best homozygote genotypes in the treated and nontreated area, respectively. With such a fitness scheme, each duplication could be maintained independently even if they have identical fitness effects. In addition, the presence of both allows them to reach a higher total frequency. More precisely, the overall D frequency when both are present should be intermediate between the frequency they could reach alone and twice this frequency, due to the overdominance of the (D₂/D₃) genotype.

The last step in the analysis was to confirm experimentally this possible pattern of double overdominance. To do so, we followed the relative success, over a single generation in the laboratory, of different pairs of genotypes in the absence of insecticide. In a first set of experiments, we showed that (D₂/S) and (D₃/S) largely outperformed (D₂/D₂) and (D₃/D₃) in terms of survival, development time, and fertility (Figure 3; Tables 5, 6 and 7). Homozygotes for both duplications indeed show an extremely low fitness with high mortality at pupation and emergence along with low fertility, such that it was nearly impossible to fix a strain for these duplications. By contrast, we showed that heterozygotes involving the two duplications (D₂/D₃) were as fit as (D₃/S) (Figure 3; Table 5). Finally, we compared the fitness of (D₃/S) and (D₃/R) with R homozygotes (Figure 4; Table 5). Here again, we found that the heterozygotes outperformed the homozygotes, although less strikingly. We were not able to perform this last comparison with the D₂ duplication because the survival and fertility of (D₃/D₂) homozygotes are too low to maintain a laboratory strain fixed for this duplication. Overall, these experiments directly confirm the fitness relationship deduced from the clinal pattern observed in the field. Homozygotes for either D₂ or D₃ have a severely reduced fitness, but heterozygotes for the two duplications, (D₂/S), (D₃/S), (D₂/D₃), and (D₃/R) perform well.

## Discussion

This study shows how *ace-1*, one of the two major genes involved in resistance to OP insecticide in the mosquito *Culex pipiens*, evolved in the last 40 y of control using insecticides. The evolution of resistance to insecticide in *C. pipiens* does not follow a classical scenario whereby a beneficial mutation

<table><caption>**Table 3.** Identifying Duplications in the Field</caption><thead><tr><th rowspan="2">Sample Site</th><th colspan="2">Susceptible Genotype</th><th colspan="6">Resistance Genotype</th></tr><tr><th>N</th><th>(S/S) frequency</th><th>N<sub>females</sub></th><th>N<sub>egg rafts</sub></th><th>[D<sub>3</sub>]</th><th>[D<sub>2</sub>]</th><th>(D<sub>2</sub>/D<sub>3</sub>)</th><th>(R/R)</th></tr></thead><tbody><tr><td>Maurin</td><td>200</td><td>0.01</td><td>700</td><td>17</td><td>1 (0.06)</td><td>8 (0.46)</td><td>1 (0.06)</td><td>7 (0.41)</td></tr><tr><td>Lattes</td><td>240</td><td>0.02</td><td>700</td><td>4</td><td>0 (0)</td><td>1 (0.24)</td><td>0 (0)</td><td>3 (0.73)</td></tr><tr><td>Distill</td><td>211</td><td>0.17</td><td>700</td><td>18</td><td>0 (0)</td><td>8 (0.37)</td><td>1 (0.05)</td><td>9 (0.42)</td></tr><tr><td>Viols</td><td>200</td><td>0.55</td><td>700</td><td>6</td><td>0 (0)</td><td>3 (0.22)</td><td>0 (0)</td><td>3 (0.22)</td></tr><tr><td>Ganges</td><td>200</td><td>0.66</td><td>700</td><td>16</td><td>0 (0)</td><td>9 (0.19)</td><td>1 (0.02)</td><td>6 (0.13)</td></tr></tbody></table>

The presence of the two duplications across the transect was demonstrated through crosses involving field samples (see Materials and Methods). Resistant females (N<sub>females</sub>) were mated with SLABTC males and the total number of egg rafts collected recorded (N<sub>egg rafts</sub>). Among them, the presence of D₂ and D₃ was detected using molecular tests. Number of individuals of each class is indicated with, in brackets, an estimate of its frequency (in addition to the very small number of egg rafts collected, a single female may have produced several egg rafts, preventing computation of accurate frequency). [D₂] represents (D₂/D₂), (D₂/S) or (D₂/R), and [D₃] represents (D₃/D₃), (D₃/S) or (D₃/R). doi:10.1371/journal.pgen.0030205.t003

<figure><img src="image_7.png" /></figure>
<table><caption>**Table 4.** Relative Fitness Estimates Fitted Using Summer Clinal Patterns Observed over the Period 1999–2002</caption><thead><tr><th>Genotype</th><th>Relative Fitness Treated Area</th><th>Support Limits</th><th>Relative Fitness Nontreated Area</th><th>Support Limits</th></tr></thead><tbody><tr><td>(S/S)</td><td>0.68</td><td>(0.58 – 0.72)</td><td>1</td><td>Reference</td></tr><tr><td>(S/R)</td><td>0.68</td><td>(0.66 – 0.83)</td><td>1</td><td>(0.99 – 1)</td></tr><tr><td>(R/R)</td><td>0.91</td><td>(0.80 – 1.02)</td><td>0.63</td><td>(0.63 – 0.72)</td></tr><tr><td>(D<sub>2</sub>/S) or (D<sub>3</sub>/S)</td><td>0.68</td><td>(0.60 – 0.80)</td><td>1</td><td>(0.99 – 1)</td></tr><tr><td>(D<sub>2</sub>/R) or (D<sub>3</sub>/R)</td><td>1</td><td>Reference</td><td>0.94</td><td>(0.86 – 0.98)</td></tr><tr><td>(D<sub>2</sub>/D<sub>2</sub>) or (D<sub>3</sub>/D<sub>3</sub>)</td><td>0.29</td><td>(0.26 – 0.69)</td><td>0.43</td><td>(0.37 – 0.53)</td></tr><tr><td>(D<sub>2</sub>/D<sub>3</sub>)</td><td>0.68</td><td>(0.39 – 1.46)</td><td>1</td><td>(0.65 – 1)</td></tr></tbody></table>

See Materials and Methods for details. “Reference” indicates which genotype is taken as the reference to compute relative fitnesses.
doi:10.1371/journal.pgen.0030205.t004

with deleterious side effects spreads and is followed by a steady “amelioration” process correcting for these side effects.

About 10 y after the beginning of OP treatments (1977), the major resistance allele *ace-I*<sup>R</sup>, which is beneficial in the treated area but has strong deleterious pleiotropic effects in absence of insecticide, appeared and spread [17,18]. Because of this fitness cost and incoming gene flow from the nontreated area favored by the absence of insecticide treatment in winter, it did not fix but remained polymorphic with a clinal pattern across treated and nontreated areas [26,27,35]. In the early nineties, at least one and probably two duplications involving a resistance and susceptible *ace-1* copy, started to spread and replace *ace-I*<sup>R</sup> [29,30]. Our lab experiments indicate that these two duplications are both severely deleterious when homozygous, but that heterozygotes involving either *ace-I*<sup>S</sup> or *ace-I*<sup>R</sup> alleles do not exhibit deleterious side effects. After ~1999, the pooled frequency of the duplications does not vary significantly, suggesting that the

four “alleles” (*ace-I*<sup>S</sup>, *ace-I*<sup>R</sup>, *ace-I*<sup>D2</sup>, and *ace-I*<sup>D3</sup>) reached a stable equilibrium. The stability of the duplicated haplotype frequencies since 1996 suggests that both duplications occurred and spread before 1996. This latter conclusion is tentative, since our data may not be accurate enough to detect the small perturbation of the frequency equilibrium caused by the spread of one of the duplications after 1996, due to the indirect estimation of duplicated haplotype frequency.

Clearly, however, considerable polymorphism in this system is maintained by both overdominance and migration. What is the mechanism by which overdominance operates? First, a duplicated haplotype restores AChE1 activity while maintaining resistance [30]. It therefore combines resistance with no deleterious pleiotropic effect caused by a deficit of AChE1 activity. An excess of AChE1 activity may be deleterious as well, but the one caused by the duplication is mild. For this reason, duplications certainly exhibit at least marginal overdominance over treated and nontreated areas.

<figure><img src="image_8.png" /><figcaption>**Figure 3.** Frequency of Homozygotes (D/D) Compared to Heterozygotes (D/S) during Development</figcaption></figure>

The frequency of (D/D) individuals is indicated for second instar and emerging adults on the left. On the right, adults are divided into early and late emerging adults. (i) MAURIN-D ((D<sub>2</sub>/D<sub>2</sub>) versus (D<sub>2</sub>/S), green), (ii) BIFACE-D ((D<sub>3</sub>/D<sub>3</sub>) versus (D<sub>3</sub>/S), orange), and (iii) BIFACE-DFix x MAURIN-D ((D<sub>2</sub>/D<sub>3</sub>) versus (D<sub>3</sub>/S), red). The bars indicate one error-type. Significant frequency differences between the different stages are indicated (n.s., *p*-value > 0.05; *, *p*-value < 0.05; **, *p*-value < 0.01; and ***, *p*-value < 0.001). doi:10.1371/journal.pgen.0030205.g003

<figure><img src="image_9.png" /></figure>
<table><caption>**Table 5.** Statistics for Larval Mortality and Development Time</caption><thead><tr><th rowspan="2"></th><th rowspan="2">Genotypes</th><th colspan="3">F (p-Value)</th></tr><tr><th>Time</th><th>Replicate</th><th>Time.Replicate</th></tr></thead><tbody><tr><td rowspan="4">Larval mortality</td><td>(D<sub>2</sub>/D<sub>2</sub>) versus (D<sub>2</sub>/S)</td><td><strong>13.84</strong> (&lt;0.001)</td><td>1.47 (0.23)</td><td>0.4 (0.82)</td></tr><tr><td>(D<sub>3</sub>/D<sub>3</sub>) versus (D<sub>3</sub>/S)</td><td><strong>14.5</strong> (&lt;0.001)</td><td>0.4 (0.67)</td><td>0.56 (0.57)</td></tr><tr><td>(D<sub>2</sub>/D<sub>3</sub>) versus (D<sub>2</sub>/S)</td><td>0.03 (0.86)</td><td>1.80 (0.17)</td><td>2.01 (0.13)</td></tr><tr><td>(D<sub>3</sub>/R) versus (R/R)</td><td>0.26 (0.61)</td><td>0.13 (0.87)</td><td>2.26 (0.10)</td></tr><tr><td rowspan="4">Development time</td><td>(D<sub>3</sub>/S) versus (R/R)</td><td><strong>15.05</strong> (&lt;0.001)</td><td>0.02 (0.98)</td><td>0.42 (0.65)</td></tr><tr><td>(D<sub>2</sub>/D<sub>2</sub>) versus (D<sub>2</sub>/S)</td><td><strong>3.91</strong> (&lt;0.05)</td><td>1.28 (0.26)</td><td>0.03 (0.87)</td></tr><tr><td>(D<sub>3</sub>/D<sub>3</sub>) versus (D<sub>3</sub>/S)</td><td><strong>4.58</strong> (&lt;0.05)</td><td>0.48 (0.62)</td><td>0.59 (0.56)</td></tr><tr><td>(D<sub>2</sub>/D<sub>3</sub>) versus (D<sub>2</sub>/S)</td><td>1.42 (0.23)</td><td>2.44 (0.09)</td><td>1.43 (0.24)</td></tr><tr><td rowspan="2"></td><td>(D<sub>3</sub>/R) versus (R/R)</td><td><strong>52.78</strong> (&lt;0.001)</td><td><strong>3.38</strong> (0.034)</td><td>2.90 (0.06)</td></tr><tr><td>(D<sub>3</sub>/S) versus (R/R)</td><td><strong>49.31</strong> (&lt;0.001)</td><td>0.53 (0.59)</td><td>1.36 (0.26)</td></tr></tbody></table>

Differential larval mortality is measured by following the frequency of two competing genotypes (Genotypes) between young larvae and emerging adults stages. Difference in development time (Time) is measured by the frequency change between early and late emerging adults. Difference among replicates is also tested (Replicate). The F-test statistics and the corresponding p-value of each variable and their interaction are indicated and bolded when significant.
doi:10.1371/journal.pgen.0030205.t005

Duplications could be described as “generalist” haplotypes that perform well in both habitats, compared to the “specialist” alleles *ace-I*<sup>S</sup> and *ace-I*<sup>R</sup>. Second, through lab experiments, we found, surprisingly, that both *ace-I*<sup>D2</sup> and *ace-I*<sup>D3</sup> duplications cause very strong deleterious effects when homozygous. However, these effects disappear in *ace-I*<sup>D2</sup> *lace-I*<sup>D3</sup> heterozygotes. The simplest explanation for this observation is that each duplication occurred independently rather than being generated by recombination (which may not be surprising, given the high duplication rate at this locus [30]), and that they carry distinct recessive sublethal mutations not necessarily related to *ace-I*-mediated resistance. Distinct unequal crossing over could have generated each duplication and in the same time disrupted different genes close to *ace-1*. Alternatively, different sublethal recessive mutations could have hitchhiked with the initial spread of each duplication. In all cases, these explanations require that the duplications cannot get rid of these deleterious mutations by recombining. It therefore suggests that recombination is very low around or within the duplicated haplotypes or that they are situated on chromosomal inversions (and thus behave like a “balancer region” in analogy with the balancer chromosomes used in laboratory *Drosophila*). The fact that we never observed in laboratory crosses a recombinant separating the two *ace-1* duplicated

copies [30] is consistent with this last hypothesis. We consider this explanation plausible since duplication events involve an important modification of the genome, which can easily disrupt other genes or regulatory regions [31,33,36].

This example of adaptation involves three successful steps (the formation of *ace-I*<sup>R</sup>, *ace-I*<sup>D2</sup> and *ace-I*<sup>D3</sup>), each of them being driven by natural selection. However, each of these three steps presents severe deleterious pleiotropic effects. The deleterious pleiotropic effect of *ace-I*<sup>R</sup> may be unavoidable if changing the AChE1 active site to decrease affinity towards OP insecticide necessarily also leads to a lower affinity towards ACh. However, the occurrence and spread of duplications that are sublethal when homozygotic is more perplexing. For instance, another *ace-1* duplication, *ace-I*<sup>D1</sup> seems to have spread and be almost fixed on the island of Martinique [37], indicating that *ace-1* duplications do not necessarily involve strong deleterious pleiotropic effects when homozygous (confirmed by laboratory analyses, P. Labbé and M. Weill, unpublished data). However, once a duplication with a recessive cost has spread, even partially, in a population, selection is likely to be less effective at replacing it with a better one because the fate of a beneficial mutation is mostly determined when it is rare and therefore

<table><caption>**Table 6.** Cross Results for Fertility</caption><thead><tr><th rowspan="2">Pairs</th><th colspan="3">D<sub>3</sub></th><th colspan="3">D<sub>2</sub></th></tr><tr><th>N</th><th>N<sub>E</sub></th><th>N<sub>H</sub></th><th>N</th><th>N<sub>E</sub></th><th>N<sub>H</sub></th></tr></thead><tbody><tr><td>♂ (D/S) × ♀ (D/S)</td><td>16</td><td>11</td><td>9</td><td>54</td><td>45</td><td>23</td></tr><tr><td>♂ (D/D) × ♀ (D/S)</td><td>24</td><td>15</td><td>10</td><td>40</td><td>30</td><td>12</td></tr><tr><td>♂ (D/S) × ♀ (D/D)</td><td>18</td><td>16</td><td>12</td><td>44</td><td>35</td><td>11</td></tr><tr><td>♂ (D/D) × ♀ (D/D)</td><td>24</td><td>13</td><td>1</td><td>24</td><td>17</td><td>1</td></tr></tbody></table>

For *ace-I*<sup>D3</sup> (D<sub>3</sub>) and *ace-I*<sup>D2</sup> (D<sub>2</sub>), the total number of each type of cross (N) is indicated with the number of crosses where the female laid eggs (N<sub>E</sub>) and with the number of crosses where the egg rafts hatched (N<sub>H</sub>). Males and females are either heterozygotes (D/S) or homozygotes (D/D).
doi:10.1371/journal.pgen.0030205.t006

<table><caption>**Table 7.** Statistics for Fertility Analysis</caption><thead><tr><th rowspan="2">Trait</th><th rowspan="2">Haplotype</th><th colspan="3">F (p-Value)</th></tr><tr><th>Male</th><th>Female</th><th>Male.Female</th></tr></thead><tbody><tr><td rowspan="2">Ne</td><td>D<sub>3</sub></td><td><strong>4.11</strong> (&lt;0.05)</td><td>0.12 (0.73)</td><td>2.37 (0.12)</td></tr><tr><td>D<sub>2</sub></td><td>1.62 (0.2)</td><td>0.36 (0.55)</td><td>0.002 (0.96)</td></tr><tr><td rowspan="2">Nh</td><td>D<sub>3</sub></td><td><strong>8.04</strong> (&lt;0.01)</td><td><strong>12.1</strong> (&lt;0.001)</td><td>3.51 (0.061)</td></tr><tr><td>D<sub>2</sub></td><td><strong>3.82</strong> (0.05)</td><td><strong>8.53</strong> (&lt;0.01)</td><td>2.08 (0.15)</td></tr></tbody></table>

The effects of male and female genotypes (Male and Female variables, respectively) on the proportion of female laying eggs (N<sub>E</sub>) and on the proportion of hatching egg rafts (N<sub>H</sub>) (see Table 6) are analyzed. Males and females are either heterozygotes (D<sub>3</sub>/S) or (D<sub>2</sub>/S) or homozygotes (D<sub>3</sub>/D<sub>3</sub>) or (D<sub>2</sub>/D<sub>2</sub>). The F-test statistics and the corresponding p-value of each variable and their interaction are indicated and bolded when significant (and put in italics when close to significance).
doi:10.1371/journal.pgen.0030205.t007

<figure><img src="image_10.png" /></figure>
<figure><img src="image_11.png" /><figcaption>**Figure 4.** Frequency of Heterozygotes (D/S) and (D/R) Compared to Homozygotes (R/R) during Development</figcaption></figure>

The frequency of (D/S) and (D/R) individuals is indicated on the left for two developmental stages, first instar and emerging adults. On the right, emerging adults are divided into early and late emerging adults. (i) (D₃/S) versus (R/R) (green), (ii) (D₃/R) versus (R/R) (orange). The bars indicate one error-type. Significant frequency differences between the different stages are indicated (n.s., *p*-value > 0.05; *, *p*-value < 0.05; **, *p*-value < 0.01; and ***, *p*-value < 0.001).

doi:10.1371/journal.pgen.0030205.g004

by its heterozygous effect [38]. Since *ace-1*<sup>D2</sup> or *ace-1*<sup>D3</sup> duplications enjoy almost no fitness cost when heterozygotic, any new duplication has a low chance to spread in these populations (although migration bringing fitter haplotypes could help escaping this apparent dead end). Thus, selection has taken the mosquito populations on a difficult path indeed.

More generally, we may wonder if such a tortuous adaptive trajectory is frequent in nature. Clearly, gene duplications may solve many genetic trade-offs and chromosome rearrangements such as inversions may strongly perturb the genetic architecture. This type of situation may be more common with strong selective pressure and strong pleiotropy, whereas less intense selection may select for more subtle variation (e.g., see [39]). Nevertheless, insecticide resistance may not be an exception, since other selective pressures appear to be intense as well (e.g., parasitism [40,41]) and thus could favor similar complex genetic responses. The new molecular tools available will allow deeper investigation of adaptation genetics evolution, and thus will help to settle the issue of the frequency of the kind of complex patterns uncovered by our study.

In the long term, selection may produce exquisite adaptations, but this study lays bare that even impressive adaptations are likely to have begun with a process of trial and error that seems to be anything but optimal. It appears that natural selection is forced to tinker with available variability, despite the costs, rather than build impressive and cost-free adaptations that are wholly novel.

## Methods

**Data collection.** Data used in cline analysis were collected in the Montpellier area since 1986 along a transect across the treated and the nontreated areas (Figure 2). Published data from the summers of 1986, 1991, 1995, and 1996, spring of 1993, and winters of 1995 and 1996 [26,29,35] were used to perform the overall analysis. This was complemented with unpublished samples: summer of 1999, 2001, and 2002, spring of 1996 and 2000, and winter of 1999 and 2002 (complete dataset in Table S1). For illustration, we also indicate the frequency observed in a single population in 1984 near the coast, but this population was not included in the cline analysis. Five populations were also sampled on the same transect in 2005 (Maurin, Lattes, Distill, Viols, and Ganges) to assess the distribution of the two duplicated haplotypes in the Montpellier area.

**Mosquito strains.** In laboratory experiments, different strains were used to compare life history traits. Two reference strains were used: SLAB, the susceptible reference strain (homozygote for the allele $S_{SLAB}$ [34]), and SR, homozygote for *ace-1*<sup>R</sup> (R), but with the same genetic background as SLAB [23]. Two other strains, named MAURIN-D and BIFACE-D and respectively harboring *ace-1*<sup>D2</sup> ($D_2$) and *ace-1*<sup>D3</sup> ($D_3$), were also used. The duplicated strains originate from MAURIN and BIFACE strains [30] backcrossed for more than 15 generations with SLAB (method in [23]). These strains are not homozygous but contain three genotypes: ($S_{SLAB}/S_{SLAB}$), ($D/S_{SLAB}$), and (D/D) ($D_2$ and $D_3$ for MAURIN-D and BIFACE-D, respectively).

Egg rafts resulting from the cross of homozygous males and females originating from BIFACE-D were used to constitute a strain homozygous for each duplicated haplotype, BIFACE-DFix. SLAB-TC strain (SLAB strain cured from *Wolbachia* bacteria [25]) was used for crosses with field samples to avoid incompatibility phenomena.

**Cline analysis.** *Identification of ace-1 phenotype.* For each mosquito, the head was used to establish the phenotype at the *ace-1* locus, using the TPP test [42], based on enzymatic activity of AChE1 in the presence or absence of insecticide. Single copy allele homozygotes (S/S) and (R/R) are easily detectable using this test. However, the heterozygotes (R/S) and the genotypes involving a duplicated haplotype (heterozygotes: (D/S) and (D/R) and homozygotes: (D/D)) can not be distinguished

<figure><img src="image_12.png" /></figure>
(they all display a [RS] phenotype). Moreover, this test does not allow distinction between the two duplicated haplotypes D₂ and D₃. The pooled frequency of D₂ and D₃ can be estimated from the apparent excess of [RS] phenotypes caused by the presence of the duplications [26]. This method assumes Hardy-Weinberg proportions for the different genotypes and is therefore not as accurate as direct observation of the different genotypes. However, it has been shown to correctly estimate D frequency in field samples where duplication frequency was independently estimated using crosses [26]. Note that no molecular test is currently available to directly detect the duplicated haplotypes in natural populations, as their D(R) and D(S) copies are not different from single copy alleles (R and S) present in the same sampling sites.

*Clines description.* We first analyzed the spatial variation in allele frequency across the treated and nontreated areas for each sample independently (i.e., for a given year and season) using a purely descriptive model. We estimated the pooled frequency of both duplications. More specifically, we assumed that the frequency (denoted *p*) of each resistance allele (indicated by the subscript *i* = R or D) at time *j* followed a scaled negative exponential as follows:

$$p_{ij} = h_i \exp[-(a_{ij}x^2 + b_{ij}x)] \quad (1)$$

where *x* is the distance from the coast and *h*<sub>ij</sub>, *b*<sub>ij</sub> and *a*<sub>ij</sub> are the estimated parameters. *h*<sub>ij</sub> measures the frequency of resistance allele *i* on the coast (i.e., at *x*=0) at time *j*. The parameters *b*<sub>ij</sub> and *a*<sub>ij</sub> describe rates of decline of the frequency of allele *i* (at time *j*) with distance and with the square of distance from the coast, respectively. We allow for a flexible clinal shape because it tends to vary with season [35].

The second step in the analysis was to compare frequency patterns across years and seasons. For this purpose, we fitted all samples simultaneously to measure the variation of clines trough time using the method developed in Labbé et al. [16]. More specifically, we assumed that *h*<sub>ij</sub> values changed smoothly as a logistic function of time (measured in months)

$$h_{ij} = \frac{\text{Exp}(\alpha_i t_{1j} + \beta_i t_{2j} + \gamma_i)}{1 + \text{Exp}(\alpha_i t_{1j} + \beta_i t_{2j} + \gamma_i)}, \quad (2)$$

where *t*<sub>1j</sub> is the number of months after January 1986 when the date of sampling is before January 1986 + *t*<sup>*</sup> and is *t*<sup>*</sup> otherwise. *t*<sub>2j</sub> is the number of months after January 1986 + *t*<sup>*</sup>. *α*<sub>i</sub>, *β*<sub>i</sub>, *γ*<sub>i</sub> and *t*<sup>*</sup> are estimated parameters. The overall change in frequency over the 1986–2002 period is measured for each resistance allele *i* by *α*<sub>i</sub> and *β*<sub>i</sub>, which measure the rate of frequency change between 1986 and 1986 + *t*<sup>*</sup> and between 1986 + *t*<sup>*</sup> and 2002, respectively. *t*<sup>*</sup> was introduced to allow for changes in the rate of allele replacement. Parameter *γ*<sub>i</sub> is related to the initial frequency *h*<sub>i0</sub> of each allele *i* (*h*<sub>i0</sub> = Exp(γ<sub>i</sub>))[1 + Exp(γ<sub>i</sub>)]]. Expected phenotypic distributions were computed using allelic frequency and assuming Hardy-Weinberg proportions in each location (see [16]). The phenotype was considered to be a three-state random variable ([RR], [RS], and [SS]). The log-likelihood of a sample was computed from the phenotypic multinomial distribution and maximized using the Metropolis algorithm (see [26,27,35]). Models were compared using *F*-tests in order to correct for overdispersion. Deviance was also corrected for overdispersion to find the support limits of each parameter [43,44].

*Fitness estimation of the different genotypes.* In order to determine the fitness of the different genotypes that would yield the stable clines observed over the period 1999–2002, we also used a maximum-likelihood approach. For a given set of fitness values, we obtained by simulation the distribution of genotypes at equilibrium at different distances from the coast. From this distribution, we computed the expected frequency of [RR], [RS], and [SS] for each population in our dataset. The likelihood was then computed and maximized in the same way as with the descriptive models above. The simulation was performed using a stepping stone across treated and nontreated areas using the dispersion kernel that has been previously estimated for *C. pipiens* as described in Lenormand et al. [26] but with a treated area of 16 km that reflects the treatment practices over this period (P. Labbé, unpublished data). Because we do not have estimates for the relative frequency of the two duplications, we could not estimate their relative fitness. We therefore report estimates assuming that the two duplications have the same fitness effects. The fitness of a given genotype was modelled with two components: a fitness cost *c* (expressed in both treated and non treated areas) and a reduction in survival due to the presence of insecticide in the treated area *s*. The full parameter range was constrained to reflect that *s* should be lower (or equal) for genotypes with an increasing number *R* copies.

At first it may be surprising that polymorphism with more than two alleles could be maintained at migration—selection equilibrium

in a situation with only two habitats (treated and nontreated areas). First, even with haploid selection, the number of alleles that can be maintained is larger than the number of habitats if dispersal is localized (as confirmed by our simulations, see also [45]). Second, particular overdominance relationships among alleles (such as the ones we find) can also increase the number of alleles that can be maintained at equilibrium.

*Distribution of the two duplicated haplotypes.* Crosses. Larvae were sampled in five populations along the sampling transect (Maurin, Lattes, Distill, Viols, and Ganges; Figure 2) and reared in the laboratory. They were exposed to a dose of insecticide that kills all [SS] individuals (25 × 10⁻⁶ M of Propoxur). About 700 resistant females from each population were crossed with about 800 males of SLAB-TC strain (S<sub>SLAB</sub>/S<sub>SLAB</sub>). They were repeatedly blood fed each week until they died. Egg rafts were collected every day, isolated, and reared to the third instar. They were then exposed to a dose of insecticide that kills all [SS] individuals in order to discard all susceptible field alleles resulting from (D/S) or (R/S) mothers. DNA was extracted from a pool of ~20 survivors of each egg raft. Molecular tests described below were used to detect the presence of each resistance haplotype (R, D₂, and D₃). Four types could be identified: genotypes (D₂/D₃) and (R/R) and phenotypes [D₂] (i.e., (D₂/D₃), (D₂/S), or (D₂/R)), and [D₃] (i.e., (D₃/D₃), (D₃/S), or (D₃/R)).

**Molecular tests.**

Using partial *ace-1* sequences of each haplotype [30], we designed RFLP tests to discriminate D₂, D₃, R, and S<sub>SLAB</sub>. PCR amplification of a 458 bp fragment of exon 3 using the primers CxEx3dir 5'-CGA CTC GGA CCC ACT CGT-3' and CpEx3rev 5'-GAC TTG CGA CAC GGT ACT GCA-3' was performed (30 cycles, 93 °C for 30s, 55 °C for 30s, and 72 °C for 1min). The amplified fragment was then digested in parallel by different restriction enzymes. First, the fragment was cut twice by the enzyme BsrBI only when S<sub>SLAB</sub> is present, generating three fragments (127 bp, 141 bp, and 190 bp; Figures S1 and S2), all the other alleles being cut only once (two fragments of 127 bp and 331 bp; Figures S1 and S2). Second, the fragment is cut by the enzyme EagI only when D₂(S) is present, generating two fragments (150 bp and 308 bp; Figures S1 and S2). Third, the fragment is cut by the enzyme Hinfl only when D₃(S) is present, generating two fragments (102 bp and 354 bp; Figures S1 and S2; note that there is a Hinfl site in the primer CxEx3dir, subtracting 2 bp in each fragment; Figure S1). Discriminating between the resistance and susceptible copies is possible using the test provided by Weill et al. [20]: the G119S mutation providing resistance creates a site for the enzyme AluI.

*Laboratory experiments.* Larval mortality: (D/D) versus (D/S). Trials between (D/D) and (D/S) individuals were performed in triplicate. Larvae of different genotypes were reared in competition under the same environmental conditions (food, temperature, etc.). They were selected at the first instar stage using Propoxur at a concentration of 25 × 10⁻⁶ M, which eliminates only (S/S) individuals. Adults were collected during the first and the second wk after the first adult emergence (~30 individuals each wk). They will be respectively referred as early (first wk) and late (second wk) emerging adults. Genotype frequency was measured at second larvae instar and both adulthood stages. Three trials were conducted: (i) (D₃/S) versus (D₃/D₃), (ii) (D₂/S) versus (D₂/D₂), and (iii) (D₃/S) versus (D₂/D₃). The different genotypes were obtained respectively from (i) a cross between males and females from BIFACE-D (progeny genotypes (D₃/D₃), (D₃/S<sub>SLAB</sub>), and (S<sub>SLAB</sub>/S<sub>SLAB</sub>)), (ii) a cross between males and females from MAURIN-D (progeny genotypes (D₂/D₂), (D₂/S<sub>SLAB</sub>), and (S<sub>SLAB</sub>/S<sub>SLAB</sub>)) and (iii) a cross between males from BIFACE-DFix and females from MAURIN-D and the reverse cross (progeny genotypes (D₂/D₃) and (D₃/S<sub>SLAB</sub>)). Each sample was analyzed using the BsrBI-based RFLP test to determine the proportion of individuals of genotype (D₃/S<sub>SLAB</sub>) (in the first and third trials) or (D₂/S<sub>SLAB</sub>) (in the second trial).

*Larval mortality: (D/S) or (D/R) versus (R/R).* Trials were conducted between (D/R) or (D/S) individuals and individuals homozygote for the single resistance allele (R/R). Trials were performed in triplicates with 500 first instar larvae of each genotype reared under the same environmental conditions (food, temperature, etc.). In each replicate, early and late emerging adults were collected as indicated above. Two trials were carried out: (i) (D₃/R) versus (R/R) and (ii) (D₃/S) versus (R/R). The different genotypes were obtained from (i) a cross between females from BIFACE-DFix and males from strain SR to obtain (D₃/R) individuals, and (ii) a cross between females from BIFACE-DFix and males from SLAB to obtain (D₃/S<sub>SLAB</sub>) individuals. (R/R) individuals were directly obtained from strain SR. Each sample was analyzed using the AluI-based RFLP test to determine the proportion of individuals of genotype (D₃/R) (first trial) or (D₃/S<sub>SLAB</sub>) (second trial).

*Fertility.* Different crosses were realized in order to determine the

<figure><img src="image_13.png" /></figure>
fertility of individuals carrying a duplicated allele at the homozygous or heterozygous state (D(D)/D or (D/S)). In each case, the proportion of females laying eggs and the proportion of hatching eggs rafts was recorded. For each series of cross, about 50 males were mated independently with five females each, all from the same strain. Five days after mating, the genotype of the males was determined using the BsrBI-based RFLP test. Females were then grouped according to the male genotype, blood fed, and kept without access to laying substrate. Six days later, they were allowed to lay eggs individually. The females were genotyped either after they laid eggs or less than 6 h after their death. Two series of crosses were realized (i) between males and females from BIFACE-D (genotype (D<sub>3</sub>/D<sub>3</sub>) or (D<sub>3</sub>⁄S<sub>SLAB</sub>) after selection), and (ii) between males and females from MAURIN-D (genotype (D<sub>2</sub>/D<sub>2</sub>) or (D<sub>2</sub>⁄S<sub>SLAB</sub>) after selection), to assess the fertility of individuals carrying *ace*-I<sup>D3</sup> and *ace*-I<sup>D2</sup>, respectively.

**Statistical analysis.** For the larval mortality analysis, the following generalized linear model (GLM) was fitted, with binomial error: DD = Time + Replicate + Time.Replicate, where DD represents the proportion of the (D/D) genotype in the population (i.e., (D<sub>3</sub>/D<sub>3</sub>), (D<sub>2</sub>/D<sub>2</sub>), and (D<sub>2</sub>/D<sub>3</sub>) for the first, second, and third crosses, respectively), Time is a factor indicating when the sample was taken, Replicate is a factor indicating the three containers in which the experiment was replicated, and Time.Replicate is the interaction between the two factors. Two analyses were performed: (i) one to test for a difference in mortality between (D/D) and (D/S) individuals (in that case, Time was 2nd instar or emerging adult), and (ii) the second to test for a difference in development time between (D/D) and (D/S) individuals (in that case, Time was early emerging or late emerging adults). A similar model was used to analyze the proportion of (D<sub>3</sub>/S) and (D<sub>3</sub>/R) in trials versus (R/R). These models were simplified according to Crawley [46]: significance of the different terms was tested starting from the higher-order terms using *F*-test. Non-significant terms (*p* > 0.05) were removed. Factor levels of qualitative variables that were not different in their estimates (using *F*-test) were grouped as described by Crawley [46]. This process yielded a minimal adequate model.

The fertility of males and females of different genotypes was analyzed by comparing the proportion of females laying eggs and the proportion of hatching egg rafts among the different types of cross. The number of females laying eggs (Ne) and the number of hatching egg rafts (Nh) were analyzed using GLM with binomial error: Male + Female + Male.Female, where Male and Female are factors indicating male (or female) genotype. These models were simplified as above. All analyses were performed using R software (v 2.0.1., http://www.r-project.org).

**Supporting Information**

**Figure S1. Susceptible and Resistance Copy Sequences**

The sequence of R, S<sub>SLAB</sub>, D<sub>2</sub>(S), and D<sub>3</sub>(S) are presented. The restriction enzyme sites specific to each copy are also indicated (AluI

**References**

1. Fisher RA (1930) The genetical theory of natural selection. Oxford: Oxford University Press.

2. Hartl D, Taubes C (1996) Compensatory nearly neutral mutations: selection without adaptation. J Theor Biol 182: 303–309.

3. Duboule D, Wilkins AS (1998) The evolution of ‘bricolage’. Trends Genet 14: 54–59.

4. Jacob F (1977) Evolution and Tinkering. Science 196: 1161–1166.

5. Caspari E (1952) Pleiotropic gene action. Evolution 6: 1–18.

6. Lande R (1983) The response to selection on major and minor mutations affecting a metrical trait. Heredity 50: 47–65.

7. Wright S (1969) Evolution and the genetics of populations. Chicago: University of Chicago Press.

8. Carrière Y, Deland J-P, Roff DA, Vincent C (1994) Life-history cost associated with the evolution of insecticide resistance. Proc R Soc Lond B Biol Sci 258: 35–40.

9. Otto SP (2004) Two steps forward, one step back: the pleiotropic effects of favoured alleles. Proc R Soc Lond B Biol Sci 271: 705–714.

10. Cohan FM, King EC, Zawadzki P (1994) Amelioration of the deleterious pleiotropic effects of an adaptive mutation in *Bacillus subtilis*. Evolution 48: 81–95.

11. Fisher RA (1928) The possible modification of the response of the wild type to recurrent mutations. Am Nat 62: 115–126.

12. Haldane JBS (1932) The Causes of Evolution. New York: Harper.

13. Lenski RE (1988) Experimental studies of pleiotropy and epistasis in *Escherichia coli*. *II*. Compensation for maladaptive effects associated with resistance to virus T4. Evolution 42: 433–440.

[blue], BsrBI [yellow], Hinfl [green], and EagI [red] for R, S<sub>SLAB</sub>, D<sub>2</sub>(S), and D<sub>3</sub>(S), respectively; see Materials and Methods and Figure S2].

Found at doi:10.1371/journal.pgen.0030205.sg001 (34 KB DOC).

**Figure S2. RFLP Profiles**

The RFLP profiles allowing the identification of the different genotypes are presented for homozygotes S<sub>SLAB</sub>, R, D<sub>2</sub>, and D<sub>3</sub>. For each genotype, three restriction enzymes are used (BsrBI, Hinfl, and EagI; Figure S1) on the CxEx3dir-CpEx3rev PCR product (see Methods) and an undigested control is also presented. Note that for D<sub>2</sub> and D<sub>3</sub> homozygotes, both the susceptible and the resistant copies are present in the PCR product, thus leading to an apparent heterozygote (e.g., for the (D<sub>2</sub>/D<sub>2</sub>) genotype only the D<sub>2</sub>(S) copy is digested by EagI, not the D<sub>2</sub>(R) copy).

Found at doi:10.1371/journal.pgen.0030205.sg002 (285 KB DOC).

**Table S1. Sampling Data**

For each year, the numbers of individuals carrying each phenotype are indicated: (A) in spring, (B) in summer, and (C) in winter. For nomenclature, [SS] and [RR] phenotypes stand for (S/S) and (R/R) genotypes, whereas [RS] phenotype covers all (R/S), (D/S), (D/R), and (D/D) genotypes. (–) indicates that the site was not sampled in the corresponding year. See Figure 2 for details of sampling location.

Found at doi:10.1371/journal.pgen.0030205.st001 (47 KB DOC).

**Accession Numbers**

The National Center for Biotechnology Information (NCBI) GenBank database (http://www.ncbi.nlm.nih.gov/sites/entrez?db=Nucleotide) accession numbers for *ace*-I are AJ489456 and AJ515147.

**Acknowledgments**

We are very grateful to T. Vines, T. Little, N. Johnson, and two anonymous reviewers for helpful comments on the manuscript. M. Marquine for technical assistance, and V. Durand for the literature search.

**Author contributions.** P. Labbé, M. Weill, and T. Lenormand conceived and designed the experiments, contributed reagents/materials/analysis tools, and wrote the paper. P. Labbé, C. Berticat, A. Berthomieu, S. Unal, C. Bernard, and M. Weill performed the experiments. P. Labbé, C. Berticat, M. Weill, and T. Lenormand analyzed the data.

Funding. This work was financed by the ANR Morevol Sante-Environnement (Ministère délégué à la Recherche). Contribution 2007.107 of the Institut des Sciences de l'Evolution de Montpellier (UMR CNRS-UM2 5554).

**Competing interests.** The authors have declared that no competing interests exist.

14. McKenzie JA (1996) Ecological and evolutionary aspects of insecticide resistance. Austin (Texas): Academic Press.

15. Guillemaud T, Lenormand T, Bourguet D, Chevillon C, Pasteur N, et al. (1998) Evolution of resistance in *Culex pipiens*: allele replacement and changing environment. Evolution 52: 443–453.

16. Labbé P, Lenormand T, Raymond M (2005) On the worldwide spread of an insecticide resistance gene: a role for local selection. J Evol Biol 18: 1471–1484.

17. Raymond M, Berticat C, Weill M, Pasteur N, Chevillon C (2001) Insecticide resistance in the mosquito *Culex pipiens*: what have we learned about adaptation? *Genetica*: 1–10.

18. Weill M, Labbé P, Duron O, Pasteur N, Fort P, et al. (2005) Insecticide resistance in the mosquito *Culex pipiens*: towards an understanding of the evolution of *ace* genes. In: Fellows MDE, Holloway GJ, Rolff J, editors. Insect evolutionary ecology. Oxon (United Kingdom): CABI publishing.

19. Weill M, Lutfalla G, Mogensen K, Chandre F, Berthomieu A, et al. (2003) Insecticide resistance in mosquito vectors. Nature 425: 366–366.

20. Weill M, Malcolm C, Chandre F, Mogensen K, Berthomieu A, et al. (2004) The unique mutation in *ace*-I giving high insecticide resistance is easily detectable in mosquito vectors. Insect Mol Biol 13: 1–7.

21. Alout H, Berthomieu A, Cui F, Tan Y, Berticat C, et al. (2007) Different amino-acid substitutions confer insecticide resistance through acetylcholinesterase 1 insensitivity in *Culex vishnui* and *Culex tritaeniorhynchus* (Diptera: Culicidae) mosquitoes from China. J Med Entomol 44: 463–469.

22. Bourguet D, Lenormand T, Guillemaud T, Marcel V, Fournier D, et al. (1997) Variation of dominance of newly arisen adaptive genes. Genetics 147: 1225–1234.

<figure><img src="image_14.png" /></figure>
23. Berticat C, Boquien G, Raymond M, Chevillon C (2002) Insecticide resistance genes induce a mating competition cost in *Culex pipiens* mosquitoes. *Genet Res* 79: 41-47.

24. Bourguet D, Guillemaud T, Chevillon C, Raymond M (2004) Fitness costs of insecticide resistance in natural breeding sites of the mosquito *Culex pipiens*. *Evolution* 58: 128-135.

25. Duron O, Labbe P, Berticat C, Rousset F, Guillot S, et al. (2006) High *Wolbachia* density correlates with cost of infection for insecticide resistant *Culex pipiens* mosquitoes. *Evolution* 60: 303-314.

26. Lenormand T, Guillemaud T, Bourguet D, Raymond M (1998) Evaluating gene flow using selected markers: a case study. *Genetics* 149: 1383-1392.

27. Lenormand T, Bourguet D, Guillemaud T, Raymond M (1999) Tracking the evolution of insecticide resistance in the mosquito *Culex pipiens*. *Nature* 400: 861-864.

28. Bourguet D, Raymond M, Bisset J, Pasteur N, Arpagaus M (1996) Duplication of the *Ace.1* locus in *Culex pipiens* from the Caribbean. *Biochem Genet* 34: 351-362.

29. Lenormand T, Guillemaud T, Bourguet D, Raymond M (1998) Appearance and sweep of a gene duplication: adaptive response and potential for new functions in the mosquito *Culex pipiens*. *Evolution* 52: 1705-1712.

30. Labbé P, Berthomieu A, Berticat C, Alout H, Raymond M, et al. (2007) Independent duplications of the acetylcholinesterase gene conferring insecticide resistance in the mosquito *Culex pipiens* Mol Biol Evol 24: 1056-1067.

31. Ohno S (1970) Evolution by gene duplication. Berlin: Springer-Verlag.

32. Lynch M, Conery JS (2000) The evolutionary fate and consequences of duplicate genes. *Science* 290: 1151-1155.

33. Otto S, Yong P (2002) The evolution of gene duplicates. *Adv Genet* 46: 451-483.

34. Georgiou GP, Metcalf RL, Gidden FE (1966) Carbamate-resistance in mosquitoes: selection of *Culex pipiens fatigans* Wied (= *Culex quinquefasciatus*) for resistance to Baygon. *Bull World Health Organ* 35: 691-708.

35. Lenormand T, Raymond M (2000) Analysis of clines with variable selection and variable migration. *Am Nat* 155: 70-82.

36. Kirkpatrick M, Barton N (2006) Chromosome inversions, local adaptation and speciation. *Genetics* 173: 419-434.

37. Yebakima A, Marquine M, Rosine J, Yp-Tcha MM, Pasteur N (2004) Evolution of resistance under insecticide selection pressure in *Culex pipiens quinquefasciatus* (Diptera, Culicidae) from martinique. *J Med Entomo* 41: 718-725.

38. Roze D, Rousset F (2003) Selection and drift in subdivided populations: a straightforward method for deriving diffusion approximations and applications involving dominance, selfing and local extinctions. *Genetics* 165: 2153-2166.

39. McKenzie JA, Parker AG, Yen JL (1992) Polygenic and single gene responses to selection for resistance to diazinon in *Lucilia cuprina*. *Genetics* 130: 613-620.

40. Little TJ, Watt K, Ebert D (2006) Parasite-host specificity: experimental studies on the basis of parasite adaptation. *Evolution* 60: 31-38.

41. Carton Y, Nappi AJ, Poirie M (2005) Genetics of anti-parasite resistance in invertebrates. *Dev Comp Immunol* 29: 9-32.

42. Bourguet D, Pasteur N, Bisset J, Raymond M (1996) Determination of Ace.1 genotypes in single mosquitoes: toward an ecumenical biochemical test. *Pestic Biochem Physiol* 55: 122-128.

43. Lebreton J-D, Burnham KP, Clobert J, Anderson DR (1992) Modeling survival and testing biological hypotheses using marked animals: a unified approach with case studies. *Ecol Monogr* 62: 67-118.

44. Anderson DR, Burnham KP, White GC (1994) AIC model selection in overdispersed capture-recapture data. *Ecology* 75: 1780-1793.

45. Nagylaki T, Lou Y (2001) Patterns of multiallelic polymorphism maintained by migration and selection. *Theor Popul Biol* 59: 297-313.

46. Crawley MJ (2003) Statistical computing. An introduction to data analysis using S-plus;JWaS Ltd, editor. Chichester (United Kingdom): Wiley.

<figure><img src="image_15.png" /></figure>