Use of Mutagenesis, Genetic Mapping and Next Generation Transcriptomics to Investigate Insecticide Resistance Mechanisms

Insecticide resistance is a worldwide problem with major impact on agriculture and human health. Understanding the underlying molecular mechanisms is crucial for the management of the phenomenon; however, this information often comes late with respect to the implementation of efficient counter-measures, particularly in the case of metabolism-based resistance mechanisms. We employed a genome-wide insertional mutagenesis screen to Drosophila melanogaster, using a Minos-based construct, and retrieved a line (MiT[w−]3R2) resistant to the neonicotinoid insecticide Imidacloprid. Biochemical and bioassay data indicated that resistance was due to increased P450 detoxification. Deep sequencing transcriptomic analysis revealed substantial over- and under-representation of 357 transcripts in the resistant line, including statistically significant changes in mixed function oxidases, peptidases and cuticular proteins. Three P450 genes (Cyp4p2, Cyp6a2 and Cyp6g1) located on the 2R chromosome, are highly up-regulated in mutant flies compared to susceptible Drosophila. One of them (Cyp6g1) has been already described as a major factor for Imidacloprid resistance, which validated the approach. Elevated expression of the Cyp4p2 was not previously documented in Drosophila lines resistant to neonicotinoids. In silico analysis using the Drosophila reference genome failed to detect transcription binding factors or microRNAs associated with the over-expressed Cyp genes. The resistant line did not contain a Minos insertion in its chromosomes, suggesting a hit-and-run event, i.e. an insertion of the transposable element, followed by an excision which caused the mutation. Genetic mapping placed the resistance locus to the right arm of the second chromosome, within a ∼1 Mb region, where the highly up-regulated Cyp6g1 gene is located. The nature of the unknown mutation that causes resistance is discussed on the basis of these results.


Introduction
Insecticide resistance is an increasing problem that compromises the control of insect pests of medical, veterinary and agricultural impact. An understanding of insecticide resistance mechanisms is essential for the subsequent development of tools and practices that can improve pest control interventions. During the last decades, extensive biochemical, genetic and molecular studies have been conducted to elucidate insecticide resistance mechanisms [1,2,3]. Knowledge of the mechanisms underlying target site resistance in major pests to some commonly used insecticides has been established to some extent [4,5,6]. The understanding of detoxification/metabolism-based insecticide resistance mechanisms has not kept similar pace, due to the complexity of the involved multi-gene systems and the lack of genome sequence data. However, in a few cases, the molecular basis of metabolism-based insecticide resistance mechanisms was identified. A single P450, CYP6P3, was over-expressed in pyrethroid resistant Anopheles gambiae mosquitoes, and it was capable of metabolizing pyrethroids [7]. Karunker et al. [8] showed that the B. tabaci cytochrome P450 BTCYP6CM1 is capable of metabolizing the neonicotinoid Imidacloprid, one of the most important insecticides worldwide, and to confer neonicotinoid resistance.
These studies have shed light on cases of metabolism-based insecticide resistance mechanisms. However, there is a number of issues which remain unsolved, such as the underlying molecular mechanisms that are responsible for over-expression of detoxification enzymes. In addition, the information on molecular changes responsible for resistance often comes too late, i.e. when resistance has been irreversibly established in pest populations and/or when the active ingredient has already been replaced by others. The use of modern molecular approaches and models for the early identification or even prediction of insecticide resistance mechanisms could improve the management of the phenomenon.
Drosophila melanogaster, although not a pest species, has been used extensively for insecticide resistance research [9,10,11]. Resistance associated with the over-expression of a single P450 gene (Cyp6g1) has been documented for field-derived Drosophila lines resistant to Imidacloprid and DDT [12]. Over-expression correlated with the presence of a single insertion of an Accord transposable element into the 59 end of the Cyp6g1 gene has also been reported [12]. A recent study of Cyp6g1 induction in transgenic Drosophila showed tissue-specific expression of this gene controlled by two distinct specific enhancers, suggesting that a single mutation event can modulate Cyp6g1 expression [13].
In contrast to field pest populations, which often possess a highly heterogeneous genetic background, the possibility for the generation of single mutations in a known and characterized background would substantially facilitate the identification of resistance-associated changes. Insertional mutagenesis using transposable elements has been an exceptionally efficient method to create mutants in phylogenetically very distant species, including Drosophila melanogaster [14]. The transposon Minos, a member of the Tc1/mariner superfamily, produces stable transformants with high efficiency in different insect species [15]. This allows genome-wide mutagenesis in insects [16] making Minos a promising genomewide transgenesis tool.
High-throughput deep sequencing transcription profiling is a powerful approach to provide genome-wide information in a very short time and a cost effective way [17]. This method is classified as an ''open'' technology [18], which in contrast to ''closed'' technologies like microarrays, does not require biological or sequence information of the analyzed organism.
Here, by combining a genome-wide insertional mutagenesis screen and next generation transcriptomics, we were able to identify genes involved in Imidacloprid resistance in Drosophila melanogaster within a reasonable time frame and at moderate cost. Gene ontology analysis identified several overrepresented functional gene groups that are differentially expressed in the resistant Drosophila line. The results of our novel approach were in line with previous findings that showed that the Cyp6g1 gene is mainly responsible for resistance. The deep sequencing information was further explored to identify transcription binding factors or microRNAs possibly associated with the over-expression of Cyp genes, which are implicated in resistance. Genetic mapping placed the resistance locus to the right arm of the second chromosome, within a ,1 Mb region in which the Cyp6g1 gene is located.

Results
During the genome-wide insertional mutagenesis, about 12,900 new TREP insertions were generated. The D. melanogaster genome is estimated to contain approximately 13,000 known or predicted genes [19]. Using the information in Metaxakis et al. [16], it can be estimated that in our screen approximately 22% of known or predicted genes of Drosophila genome were hit at least once, directly or within 2 Kb upstream and downstream, excluding introns. Flies with new TREP insertions were selected on medium with 3 mg/ml of Imidacloprid (3 times higher than the LC99 of the susceptible line).
One female carrying a novel TREP insertion on the Xchromosome and exhibiting high resistance to Imidacloprid, was retrieved during the mutagenesis. Line MiT[w 2 ]3X, originating from the retrieved resistant female, was established. Genetic analysis of the MiT[w 2 ]3X line showed that the resistance trait mapped to the second chromosome (see below) but showed no linkage between the TREP insertion and the Imidacloprid resistance.

Resistance to Imidacloprid and cross-resistance to DDT
The resistance was characterized by determining the susceptibility of the resistant mutant MiT[w 2 ]3R2 (derived from MiT[w 2 ]3X line) homozygous for the second chromosome and the control line iso31 to Imidacloprid. The Imidacloprid resistance of line MiT[w 2 ]3R2 was found to be about 18-fold higher than that of the wild-type line iso31 (table 1). Mutant MiT[w 2 ]3R2 showed resistance to Imidacloprid also at the adult stage, as well as cross resistance to DDT (table 1).

Biochemical assays
In order to assess if there is a contribution of known resistance pathways in resistance of mutant MiT[w 2 ]3R2, the activities of cytochrome P450 monooxygenase, esterases and GSTs were analyzed (table 2). Esterase activity was measured using aand bnaphthol, and GST activity was measured using 1-chloro-2,4dinitrobenzene [20]. For both enzymes, no significant difference in activity was detected in the resistant line compared to line iso31. Cytochrome P450-dependent monooxygenase activity was determined by analyzing living third instar larvae, following the protocol of Inceoglu et al. [21]. The activity of cytochrome P450 was 3-fold higher in resistant MiT[w 2 ]3R2 third instar larvae compared to third instar susceptible larvae (table 2). Gene functional classification analysis by grouping genes based on functional similarities identified three functional groups in the up-regulated genes (table 3) and two functional groups in the down-regulated genes (table 4). The cytochrome P450 genes, proteolytic genes and genes showing peptidase activity were overrepresented in the up-regulated genes (table 3). Cuticular protein genes and genes showing peptidase activity were overrepresented in the down-regulated genes (table 4).

Deep sequencing analysis
Functional annotation clustering, which groups genes with similar predicted biological functions, identified 10 overrepresented groups in the up-regulated genes (table 5) and 13 overrepresented groups in the down-regulated genes (table 6). Among the functional groups overrepresented in the up-regulated genes, four clusters are connected to peptidase activity and three functional clusters are connected to P450 gene family activity (table 5). There were also other functional groups with significantly overrepresent-ed members in the over-expressed genes, like oxidoreductase activity, mitotic sister chromatid segregation, electron carrier activity and response to DNA damage. In the down-regulated genes, groups like nutrient reservoir activity, chitin and aminoglycan metabolic processes, response to bacteria and immune response activity were identified (table 6).
We performed quantitative real time PCR to validate the expression difference of two representative cytochrome P450 genes (Cyp6g1 and Cyp6a2), already known to play important role in insecticide resistance. Cyp6g1 showed an expression difference between the resistant and susceptible lines of 8.4 (60.7), while gene Cyp6a2 showed an expression difference of 10.3 (6 1.7) fold. The expression difference of Cyp4p2 was also validated with quantitative PCR, and showed 4.9 (6 0.3) fold elevated expression in the resistant line.

Chromosomal mapping of the Imidacloprid-resistance locus
Mapping of the resistance locus of MiT[w 2 ]3X flies to a chromosome was done with standard genetic tools. Males from Imidacloprid-resistant line MiT[w 2 ]3X were individually crossed with w1118iso/Dp(1;Y)y+; nocSco/SM6a females, who carry chromosome 2 balancer SM6a. Progeny from this cross, heterozygous for the second chromosome, was selected on Imidacloprid as described. Both resistant male and female progeny emerged, implying that the resistance locus does not map to the sex chromosome. Resistant heterozygous male progeny carrying SM6a was individually crossed with iso31 (susceptible line) females. Progeny from this cross was also selected on Imidacloprid. None of the progeny carried the chromosome 2 balancer, which places the resistance on the second chromosome. Equivalent crosses were performed between resistant MiT[w 2 ]3X males and w1118/Dp(1;Y)y+; TM2/TM6C, Sb1 female flies, carrying two balancers of chromosome 3. This cross showed that there is no correlation between the resistance locus and the third chromosome, confirming that the resistance locus is located on the second chromosome. Line MiT[w 2 ]3R2/CyO, which lacked the Minos insertion but carried the resistance trait and a lethal locus was derived from MiT[w]3X. The lethality locus was mapped to the right arm of the second chromosome, using the Bloomington Stock Center Drosophila deletion kit, between position 49C1-4; 50C23-D2 (8.5 Mb -9.9 Mb; figure 1). Recombination mapping placed the lethality locus on the same chromosome arm that harbors the resistance locus. Line MiT[w 2 ]3R2, homozygous for the second chromosome resistance locus, was established during the recombination analysis. MiT

P element mapping
In order to narrow down the position of the resistance locus in the MiT[w 2 ]3R2 line on the 2R chromosome, genetic mapping relative to P element insertions was used. The distance between a P element located at ,0.5 Mb and the resistance locus was determined to be 8.2 cM. The distance between a P element located at ,6.1 Mb and the resistance locus is 3.5 cM. The distance between a P element located at , 6.5 Mb and the resistance locus was determined to be 2.8 cM. The distance between a P element located at , 11.2 Mb and the resistance locus was determined to be 3.0 cM. These genetic distances were converted into Mb using estimates of local recombination rates according to Fiston-Lavier et al. [22] and Singh et al. [23] (figure 1). The genetic mapping relative to the P element insertions places the resistance locus roughly between 8 Mb and 9.7 Mb on  under selection with 3 mg/ml of Imidacloprid in order to homogenize the genetic background. The SNP comparison indicates a hybrid origin of the 2R chromosome, where the right half comes from iso31, while the left half comes from a different line, most likely yw (figure 1). This result indicates a recombination event on 2R, close to the region between 8.5 Mb and 9.9 Mb, to which the lethality was mapped (figure 1). The position of the lethality, resistance and recombination break point shows that the recombination event occurred between the resistance and lethality loci (figure 1). The highly over-expressed Cyp6g1 gene (16.3-folddeep sequencing analysis; 8.4-fold-real time PCR analysis) lies close to the recombination break point to the region into which the resistance locus has been placed (figure 1).
Comparison of the sequences of the Cyp genes differently expressed in the resistant versus the susceptible line showed no sequence changes of the P450 proteins. We also analyzed the flanking sequences of differentially expressed genes for possible common transcription factor binding sites. In silico analysis, using the JASPAR database [24], did not detect common transcription binding factors either for the subgroup of Cyp genes or for all over-expressed genes.
Similarly, a search for predicted targets of microRNAs, performed with DIANA-microT version 3.0 [25], in the 39UTRs of all up-regulated and down-regulated genes, or of just the upregulated and down-regulated Cyp genes, did not identify any significantly overrepresented common target sites.

Discussion
We tested a combined approach of mutagenesis and next generation transcriptomics to study insecticide resistance in the model organism Drosophila melanogaster.
A Minos-based construct was used for genome-wide insertional gene activation mutagenesis. During this screen, an Imidacloprid- Table 3. Gene functional groups in the up-regulated genes (analyzed with the DAVID 6.7 BETA bioinformatics resource).  Table 4. Gene functional groups in the down-regulated genes (analyzed with DAVID 6.7 BETA bioinformatics resource).  Table 5. Functional annotation clusters in the up-regulated genes (analyzed with the DAVID 6.7 BETA bioinformatics resource).
No. clusters 10 Annotation Cluster 1  Table 6. Functional annotation clusters in the down-regulated genes (analyzed with the DAVID 6.7 BETA bioinformatics resource).

No clusters 13
Annotation Cluster 1  (table 1) suggests metabolic resistance as the mechanism of resistance in this line. Furthermore, biochemical analysis showed increased P450 activity in the resistant line compared to the susceptible line (table 2).
The Illumina parallel short-sequencing technology was used to obtain total cDNA sequences of the resistant line and of the nonresistant isogenic line iso31 (w 1118 iso ; 2 iso ; 3 iso ) [26]. This approach was used in order to identify and quantify differences in expression between mutant line MiT[w 2 ]3R2 and susceptible line iso31, covering nearly all D. melangaster genes. Out of 357 genes differently expressed in the resistance line, 150 were up-regulated, and 207 genes were down-regulated in comparison to the susceptible line.
Gene ontology functional classification of the sequenced transcripts identified a significantly up-regulated P450 family group and two groups of genes coding for peptidase activity in the resistant line. Significantly overrepresented down-regulated groups of genes were cuticular protein genes and other peptidase genes.
Deep sequencing analysis detected eight members of the P450 family, Cyp4p2, Cyp6a2, Cyp6g1, Cyp6w1, Cyp4e3, Cyp309a2, Cyp6g2 and Cyp4d14, with elevated expression in the resistant line. Genes encoding glutathione-S-transferases, as well as esterases did not show elevated expression in the resistant line. The most highly over-expressed P450 genes as detected with deep sequencing and confirmed by real time PCR, are Cyp4p2, Cyp6a2 and Cyp6g1 (GSE28560_resistant_vs._susceptible_UPREGU-LATED_GENES.txt.gz). The cytochrome P450 genes play an important role in insecticide resistance, because of their variety and the broad substrate specificity of some P450 genes [27]. We report for the first time elevated expression of the Cyp4p2 gene in a D. melanogaster line resistant to Imidacloprid and DDT, although its role in resistance (if any) remains to be elucidated. The detoxification function of Cyp6a2 and Cyp6g1 in Drosophila is well documented. Over-expression of Cyp6g1 in Drosophila confers resistance to DDT and neonicotinoids [12,28,29] which provides a certain degree of validation in our approach for detecting genes conferring insecticide resistance. The Cyp6a2 is also highly expressed in different insecticide resistant Drosophila strains [30,31,32,33] and the CYP6A2 encoded enzyme can metabolize insecticides [34,35].
Five other P450 genes (Cyp6w1, Cyp4e3, Cyp309a2, Cyp6g2 and Cyp4d14) detected in the resistant line are over-expressed up to 6-fold. Microarray analysis showed that expression of Cyp6w1 is higher in DDT resistant Drosophila strain compared to a susceptible line [33]. Over-expression of the Cyp6g2 gene confers resistance to diazonin and nitenpyram in transgenic Drosophila [28]. To date, Cyp4e3, Cyp309a2 and Cyp4d14 have not been implied in insecticide resistance.
A number of cuticular protein genes were down-regulated in the resistant mutant compared to the susceptible line (table 4). This could occur as a result of the general stress response induced by the up-regulated detoxification system. It is not likely that the down-regulation of cuticular protein genes plays a role in the insecticide resistance mechanism. It would be in disaccord with the fact that reduced cuticular penetration of insecticides can contribute to resistance in some insect species [36,37,38]. The identification of a group of 21 up-regulated genes involved in peptidase activity is consistent with the finding that genes coding for peptidase activity are also significantly over-expressed in DDT resistant Drosophila [33]. The role of the proteolytic genes and genes showing peptidase activity in insecticide resistance is still poorly understood. There is increasing evidence of involvement of protein metabolism in insecticide resistances of different insect species [39,40,41,42,43]. Proteases may be involved in modification of enzyme conformation and protein biosynthesis, in order to meet energy requirements during xenobiotic stress [39].
Other groups of overrepresented members among the up-or down-regulated genes belong to the following categories: oxidoreductase activity, chromosome establishment, organelle localization and cellular response to DNA damage (up-regulated; table 5) and nutrient reservoir activity, response to bacteria, biotic stimulus and immune response (down-regulated; table 6). Down-regulation of genes involved in immune response was not seen in other DDTresistant Drosophila lines [33]. Oxidoreductase activity plays a role in detoxification, while the other biological processes could be an indication of general stress response.
The line MiT[w 2 ]3R2, homozygous for the resistance chromosome, derives from the mutant line MiT[w 2 ]3R2/CyO heterozygous for the second chromosome carrying both resistance and lethality. Genetic analysis of the mutant line MiT[w 2 ]3R2/ CyO line placed the lethality locus to the region between 8.5 and 9.9 Mb on the right arm of the second chromosome of D. melanogaster. Single nucleotide polymorphism analysis between the resistant line homozygous for resistant chromosome MiT[w 2 ]3R2 and the susceptible line iso31 indicates that the recombination event that separated the lethality from the resistance locus occurred in close vicinity to the lethality locus (figure 1). The three highest up-regulated P450 genes Cyp4p2, Cyp6a2 and Cyp6g1 are also located on the right arm of the second chromosome, but they are not closely linked (figure 1). Mapping against P element insertions confirmed that the resistance locus lies on the right arm of the second chromosome between 8 Mb and 9.7 Mb. Chromosomal mapping of the resistance in DDT and Imidacloprid resistant Drosophila lines [44] placed the DDT resistance locus (Rst (2) DDT) in an area that overlaps the interval in which the resistance locus of MiT[w 2 ]3R2 is located. The position of the lethality locus (between 8.5 and 9.9 Mb), together with the SNP analysis and P element mapping, suggests that the resistance locus in MiT[w 2 ]3R2 lies within an interval of less than ,1 Mb (figure 1). Interestingly, the highly up-regulated Cyp6g1 (16.3-fold -deep sequencing analysis; 8.4-fold -real time PCR analysis) gene is located within this range. In the mentioned study of Daborn and colleagues [44] the Cyp6g1 is strongly suggested as the main candidate gene responsible for the resistance in DDT and Imidacloprid resistant Drosophila lines.
The mutation event which causes the resistance in MiT[w 2 ]3R2 remains to be identified. The resistance locus is not linked to an insertion of the transposon used in the screen. It is conceivable that a ''hit and run'' Minos insertion effect might be responsible for the mutation, where the transposon first integrated and then re-excised. In Drosophila, Minos often leaves behind upon excision either a characteristic six bp ''footprint'' or a deletion around the site of insertion [45], both of which can be mutagenic. It has been suggested that mutations of trans-regulating factor/s, or of cis-acting elements of some of the Cyp genes are responsible for insecticide resistance in Drosophila [46,47,48]. A recent report suggests that a single mutation event in a specific enhancer can modulate Cyp6g1 tissue-specific induction in Drosophila flies [13]. One might thus speculate that a single mutation event occurred in a cis-acting element of the Cyp6g1 gene, increasing the expression of this gene. This in turn could activate a resistance cascade, affecting the expression of other Cyp genes involved in resistance. Alternatively, the mutation might involve a gene encoding a transcription factor or a microRNA which regulates in trans the Cyp genes involved. We have so far no evidence for the latter assumption, since an in silico search failed to identify common transcription factor motifs regulating the over-expressed P450 genes. The same is true for common predicted microRNA targets in the 39UTRs. Additional analyses are required in order to pinpoint the exact cause of resistance in the MiT[w 2 ]3R2 mutant.

Drosophila lines
D. melanogaster stocks were maintained on standard cornmealagar-yeast medium at 24uC with a 12-hour light/12-hour dark cycle. We analyzed a Drosophila melanogaster line (MiT[w 2 ]3R2) resistant to the neonicotinoid Imidacloprid, retrieved during a transposon Minos-based insertional mutagenesis screen. Three transgenic Drosophila lines were used for Minos-based insertional mutagenesis: TREP 2.30, BOEtTA [49] and an iso31 derivative [SM6a, MiT 2.4]/Sco] [16]. Line 2.30 carries a single insertion of the Minos transposon TREP, which contains a minimal promoter driven by the tetO operator. When inserted next to or into a gene, TREP can cause over-expression of the gene in the presence of the tetracycline trans-activator tTA [50, Kiupakis, Oehler and Savakis, manuscript in preparation]. Line BOEtTA carries a single insertion of a Minos transposon which produces tTA [49]. The mobilization of the TREP construct and generation of flies with new insertions was performed with a standard ''jumpstarter'' system [51] (figure 2). Flies with new TREP insertions were selected for insecticide resistance during egg to adult development on medium with Imidacloprid. We aimed to generate highly resistant flies, thus mutagenized flies were selected on a concentration of Imidacloprid 3 times higher than the LC99 of susceptible line iso31 (3mg/ml). Female individuals carrying both TREP 2.30 and BOEtTA construct were scored for resistance (figure 2). The insertional site distribution of new TREP insertions was estimated according to Metaxakis et al. [16]. The correction for multiple insertions into the same genes was done using the Poisson distribution, assuming the same probability of recovering insertions for all loci [52]. One resistant female with a Minos insertion located on the X chromosome was retrieved from the screen and further analyzed.

Bioassays
The insecticide Imidacloprid (98.7%, Bayer CropScience GmbH) was added directly to the standard medium for the screening experiments. Per vial, 50 eggs were added, and egg-toadult viability was determined. Six different concentrations were used with 8 replicas per concentration.
A contact assay was used for Imidacloprid and DDT assay. 35 ml glass vials were coated with DTT (DDT 4,49-DDT PESTANALH, SIGMA-ALDRICH) on the inside by evaporating 200 ml of acetone containing the required amounts of DDT. Coated vials were plugged with cotton wool soaked with 5% sucrose. The mortality of 25 flies (1-3 days old) per vial was scored after 24 hours. Each DTT amount was assayed for sex different concentration with 4 replicas per concentration.

Biochemical assays
The activities of esterases, GSTs and cytochrome P450 dependent monooxygenases were determined as previously described [20,21]. The activity of cytochrome P450 dependent monooxygenases was measured in live third instar larvae. For all assays, activity was measured at 25uC on a SpectraMax M2 microplate reader and quantified with the integrated software SoftMax pro v5.

Deep sequencing analysis
Deep sequencing data were deposited to the GEO site (series record GSE28560). The following URL allows review of the record.
Total RNA from 3 day old Drosophila melanogaster adults was extracted using standard Trizol RNA isolation protocol (http:// quantgen.med.yale.edu/). Preparation of cDNA for sequencing was done with the Illumina mRNA Seq V2 kit (Illumina, Inc). Formation of single molecule arrays, cluster growth and sequencing were done according to the standard protocols from Illumina, Inc. Sequencing was performed with a 2008 Illumina Genome Analyzer, version 2 (GA2). Mapping of the 51 nucleotides (nt) long sequencing reads of both lines, MiT[w 2 ]3R2 and iso31, to the reference genome (Drosophila release 5 sequence assembly Flybase) was performed with software RMAP, version 2.05 [55]. Genes with 10 or less reads for one line and 50 reads or less for the other line were excluded from further analysis. The minimum difference threshold between lines was set to 2-fold. Analysis of up-regulated and down-regulated genes was performed with the DAVID 6.7 BETA release bioinformatics resources [56].

Quantitative Real Time PCRs
Total RNA from D. melanogaster 1-3 day old adult flies was extracted using same protocol as for deep sequencing analysis. Synthesis of first-strand cDNA and PCR reactions were performed on a Techne TC-412 PCR machine. Synthesis of the First-Strand cDNA was done with the AccuScriptH High Fidelity RT-PCR System. Expression of Cyp6g1 and Cyp6a2 was measured relative to the housekeeping ribosomal protein gene Rp49. For this purpose, three sets of primers were designed. In order to obtain products specific for the cDNA templates, primers were designed to span exon-intron junctions. The forward and reverse primer sequences were as follows: Cyp6g1 -59ACCCTTATGCAGGA-GATTG39 and 59TAGGCTGTTAGCACGAATG39; Cyp6a2 -59GTTACTGCCTGTATGAGTTGG39 and 59TAGAGCCT-CAGGGTTTCTG39; Rp49 -59CGGTTACGGATCGAA-CAAGCG39 and 59TTGGCGCGCTCGACAATCT39. Quantitative real time PCR was performed using the QIAGEN SYBR green kit with the DNA Engine Opticon TM MJ Research analyzer. Three technical replicates were performed on each of three biological samples. The efficiency of PCR amplification with each gene-specific primer pair was analyzed with five serial dilutions in three technical replicates. Cycling conditions were: 94uC for 5 min, then 37 cycles of 94uC for 30 sec, 52uC for 30 sec and 72uC for 30 sec (plate reading at 78uC, 80uC and 82uC). Data were analyzed with the MJ Opticon Monitor 3.1 analysis software. Calculations were done with software REST-MCS [57]. Additionally, relative expression of the Cyp4p2 in the resistant line was analyzed. Quantitative real time PCR for Cyp4p2 was performed using same protocols as for Cyp6g1 and Cyp6a2 expression analysis. Flies maintained for more than 25 generations on standard medium after deep sequencing and Cyp6g1 and Cyp6a2 expression analysis, were used for Cyp4p2 expression analysis. The forward and reverse primer sequences were as follows: Cyp4p2 -59 CTGAAAAGGCATCCTTACGC 39 and 59 TTGGGATC-GATAACAGGCAG 39. Quantitative real time PCR was performed on the Bio-Rad CFX analyzer with cycling conditions: 95uC for 2 min, then 35 cycles of 95uC for 15 sec, 55uC for 30 sec and 60uC for 30 sec (melt curve 60 to 95 C, increment 1.0 C).

Mapping of the lethality
We used the Bloomington Stock Center Drosophila deletion kit (111 lines) for the second chromosome to map the lethality locus (text S1). Individual crosses were set up between MiT[w 2 ]3R2/ CyO flies carrying a balancer of the second chromosome (SM6a/ lethality) and the lines from the deletion kit (also carrying a balancer of the second chromosome). The progeny was scored for the presence of the balancer chromosome.

P element mapping of the resistance locus
To narrow down the resistance locus on the right arm of the second chromosome, four lines with P elements insertions were employed (Bloomington Drosophila Stock Center, IU; text S2). The resistant MiT[w 2 ]3R2 flies do not carry any visible marker gene, while the flies carrying P element insertions have w + as phenotypic marker. Resistant flies were mass crossed with flies carrying the P element insertion. Virgin female progeny with red eyes (one chromosome deriving from the resistant line and the other from the P element line) were collected and crossed with iso31 males. For each experiment, 50 female flies, heterozygous for the resistance chromosome were crossed with 25 iso31 males, per replica. Each experiment had eight replicas with a total number of 250 females for each P element line. After 2-3 days, crossed flies were transferred to medium with 3 mg/ml of Imidacloprid. Progeny was scored for recombination events. At least 1000 emerged flies were analyzed per replica. Recombination rates were calculated as the ratios of the total number of recombinant flies over the total number of emerged flies. The distance between the P element insertions and resistance was calculated in centimorgans (cM), from which the physical distance was calculated using estimates of the local recombination rates at the sites of the P element insertions after Fiston-Lavier et al. [22] and Singh et al. [23]. This estimate was not possible for one of the P elements, which is too close (about 0.5 Mb) to the centromere. Here, the recombination rate for the interval between the P element and the average position of the resistance locus, as determined relative to the other three P elements, was calculated.

In silico analysis
A comparison of single nucleotide polymorphisms (SNP) of the deep sequence data between the resistant line MiT[w 2 ]3R2 and the susceptible line iso31 was done for the right arm of the second chromosome. Genomic SNP analysis of the pooled assembly of the resistant and the susceptible strains reads were done with the Gigabayes SNP discovery algorithm (improved PolyBayes algorithm [58]) and MOSAIC algorithm [59], using all Refseq mRNA transcripts of the dm3 assembly [60] as a reference. A polymorphism probability threshold of 0.9 is used, with alleles requiring a minimal overall coverage of 10 and of 5 for the minor allele. A SNP density track with the number of SNPs in 1 Kb tiling windows was created. The SNP density was visualized with the UCSC Genome Browser on D. melanogaster release 5 (http:// genome.ucsc.edu/) [61] and is presented for chromosome 2R.
The in silico search for overrepresented transcription factor binding sites was conducted using the JASPAR database [24]. All up-regulated and down-regulated genes, as well as the subset of Cyp genes (up-regulated, down-regulated and all), were analyzed for the presence of transcription factor binding sites. The sequence of all genes was retrieved from Flybase (Drosophila release 5 sequence assembly) [62]. Regions from of 3Kb upstream to1Kb downstream of the gene start and the 39UTR sequences were analyzed.
A survey of predicted targets of microRNAs in the 39UTR of all up-regulated and down-regulated genes, as well as the subsets of up-regulated and down-regulated Cyp genes, was performed with DIANA-microT (version 3.0) [25].
We also compared the sequences of Cyp genes differently expressed in the resistant versus the susceptible line for nucleotide differences. Comparison of the DNA sequences and translation to amino acids were done with the APE software (http://biologylabs. utah.edu/jorgensen/wayned/ape/).

Supporting Information
Text S1 Deletion kit stocks (second chromosome).