Investigation of DDT resistance mechanisms in Anopheles funestus populations from northern and southern Benin reveals a key role of the GSTe2 gene

Understanding the molecular basis of insecticide resistance in mosquito, such as Anopheles funestus, is an important step in developing strategies to mitigate the resistance problem. This study aims to assess the role of the GSTe2 gene in DDT resistance and determine the genetic diversity of this gene in An. funestus. Gene expression analysis was performed using microarrays and PCR while the potential mutation associated with resistance was determined using sequencing. Low expression level of GSTe2 gene was recorded in Burkina-Faso samples with a fold change of 3.3 while high expression (FC 35.6) was recorded in southern Benin in Pahou (FC 35.6) and Kpome (FC 13.3). The sequencing of GSTe2 gene in six localities showed that L119F-GSTe2 mutation is almost getting fixed in highly DDT-resistant Benin (Pahou, Kpome, Doukonta) and Nigeria (Akaka Remo) mosquitoes with a low mutation rate observed in Tanongou (Benin) and Burkina-Faso mosquitoes. This study shows the key role of the GSTe2 gene in DDT resistant An. funestus in Benin. Polymorphism analysis of this gene across Benin revealed possible barriers to gene flow, which could impact the design and implementation of resistance management strategies in the country.


Background
Malaria remains the most severe infectious disease and a major public health challenge in sub-Saharan Africa [1]. The mortality and the loss of productivity due to the illness, has devastating effects on cognitive development in children surviving the disease, leaving many disabled for life [2]. Since the discovery of the connection between Anopheles vectors and malarial transmission in 1897, vector control strategies have been the most widely used malarial control measures [3]. These measures (based on insecticide use) are insecticide treated bed-nets (ITN) and indoor residual spraying (IRS), both of which have been shown to be effective for reducing malaria prevalence in Africa [4]. One of the insecticides of choice for IRS is DDT (dichloro-diphenyl-trichloro-ethane) because of its high insecticidal activity, low acute mammalian toxicity, wide spectrum use, low price, and long duration of activity.
The availability of dichlorodiphenyltrichloroethane (DDT) and other insecticides in the 1940s marked a new era for malarial control in the world. The effectiveness of DDT against indoor resting mosquitoes led to the adoption of the Global Eradication Programme of Malaria in 1955, coordinated and supported by the World Health Organization (WHO). Although the use of DDT raises concerns of potential harm to the environment and human health, mainly because of the persistent and bioaccumulative nature of DDT and its potential to magnify through the food chain, it continued to be used for pest control, for which exemptions were granted by the federal government and it is still available for public health use today [5].
Inevitably, the major malaria vectors, Anopheles gambiae and Anopheles funestus, have developed resistance to this insecticide. The basic mechanisms underlying insecticide resistance include insecticide target-site mutations, and increased metabolic detoxification of the insecticide through overproduction or elevated enzymatic activity [6]. Three enzyme families are primarily involved in insecticide detoxification: the carboxylesterases (COEs), glutathione-S-transferases (GSTs) and cytochrome P450s (P450s). DDT resistance in An. gambiae can be due either to a specific detoxification mechanism (glutathione-Stransferase) or to a nerve insensitivity resulting from a modification of the target site (sodium channel). The latter, governed by the kdr gene, reduces both the knockdown and lethal effects of DDT [7]. In West Africa, it induces a cross-resistance to pyrethroids, which also depends on kdr mutation [7,8]. In contrary, no kdr mutation has been detected in An. funestus so far [9][10][11]. Indeed, a single amino acid change in the binding pocket of the glutathione-S-transferase epsilon 2 (GSTe2) gene, coupled with increased transcription of this gene, confers a high level of DDT resistance and also cross-resistance to pyrethroids in An. funestus. Furthermore, analysis of GSTe2 polymorphism established that the L119F-Gste2 mutation is tightly associated with metabolic resistance to DDT and its geographical distribution strongly correlates with DDT resistance patterns across Africa [12]. Nevertheless, the strong contrast in the allele frequencies of the L119F-GSTe2 frequencies despite the similar resistance profile recorded in An. funestus populations from two localities in Ghana [13] suggest that possible barriers to gene flow could exist between populations of the same country. Such differences in the underlying resistance mechanisms should be taken into account when designing suitable insecticide resistance management strategies. In southern Benin (Kpome and Pahou), An. funestus was found to be highly resistant to DDT [14] [15] while the population from Tanongou was moderately resistant to DDT with 90% mortality [16]. Also, as GSTe2 gene has been associated with DDT resistance patterns across Africa, this study aimed to investigate the role of the GSTe2 gene in DDT resistance across Benin to fill the knowledge gap by checking if this resistance is driven by the same mechanism.

Samples description
In this study, mosquitoes from the previously published research results were used to further describe the molecular basis of DDT resistance in An. funestus population from different localities [12,[14][15][16]. Mosquito samples generated from the previous investigation were used for genetic analysis in this work.

Study area and mosquito collection
Adult anopheles mosquitoes were collected from three (3) locations in Benin: Kpome in South-Est (6° 55′ N, 2° 19′E), Pahou (6° 23′ N, 2° 13′E) in South-West and Tanongou in North West (10° 48′ N, 1° 26′ E). Mosquitoes were also collected in South-West Burkina Faso (11° 23′ N, 4° 24′ E) (Fig. 1). The selected sites are located in close proximity with rivers, swamps as these permanent water bodies are suitable breeding sites for An. funestus. After obtaining consent from village chiefs and house owners, indoor resting mosquitoes were collected from December 2013 to March 2014 inside households using electric aspirator. Blood fed mosquitoes collected were kept in cups until fully gravid before being subjected to the forced-egg laying technique [17]. The eggs obtained were pooled and reared in a mineral water. Larvae were reared under standard insectary conditions (26 ± 2 °C with a relative humidity of 80%) and were fed daily with Tetramin ™ baby fish food. The water of each larvae bowl was changed every two days to reduce mortality. F1 adult generated were pooled in cages for subsequent analyses.

Microarrays
A custom microarray chip containing 44,000 probes (4 × 44 k) [18] was used to identify the set of genes associated with DDT resistance in Pahou and Burkina-Faso. The 8 × 60 k (60 mer) Agilent An. funestus chip was used to screen for the genes involved in resistance of An. funestus from Kpome. This Agilent microarray chip was designed using the eArray program (Agilent, Santa Clara, CA, USA) (A-MEXP-2374) by adding the 15,527 expressed sequence tags (ESTs) generated from another transcriptome sequencing of An. funestus [19] to the previous 4 × 44 k array (A-MEXP-2245) [18]. Labelled cRNA was obtained from three biological replicates (10 mosquitoes per replicate) for the following samples: (i) resistant (R) (mosquitoes alive after a 1-h exposure to 4% DDT); (ii) control (C) (mosquitoes unexposed to insecticide and thus representative of the wild-type population); and (iii) susceptible (S) (unexposed mosquitoes from the fully susceptible laboratory strain of An. funestus: FANG) making a total of 60 mosquitoes per locality (Pahou, Kpome and Burkina-Faso). Complementary RNA (cRNA) was amplified from each sample using the Agilent Quick Amp Labelling Kit (two-colour) following the manufacturer's protocol. These cRNA were reciprocally hybridized against each other comparing R-S for resistant vs. susceptible, C-S for control vs. susceptible. Microarray data were analyzed using Genespring GX 13.0 software. To identify differentially expressed genes, a cut-off of twofold-change (FC) and a statistical significance of P < 0.05 using Storey with boostrapping correction for multiple testing were applied. These results were compared to those obtained from Kpome [20].

Quantitative reverse transcriptase PCR
Three genes (GSTe2, CYP6P9a and CYP6P9b) ( Table 1) up-regulated from the microarray analysis and mostly associated with DDT and pyrethroids resistance [12,14,18,21] were assessed by qRT-PCR to validate their expression pattern using the three biological replicates for resistant, control, and FANG. cDNA from the Resistant (R), Control (C) and FANG (S) populations were synthesized using one microgram of total RNA from each of the three biological replicates. The relative expression level and FC of each target gene in R and C relative to S were calculated according to the 2-ΔΔCT method incorporating the PCR efficiency [22] after normalization with the housekeeping genes ribosomal protein S7 (RSP7; AFUN007153-RA), and actin (Actin; AFUN006819) ( Table 1). The results were compared to those obtained in Kpome.

Genotyping of L119F-GSTe2 resistance
The role of the L119F-Gste2 mutation recently shown to play a major role in the DDT resistance was assessed. Field-collected An. funestus sensu stricto (s.s.) females from each selected location were genotyped using a Taqman assay [12]. The reaction was performed in a 10-μl final volume containing 1 × SensiMix (Bioline, London, UK), 800 nM of each primer and 200 nM of each probe using an Agilent MX3005P machine. The following cycling conditions were used: 10 min at 95 °C, 40 cycles of 15 s at 92 °C and 1 min at 60 °C. Two probes labelled

Genetic diversity of GSTe2 across Benin
A full-length GSTe2 (exons and introns) was amplified from 10 field-collected female mosquitoes from each location using Phusion High-Fidelity DNA Polymerase (Fermentas, Burlington, Ontario, Canada) and the following conditions: 1 cycle at 95 °C for 5 min; 35 cycles of 94 °C for 20 s, 57 °C for 30 s and 72 °C for 60 s; and 1 cycle at 72 °C for 5 min. The PCR products were purified using the QIAquick PCR Purification Kit (Qiagen, Valencia, CA, USA) and subsequently sequenced. The GSTe2-L119F polymorphic position was detected through a manual analysis of sequence traces and sequences alignments were done using BioEdit. Data were exported to the software DnaSp-version 5.10.01 to detect genetic variability of the GSTe2 gene among the different populations. A maximum likelihood phylogenetic tree for the coding sequences of GSTe2 in the five localities was constructed using MEGA 5.2 [23]. The best model was firsly assessed and this indicated that the Jukes-Cantor model best describes the GSTe2 haplotypes. This was then used to generate the maximum likehood tree using MEGA 5.2.
In addition, the level of pairwise genetic differentiation between the populations were determined in dnasp 5. 10 using the Kst statistic [24] and the neighbour-joining tree was built using Mega 6.06 [24].

Susceptibility profiles to insecticides
The Pahou population (Benin) had previously been described as highly DDT resistant [14] with no mortality 24 h after 1 h of exposure. The WHO bioassays conducted in Kpome [15] indicated that this An. funestus population, which is located approximately 100 km from Pahou, was also resistant to DDT, with 9.1 ± 2.5% mortality 24 h after 1 h of exposure to 4% DDT for females. The population from Tanongou was moderately resistant to DDT with 90% mortality [16].

Genome-wide transcription microarray analysis
A genome-wide transcription analysis enabled us to identify the set of genes associated with DDT resistance in Pahou (Benin) and Burkina-Faso (Table 2; Fig. 2). These results were compared to Kpome (Benin) ( Table 3; Fig. 3) results where high level of DDT resistance were recorded recently [15]. A total of 6610 probes were differentially expressed (FC ≥ 2 at P < 0.05) between the DDT-resistant samples from Pahou and the susceptible strain FANG with 4637 up regulated and 1973 down regulated. The comparison between the control wild type samples (Control) from Pahou and the susceptible strain FANG showed 9756 probes differentially expressed with 7489 up regulated and 2267 down regulated. In Burkina-Faso, a total of 3602 probes were differentially expressed between the DDT-resistant samples and the susceptible FANG. When comparison was made between samples from Pahou and Burkina-Faso, 1007 probes were differentially expressed with 779 over expressed and 228 down expressed as presented in (Table 2, Fig. 2). On the other hand, samples from Pahou were also compared to those from Kpome (20) and 852 common probes were differentially expressed with 326 up regulated and 526 down regulated (Table 3; Fig. 2  analysis in Tanongou (Benin). Several gene families among which the most preeminent were the cytochrome P450 genes were also over expressed. Beside cytochrome P450s, other genes belonging to multiple gene families included alcohol and aldehyde dehydrogenases were up-regulated.

Validation of the microarray upregulation with qRT-PCR
Transcription analysis of the candidate resistance genes GSTe2, CYP6P9a and CYP6P9b revealed that these genes are significantly upregulated in An. funestus from Pahou and Kpome. Indeed the GSTe2 was the most upregulated gene with a fold-change FC of 44.8 [12] in Pahou    [20] and this expression pattern goes with DDT resistance observed in both localities. The two P450 duplicated genes CYP6P9a and CYP6P9b were also upregulated with a FC of 2.9; 7.1 and 3.7; 3.4, respectively in Pahou and Kpome [12,20].

Correlation between the L119F mutation and DDT resistance
The genotyping of the GSTe2-L119F mutation in Kpome, and Doukonta in the southern Benin where high resistance was observed against DDT revealed a high frequency of 96% and 93% of the 119F in these locations compared to Tanongou in the North Benin (35%) where moderate resistance was observed to DDT as reported by Djouaka et al. [16]. Also, similar results were reported in Burkina Faso with 25% of the 119F mutation in correlation with the prevalence of DDT resistance [12].

Role of the GSTe2-L119F mutation in DDT resistance and genetic diversity of Gste2 gene
Full length GSTe2 (exons and introns) was successfully amplified and directly sequenced in ten samples from each locality. These localities are Kpome, Pahou, Doukonta in the southern Benin, Tanongou in the northern Benin, Akaka-Remo in the southern Nigeria and Burkina-Faso. The L119F-GSTe2 mutation is the replacement of leucine (CTT) with phenylalanine (TTT) at the position 119. The C/C is the homozygote susceptible wild type, the T/T is the homozygote mutant genotype while the C/T is a codominant genotype. Interestingly, no T/T genotype (the homozygous resistant allele) was detected in Tanongou (North Benin) and Burkina Faso population where moderate resistance was recorded against DDT while in others populations highly resistant to DDT were almost all homozygote T/T (Fig. 4). The alignment of 739 bp of the sequenced samples showed a heterogeneity between the An. funestus population analysed as reported in Table 4. The analysis of maximum likehood phylogenetic tree of GSTe2 indicated that An. funestus populations are structured according to their pattern of DDT resistance. The ML tree shows that sequences from southern Benin cluster closer to those from southern Nigeria where high resistance level was recorded and sequences from Tanongou cluster with those from Burkina-Faso where moderate resistance level was observed (Fig. 5a). This pattern is also supported by the neighbour-joining tree with genetic distances based on Fst estimates (proportion of the total genetic variance contained in a subpopulation) (Fig. 5b). This result suggest the presence of barriers to gene flow that are affecting the spread of resistance genes. In addition, the presence of a large indel in the GSTe2 gene was

Discussion
Insecticide resistance is a complex trait and factors involved vary depending on species, insecticide and population. This research was designed to assess the underlying molecular basis driving DDT resistance in the South-North transect of Benin to improve performances of malaria controls tools. Physiological resistance to insecticides often involves either mutations in the insecticide target site (target-site resistance), or elevated activity of detoxifying enzymes that metabolise and/or sequester insecticides (metabolic resistance). In the absence of knockdown resistance mutations in the voltage-gate sodium channel in An. funestus [10,11,14], this study showed how DDT resistance in this mosquito species is a result of both target-site resistance and up-regulation of a DDT detoxifying enzyme. Overall, this study has revealed GSTe2 gene as a key gene implicated in DDT as a result of elevated expression rather than allelic variation through the GSTe2-L119F mutation. This is in line with previous studies [2-7, 11, 12] showing that over transcription of the GSTe2 gene confers DDT resistance and cross-resistance to permethrin. It has also been shown that the overexpression of GSTe2 gene in DDT resistant strain of An. gambiae [8,25] and the overproduction of this gene is very efficient at metabolizing DDT [26]. Also, the GSTe2 gene has been implicated in DDT resistance in Aedes aegypti species from Thailand [27] showing the important role of the up-regulation of this gene in DDT resistance. Such observations are in accordance with the resistance profile of Pahou and Kpome An. funestus populations, which are highly resistant to DDT compared to Tanongou closer to Burkina-Faso where moderate resistance was recorded. Beside the overexpression of the GSTe2 gene, it is acknowledged that the presence of L119F-GSTe2 mutation confers DDT resistance in An. funestus s.s. populations in West/Central and East Africa [12] and this is in line with the high allelic frequency of this mutation recorded in Kpome and Doukonta. High allelic frequency of the GSTe2-L119F mutation and upregulation of the GSTe2 gene were observed in Kpome, Pahou and Doukonta while low regulation and low allelic frequency of the mutation were observed in Tanongou and Burkina-Faso in correlation with DDT profile observed. The consistent difference observed for this gene between the population of southern Benin (Kpome, Pahou and Doukonta) and that of Tanongou (North Benin) and Burkina-Faso suggest that possible barriers to gene flow exist between these populations. These barriers to gene flow could be due to geographic distance, because gene flow can be restricted by physical barriers separating the populations. Clarke [28] and Duke et al. [29]. also reported that habitat discontinuities may present barriers to gene flow. Furthermore, the genetic and behavioural divergence may be related to differences in the scale of vector control interventions between the regions or an effect of climate change could also explain this phenomenon. Such anti-vector interventions have been found to impact population size of vector populations [30]. Insecticide resistance in vector populations has been widespread with large scale exposure resulting in altered abundance, behavioural shifts and general ecology of major vector populations (e.g. An. funestus, An. gambiae) [31]. However, this observation needs to be confirmed in future studies.
Analysis of the full-length GSTe2 gene shows a possible association between the GSTe2 polymorphism and observed DDT resistance in the 6 localities. The 119F resistant allele is fixed in highly DDT-resistant Benin mosquitoes especially in Pahou, Kpome and Doukonta and in Nigeria (Akaka Remo), but very low in moderate resistant mosquitoes in Tanongou and Burkina-Faso showing the key role of this mutation in DDT resistance as reported by [12]. This study revealed that southern Benin and Nigeria populations of An. funestus are more genetically differentiated as they form a unique cluster compared to North populations. This pattern of genetic diversity of the GSTe2 gene observed in this study support the contrast in resistance patterns between populations of An. funestus. In addition, a significant shift in the over-expression profile of this gene was detected across a South/North transect of Benin in line with the DDT resistance profile observed, showing that the L119F-GSTe2 mutation coupled with up-regulation of this gene confer a high level DDT resistance in An. funestus [12]. The consistent differences between the An. funestus population across Benin is likely to impact the design and implementation of resistance management strategies in Benin.

Conclusion
Effective management of resistance requires an understanding of the dynamics and mechanisms driving resistance. This study shows that molecular basis of DDT resistance in southern Benin An. funestus is associated with L119F-GSTe2 mutation and over-expression of this gene. The variations observed between southern and northern populations of An. funestus could suggest the presence of barriers to gene flow that are affecting the spread of resistance and associated genes.