Prerpints.org logo
Preprint
Article

Comparative Transcriptomic Analysis of Insecticide Resistant Aedes aegypti from Puerto Rico Reveals Insecticide-Specific Patterns of Gene Expression

Altmetrics

Downloads

182

Views

97

Comments

0

A peer-reviewed article of this preprint also exists.

supplementary.zip (25.03MB )

Submitted:

10 July 2023

Posted:

12 July 2023

You are already at the latest version

Alerts
Abstract
Aedes aegypti transmits major arboviruses of public health importance, including dengue, chikungunya, Zika, and yellow fever. The use of insecticides represents the cornerstone of vector control, however insecticide resistance in Ae. aegypti has become widespread. Understanding the molecular basis of insecticide resistance in this species is crucial to design effective resistance management strategies. Here, we applied Illumina RNA-Seq to study the gene expression patterns associated with resistance to three widely used insecticides (malathion, alpha-cypermethrin and lambda-cyhalothrin) in Ae. aegypti populations from 2 sites (Manatí and Isabela) in Puerto Rico (PR). Cytochrome P450s were the most over-expressed detoxification genes across all resistant phenotypes. Some detoxification genes (CYP6Z7, CYP28A5, CYP9J2, CYP6Z6, CYP6BB2, CYP6M9, and two CYP9F2 orthologs) were commonly over-expressed in mosquitoes that survived exposure to all three insecticides (independent of geographical origin) while others including CYP6BY1 (malathion), GSTD1 (alpha-cypermethrin), CYP4H29 and GSTE6 (lambda-cyhalothrin) were uniquely overexpressed in mosquitoes that survived exposure to specific insecticides. The gene ontology (GO) terms associated with monooxygenase, iron binding, and passive transmembrane transporter activities were significantly enriched in 4 out of 6 resistant vs susceptible comparisons while serine protease activity was elevated in all insecticide resistant groups relative to the susceptible strain. Interestingly, cuticular-related protein genes (chinase and chitin) were predominantly downregulated, which was also confirmed in the functional enrichment analysis. This RNA-Seq analysis presents a detailed picture of the candidate detoxification genes and other pathways that are potentially associated with pyrethroid and organophosphate resistance in Ae. aegypti populations from PR. These results could inform development of novel molecular tools for detection of resistance-associated gene expression in this important arbovirus vector and guide the design and implementation of resistance management strategies.
Keywords: 
Subject: Biology and Life Sciences  -   Insect Science

1. Introduction

Chemical insecticides are widely used to control insects of public health significance. However, the intensive and widespread use of these insecticides has led to rapid and widespread evolution of insecticide resistance in major insect disease vectors, raising concerns about their sustained use in public health. Developing strategies to maintain efficacy of the limited number of chemical insecticides approved for public health use in the face of growing insecticide resistance is a major priority for the control of malaria. This effort could benefit from an in-depth understanding of the molecular basis of resistance to key insecticides used in public health vector control.
Mosquitoes are among the major disease vectors that have become increasingly resistant to the principal classes of chemical insecticides approved for vector control (Ranson and Lissenden 2016). Target site insensitivity and increased metabolic detoxification are the main mechanisms responsible for mosquito adaptation to insecticides, but other mechanisms such as behavioral avoidance of insecticide-treated surfaces and cuticular alterations have also been reported (Hemingway et al. 2004, Chareonviriyaphap et al. 2013, Liu 2015, Simma et al. 2019). Target site insensitivity occurs when a mutation of the protein targeted by the insecticide prevents the insecticide from binding to the molecule. These mutations are well characterized in mosquitoes and molecular tests for their detection are available for many major vectors (Weill et al. 2004, Bass et al. 2007, Mathias et al. 2011). In contrast, metabolic resistance is much less understood and is primarily conferred via gene amplifications, transcriptional enhancements, and coding mutations in detoxification enzymes, particularly cytochrome P450 monooxygenases, carboxylesterases and glutathione-S-transferases (Li et al. 2007). Because metabolic resistance can arise through multiple mechanisms and involves many gene families, understanding its molecular basis remains a complex public health challenge. Moreover, the specific detoxification genes involved in metabolic resistance can vary greatly between insecticides and across mosquito populations and species (David et al. 2013, Ishak et al. 2017, Riveron et al. 2017). These complexities hamper the development of universal molecular diagnostic markers for metabolic resistance and highlight the need for site-specific understanding of the molecular mechanisms responsible for insecticide resistance to facilitate development of evidence-based resistance management strategies.
Here, we use RNA-Seq to examine the gene expression profiles of malathion, alpha-cypermethrin, and lambda-cyhalothrin resistant populations of Aedes aegypti from two locations in Puerto Rico. This mosquito species lives near humans in urban settings and is the primary vector of several globally important arthropod-borne viruses (arboviruses) including dengue, chikungunya, Zika and yellow fever. Resistance of Ae. aegypti to multiple insecticide classes has been reported worldwide posing serious operational challenges towards its control (Braga et al. 2004, Vontas et al. 2012, Amelia-Yap et al. 2018, Schluep and Buckner 2021). Both target site insensitivity and metabolic resistance have been identified as important mechanisms of insecticide resistance in Ae. aegypti and a few previous studies have investigated the specific genes involved in metabolic detoxification of insecticides in this mosquito species (Strode et al. 2008, Lertkiatmongkol et al. 2010, Lima et al. 2011, David et al. 2014, Faucon et al. 2015, Aponte et al. 2019, Badolo et al. 2019, Lien et al. 2019, Saavedra-Rodriguez et al. 2021, Sene et al. 2021, Samal et al. 2022). This study seeks to disentangle the transcriptional basis of resistance to multiple insecticides in Ae. aegypti from two geographically distinct populations, with the aim of better understanding the specificity of metabolic resistance with respect to insecticide and geography.

2. Materials and Methods

2.1. Study area and mosquito collections

Three populations of Ae. aegypti were used in this study: the Rockefeller (ROCK) reference susceptible strain (obtained from the Malaria Research and Reference Reagent Resource Center, MR4) and two Ae. aegypti populations collected from the towns of Isabela (18.4704° N, 67.0242° W) and Manatí (18.4330° N, 66.4759° W) in Puerto Rico in spring of 2017. Eggs were collected from the adults and shipped to the U.S. Centers for Disease Control and Prevention (CDC) in Atlanta, USA. At CDC, the mosquitoes were reared at 27 ± 2°C and 70 ± 10% relative humidity, and a 14L:10D photoperiod. Adults were provided 10% sucrose ad libitum. All three strains were reared in parallel for each of the biological replicates for the RNA-Seq libraries.

2.2. Insecticide bioassays

Three-to-four-day old non-blood fed and unmated Ae. aegypti females from Isabela (ISA) and Manatí (MAN) were tested for resistance to the diagnostic doses of alpha-cypermethrin (10 μg/bottle), malathion (50 μg/bottle), and lambda-cyhalothrin (10 μg/bottle) using the CDC bottle bioassay (Brogdon and Chan 2010). Assays were carried out in four replicates, each containing approximately 25 individuals per bottle. For each test, a control bottle coated with only 1 ml of acetone was included. Mosquitoes alive after 30 minutes of insecticide exposure were considered resistant. Assays were repeated to obtain sufficient numbers of phenotyped mosquitoes for RNA extraction. For the susceptible ROCK strain, three-to-four-day-old adult non-blood fed and unmated female mosquitoes were killed by freezing and stored at -80˚C until RNA extraction was conducted. For ISA and MAN, mosquitoes that were alive (survivors) at the end of the 30-minute exposure period were separated and tallied, then frozen at -80 °C. In addition, individuals from ISA and MAN that were used in the control bottle and not exposed to insecticide were killed by freezing and stored at -80˚C for RNA extraction and sequencing.

2.3. RNA and library preparation methods

Three biological replicates with pools of 5 mosquitoes each were prepared for three groups of mosquitoes: the susceptible Ae. aegypti Rockefeller strain (ROCK), field-collected mosquitoes that were not exposed to insecticides during the bottle bioassay, and field-collected mosquitoes that survived exposure to insecticides at the diagnostic dose and diagnostic time during the bottle bioassay. RNA was extracted using the Applied Biosystems Arcturus PicoPure RNA isolation kit (Arcturus, Applied Biosystems, USA) and quantified using the Agilent RNA ScreenTape 4200 assay, according to the manufacturer’s protocols. Two micrograms of starting material were treated with Baseline-ZERO™ DNase (Lucigen, USA) and ribosomal RNA was removed using the Ribo-Zero™ Magnetic Core Kit and Ribo-Zero™ rRNA Removal kit (Illumina, USA), according to the manufacturer’s protocols. RNA-Seq libraries were prepared from each pool of extracted RNA using the ScriptSeq™ v2 RNA-Seq library preparation kit (Illumina, USA), using 12 cycles of PCR amplification, according to the manufacturer’s protocol. Libraries were purified using Agencourt AMPure XP beads (Beckman Coulter, USA) and assessed for quantity and size distribution using the Agilent DNA ScreenTape D5000 assay. Each library was barcoded, and equal amounts of each library were pooled and sequenced (2x125bp paired end) on an Illumina HiSeq 2500 sequencer, using v2 chemistry. Sequencing was done at the Biotechnology Core Facility at CDC, Atlanta, USA.

2.4. Measurement of gene expression by quantitative PCR

Five genes that were significantly differentially expressed in the RNA-Seq analysis were validated using quantitative PCR (qPCR). One microgram of RNA from 3 replicates of samples resistant to lambda-cyhalothrin from Manatí and Isabela, and samples resistant to alphacypermethrin from Isabela were used to synthesize cDNA using the High-Capacity cDNA reverse transcription kit (Applied Biosystems) with oligo-dT20 (NEB), according to the manufacturer’s instructions. RNA and cDNA synthesis was also prepared for the susceptible laboratory strain ROCK. There was insufficient material remaining from the unexposed samples from both locations hence this was not included in the qPCR validation. The primers used are listed in Table S1. Standard curves of Ct values for each gene were generated using a serial dilution of cDNA, allowing assessment of PCR efficiency. The qPCR amplification was carried out on a QuantStudio 6 Flex Real-Time PCR system (Applied Biosystems) using PowerUp SYBR Green Master Mix (Applied Biosystems). cDNA from each sample was used as a template in a three-step program as follows. Uracil-DNA glycosylase (UDG) activation at 50°C for 2 minutes, DNA polymerase activation at 95°C for 10 minutes, followed by 40 cycles of DNA denaturation for 15 seconds at 95°C, DNA annealing and extension for 1 minute at 60°C, and a last DNA extension step of 15 seconds at 95°C. The relative expression level and Fold Change (FC) of each target gene from resistant field samples relative to the susceptible lab strain were calculated using the 2−ΔΔCT method (Schmittgen et al 2008) incorporating the PCR efficiency. Housekeeping genes RPS3 (Smith et al 2018) and Actin (Mitra et al 2021) were used for normalization.

2.5. Data analysis

2.5.1. Quality control filtering and mapping

De-multiplexed Paired-End (PE) RNA-Seq reads for each sample were evaluated for quality using FastQC v0.11.5 (Andrews 2010). Subsequently, the PE reads were trimmed and filtered using fastp v0.21.0 (Chen et al. 2018) to remove adapter and low-quality reads according to the following criteria: minimum base quality score = 20, minimum length required = 25, polyG and poly tail trimming = True. The resulting trimmed and filtered read pairs (R1/R2) were aligned against the reference genome of Ae. aegypti (genome assembly version = AaegL5.1, GeneBank assembly identifier = GCA_002204515.1) (Giraldo-Calderón et al. 2015), using ‘subjunc’ v2.0.1, part of the subread aligner v2.0.1 (Liao et al. 2013), with default parameters. The resulting alignment was filtered to remove reads with low mapping quality (q < 10) and sorted successively using Samtools v1.10 (Li et al. 2009). Summary statistics for the QC filtering and sequencing alignments corresponding to malathion, alpha-cypermethrin, and lambda-cyhalothrin are reported in Table S2. Tags (a read pair or single, unpaired read) mapped to the sense orientation of the annotated Ae. aegypti gene set AaegL5.1 were quantified with ‘FeatureCounts’, as part of the subread-aligner package v2.0.1 (Liao et al. 2013), by using the following criteria: 1) count only read pairs that have both ends aligned; 2) count fragment instead of reads; 3) minimum number of overlaps required = 1; 4) feature type = exon; 5) attribute type = gene_id; and 6) strandness = sense. The FeatureCount analysis generated a tag count matrix table which was used as input to edgeR (Robinson et al. 2010) for differential gene expression analysis.

2.5.2. Differential Gene Expression Analysis

For each group, three pairwise comparisons of differential gene expression (DGE) analysis were conducted to identify genes associated with resistance to each insecticide: 1) resistant vs unexposed control (R-C), 2) resistant vs ROCK susceptible (R-S), and 3) unexposed control vs ROCK susceptible (C-S). The R-C comparison was used to detect transcription induced by insecticide exposure, while the R-S and C-S comparisons were used to distinguish genes that are constitutively associated with insecticide resistance as well as those induced by insecticide exposure.
The DGE analysis was conducted using the R package edgeR (Robinson et al. 2010). To remove the effect of noise and low-expressed genes, for each pairwise comparison, genes with a total tag count less than 50 across all libraries (any pairwise comparison) were filtered out before further analysis as previously suggested (Robinson et al. 2010, Chen et al. 2020). The function calcNormFactors was used to normalize tag counts among samples, by using the TMM (Trimmed Mean M-values) method. This approach normalizes the dataset for sequencing depth and library size, the two most important technical factors that influence DGE. The tag count was not normalized for gene length and GC content, as these values do not vary from sample to sample and would be expected to have little effect on differentially expressed genes (DEGs). The tag-wise and common dispersion of the read count distribution were estimated using the estimateDisp function from the edgeR package. The DEGs were selected after multiple testing using the decideTests function in the limma package (Ritchie et al., 2015). A critical value absolute log2 fold-change |log2FC| = 1 (|FC| =2), and a (False Discovery Rate) FDR adjusted P value ≤ 0.01 were used to tag a gene as significantly differentially expressed.

2.5.3. Gene ontology annotation and functional enrichment analysis

The Ae. aegypti reference genome corresponding to genome assembly version ‘AaegL5.1’ and GeneBank assembly identifier ‘GCA_002204515.1’ (directly downloaded from Vectorbase, Release 51) was used. AaegL5.1 contains 14,718 protein coding genes, 4,704 ncRNA genes and 382 pseudogenes (corresponding to 19804 genes) (Matthews et al. 2018). However, functional annotation is provided for only 6,397 (36.68%), and gene ontology annotation is provided for only 10,910 (74.13%) genes in VectorBase. To help with the interpretation of our results, predicted genes of AaegL5.1 were functionally annotated using Blast2GO as follows. First, a local BLASTp (v2.9) search of the predicted protein coding sequences was conducted against the Arthropoda (taxid=6665) category of the nr protein NCBI database with maximum e-value cut-off 10-3. Second, the protein sequences were searched against the InterPro database (Hunter et al. 2008), using InterProScan v5 (Jones et al. 2014). The Blastp and InterProScan outputs were simultaneously provided to Blast2GO command line as input, which map the RefSeq and InteProScan identifiers to the GO database as curated and updated in the Blast2GO database (July 2021). This analysis assigned putative descriptions to 13,302 (90.37%) protein coding genes and gene ontology annotation to 11,672 (79.30%), which could be considered a significant improvement. All annotation results (VectorBase and Blast2GO) are provided in Table S3. GO term enrichment analysis of up- and down-regulated genes was carried out using Goatools (Klopfenstein et al. 2018), based on the go-basic database (release 2021-02-01). The list of 11,672 blast2GO annotated genes of Ae. aegypti with their associated GO terms was used as the background reference set. The P-values used to evaluate significantly enriched GO terms were calculated based on Fisher’s exact test and corrected with the Benjamini-Hochberg multiple test correction method(Klopfenstein et al. 2018). Finally, we used an FDR adjusted P-value <0.05 to define statistically significantly overrepresented (enriched) and underrepresented (purified) GO terms associated with the list of DEGs.

3. Results

3.1. Bioassay results

The mosquitoes from both sites showed different patterns of resistance to the three insecticides tested. After 30 min of exposure to alpha-cypermethrin, lambda-cyhalothrin, or malathion the mortality rate of mosquitoes from Manatí was 46%, 85% and 56%, respectively, compared to 81%, 48%, and 25%, respectively, for mosquitoes from Isabela (Figure 1, Table S4).

3.2. RNA sequencing, quality control filtering and alignment rate

The Table S2 shows the statistics of the RNA-Seq sequencing results (before and after quality filtering), the alignment of the filtered reads to the reference genome and the read quantification. The Illumina sequencing generated approximatively 27.6-69.4 million of reads per sample. After adapter trimming and read quality filtering, an average of 41.52 ± 10.43 million reads per sample were retained for subsequent analysis. An average of 24.12 ± 9.43 million quality filtered PE reads per sample were mapped to the reference genome of Ae. aegypti, while an average of 7.2 ± 3.03 million tags (read pairs) per sample, representing 50% to 74% of the total alignment, were successfully assigned to the exonic regions of the gene set AaegL5.1 (Figure S1).

3.3. Differential gene expression analysis

3.3.1. Differential gene expression associated with malathion resistance

For ISA samples, 1132 genes (441 upregulated and 691 downregulated), 563 genes (256 upregulated and 307 downregulated), and 119 genes (33 up regulated and 86 down regulated) were significantly differentially expressed in R-S, C-S and R-C comparisons, respectively. For MAN samples, a total of 1,396 genes (739 upregulated and 657 downregulated), 653 genes (197 up regulated and 456 downregulated), and 215 genes (174 upregulated and 41) were significantly differentially expressed in R-S, C-S, and R-C comparisons, respectively (Table 1, Figure 2, Table S5).
The DEGs that overlapped 2 or more comparisons for each experiment (4 intersections in total) were extracted and inspected for important gene expression patterns (Figure S2). For the ISA samples, there were 69 DEGs (24 upregulated and 45 downregulated) that were shared between R-S and R-C groups (Figure S2A). The top 5 shared over-expressed genes between R-S and R-C groups included 2 odorant binding proteins (OBPs) (the duplicated OBP66), a kDa salivary protein, a zinc finger 512B-like, and an uncharacterized protein (Table S6). However, for MAN samples, there were 148 differentially expressed genes (126 upregulated and 22 downregulated) that overlapped between R-S and R-C groups (Figure S2B); the top 10 shared upregulated genes included 7 OBPs (including duplicated OBP66, and OBP67), serine protease snake-like, and two uncharacterized proteins (Table S6).
Focusing on R-S comparisons, we identified 753 genes (309 upregulated and 444 downregulated) that were differentially expressed in malathion survivors from both ISA and MAN relative to the susceptible ROCK strain (Figure 3A). The top 15 upregulated genes with retrievable annotations included 2 serine proteases, 3 uncharacterized proteins, 3 odorant binding proteins, a 39S ribosomal mitochondrial protein, a digestive organ expansion factor homolog, a kDa secreted -1, a nuclear valosin-containing -like, nucleolar dao-5, peritrophin-1-like, pre-mRNA-processing factor 39, and peritrophin-1-like protein coding genes (Table S7).
The detoxification genes annotated as cytochrome P450 monooxygenases were predominantly over-expressed, while the cuticular-related protein genes were predominantly downregulated in malathion resistant samples from both ISA and MAN (Figure 2A-B). 83% (20/24) and 82% (27/33) of differentially expressed CYP450s were upregulated in the R-S comparisons from ISA and MAN, respectively, while 78% (18/23) and 91% (39/33) of cuticular-related DEGs were significantly down regulated in R-S comparisons from ISA and MAN, respectively (Figure S3).
From the list of DEGs shared between the R-S groups from the two sites, we have identified a list of core detoxification genes that are associated with malathion resistance in Ae. aegypti in both populations. These core candidate detoxification genes include 22 CYP450s (19 upregulated: CYP6Z7, CYP6D5, CYP6BY1, CYP325Q1, CYP9J2, CYP6N9, CYP6Z6, CYP325V1, CYP6N12, CYP6BB2, CYP6M9, CYP9J23, CYP9J22, CYP9J9, CYP4G36 and four CYP9F2 orthologs; and 3 downregulated: CYP12F5, CYP325N2, and CYP6N15), 1 upregulated carboxylesterase and 3 glutathione-s-transferases (including 1 upregulated: GSTD6, and 2 downregulated: GSTS1 and GSTX1) (Figure 4, Table 2). Although ISA Ae. aegypti showed higher resistance to malathion than MAN (Figure 1), we did not detect any detoxification genes that were uniquely overexpressed in malathion-resistant ISA. However, several CYP450s (including CYP4D39, CYP28A5, CYP4J16, and CYP9J23) were overexpressed in malathion-resistant MAN, but not in malathion-resistant ISA. A detailed summary of all the DEGs shared between ISA and MAN is shown in Table S7.

3.3.2. Differential gene expression associated with alpha-cypermethrin resistance

For ISA, 1171 genes (642 upregulated and 529 downregulated), 478 genes (158 upregulated and 322 downregulated), and 433 genes (87 up regulated and 346 down regulated) were significantly differentially expressed in R-S, C-S and R-C comparisons, respectively. For MAN, a total of 1175 genes (605 upregulated and 570 downregulated), 811 genes (222 upregulated and 589 downregulated), and 602 genes (154 upregulated and 448 downregulated) were significantly differentially expressed in R-S, C-S, and R-C comparisons, respectively (Table 1, Table S5). The list of DEGs that overlapped 2 or more comparisons for each experiment (4 intersections) (Figure S2B–C), were manually extracted and are reported in Table S6. Unlike the malathion experiment, no DEGs were shared between the R-S and R-C groups.
Focusing on the R-S comparisons, a total of 801 differentially expressed genes (421 upregulated and 380 downregulated) were shared between the MAN and ISA alpha-cypermethrin survivors, relative to the susceptible strain (Figure 3B). The top 10 overlapping DEGs with retrievable annotation (FC=15.1-2101.6) included a serine protease, a peritrophin-1-like, a kDa secreted salivary protein, an uncharacterized protein LOC5569546 isoform, an uncharacterized protein LOC5569546 isoform X1, NPC intracellular cholesterol transporter, a nuclear valosin-containing -like, a V-type proton ATPase subunit d1, and 2 cytochrome P450s (CYP67 and CYP9F2) (Table S7).
Similar to the malathion analysis, the detoxification genes annotated as cytochrome P450 monooxygenases were predominantly over-expressed, while the cuticular protein genes were predominantly downregulated in both alpha-cypermethrin survivors (R-S) and unexposed groups (C-S) relative to the susceptible strain for both study sites (Figure 2C-D). In particular, 76% (22/29) and 74% (20/27) of differentially expressed CYP450s were upregulated in the R-S comparisons for ISA and MAN, respectively, while 97% (34/35) and 97% (33/34) of cuticular related DEGs were significantly down regulated in R-S comparisons for ISA and MAN, respectively (Figure 2, Figure S3). From the list of DEGs that overlapped the two R-S groups, a list of core detoxification genes associated with alpha-cypermethrin resistance was identified due to their consistent differential expression in both MAN and ISA survivors (Figure 4, Table 2). These candidate alpha-cypermethrin detoxification genes included 21 CYP450s, of which 15 were significantly upregulated (CYP6Z7, CYP9J2, CYP6N12, CYP6BB2, CYP9J22, CYP6N9, CYP6M9, CYP6D5, CYP325V1, CYP6Z6, CYP9J9, CYP4G36 and three CYP9F2 orthologs) and 6 were downregulated (including CYP12F5, CYP6Y3, CYP4AC1, CYP325N, CYP6N15, and CYP4H32) (Figure 4A). Additionally, we identified 4 GST genes (GSTS1, GSTD1, GSTD11 and GSTE4) that were downregulated in both ISA and MAN alpha-cypermethrin survivors as compared to the susceptible strain (Figure 4B).

3.3.3. Differential gene expression associated with lambda-cyhalothrin resistance

For ISA samples, 1261 genes (341 upregulated and 920 downregulated), 480 genes (158 upregulated and 322 downregulated), and 286 genes (30 upregulated and 256 downregulated) were significantly differentially expressed in R-S, C-S and R-C comparisons, respectively. For MAN samples, a total of 1413 genes (362 upregulated and 1051 downregulated), 811genes (222 up regulated and 589 downregulated), and 254 genes (67 upregulated and 197) were significantly differentially expressed in R-S, C-S, and R-C comparisons, respectively (Table 1, Table S5). The DEGs that overlapped 2 or more comparisons for each experiment were manually extracted and are shown in Table S6. Although 174 DEGs (8 upregulated and 166 downregulated) and 139 DEGs (12 upregulated and 127 downregulated) were identified in the intersection R-S/R-C, for ISA and MAN, respectively, no upregulated genes were shared between the resistant, susceptible and controls from each study site (Figure S2B-C). Among the over-expressed genes in the R-S/R-C intersection for MAN was a glutathione-s-transferase gene (GSTD6) with fold change expression of 10.56 and 2.73 in R-S and R-C, respectively. However, of the 8 upregulated genes that overlapped the R-S and C-S groups for ISA, none were previously reported to be associated with insecticide detoxification (Table S6).
Focusing on only the R-S groups, a total of 943 differently expressed genes (203 upregulated and 743 downregulated) were shared between MAN and ISA lambda-cyhalothrin survivors relative to the susceptible strain (Figure 3C). The top 10 shared upregulated genes (FC = 14.2-1901.7) included a serine protease, a peritrophin-1-like, a digestive organ expansion factor homolog, an uncharacterized protein, nucleolar dao-5, pre-mRNA-processing factor 39, kDa secreted salivary protein, nuclear valosin-containing -like, NPC intracellular cholesterol transporter 2 homolog a, and a CYP450 (CYP6Z7) (Table S7). Similar to the malathion and alpha-cypermethrin experiments, the detoxification genes annotated as cytochrome P450 monooxygenases were predominantly over-expressed, while the annotated cuticular protein genes were predominantly downregulated in both R-S and C-S groups (Figure 2). In particular, 69% (18/26) and 55% (12/22) of differentially expressed CYP450s were significantly overexpressed in the R-S comparisons from ISA and MAN, respectively, while 97% (34/35) and 97 % (33/34) of cuticular-related DEGs were significantly downregulated in R-S comparisons of ISA and MAN lambda-cyhalothrin survivors, respectively (Figure 2, Figure S3).
From the DEGs that overlapped the R-S comparisons of both study sites, several core insecticide detoxification genes were identified. These included 17 cytochrome P450s (10 upregulated: CYP6Z7, CYP9J2, CYP28A5, CYP6Z6, two CYP9F2 orthologs, CYP6M9, CYP6BB2, CYP307A1, CYP4C2-like; and 7 downregulated: CYP6P12, CYP6Y3, CYP4H29, CYP6S3, CYP6N15, CYP325N2, CYP4C1-like). Additionally, 5 GSTs (1 upregulated: GSTD6; and 4 downregulated: GSTX1, GSTE6, GSTD11 and GSTS1) were significantly differentially expressed in mosquitoes that survived lambda-cyhalothrin exposure relative to the susceptible strain for both ISA and MAN (Figure 4, Table 2).

3.3.4. Genes associated with resistance to multiple insecticides

A total of 234 genes (43 upregulated and 191 downregulated) were significantly differentially expressed in mosquitoes that survived either malathion, alpha-cypermethrin or lambda cyhalothrin exposure compared to the susceptible laboratory strain (Figure 5A, Table S8). Two of the top 10 upregulated DEGs that overlapped all the R-S groups included a putative serine protease (AAEL029072; FC=1197-6081) and a peritrophin-1-like (AAEL021372; FC=189-1078), reflecting their potential association with the multiple-insecticide resistance phenotypes (Figure 5B). Eight key genes coding for CYP450 enzymes including CYP6Z7, CYP28A5, CYP9J2, CYP6Z6, CYP6BB2, CYP6M9, and two CYP9F2 orthologs were upregulated across the six R-S comparisons. Two additional CYP450s (CYP325N2, CYP6N15) and one GST (GSTS1) were downregulated in all six resistant groups (Figure 5C). A total of 131, 325, and 240 genes were uniquely differentially expressed in malathion, alpha-cypermethrin, lambda cyhalothrin survivors, respectively, suggesting their association with a specific insecticide (Figure 5A). Two detoxification DEGs, a CYP450 (CYP6BY1) and a carboxylesterase (AAEL004022), were unique to malathion resistant samples (Figure 4A). Additionally, one downregulated detoxification gene (GSTD1) was only detected in alpha-cypermethrin-resistant samples. Unique downregulated detoxification genes associated with lambda-cyhalothrin resistance were CYP4H29, CYP307A1, CYP6P12, and GSTE6 (Table 2).

3.4. Gene Ontology enrichment analysis

Gene ontology enrichment analysis (GOEA) was performed on the list of DEGs associated with all the 6 R-S comparisons (Table S7). Not surprisingly, the differences seen in the transcriptomic profiles between the different comparisons was also reflected in the total number of enriched gene ontology (GO) terms associated with DE genes (Table 1, Table S9). GO analysis of the upregulated genes for the R-S comparisons, identified the enrichment of some relevant molecular functions that are strongly associated with multiple-insecticide resistance phenotypes, including “monooxygenase activity” (GO:0004497), “oxidoreductase activity” (GO:0016705), “heme binding” (GO:0020037), and “tetrapyrrole binding”(GO:0046906) (Figure 6). Additionally, the molecular functions “catalytic activity” (GO:0003824) and “iron ion binding” (GO:0005506) were enriched in 4 out of 6 upregulated gene lists analyzed, reflecting their association with multiple-insecticide resistance. Interestingly, the enrichment of some relevant GO terms was found to be specific to malathion resistant mosquitoes, including several GO terms associated with postsynaptic neurotransmitter activities (GO:1904315, GO:0099529, GO:0098960, GO:0030594), transmembrane transporter activities (GO:0022803, GO:002285, GO:0022857), ion channel activities (GO:0005216 , GO:0005261, GO:0015267, GO:0022848, GO:0099094, GO:0005231, GO:0022848) and protein receptor activities (GO:0004930, GO:0005102) (Figure 6). For the list of downregulated genes, several GO terms were overrepresented in 4 out of the 6 R-S comparisons. These included “catalytic activity” (GO:0003824), “carboxypeptidase activity” (GO:0004180), “metallocarboxypeptidase activity” (GO:0004181), “structural constituent of chitin-based cuticle” (GO:0005214), “chitin binding” (GO:0008061) and “exopeptidase activity” (GO:0008238) (Figure S4). GO terms associated with chitin binding activities were overrepresented in the list of downregulated genes of all the six R-S DEG comparisons, which was also clearly reflected in the functional annotation of the DEGs displayed in Figure 2, where genes associated with chitin and chitinase activities were mostly downregulated (Figure S4).

3.5. Validation of relative expression levels estimated by RNA-Seq with quantitative RT-PCR

Five transcripts differentially expressed in resistant samples (I-ACYP, I-LCT and M-LCT) including GSTD6, GSTE4, CYPBB2, CYP6M9 and CYP6N9, were used to validate the RNA-Seq results by quantitative RT-PCR. The threshold cycle (ct) value was not detected for 2 of the 15 qPCR assays, likely associated to the low expression level of these transcripts in the mosquitoes. This is evident in the low read count from the RNA-seq. The qRT-PCR results broadly supported the directionality of the changes in expression levels (75% of the essays with ct values), although for several genes, the magnitude of the expression difference was smaller than estimated by RNA-Seq (Figure S5).

4. Discussion

Insecticide resistance in Ae. aegypti continues to expand globally due the extensive use of insecticides for its control. In this study, we used whole transcriptomic approach to investigate the molecular basis of resistance in Ae. aegypti populations from two locations in Puerto Rico, Isabela and Manatí, that were resistant to an organophosphate (malathion) and two pyrethroids (lambda-cyhalothrin and alpha-cypermethrin).
Differential gene expression analysis of genes associated with resistance to all the three insecticides from the two locations showed similarities in overexpression of the detoxification genes belonging mainly to the cytochrome P450 gene family and downregulation in cuticular protein genes. The insect cuticle is comprised mainly of chitin and cuticular proteins, and modifications to the cuticle can lead to thickening of the cuticle hence slowing the penetration of insecticides (Balabanidou et al. 2018). Cuticular thickening has been associated with resistance to pyrethroids in Anopheles funestus (Wood et al. 2010) and with reduced penetration of deltamethrin in highly resistant strain of Anopheles gambiae (Balabanidou et al. 2016). The downregulation of cuticular proteins in this study suggests that resistance in these populations is predominantly driven by cytochrome P450-mediated detoxification of all three insecticides. A study on the transcriptomic profile of a resistant strain of Ae. aegypti from Vietnam similarly reported downregulation of cuticular proteins in resistant versus susceptible strains (Lien et al. 2019).
In all the resistant groups from all three insecticides, a serine protease (AAEL029072) was significantly upregulated with a high fold change (Figure 5B). Serine protease has been reported to be highly upregulated in An. gambiae resistant to DDT (Vontas et al. 2005) as well as in Culex pipiens pallens resistant to deltamethrin (Gong et al. 2005), where it was suggested to play a crucial role in the innate immunity in this species and other insects (El Moussawi et al. 2019, Wang et al. 2021, Shan et al. 2022). This finding suggests that the increase in protease activity may be important in modulating Ae. aegypti resistance to multiple insecticides, however this remains speculative as this has not yet been functionally proven.
Several CYP450s were also commonly differentially expressed in the three resistant groups. These included CYP6Z7, CYP28A5, CYP9J2, CYP6Z6, CYP6BB2, CYP6M9, and two CYP9F2 orthologs, suggesting that these CYPs may confer cross resistance to pyrethroids and organophosphates. The majority of these CYPs have been associated with resistance to insecticides in previous studies. CYP6BB2 has consistently been reported to be overexpressed in Ae. aegypti resistant to pyrethroids in Asia, Europe, and the Americas (Bariami et al. 2012, Kasai et al. 2014, Reid et al. 2014, Dusfour et al. 2015, Moyes et al. 2017). The overexpression of CYP6BB2 has also shown strong metabolic activity for permethrin in in vitro metabolic studies (Kasai et al. 2014), highlighting its importance as a candidate metabolic resistance marker for pyrethroids in Ae. aegypti populations. In the same study, CYP6Z7 and CYP6Z6 were also overexpressed in permethrin-resistant Ae. aegypti. In another study, CYP6Z7 and CYP6BB2 were among four of the only CYPs upregulated in three tested strains of Ae. aegypti resistant to deltamethrin while CYP9J2 also appeared in one of the resistant strains in the same study (Epelboin et al. 2021). Multiple studies have also reported the upregulation of CYP6Z6 in pyrethroid resistant Ae. aegypti populations around the world (Marcombe et al. 2009, Bariami et al. 2012, Saavedra-Rodriguez et al. 2012, Reid et al. 2014, Dusfour et al. 2015).
Despite some shared gene expression patterns across the three insecticides, some of the genes were only found to be differentially regulated in survivors of specific insecticides and not others. In the malathion resistant samples from both locations, CYP6BY1 was found to be upregulated but was not upregulated in the lambda-cyhalothrin or alpha-cypermethrin resistant samples. CYP6BY1 has previously been reported to be among only two overexpressed P450s in Colombian Ae. aegypti larvae resistant to temephos, which is also an organophosphate insecticide (Morgan et al. 2022). In samples resistant to alpha-cypermethrin, the glutathione S-transferase gene GSTD1 was found to be uniquely overexpressed in this group compared to survivors of the other two insecticides. GSTD1 has been implicated in resistance to DDT in both An. gambiae (Ranson et al. 2001) and D. melanogaster (Tang and Tu 1994). Two genes, CYP4H29 and GSTE6, were uniquely overexpressed in samples resistant to lambda-cyhalothrin. In a targeted RNA-seq study, GSTE6 was the only GST reported to be associated with deltamethrin resistance across eight Ae. aegypti strains from distinct geographical origins (Faucon et al. 2017). CYP4H29 has previously been associated with resistance to pyrethroids in Puerto Rican Ae. aegypti (Reid et al. 2014) The same study reported the overexpression of four CYP450s belonging to family four. Studies on the role of cytochrome P450 4 (CYP4) family in resistance are scarce, however, Reid et al showed increased survival of D. melanogaster expressing CYP4H29 when exposed to permethrin.
Another difference that stands out is the enrichment of GO terms belonging to the postsynaptic neurotransmitter receptor activity mainly in the alpha-cypermethrin resistant samples from both Isabela and Manatí. This supports the evidence that type II pyrethroids such as alpha-cypermethrin have been shown to mainly inhibit GABAA receptor function, hence leading to higher excitation in the insect’s central nervous system (Eldefrawi and Eldefrawi 1989). However, the same GO terms did not appear enriched in the lambda-cyhalothrin resistant group which is also a type II pyrethroid. According to the WHO report 1998 (World Health 1998), alpha-cypermethrin was more toxic than other pyrethroids and was significantly more effective at killing anopheline mosquitoes than permethrin or lambda-cyhalothrin and this may explain the high activity in the enriched GO terms. In addition, the bottle bioassay results showed a higher frequency of resistance to alpha-cypermethrin in the Manatí population than in the Isabela population, potentially also explaining the higher enrichment score in this group.

6. Conclusions

This study investigated the gene expression profiles associated with pyrethroid and organophosphate resistance in Ae. aegypti populations from two locations in Puerto Rico. We identified several detoxification genes which exhibited similar differential over-expression patterns across all the insecticide resistant phenotypes from both populations, suggesting that these genes could be used as expression-based markers for multiple-insecticide resistance in Ae. aegypti. We also identified some detoxification genes that were uniquely overexpressed in mosquitoes that survived exposure to malathion, alpha-cypermethrin, or lambda-cyhalothrin, indicating that their overexpression is more closely associated with resistance to specific insecticides. Additionally, we report the association of significant overexpression of serine protease genes and downregulation of cuticular-related protein genes (chinase and chitin) with organophosphate and pyrethroid resistance in both Ae. aegypti populations. The genes identified in this study, particularly in the cytochrome P450 family, should be functionally validated to confirm their importance as molecular markers for resistance detection.

Supplementary Materials

The following supporting information can be downloaded at the website of this paper posted on Preprints.org. Figure S1: Bar plot of the feature count results. Percentage shows the number of tags successfully mapped to exonic regions (pink) and non exonic region (orange) of the reference genome, as well as those that were not assigned to exons because they were tagged as singleton (green), chimera (blue) and ambiguous (black) alignments. Figure S2: Venn diagrams summarizing the number of differentially expressed genes (DE) and the possible sharing of DE genes in the resistant samples vs susceptible (R-S): control vs susceptible (C-S) and resistant vs control (R-C) for each experiment including Isabela malathion (A), Manatí malathion (B), Isabela alpha-cypermethrin (C), Manatí alpha-cypermethrin (D), Isabela lambda-cyhalothrin (E) and Manatí lambda-cyhalothrin (F). Figure S3: Bar plot showing the percentage of significant up- and downregulated cytochrome P450s (left panel) and the cuticular related protein genes (left panel) in the resistant samples: compared the susceptible strain. Each bar plot is divided into three groups, representing the three insecticides (ACYP, LCT and MAL). No consistent patterns were found between the number of differentially expressed CYP450s, cuticular protein genes, and the survival rate of the mosquitoes to the insecticides. Figure S4: Heatmap of the gene ontology (GO) enrichment analysis of the significantly down-regulated genes in the malathion (I-MAL and M-MAL): alpha-cypermethrin (I-ACYP and -M-CYP) and lambda-cyhalothrin (I-LCT and M-MLCT) resistant samples compared to the susceptible strain of Ae. aegypti. The white-red gradient indicates the enrichment score of each GO term, while the gray color indicates that this GO term was not significantly enriched. Only the core enriched GO-terms (GO terms overlapping across both MAN and ISA) for each insecticide are shown. Figure S5: qPCR validation. Comparison of the log2FC expression of 5 detoxification genes measured by qPCR and RNA-seq for the insecticide resistant samples compared to the susceptible. The colors of the bars indicate the type of the pairwise comparison (insecticide resistant vs susceptible) and the patterns of the bars indicate the gene expression measurement approach (qPCR and RNA-seq). The red and the black asterisks (*) indicate genes for which the expression was not detected for the qPCR and the RNA-seq, respectively. Table S1: qRT-PCR primers used in this study. *Candidate internal reference genes. Table S2: Summary statistics of the RNA-Seq reads, alignment, and tag counting. Table shows the key statistics describing the RNA-Seq Reads before and after quality filtering, the alignment of the filtered reads to the reference genome and the read quantification. Table S3: Functional annotation of AaegL5.1 using Blast2Go and VectorBase. Table S4: Bioassay results. Table S5: Results of differential gene expression analyses from edgeR. Full results of pairwise differential gene expression analyses of RNA-Seq datasets. Table S6: Identifying the overlapped DEGs between R-S , C-S, or R-C comparisons. All overlapping DEGs belonging to different intersections, for all six comparisons, were mapped to their functional description and gene name. Table S7: Identifying overlapping DEGs between 2 R-S comparisons (2 geographical locations) for each insecticide. The core DEGs for each insecticide mapped to their functional description and gene symbol. Table S8: Identifying genes that are associated with multiple-insecticide resistance. Differentially expressed genes that overlapped across all three insecticides (independently of geographical location) identified and mapped to their fold change expression and their functional description. Table S9: Gene ontology enrichment analysis results. Full results of gene ontology enrichment analysis for differentially expressed gene sets of Aedes aegypti (AaegL5.1).

Author Contributions

Conceptualization, LMI, LK and AL; Data curation, DD; Formal analysis, DD; Funding acquisition, AL; Investigation, DD, LMI, EJM, and AL; Methodology, DD, LMI, JM, JK, HRR, LK and AL; Project administration, AL; Resources, LMI and AL; Software, DD and SEM; Supervision, AL; Validation, DD and LMI; Visualization, DD; Writing – original draft, DD, LMI, EJM and AL; Writing – review & editing, DD, LMI, EJM, SEM and AL. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by CDC’s Advanced Molecular Detection (AMD) program. DD was supported by the Association of Public Health Laboratories (APHL) through the Public Health Bioinformatics Fellowship.

Code availability

Custom scripts used for all the analysis are available from the authors on request.

Informed Consent Statement

Not Applicable.

Data Availability Statement

Sequence data generated by this study is available at Sequence Read Archive (SRA) under the Bio Project accession number PRJNA801226.

Acknowledgments

We acknowledge Jonathan Gerhart from Division of Bacterial Diseases, Centers for Disease Control and Prevention, Atlanta, Georgia, for performing bottle bioassays.

Conflicts of Interest

The authors declare no conflict of interest.

Disclaimer

The views expressed in this manuscript are those of the authors and do not necessarily reflect the official policy or position of the Centers for Disease Control and Prevention.

References

  1. Amelia-Yap, Z. H.; Chen, C. D.; Sofian-Azirun, M.; Low, V. L. Pyrethroid resistance in the dengue vector Aedes aegypti in Southeast Asia: present situation and prospects for management. Parasit Vectors 2018, 11, 332. [Google Scholar] [CrossRef] [PubMed]
  2. Andrews, S. FastQC: a quality control tool for high throughput sequence data. Version 0.11. 2. 2010. Available online: http://www.bioinforma tics.babraham.ac.uk/projects/fastqc.
  3. Aponte, A.; Penilla, R. P.; Rodríguez, A. D.; Ocampo, C. B. Mechanisms of pyrethroid resistance in Aedes (Stegomyia) aegypti from Colombia. Acta Tropica 2019, 191, 146–154. [Google Scholar] [CrossRef] [PubMed]
  4. Badolo, A.; Sombié, A.; Pignatelli, P. M.; Sanon, A.; Yaméogo, F.; Wangrawa, D. W.; Sanon, A.; Kanuka, H.; McCall, P. J.; Weetman, D. Insecticide resistance levels and mechanisms in Aedes aegypti populations in and around Ouagadougou, Burkina Faso. PLoS Negl Trop Dis 2019, 13, e0007439. [Google Scholar] [CrossRef]
  5. Balabanidou, V.; Grigoraki, L.; Vontas, J. Insect cuticle: a critical determinant of insecticide resistance. Curr Opin Insect Sci 2018, 27, 68–74. [Google Scholar] [CrossRef] [PubMed]
  6. Balabanidou, V.; Kampouraki, A.; MacLean, M.; Blomquist, G. J.; Tittiger, C.; Juarez, M. P.; Mijailovsky, S. J.; Chalepakis, G.; Anthousi, A.; Lynd, A.; Antoine, S.; Hemingway, J.; Ranson, H.; Lycett, G. J.; Vontas, J. Cytochrome P450 associated with insecticide resistance catalyzes cuticular hydrocarbon production in Anopheles gambiae. Proc Natl Acad Sci U S A 2016, 113, 9268–9273. [Google Scholar] [CrossRef]
  7. Bariami, V.; Jones, C. M.; Poupardin, R.; Vontas, J.; Ranson, H. Gene amplification, ABC transporters and cytochrome P450s: unraveling the molecular basis of pyrethroid resistance in the dengue vector, Aedes aegypti. PLoS Negl Trop Dis 2012, 6, e1692. [Google Scholar] [CrossRef]
  8. Bass, C.; Nikou, D.; Donnelly, M. J.; Williamson, M. S.; Ranson, H.; Ball, A.; Vontas, J.; Field, L. M. Detection of knockdown resistance (kdr) mutations in Anopheles gambiae: a comparison of two new high-throughput assays with existing methods. Malar J 2007, 6, 111. [Google Scholar] [CrossRef]
  9. Braga, I. A.; Lima, J. B.; Sda, S. Soares; Valle, D. Aedes aegypti resistance to temephos during 2001 in several municipalities in the states of Rio de Janeiro, Sergipe, and Alagoas, Brazil. Mem Inst Oswaldo Cruz 2004, 99, 199–203. [Google Scholar] [CrossRef]
  10. Brogdon, W.; Chan, A. Guideline for evaluating insecticide resistance in vectors using the CDC bottle bioassay; CDC Atlanta: USA, 2010. [Google Scholar]
  11. Chareonviriyaphap, T.; Bangs, M. J.; Suwonkerd, W.; Kongmee, M.; Corbel, V.; Ngoen-Klan, R. Review of insecticide resistance and behavioral avoidance of vectors of human diseases in Thailand. Parasit Vectors 2013, 6, 280. [Google Scholar] [CrossRef]
  12. Chen, S.; Zhou, Y.; Chen, Y.; Gu, J. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 2018, 34, i884–i890. [Google Scholar] [CrossRef]
  13. Chen, Y.; McCarthy, D.; Ritchie, M.; Robinson, M.; Smyth, G.; Hall, E. edgeR: differential analysis of sequence read count data User’s Guide. 2020. [Google Scholar]
  14. David, J. P.; Ismail, H. M.; Chandor-Proust, A.; Paine, M. J. Role of cytochrome P450s in insecticide resistance: impact on the control of mosquito-borne diseases and use of insecticides on Earth. Philos Trans R Soc Lond B Biol Sci 2013, 368, 20120429. [Google Scholar] [CrossRef] [PubMed]
  15. David, J. P.; Faucon, F.; Chandor-Proust, A.; Poupardin, R.; Riaz, M. A.; Bonin, A.; Navratil, V.; Reynaud, S. Comparative analysis of response to selection with three insecticides in the dengue mosquito Aedes aegypti using mRNA sequencing. BMC Genomics 2014, 15, 174. [Google Scholar] [CrossRef] [PubMed]
  16. Dusfour, I.; Zorrilla, P.; Guidez, A.; Issaly, J.; Girod, R.; Guillaumot, L.; Robello, C.; Strode, C. Deltamethrin resistancemMechanisms in Aedes aegypti populations from three French overseas territories worldwide. PLoS Negl Trop Dis 2015, 9, e0004226. [Google Scholar] [CrossRef] [PubMed]
  17. El Moussawi, L.; Nakhleh, J.; Kamareddine, L.; Osta, M. A. The mosquito melanization response requires hierarchical activation of non-catalytic clip domain serine protease homologs. PLoS Pathog 2019, 15, e1008194. [Google Scholar] [CrossRef] [PubMed]
  18. Eldefrawi, M.E.; Eldefrawi, A.T. Insecticide Actions on Gaba Receptors and Voltage-Dependent Chloride Channels. 1989. [Google Scholar]
  19. Epelboin, Y.; Wang, L.; Gianetto, Q. Giai; Choumet, V.; Gaborit, P.; Issaly, J.; Guidez, A.; Douche, T.; Chaze, T.; Matondo, M.; Dusfour, I. CYP450 core involvement in multiple resistance strains of Aedes aegypti from French Guiana highlighted by proteomics, molecular and biochemical studies. PLoS One 2021, 16, e0243992. [Google Scholar] [CrossRef] [PubMed]
  20. Faucon, F.; Gaude, T.; Dusfour, I.; Navratil, V.; Corbel, V.; Juntarajumnong, W.; Girod, R.; Poupardin, R.; Boyer, F.; Reynaud, S.; David, J. P. In the hunt for genomic markers of metabolic resistance to pyrethroids in the mosquito Aedes aegypti: An integrated next-generation sequencing approach. PLoS Negl Trop Dis 2017, 11, e0005526. [Google Scholar] [CrossRef] [PubMed]
  21. Faucon, F.; Dusfour, I.; Gaude, T.; Navratil, V.; Boyer, F.; Chandre, F.; Sirisopa, P.; Thanispong, K.; Juntarajumnong, W.; Poupardin, R.; Chareonviriyaphap, T.; Girod, R.; Corbel, V.; Reynaud, S.; David, J. P. Identifying genomic changes associated with insecticide resistance in the dengue mosquito Aedes aegypti by deep targeted sequencing. Genome Res 2015, 25, 1347–1359. [Google Scholar] [CrossRef]
  22. Giraldo-Calderón, G. I.; Emrich, S. J.; MacCallum, R. M.; Maslen, G.; Dialynas, E.; Topalis, P.; Ho, N.; Gesing, S.; Madey, G.; Collins, F. H.; Lawson, D. VectorBase: an updated bioinformatics resource for invertebrate vectors and other organisms related with human diseases. Nucleic Acids Res 2015, 43, D707–D713. [Google Scholar] [CrossRef]
  23. Gong, M.; Shen, B.; Gu, Y.; Tian, H.; Ma, L.; Li, X.; Yang, M.; Hu, Y.; Sun, Y.; Hu, X.; Li, J.; Zhu, C. Serine proteinase over-expression in relation to deltamethrin resistance in Culex pipiens pallens. Arch Biochem Biophys 2005, 438, 53–62. [Google Scholar] [CrossRef]
  24. Hemingway, J.; Hawkes, N. J.; McCarroll, L.; Ranson, H. The molecular basis of insecticide resistance in mosquitoes. Insect Biochem Mol Biol 2004, 34, 653–665. [Google Scholar] [CrossRef]
  25. Hunter, S.; Apweiler, R.; Attwood, T. K.; Bairoch, A.; Bateman, A.; Binns, D.; Bork, P.; Das, U.; Daugherty, L.; Duquenne, L.; Finn, R. D.; Gough, J.; Haft, D.; Hulo, N.; Kahn, D.; Kelly, E.; Laugraud, A.; Letunic, I.; Lonsdale, D.; Lopez, R.; Madera, M.; Maslen, J.; McAnulla, C.; McDowall, J.; Mistry, J.; Mitchell, A.; Mulder, N.; Natale, D.; Orengo, C.; Quinn, A. F.; Selengut, J. D.; Sigrist, C. J. A.; Thimma, M.; Thomas, P. D.; Valentin, F.; Wilson, D.; Wu, C. H.; Yeats, C. InterPro: the integrative protein signature database. Nucleic Acids Research 2008, 37, D211–D215. [Google Scholar] [CrossRef]
  26. Ishak, I. H.; Kamgang, B.; Ibrahim, S. S.; Riveron, J. M.; Irving, H.; Wondji, C. S. Pyrethroid resistance in Malaysian populations of Dengue vector Aedes aegypti is mediated by CYP9 Family of Cytochrome P450 Genes. PLoS Negl Trop Dis 2017, 11, e0005302. [Google Scholar] [CrossRef] [PubMed]
  27. Jones, P.; Binns, D.; Chang, H.-Y.; Fraser, M.; Li, W.; McAnulla, C.; McWilliam, H.; Maslen, J.; Mitchell, A.; Nuka, G.; Pesseat, S.; Quinn, A. F.; Sangrador-Vegas, A.; Scheremetjew, M.; Yong, S.-Y.; Lopez, R.; Hunter, S. InterProScan 5: genome-scale protein function classification. Bioinformatics (Oxford, England) 2014, 30, 1236–1240. [Google Scholar] [CrossRef] [PubMed]
  28. Kasai, S.; Komagata, O.; Itokawa, K.; Shono, T.; Ng, L. C.; Kobayashi, M.; Tomita, T. Mechanisms of pyrethroid resistance in the dengue mosquito vector, Aedes aegypti: target site insensitivity, penetration, and metabolism. PLoS Negl Trop Dis 2014, 8, e2948. [Google Scholar] [CrossRef] [PubMed]
  29. Klopfenstein, D. V.; Zhang, L.; Pedersen, B. S.; Ramírez, F.; Vesztrocy, A. Warwick; Naldi, A.; Mungall, C. J.; Yunes, J. M.; Botvinnik, O.; Weigel, M.; Dampier, W.; Dessimoz, C.; Flick, P.; Tang, H. GOATOOLS: A Python library for Gene Ontology analyses. Scientific Reports 2018, 8, 10872. [Google Scholar] [CrossRef]
  30. Lertkiatmongkol, P.; Pethuan, S.; Jirakanjanakit, N.; Rongnoparut, P. Transcription analysis of differentially expressed genes in insecticide-resistant Aedes aegypti mosquitoes after deltamethrin exposure. J Vector Ecol 2010, 35, 197–203. [Google Scholar] [CrossRef] [PubMed]
  31. Li, H.; Handsaker, B.; Wysoker, A.; Fennell, T.; Ruan, J.; Homer, N.; Marth, G.; Abecasis, G.; Durbin, R.; Subgroup, G. P. D. P. The Sequence Alignment/Map format and SAMtools. Bioinformatics 2009, 25, 2078–2079. [Google Scholar] [CrossRef]
  32. Li, X.; Schuler, M. A.; Berenbaum, M. R. Molecular mechanisms of metabolic resistance to synthetic and natural xenobiotics. Annu Rev Entomol 2007, 52, 231–253. [Google Scholar] [CrossRef]
  33. Liao, Y.; Smyth, G. K.; Shi, W. The Subread aligner: fast, accurate and scalable read mapping by seed-and-vote. Nucleic acids research 2013, 41, e108. [Google Scholar] [CrossRef]
  34. Lien, N. T. K.; Ngoc, N. T. H.; Lan, N. N.; Hien, N. T.; Tung, N. V.; Ngan, N. T. T.; Hoang, N. H.; Binh, N. T. H. Transcriptome sequencing and analysis of changes associated with Insecticide resistance in the dengue mosquito (Aedes aegypti) in Vietnam. Am J Trop Med Hyg 2019, 100, 1240–1248. [Google Scholar] [CrossRef]
  35. Lima, E. P.; Paiva, M. H.; de Araújo, A. P.; da Silva, E. V.; da Silva, U. M.; de Oliveira, L. N.; Santana, A. E.; Barbosa, C. N.; de Paiva Neto, C. C.; Goulart, M. O.; Wilding, C. S.; Ayres, C. F.; de Melo Santos, M. A. Insecticide resistance in Aedes aegypti populations from Ceará, Brazil. Parasit Vectors 2011, 4, 5. [Google Scholar] [CrossRef] [PubMed]
  36. Liu, N. Insecticide resistance in mosquitoes: impact, mechanisms and research directions. Annu Rev Entomol 2015, 60, 537–559. [Google Scholar] [CrossRef] [PubMed]
  37. Marcombe, S.; Poupardin, R.; Darriet, F.; Reynaud, S.; Bonnet, J.; Strode, C.; Brengues, C.; Yébakima, A.; Ranson, H.; Corbel, V.; David, J. P. Exploring the molecular basis of insecticide resistance in the dengue vector Aedes aegypti: a case study in Martinique Island (French West Indies). BMC Genomics 2009, 10, 494. [Google Scholar] [CrossRef] [PubMed]
  38. Mathias, D. K.; Ochomo, E.; Atieli, F.; Ombok, M.; Bayoh, M. N.; Olang, G.; Muhia, D.; Kamau, L.; Vulule, J. M.; Hamel, M. J.; Hawley, W. A.; Walker, E. D.; Gimnig, J. E. Spatial and temporal variation in the kdr allele L1014S in Anopheles gambiae s.s. and phenotypic variability in susceptibility to insecticides in Western Kenya. Malar J 2011, 10, 10. [Google Scholar] [CrossRef]
  39. Matthews, B. J.; Dudchenko, O.; Kingan, S. B.; Koren, S.; Antoshechkin, I.; Crawford, J. E.; Glassford, W. J.; Herre, M.; Redmond, S. N.; Rose, N. H.; Weedall, G. D.; Wu, Y.; Batra, S. S.; Brito-Sierra, C. A.; Buckingham, S. D.; Campbell, C. L.; Chan, S.; Cox, E.; Evans, B. R.; Fansiri, T.; Filipović, I.; Fontaine, A.; Gloria-Soria, A.; Hall, R.; Joardar, V. S.; Jones, A. K.; Kay, R. G. G.; Kodali, V. K.; Lee, J.; Lycett, G. J.; Mitchell, S. N.; Muehling, J.; Murphy, M. R.; Omer, A. D.; Partridge, F. A.; Peluso, P.; Aiden, A. P.; Ramasamy, V.; Rašić, G.; Roy, S.; Saavedra-Rodriguez, K.; Sharan, S.; Sharma, A.; Smith, M. L.; Turner, J.; Weakley, A. M.; Zhao, Z.; Akbari, O. S.; Black, W. C.; Cao, H.; Darby, A. C.; Hill, C. A.; Johnston, J. S.; Murphy, T. D.; Raikhel, A. S.; Sattelle, D. B.; Sharakhov, I. V.; White, B. J.; Zhao, L.; Aiden, E. L.; Mann, R. S.; Lambrechts, L.; Powell, J. R.; Sharakhova, M. V.; Tu, Z.; Robertson, H. M.; McBride, C. S.; Hastie, A. R.; Korlach, J.; Neafsey, D. E.; Phillippy, A. M.; Vosshall, L. B. Improved reference genome of Aedes aegypti informs arbovirus vector control. Nature 2018, 563, 501–507. [Google Scholar] [CrossRef]
  40. Morgan, J.; Salcedo-Sora, J. E.; Triana-Chavez, O.; Strode, C. Expansive and diverse phenotypic landscape of field Aedes aegypti (Diptera: Culicidae) larvae with differential susceptibility to Temephos: Beyond metabolic detoxification. J Med Entomol 2022, 59, 192–212. [Google Scholar] [CrossRef]
  41. Moyes, C. L.; Vontas, J.; Martins, A. J.; Ng, L. C.; Koou, S. Y.; Dusfour, I.; Raghavendra, K.; Pinto, J.; Corbel, V.; David, J. P.; Weetman, D. Contemporary status of insecticide resistance in the major Aedes vectors of arboviruses infecting humans. PLoS Negl Trop Dis 2017, 11, e0005625. [Google Scholar] [CrossRef]
  42. Ranson, H.; Lissenden, N. Insecticide Resistance in African Anopheles Mosquitoes: A Worsening Situation that Needs Urgent Action to Maintain Malaria Control. Trends Parasitol 2016, 32, 187–196. [Google Scholar] [CrossRef]
  43. Ranson, H.; Rossiter, L.; Ortelli, F.; Jensen, B.; Wang, X.; Roth, C. W.; Collins, F. H.; Hemingway, J. Identification of a novel class of insect glutathione S-transferases involved in resistance to DDT in the malaria vector Anopheles gambiae. Biochem J 2001, 359, 295–304. [Google Scholar] [CrossRef]
  44. Reid, W. R.; Thornton, A.; Pridgeon, J. W.; Becnel, J. J.; Tang, F.; Estep, A.; Clark, G. G.; Allan, S.; Liu, N. Transcriptional analysis of four family 4 P450s in a Puerto Rico strain of Aedes aegypti (Diptera: Culicidae) compared with an Orlando strain and their possible functional roles in permethrin resistance. J Med Entomol 2014, 51, 605–615. [Google Scholar] [CrossRef]
  45. Riveron, J. M.; Ibrahim, S. S.; Mulamba, C.; Djouaka, R.; Irving, H.; Wondji, M. J.; Ishak, I. H.; Wondji, C. S. Genome-Wide Transcription and Functional Analyses Reveal Heterogeneous Molecular Mechanisms Driving Pyrethroids Resistance in the Major Malaria Vector Anopheles funestus Across Africa. G3 (Bethesda) 2017, 7, 1819–1832. [Google Scholar] [CrossRef] [PubMed]
  46. Robinson, M. D.; McCarthy, D. J.; Smyth, G. K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics (Oxford, England) 2010, 26, 139–140. [Google Scholar] [CrossRef] [PubMed]
  47. Saavedra-Rodriguez, K.; Suarez, A. F.; Salas, I. F.; Strode, C.; Ranson, H.; Hemingway, J.; Black, W. C. t. Transcription of detoxification genes after permethrin selection in the mosquito Aedes aegypti. Insect Mol Biol 2012, 21, 61–77. [Google Scholar] [CrossRef] [PubMed]
  48. Saavedra-Rodriguez, K.; Campbell, C. L.; Lozano, S.; Penilla-Navarro, P.; Lopez-Solis, A.; Solis-Santoyo, F.; Rodriguez, A. D.; Perera, R.; Iv, W. C. Black. Permethrin resistance in Aedes aegypti: Genomic variants that confer knockdown resistance, recovery, and death. PLoS Genet 2021, 17, e1009606. [Google Scholar] [CrossRef] [PubMed]
  49. Samal, R. R.; Panmei, K.; Lanbiliu, P.; Kumar, S. Metabolic detoxification and ace-1 target site mutations associated with acetamiprid resistance in Aedes aegypti L. Front Physiol 2022, 13, 988907. [Google Scholar] [CrossRef]
  50. Schluep, S.M.; Buckner, E.A. Metabolic resistance in permethrin-resistant Florida Aedes aegypti (Diptera: Culicidae). Insects 2021, 12. [Google Scholar] [CrossRef]
  51. Sene, N. M.; Mavridis, K.; Ndiaye, E. H.; Diagne, C. T.; Gaye, A.; Ngom, E. H. M.; Ba, Y.; Diallo, D.; Vontas, J.; Dia, I.; Diallo, M. Insecticide resistance status and mechanisms in Aedes aegypti populations from Senegal. PLoS Negl Trop Dis 2021, 15, e0009393. [Google Scholar] [CrossRef]
  52. Shan, T.; Wang, Y.; Dittmer, N. T.; Kanost, M. R.; Jiang, H. Serine Protease Networks Mediate Immune Responses in Extra-Embryonic Tissues of Eggs in the Tobacco Hornworm, Manduca sexta. J Innate Immun 2022, 1–15. [Google Scholar] [CrossRef]
  53. Simma, E. A.; Dermauw, W.; Balabanidou, V.; Snoeck, S.; Bryon, A.; Clark, R. M.; Yewhalaw, D.; Vontas, J.; Duchateau, L.; Leeuwen, T. Van. Genome-wide gene expression profiling reveals that cuticle alterations and P450 detoxification are associated with deltamethrin and DDT resistance in Anopheles arabiensis populations from Ethiopia. Pest Manag Sci 2019, 75, 1808–1818. [Google Scholar] [CrossRef]
  54. Strode, C.; Wondji, C. S.; David, J.-P.; Hawkes, N. J.; Lumjuan, N.; Nelson, D. R.; Drane, D. R.; Karunaratne, S. H. P. P.; Hemingway, J.; Black, W. C.; Ranson, H. Genomic analysis of detoxification genes in the mosquito Aedes aegypti. Insect Biochemistry and Molecular Biology 2008, 38, 113–123. [Google Scholar] [CrossRef]
  55. Tang, A. H.; Tu, C. P. Biochemical characterization of Drosophila glutathione S-transferases D1 and D21. J Biol Chem 1994, 269, 27876–27884. [Google Scholar] [CrossRef]
  56. Vontas, J.; Kioulos, E.; Pavlidi, N.; Morou, E.; Torre, A. della; Ranson, H. Insecticide resistance in the major dengue vectors Aedes albopictus and Aedes aegypti. Pesticide Biochemistry and Physiology 2012, 104, 126–131. [Google Scholar] [CrossRef]
  57. Vontas, J.; Blass, C.; Koutsos, A. C.; David, J. P.; Kafatos, F. C.; Louis, C.; Hemingway, J.; Christophides, G. K.; Ranson, H. Gene expression in insecticide resistant and susceptible Anopheles gambiae strains constitutively or after insecticide exposure. Insect Mol Biol 2005, 14, 509–521. [Google Scholar] [CrossRef]
  58. Wang, H.-C.; Wang, Q.-H.; Bhowmick, B.; Li, Y.-X.; Han, Q. Functional characterization of two clip domain serine proteases in innate immune responses of Aedes aegypti. Parasites & Vectors 2021, 14, 584. [Google Scholar] [CrossRef]
  59. Weill, M.; Malcolm, C.; Chandre, F.; Mogensen, K.; Berthomieu, A.; Marquine, M.; Raymond, M. The unique mutation in ace-1 giving high insecticide resistance is easily detectable in mosquitoes. Insect Mol Biol 2004, 13, 1–7. [Google Scholar] [CrossRef] [PubMed]
  60. Wood, O.; Hanrahan, S.; Coetzee, M.; Koekemoer, L.; Brooke, B. Cuticle thickening associated with pyrethroid resistance in the major malaria vector Anopheles funestus. Parasit Vectors 2010, 3, 67. [Google Scholar] [CrossRef] [PubMed]
  61. World Health Organization. The World health report: 1998: Life in the 21st century: a vision for all: report of the Director-General. World Health Organization, Geneva. 1998. [Google Scholar]
Preprints 79101 g001
Figure 1. Bioassay results for alpha-cypermethrin (ACYP), lambda-cyhalothrin (LCT) and malathion (MAL) insecticides for Aedes aegypti from Isabela and Manatí in Puerto Rico. Bars show the percent mortality after 30 minutes of insecticide exposure.
Figure 1. Bioassay results for alpha-cypermethrin (ACYP), lambda-cyhalothrin (LCT) and malathion (MAL) insecticides for Aedes aegypti from Isabela and Manatí in Puerto Rico. Bars show the percent mortality after 30 minutes of insecticide exposure.
Preprints 79101 g001
Preprints 79101 g002
Figure 2. Volcano plots of gene expression profiles of mosquitoes that survived insecticide exposure compared to the susceptible laboratory strain (ROCK). The level of gene expression is plotted in the x-axis as the log2 fold-change, and the significance the statistical test is plotted in the y-axis as -log10 of the corrected p-value (-log10FDR values greater than 50 are displayed as 30). The panels include gene expression profiles compared to ROCK for Isabela malathion (I-MAL), Manatí malathion (M-MAL), Isabela alpha-cypermethrin (I-ACYP), Manatí alpha-cypermethrin (M-ACYP), Isabela lambda-cyhalothrin (I-LCT), and Manatí lambda-cyhalothrin (M-LCT). Detoxification gene families are denoted in red (COE: carboxylesterases), blue (CYP: cytochrome P450s) and black (GST: glutathione-S-transferases). Cuticular proteins (CP) are shown in green, salivary gland proteins (SGP) in cyan, odorant binding proteins (OBP) in violet, and odorant receptor genes in orange. In each plot, genes overexpressed in the resistant groups have a log2FC >0. The vertical dotted lines indicate 2-fold change expression differences, and the horizontal lines indicate an FDR of 0.01.
Figure 2. Volcano plots of gene expression profiles of mosquitoes that survived insecticide exposure compared to the susceptible laboratory strain (ROCK). The level of gene expression is plotted in the x-axis as the log2 fold-change, and the significance the statistical test is plotted in the y-axis as -log10 of the corrected p-value (-log10FDR values greater than 50 are displayed as 30). The panels include gene expression profiles compared to ROCK for Isabela malathion (I-MAL), Manatí malathion (M-MAL), Isabela alpha-cypermethrin (I-ACYP), Manatí alpha-cypermethrin (M-ACYP), Isabela lambda-cyhalothrin (I-LCT), and Manatí lambda-cyhalothrin (M-LCT). Detoxification gene families are denoted in red (COE: carboxylesterases), blue (CYP: cytochrome P450s) and black (GST: glutathione-S-transferases). Cuticular proteins (CP) are shown in green, salivary gland proteins (SGP) in cyan, odorant binding proteins (OBP) in violet, and odorant receptor genes in orange. In each plot, genes overexpressed in the resistant groups have a log2FC >0. The vertical dotted lines indicate 2-fold change expression differences, and the horizontal lines indicate an FDR of 0.01.
Preprints 79101 g002
Preprints 79101 g003
Figure 3. Venn diagram of differentially expressed genes (log2FC >1 and FDR <0.01) in insecticide resistant mosquitoes compared to the susceptible laboratory colony. A) Malathion resistant samples from Manatí and Isabela; B) Alpha-cypermethrin resistant samples from Manatí and Isabela; and C) Lambda-cyhalothrin resistant samples from Manatí and Isabela. The overlapping areas represent genes that that are differentially expressed across both sites. Upregulated (up) genes are in red and downregulated (down) genes are in blue.
Figure 3. Venn diagram of differentially expressed genes (log2FC >1 and FDR <0.01) in insecticide resistant mosquitoes compared to the susceptible laboratory colony. A) Malathion resistant samples from Manatí and Isabela; B) Alpha-cypermethrin resistant samples from Manatí and Isabela; and C) Lambda-cyhalothrin resistant samples from Manatí and Isabela. The overlapping areas represent genes that that are differentially expressed across both sites. Upregulated (up) genes are in red and downregulated (down) genes are in blue.
Preprints 79101 g003
Preprints 79101 g004
Figure 4. Heatmaps summarizing the log2 fold-change expression of detoxifications genes including cytochrome P450s (A), glutathione-S transferases (B) and carboxylesterases (C). The log2FC expression is shown in white-red scale (red =overexpressed, white= under-expressed, and gray= not detected or not differentially expressed) for the comparisons of insecticide resistant samples vs. susceptible lab strain (R-S).
Figure 4. Heatmaps summarizing the log2 fold-change expression of detoxifications genes including cytochrome P450s (A), glutathione-S transferases (B) and carboxylesterases (C). The log2FC expression is shown in white-red scale (red =overexpressed, white= under-expressed, and gray= not detected or not differentially expressed) for the comparisons of insecticide resistant samples vs. susceptible lab strain (R-S).
Preprints 79101 g004
Preprints 79101 g005
Figure 5. Identification of DEGs associated with multiple-insecticide resistance. A) The upset plot represents the intersection of differentially expressed genes from Isabela (I) and Manatí (M) between the malathion (MAL), alpha-cypermethrin (ACYP), lambda-cyhalothrin (LCT) resistant samples when compared to the insecticide susceptible strain (R-S). The left horizontal bar plot reports the total number of DEGs in each comparison (set size), the circles represent the set of comparisons associated to the intersection, while the vertical bar plot reports the number of unique and overlapping DEGs (intersection size) between the different combinations of R-S comparisons. The bar plots represent core DEGs (green), DEGs specific to MAL (blue), ACYP (red) and LCT (orange). B) Heatmap showing the log2 fold change (log2FC) expression of the top 10 DEGs in any R-S comparison (core top 10). C) Heatmap showing the log2FC expression of the detoxification genes shared across all R-S comparisons. The heatmaps are in a blue-red color gradient (red=over-expressed and blue=under-expressed).
Figure 5. Identification of DEGs associated with multiple-insecticide resistance. A) The upset plot represents the intersection of differentially expressed genes from Isabela (I) and Manatí (M) between the malathion (MAL), alpha-cypermethrin (ACYP), lambda-cyhalothrin (LCT) resistant samples when compared to the insecticide susceptible strain (R-S). The left horizontal bar plot reports the total number of DEGs in each comparison (set size), the circles represent the set of comparisons associated to the intersection, while the vertical bar plot reports the number of unique and overlapping DEGs (intersection size) between the different combinations of R-S comparisons. The bar plots represent core DEGs (green), DEGs specific to MAL (blue), ACYP (red) and LCT (orange). B) Heatmap showing the log2 fold change (log2FC) expression of the top 10 DEGs in any R-S comparison (core top 10). C) Heatmap showing the log2FC expression of the detoxification genes shared across all R-S comparisons. The heatmaps are in a blue-red color gradient (red=over-expressed and blue=under-expressed).
Preprints 79101 g005
Preprints 79101 g006
Figure 6. Heatmap of the gene ontology (GO) enrichment analysis of the significantly upregulated genes in the malathion (I-MAL and M-MAL), alpha-cypermethrin (I-ACYP and-M-ACYP) and lambda-cyhalothrin (I-LCT and M-LCT) resistant samples compared to the susceptible strain of Ae. aegypti. The white-red gradient indicates the enrichment score of each GO term, while the gray color indicates that the GO term was not significantly enriched. Only the GO terms that overlapped for both MAN and ISA for specific insecticides are shown. The complete GOEA report is shown in Table S9.
Figure 6. Heatmap of the gene ontology (GO) enrichment analysis of the significantly upregulated genes in the malathion (I-MAL and M-MAL), alpha-cypermethrin (I-ACYP and-M-ACYP) and lambda-cyhalothrin (I-LCT and M-LCT) resistant samples compared to the susceptible strain of Ae. aegypti. The white-red gradient indicates the enrichment score of each GO term, while the gray color indicates that the GO term was not significantly enriched. Only the GO terms that overlapped for both MAN and ISA for specific insecticides are shown. The complete GOEA report is shown in Table S9.
Preprints 79101 g006
Table 1. Summary of the total number of differentially expressed genes detected for malathion, alpha-cypermethrin, and lambda-cyhalothrin comparisons: R-S (Resistant vs Susceptible), C-S (Unexposed Control vs Susceptible) or R-C (Resistant vs Unexposed Control).
Table 1. Summary of the total number of differentially expressed genes detected for malathion, alpha-cypermethrin, and lambda-cyhalothrin comparisons: R-S (Resistant vs Susceptible), C-S (Unexposed Control vs Susceptible) or R-C (Resistant vs Unexposed Control).
# Of genes tested DE genes(|log2FC|>1 & FDR<0.05) DE genes(|log2FC|>1 & FDR<0.01)
Insecticide Sites Comparisons Up Down Up Down
Malathion Isabela I-MAL vs Rock (R-S) 10153 699 939 441 691
I-MAL vs -IU1 (R-C) 10161 120 186 33 86
I-U1 vs Rock (C-S) 10088 337 403 256 307
Manatí M-MAL vs Rock (R-S) 10573 961 792 739 657
M-MAL vs M-U (R-C) 10885 322 98 174 41
M-U vs Rock (C-S) 10039 333 634 197 456
Alpha-cypermethrin Isabela I-ACYP vs Rock (R-S) 9825 711 625 642 529
I-U vs Rock (R-C) 9711 230 426 158 322
I-ACYP vs I-U (R-C) 9308 127 426 87 346
Manatí M-ACYP vs Rock (R-S) 9599 686 692 605 570
M-ACYP vs M-U (R-C) 7835 237 561 154 448
M-U vs Rock (C-S) 9561 301 758 222 589
Lambda-cyhalothrin Isabela I-LCT vs R vs Rock (R-S) 9916 429 1008 341 920
I-U vs Rock (C-S) 9711 230 426 158 322
I-LCT vs I-U (R-C) 9550 51 366 30 256
Manatí M-LCT vs Rock (R-S) 9710 461 1163 362 1051
M-U2 vs Rock (C-S) 9561 301 758 222 589
M-LCT vs M-U (R-C) 8618 126 271 67 197
Table 2. Fold change expression of significantly differentially expressed genes of interest (|logFC| >=1, FDR <=0.01) in the comparison of Resistant vs Susceptible (R-S) and Unexposed Control vs Susceptible (C-S) groups for malathion (MAL), alpha-cypermethrin (ACYP) and lambda-cyhalothrin (LCT) for the two locations (I: Isabela, M: Manatí). U1 and U2 denote the Unexposed (Control) groups used for the MAL and pyrethroid (ACYP and LCT) comparisons, respectively.
Table 2. Fold change expression of significantly differentially expressed genes of interest (|logFC| >=1, FDR <=0.01) in the comparison of Resistant vs Susceptible (R-S) and Unexposed Control vs Susceptible (C-S) groups for malathion (MAL), alpha-cypermethrin (ACYP) and lambda-cyhalothrin (LCT) for the two locations (I: Isabela, M: Manatí). U1 and U2 denote the Unexposed (Control) groups used for the MAL and pyrethroid (ACYP and LCT) comparisons, respectively.
Gene ID Description I-MAL vs Rock (R-S) I-U1 vs Rock (C-S) M-MAL vs Rock (R-S) M-U1 vs Rock (C-S) I-ACYP vs Rock (R-S) I-LCT-R vs Rock (R-S) I-U2 vs Rock (C-S) M-ACYP vs Rock (R-S) M-LCT vs Rock (R-S) M-U2 vs Rock (C-S)
Cytochrome P450s monooxygenases
AAEL009130 CYP6Z7 17.19 17.45 36.13 14.77 15.5 13.63 14.92 19.44 14.69
AAEL006805 CYP9J2 9.25 8.16 8.8 11.56 10.58 8.29 7.15 23.33 8.87 8.96
AAEL009121 CYP6N9 5.57 3.17 8.42 4.5 5.16 7.28
AAEL006044 CYP325Q1 5.18 3.23 9.48 3.23 2.77 3.55
AAEL014893 CYP6BB2 4.27 7.01 4.37 6.73 5.48 5.34 6.01 4.71 3.37 3.75
AAEL009123 CYP6Z6 4.22 7.4 2.81 2.31 5.14 3.69
AAEL010154 CYP4AR2 3.62 2.84 2.54 2.24 2.35
AAEL014605 CYP9J9 3.43 2.47 2.91 3.32 2.57 2.2 6.79 2.62
AAEL017297 CYP6M9 3.35 6.91 4.36 4.67 4.84 2.68 5.39 4.15 3.41 3.38
AAEL009127 CYP6M11 3.34 5.13 3.69 2.23 2.96
AAEL012772 CYP325G3 3.19 7.72
AAEL009124 CYP6N12 3.18 5.06 4.71 4.18 5.56 2.86 3.91 3.45
AAEL014619 CYP9J22 3.13 2.44 3.09 5.27 2.39 2.27 8.75
AAEL014615 CYP9J23 3.08 3.31 9.06 0.03 0.3
AAEL014019 CYP4J16 2.86 2.54
AAEL017136 CYP325V1 2.79 2.95 6.14 2.94 2.68
AAEL007808 CYP4D39 2.38 2.91
AAEL004054 CYP4G36 2.37 2.13 2.54 2.18 2.57
AAEL009125 CYP6M10 0.46 0.4
AAEL009132 CYP6Y3 0.39 0.3 0.41 0.34 0.45 0.38 0.38 0.34
AAEL003748 CYP9AE1 0.38
AAEL012762 CYP325N2 0.28 0.29 0.22 0.31 0.3 0.34 0.1 0.22
AAEL009762 CYP307A1 2.73 2.71 2.38
AAEL014609 CYP9J26 2.73 3.04 2.23 2.77 2.02
AAEL004941 CYP6AK1 0.38 0.47
AAEL007010 CYP6AG4 0.34
AAEL009120 CYP6S3 0.28 0.31 0.33 0.27 0.22
AAEL007812 CYP4H32 0.12 0.11 0.05 0.08 0.1 0.02 0.26 0.12
AAEL006989 CYP6AG7 3.24 3.24
AAEL009126 CYP6N6 3.49 4.87 4.3 3.65
AAEL010151 CYP6N16 2.31
AAEL014891 CYP6P12 0.44 0.45
AAEL007830 CYP4H29 0.32 0.38
AAEL009131 CYP6Z8 0.3
AAEL006815 CYP9J16 3.58
AAEL009129 CYP6Z9 2.13
AAEL002633 CYP9J31 2.27
AAEL007024 CYP6AG3 2.88
AAEL026582 CYP6AA5 2.17 2.08
AAEL001960 CYP12F5 0.37 0.42 0.49 0.4 0.38 0.47
AAEL014604 CYP9f2* 10.73 12.18 7.65 14.15 16.08 8.03 12.05 36.68
AAEL009018 CYP6d5* 6.24 3.07 12.91 4.33 4.81 4.94 3.56 7.26 6.66 5.05
AAEL019504 CYP9f2* 5.12 4.14 5.81 6.98 4.16 3.89 4.04 3.43
AAEL014614 CYP9f2* 3.81 4 3.96 5.69 3.72 2.54 3.05 9.18 2.92 4.28
AAEL028635 CYP9f2* 3.16 3.98 4.18 5.02 5.34 2.3 3.15 9.43 2.45 4.78
AAEL003890 CYP28a5* 2.39 2.82 2.67
AAEL009122 CYP6a14* 0.02 0.07 0.04 0.09 0.08 0.1 0.01 0.01 0.19 0.11
AAEL021861 CYP28a5* 2.89 2.4 5.03 3.85 3.02 2.37 2.32
AAEL014830 CYP4ac1* 0.47 0.29 0.37 0.19 0.33
AAEL000340 CYP4C1-like 0.24 0.15 0.06 0.19
AAEL019603 CYP9f2* 2.05
AAEL014924 CYP6d5* 0.21
AAEL017539 CYP6BY1* 17.22 12.48
Glutathione S-transferases
AAEL010591 GSTD6 7.03 5.58 4.81 3.98 8.87 5.23 10.56
AAEL010582 GSTD11 0.43 0.33 0.41 0.37 0.41 0.3 0.38
AAEL000092 GSTX1 0.4 0.42 0.42 0.39 0.39
AAEL007962 GSTE4 0.37 0.27 0.23 0.26 0.31 0.24 0.24 0.21
AAEL011741 GSTS1 0.35 0.44 0.37 0.49 0.13 0.26 0.1 0.29
AAEL007946 GSTE6 0.24 0.41 0.28 0.18
AAEL001061 GSTD1 0.44 0.48
AAEL001054 GSTD4 5.97 8.34
Carboxylesterases
AAEL004022 CES5A 3.19 4.07
AAEL005199 CES6-like 0.36 0.33 0.3 0.48 0.24
AAEL005200 CES6-like 0.48
AAEL012886 CES-6 0.41 0.43
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
Copyright: This open access article is published under a Creative Commons CC BY 4.0 license, which permit the free download, distribution, and reuse, provided that the author and preprint are cited in any reuse.
Prerpints.org logo

Preprints.org is a free preprint server supported by MDPI in Basel, Switzerland.

facebook logotwitter logolinkedin logo
MDPI Initiatives
Important Links
Subscribe

Disclaimer

Terms of Use

Privacy Policy

© 2024 MDPI (Basel, Switzerland) unless otherwise stated