Skip to main content

Molecular basis of permethrin and DDT resistance in an Anopheles funestus population from Benin

Abstract

Background

Insecticide resistance in Anopheles mosquitoes is threatening the success of malaria control programmes. In order to implement suitable insecticide resistance management strategies, it is necessary to understand the underlying mechanisms involved. To achieve this, the molecular basis of permethrin and DDT resistance in the principal malaria vector, Anopheles funestus from inland Benin (Kpome), was investigated.

Results

Here, using a microarray-based genome-wide transcription and qRT-PCR analysis, we showed that metabolic resistance mechanisms through over-expression of cytochrome P450 and glutathione S-transferase genes (GSTs) are a major contributor to DDT and permethrin resistance in Anopheles funestus from Kpome. The GSTe2 gene was the most upregulated detoxification gene in both DDT- [fold-change (FC: 16.0)] and permethrin-resistant (FC: 18.1) mosquitoes suggesting that upregulation of this gene could contribute to DDT resistance and cross-resistance to permethrin. CYP6P9a and CYP6P9b genes that have been previously associated with pyrethroid resistance were also significantly overexpressed with FC 5.4 and 4.8, respectively, in a permethrin resistant population. Noticeably, the GSTs, GSTd1-5 and GSTd3, were more upregulated in DDT-resistant than in permethrin-resistant Anopheles funestus suggesting these genes are more implicated in DDT resistance. The absence of the L1014F or L1014S kdr mutations in the voltage-gated sodium channel gene coupled with the lack of directional selection at the gene further supported that knockdown resistance plays little role in this resistance.

Conclusions

The major role played by metabolic resistance to pyrethroids in this An. funestus population in Benin suggests that using novel control tools combining the P450 synergist piperonyl butoxide (PBO), such as PBO-based bednets, could help manage the growing pyrethroid resistance in this malaria vector in Benin.

Background

There were an estimated 216 million cases of malaria worldwide in 2016 and 445,000 deaths with 80% of all malaria deaths occurring in Africa [1]. Despite extensive control efforts, over half the world’s population remains at risk and the disease has a massive impact on health and economic development, particularly in Africa [2]. Four species of Anopheles, An. gambiae Giles, An. coluzzii Coetzee & Wilkerson, An. arabiensis Patton and An. funestus Giles, are responsible for most of the malaria transmission in this continent. The transmission role of An. funestus (sensu stricto) is further supported by observations that in many regions of Africa, its infection rate even surpasses that of An. gambiae [3]. Long-lasting insecticidal nets (LLINs) and indoor residual spraying (IRS) are the main malaria prevention interventions [4]. However, the success of these control methods is jeopardised by the development of resistance by Anopheles species to insecticides such as pyrethroids and DDT as seen in the past with loss of efficacy of dieldrin for IRS in west/central Africa [5, 6]. In Benin, due to vector resistance to pyrethroids across the country [7,8,9,10], two insecticides of two different classes were used in rotation for IRS: bendiocarb (a carbamate) and pirimiphos-methyl (an organophosphate) [11]. However, DDT is still retained for use in IRS, due to the limited number of cost-effective alternatives [12]. Field studies on insecticide susceptibility carried out as baseline surveys for malaria control programs showed An. funestus to be resistant to various insecticides at various localities [9, 10, 13]. Metabolic resistance is the main resistance mechanism recorded, and cytochrome P450 genes are playing a major role while target-site resistance like the knockdown resistance (kdr) is absent [9, 14,15,16,17]. Target site resistance was investigated in the pyrethroid/DDT-resistant population and the fifth and sixth segments of domain II in the sodium channel sequence and no polymorphism associated with resistance was identified [9, 14, 16, 18,19,20] although some mutations have been detected in Cameroon [21] and Uganda [22] but with no association yet established with resistance phenotype. A significant increase in mono-oxygenase activities confers resistance to pyrethroids in An. funestus in Mozambique and South-Africa (southern Africa) [15, 23,24,25], Uganda (east Africa) [16, 22, 26], western Kenya (east Africa) [18, 22], Cameroon (central Africa) [21] and in Senegal (west Africa) [19]. Characterisation of the molecular basis of resistance to pyrethroids and DDT in An. funestus from the coastal area of Benin [27] revealed a predominant role played by the glutathione S-transferase GSTe2. The enzyme encoded by this gene was highly upregulated in association with the presence of an L119F mutation in the substrate binding site conferring a high level of DDT and permethrin resistance [27]. Recently, An. funestus from Kpome, an inland region of Benin has also been shown to be resistant to DDT and permethrin with mortality rates of 9.1 ± 2.5% and 13.0 ± 3%, respectively [10]. However, it remains unknown if the resistant mechanism is the same as observed in Pahou, the coastal area of Benin [9] or whether different resistance mechanisms are responsible for the spread of resistance across the country. Such information is essential in designing suitable control interventions nationwide and to improve the implementation of resistance management strategies against An. funestus. In this study, using a genome-wide microarray-based transcription analysis, the underlying molecular mechanisms conferring DDT and permethrin resistance in Kpome were characterised revealing that both glutathione S-transferases, notably the GSTe2, and cytochrome P450 genes are the main drivers of resistance.

Methods

Area of study and mosquito collection

Blood-fed adult female An. funestus mosquitoes resting indoors were collected in houses between 06:00 and 10:00 h from the rural area of Kpome (6°23'N, 2°13'E) a southern inland region of Benin from December 2013 to February 2014. Kpome is a large agricultural setting with the intensive production of tomatoes which mainly serve a commercial purpose and insecticides, mainly pyrethroids (deltamethrin and lambda-cyalothrin), used in public health are also used to protect tomatoes against pest attacks. Mosquito collections and rearing were performed as described previously [10]. Adult female mosquitoes were collected indoor between 06:00 to 10:00 h using electric aspirators. F1 adults were generated from field-collected female mosquitoes using the forced-egg laying method [16] and were randomly mixed in cages for subsequent experiments. All females used for individual oviposition were morphologically and molecularly identified as An. funestus (s.s.) [10, 28].

Insecticide susceptibility assays

Two to five day-old F1 adult female mosquitoes pooled from different F0 mosquitoes collected from Kpome were used for this test. Twenty to twenty-five mosquitoes per tube with at least four replicates were exposed to DDT (4%) and permethrin (0.75%) for 1 h before transferring into clean holding tube with 10% sugar solution. Final mortality was scored 24 h post-exposure. Survivors and dead mosquitoes were conserved for further analysis. In this study, we used the same mosquitoes from the previously published research results [10] to further describe the molecular basis of the observed resistance in this malaria vector population. Mosquito samples generated from the previous investigation were used for genetic analysis in this work. The susceptible strain used in this study has the same age as the exposed mosquitoes.

Microarrays

A genome-wide transcription profiling was performed to detect the sets of genes differentially expressed in relation to observed resistance phenotypes in An. funestus populations from Kpome. The microarray hybridisation for An. funestus was carried out using the 8 × 60 k (60 mer) Agilent An. funestus chip. This Agilent microarray chip was designed using the eArray program (Agilent, Santa Clara, CA, US ) (A-MEXP-2374) by adding the 15,527 expressed sequence tags (ESTs) generated from another transcriptome sequencing of An. funestus [29] to the previous 4 × 44 k array (A-MEXP-2245) [24]. Each array was designed using 8540 ESTs generated from An. funestus transcriptome 454 sequencing [30], 2850 An. funestus cDNAs from GenBank and a set of P450 genes from the rp1 and rp2 QTL BAC sequence [31, 32], and the 13,000 transcripts of the complete An. gambiae genome. Moreover, An. gambiae detoxification genes previously present on the An. gambiae detox chip [33] were added to this. RNA was extracted from three batches of 10 An. funestus females (2–5 days-old) from the following sample sets: alive after exposure to 0.75% permethrin (Resistant permethrin: Rperm) and 4% DDT (Resistant DDT: RDDT), un-exposed to insecticides (Control: C), and from the fully susceptible laboratory strain FANG (Susceptible: S) [a strain originating from Calueque in southern Angola in 2002 which is completely susceptible to insecticides (Hunt et al. [34]) using the Picopure RNA Isolation Kit (Arcturus, Waltham, MA, USA). The FANG colony was kept in insectary conditions (temperature of 25–28 °C with a relative humidity of 80%) at the Liverpool School of Tropical Medicine where assays were performed. RNA quantity and quality were assessed using a NanoDrop ND1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) and Bioanalyzer (Agilent, Santa Clara, CA, USA), respectively. Each RNA sample was used to generate complementary RNA (cRNA) using the Agilent Quick Amp Labeling Kit (two-colour) following the manufacturer’s protocol. cRNAs from the resistant samples (R) were labelled with a cy5 dye, and cRNAs from the susceptible strain FANG (S) were labelled with the cy3 dye. cRNA quantity and quality were assessed using the NanoDrop and Bioanalyzer before labelling. Labelled cRNAs were hybridised to the arrays for 17 h at 65 °C according to the manufacturer’s protocol. Five hybridisations were performed for each comparison including two dyes swap per comparison. Microarray data were analysed using Genespring GX 13.0 software purchased from Agilent (Agilent). A normalisation of the data was performed with the GeneSpring feature extraction program using the lowess normalisation approach. To identify differentially expressed genes, a cut-off of 2-fold-change (FC) and a statistical significance of P ≤ 0.05 with Benjamini-Hochberg correction for multiple testing and Storey with bootstrapping (with a cut-off of 1.5-fold-change for the R-C comparison) were applied. A false discovery analysis was also performed using the Benjamini-Hochberg correction test for multiple testing as implemented in GeneSpring.

Quantitative reverse transcriptase PCR

A qRT-PCR was used to validate the microarrays results for seven of some upregulated detoxification genes. These genes include three cytochromes P450s (CYP6M7, CYP6P9a, CYP9K1), two glutathione transferases (GSTd1-5, GSTe2), and one aldehyde oxidase (Ald oxi). Also, the expression level of CYP6P9b previously found to be strongly associated with pyrethroid resistance [19, 31, 35] and GSTd1-3, a gene shown as possibly associated with pyrethroid resistance by 454 transcriptome profiling [19, 30] were analysed (Table 1). One microgram of the three biological replicates for resistant (permethrin: Rperm and DDT: RDDT), control (C), and FANG (S) was used for cDNA synthesis using Superscript III (Invitrogen, Carlsbad, CA, USA) with oligodT20 and RNase H according to the manufacturer’s instructions. Standard curves were established for each gene using serial dilution of cDNA. qRT-PCR amplification was performed using the MX 3005P (Agilent) system. The relative expression level and fold change (FC) of each target gene in resistant and control samples relative to susceptible were calculated according to the 2-ΔΔCT method incorporating the PCR efficiency [36] after normalisation with the housekeeping genes ribosomal protein S7 (RSP7; AFUN007153-RA), and actin (Act5C; AFUN006819).

Table 1 List of primers used for qRT-PCR

Investigation of the role of the knockdown resistance mutation in DDT and permethrin resistance

A fragment spanning a portion of the voltage-gated sodium channel gene (VGSC), including the 1014 codon (a portion of intron 19 and the entire exon 20, domain II, segment 6 of the VGSC gene) associated with resistance in An. gambiae, was amplified in eleven field-collected females of An. funestus from Kpome using the Exon 19F/KdrFunR2 primers (Exon 19F: 5'-TTT TTA AGC TCG CTA AAT CGT G-3'; KdrFunR2: 5'-CCG AAA TTT GAC AAA AGC AAA-3') [9, 15, 22]. The PCR products were purified using the QIAquick Purification Kit (Qiagen, Hilden, Germany) and directly sequenced. The polymorphic positions were detected through manual analysis of sequence traces using BioEdit. Sequences were aligned using ClustalW [37], and haplotype reconstruction and analysis were done using DnaSP v5.10 [38]. Kpome sequences were also compared to those previously obtained from Pahou in coastal Benin [9] and Gounougou in North Cameroon [21]. After the selection of the best model, a maximum likelihood phylogenetic tree was generated for the VGSC haplotypes in the different populations using Mega 6.06 [39]. The level of Kst of pairwise genetic differentiation between the populations was estimated in Dnasp v5.10, and the neighbour-Joining tree was built using Mega 6.06.

Results

Genome-wide microarray-based transcriptional profiling of permethrin resistance

A set of transcripts were differentially expressed (≥ 2-fold change, P ≤ 0.05) between permethrin resistant (Rperm), susceptible (S) and control (C) mosquitoes (Fig. 1) from microarray analysis. Overall, 3669 transcripts were differentially expressed when permethrin resistant mosquitoes were compared to susceptible strains (1460 overexpressed and 2209 underexpressed), 3617 between C-S (1211 overexpressed and 2406 underexpressed) and 210 for Rperm-C (91 overexpressed and 119 underexpressed). Overall, a total of 40 transcripts were differentially expressed in all the three comparisons.

Fig. 1
figure 1

Summary of transcripts differentially expressed in permethrin resistance. The Venn diagrams show the number of transcripts significantly (P ≤ 0.05) up- or downregulated (FC ≥ 2) in each comparison as well as the commonly expressed transcripts

Genes commonly overexpressed in Rperm-S, C-S and Rperm-C comparisons

Only the transcript, AGAP011496-PA (Anopheles gambiae str. PEST) was upregulated in Rperm-S, Rperm-C and C-S comparisons with FC of 2.7, 1.5 and 2.3, respectively. No detoxification gene was found in the comparison of the three groups of mosquitoes.

Genes commonly overexpressed in Rperm-S and C-S comparisons

Several detoxification genes were commonly over-expressed in the Rperm-S and C-S. Among these genes, the transcripts for glycogenin (Afun000500) and CYP6M7, belonging to the cytochrome P450 gene family, were the most commonly overexpressed detoxification gene in Rperm-S (FC of 36.2 for Afun000500 and 34.6 for CYP6M7) and C-S (FC of 38.1 for Afun000500 and 28.9 for CYP6M7). A transcript from CYP6AA1 was also commonly upregulated in the Rperm-S and C-S samples with FC of 5.2 and 3.8, respectively. This gene is the ortholog of CYP6AA3 in An. minimus, which was shown to metabolise pyrethroids [40]. Other upregulated detoxification genes identified in Rperm-S and C-S comparisons include CYP6Y2 located in the genomic region spanning the pyrethroid resistance rp2 QTL detected in the FUMOZ-R pyrethroid-resistant laboratory strain [32], as well as CYP6M4, CYP6S1, CYP4K2 and CYP315A1 that were expressed at lower levels (FC < 5) (Table 2). Transcripts from the glutathione S-transferase (GST) family notably from epsilon [GSTe2 (FC of 10.9 in Rperm-S and 15.2 in C-S) and GSTe4 (FC of 3.9 in Rperm-S and 4.7 in C-S)] and delta, GSTd1-5 (FC 2.0 and 2.6 in Rperm-S and C-S, respectively) classes were upregulated. Other genes such as ABC transporter (Afun009697), alcohol dehydrogenases (Afun012461, Afun013797), carboxylesterase (Afun009492), short chain dehydrogenase (combined_c738) and chymotrypsin (Afun013921) were also overexpressed (Table 2 and Additional file 1: Table S1).

Table 2 Top detoxification genes upregulated in Rperm-S, C-S and Rperm-C

Comparative transcriptional profiling of DDT and permethrin resistance

A comparative analysis of transcriptomes associated with DDT resistance in Kpome showed a total of 2158 transcripts differentially expressed (≥ 2-fold change, FC at P ≤ 0.05) in RDDT-S with 715 upregulated and 1443 downregulated. The comparison between RDDT-S, Rperm-S and C-S revealed 1311 transcripts differentially expressed. A set of 212 transcripts was commonly expressed between RDDT-S and Rperm-S, 958 transcripts and 270 transcripts were also expressed when comparing Rperm-S and RDDT-S to C-S (Fig. 2).

Fig. 2
figure 2

Summary of transcripts differentially expressed in DDT and permethrin resistance. The Venn diagrams show the number of transcripts significantly (P ≤ 0.05) up- or downregulated (FC ≥ 2) in each comparison as well as the commonly expressed transcripts permethrin

Genes commonly overexpressed in RDDT-S, Rperm-S and C-S comparisons

Many genes were commonly overexpressed in the three comparisons with glycogenin (Afun000500) and arginosuccinate lyase (Afun009227) been the most overexpressed and were consistent in all the three comparisons with fold change of (FC 36.2, 38.1 and 41.4) and (FC 22.2, 23.3 and 25.6) in RDDT-S, Rperm-S and C-S, respectively. Among the common detoxification gene families identified, the GSTe2 gene (Afun000045) was the one that remained highly upregulated with FC of 13.3, 10.9 and 15.2 for RDDT-S, Rperm-S and C-S, respectively. The CYP6M7 gene was also overexpressed in all three comparisons; FC was 34.6 in Rperm-S and FC 28.9 in C-S but with lower FC of 4 in RDDT-S. Several other genes of the family of P450 such as CYP6AA1, CYP6M4, CYP6Y2, CYP4K2 and CYP315A1 were also upregulated. Other transcripts of the GST family were also upregulated such as the Gste4 (FC 4.3, 3.9 and 4.7 in RDDT-S, Rperm-S and C-S, respectively) and the GSTd1 (FC 2.2, 2.0 and 2.6 in RDDT-S, Rperm-S and C-S, respectively (Table 3 and Additional file 1: Table S2).

Table 3 Top detoxification genes upregulated in RDDT-S, Rperm-S and C-S

Overexpressed genes in the comparison between RDDT-S and Rperm-S

A limited number of genes were commonly overexpressed in these two “resistant” comparisons, potentially as a result of induction or greater involvement in resistance than in control samples. The transcript of alcohol dehydrogenase (Afun011037) was the most overexpressed in RDDT-S, and Rperm-S comparison with FC 6.5 and 3.7, respectively. Other genes of the P450 family including CYP4H18, CYP9K1, CYP9J3, and the esterase b1 were also overexpressed (Table 3 and Additional file 1: Table S2).

Overexpressed genes in the comparison between RDDT-S and C-S

The transcript (combined_c738) belonging to the short-chain dehydrogenase gene was the most commonly upregulated in the RDDT-S, and C-S comparison with FC 10.4 and 15.9, respectively. The cytochrome P450 gene CYP307A1 was also upregulated with FC 3.3 in RDDT-S and FC 2.2 in C-S. The CYP6P9a gene, reported to be upregulated in southern Africa mosquitoes, was also upregulated with FC 2.8 and 2.3 in RDDT-S and C-S, respectively. Other genes with a known association with insecticide resistance were also upregulated in Kpome mosquitoes as shown in Table 3 and Additional file 1: Table S2.

Overexpressed genes in the comparison between Rperm-S and C-S

Genes identified from these comparisons include the alcohol dehydrogenase gene (Afun012461 ortholog of AGAP000288 in An. gambiae) which is consistently upregulated in Rperm-S (FC8.4) and C-S (FC10.9) but not in RDDT-S. Similarly, the aldehyde dehydrogenase (Afun008026) and the ABC transporter (Combined_c1762) were also upregulated in Rperm-S and C-S exclusively although with lower FC.

Overexpressed genes in the comparison between DDT resistant mosquitoes and a susceptible strain

Some detoxification genes specific for RDDT-S comparison were identified in Kpome. These genes include alcohol dehydrogenase Afun013475 (FC 3.2), the most upregulated gene (FC 3.2) as well as CYP302A1 (FC 2.6), CYP4D15 (FC 2.4), CYP4J5 (FC 2.2), gb-CYP12F3 (FC 2.1) and the cytochrome P450 (Afun010630) (FC 2.1).

Overexpressed genes in the comparison between permethrin resistant mosquitoes and a susceptible strain

For Rperm-S comparison, the ABC transporter (Afun007722) was the most upregulated gene with a fold change of 3.2. Also, the carboxylesterase (Afun012261) (FC 2.8), the cuticle protein (CD5773231) (FC 2.4) and the aldehyde oxidase (Afun00493) (FC 2.4) were upregulated.

Quantitative reverse transcriptase PCR

A qRT-PCR analysis revealed that GSTe2 is the most overexpressed genes in DDT and permethrin resistant samples (FC 16.0 and 18.1, respectively) analysed in this study. However, the expression level of the detoxification gene, CYP6M7 recorded with the qRT-PCR analysis was lower (FC 1.4 and 1.7 in permethrin and DDT resistant samples, respectively) compared to expression level from microarray analysis (FC 34.6 and 4.0 in permethrin and DDT resistant samples, respectively). The GSTd1-5 and the GSTd3 were significantly upregulated in DDT resistant samples (FC 12.5 and 6.2; P < 0.01) compared to permethrin resistant samples (FC 0.72 and 0.86; P < 0.05) (Fig. 3a, b).

Fig. 3
figure 3

Quantitative PCR results: differential expression by qRT-PCR of some candidate genes upregulated in microarray assays in An. funestus resistant to DDT (a) and permethrin (b). Error bars represent SD (n = 3). The presence of * on top of the two-fold changes for each gene indicates a statistically significant (P < 0.05) over-expression in resistant or non-exposed mosquitoes compared to the FANG susceptible strain; “ns” is shown when the difference was not significant

Role of the knockdown resistance gene in DDT and permethrin resistance profiles recorded in An. funestus from Kpome

Amplification and sequencing of a fragment (a portion of intron 19 and the entire exon 20, domain II, segment 6) of the VGSC gene showed that both L1014F (TTA-to-TTT) and L1014S (TTA-to-TCA) kdr mutation commonly found in An. gambiae are absent in this mosquito species. However, further analysis with 837 bp of common sequences generated from the eleven mosquitoes analysed revealed 12 polymorphic sites with 12 haplotypes (Table 4, Fig. 4) showing a high genetic diversity of this gene within the An. funestus mosquitoes of Kpome. No amino acid change was detected in the Benin population. Furthermore, the Tajima D and Fu and Li D* statistics (Table 4) were not statistically significant. The Neighbour-joining tree showed low genetic differentiation between Benin populations compared to Cameroon samples with respect to geographical distance (Fig. 4).

Table 4 Summary statistics for polymorphism in the voltage-gated sodium channel gene in F0 An. funestus from Kpome compared to Pahou and Cameroon population
Fig. 4
figure 4

kdr polymorphism in Anopheles funestus of Kpome. a Schematic representation of haplotypes of Exon20 fragment of the voltage-gated sodium channel gene (VGSC) observed in wild type An. funestus from Kpome. Only polymorphic sites are shown and these are numbered from the beginning of each aligned sequence. Dots indicate identity with the first sequence. A number has been given to each haplotype. The column (N) indicates the number of individuals sharing the haplotype. b Maximum-likelihood tree based on kdr data for samples from Kpome, Pahou and Cameroon. c Neighbour-joining tree based on the VGSC gene data for samples of Kpome, Pahou and Cameroon. Abbreviations: Kp, Kpome; Ph, Pahou, Cam: Cameroon

Discussion

The WHO Global Plan for Insecticide Resistance Management [12] highlights the necessity to detect and monitor the development of insecticide resistance and characterise the underlying resistance mechanisms to maintain the successes in the reduction of malaria cases across Africa. In this study, we elucidated the molecular basis driving DDT and pyrethroids resistance in an inland An. funestus population from Benin.

Permethrin and DDT resistance in An. funestus from Kpome is driven by metabolic resistance

Transcriptional profiling results and the absence of kdr mutation points out to a possible metabolic mechanism driving permethrin and DDT resistance in An. funestus from Kpome notably through overexpression of genes involved in insecticide detoxification such as cytochrome P450 genes, GSTs, aldehyde oxidases, and other gene families previously associated with resistance of An. funestus to insecticides [41].

The cytochrome P450, CYP6M7 was the most overexpressed detoxification gene in permethrin exposed mosquitoes compared to the FANG susceptible strain as well as to the unexposed mosquitoes from microarrays analysis. However, this consistent information from microarray analysis was not confirmed by qRT-PCR. This discrepancy could be associated with the high genetic redundancy of some P450 genes and a high level of sequence similarity between different genes. These factors could lead to cross amplification contributing to homogenise the expression between samples. It is more likely that CYP6M7 could also be a key pyrethroid resistance gene in Benin, but further validation is required. New approaches with RNAseq could help to confirm the expression pattern of this gene. On the other hand, only one gene (AGAP011496-PA) was commonly expressed in Rperm-S, C-S and Rperm-C comparisons. This could be due to the high resistance recorded against permethrin suggesting a less significant difference between the resistant and the control. Such a low number of differentially expressed genes between R-S and C-S comparison is commonly observed in similar studies when the resistance level is high in the population as it is the case here in Kpome [24, 27]. This is because both non-exposed (control; C) and alive mosquitoes after exposure (resistant; R) are all resistant which leads to less difference in the level of gene expression fold change between these two comparisons. The GSTe2 gene was the most upregulated gene in both DDT and permethrin resistant samples showing that this gene is likely to be involved in both DDT and permethrin resistance as previously shown in coastal Benin [27]. Hence, such genes could confer cross-resistance to both insecticides, and this represents a challenge for the success of the malaria control program. The enzyme encoded by this gene was shown to be able to metabolize both DDT and permethrin [27]. QRT-PCR confirmed the expression level of GSTe2 in permethrin, DDT resistant and unexposed mosquitoes compared to the susceptible FANG. This observation is in line with the common implication of the GSTe2 gene in permethrin and DDT resistance in the An. funestus population from Pahou (Benin) reported by Riveron et al. [27]. Furthermore, the potential role of GSTe2 in permethrin resistance observed here has been shown in previous reports that suggested that orthologs of GSTe2 in other insects are associated with pyrethroids resistance. Indeed, the elevated GSTe2 expression has been associated with pyrethroid resistance by acting as a pyrethroid-binding protein and sequestering the insecticide [42] or by protecting against oxidative stress and lipid peroxidation induced by pyrethroid exposure [43]. In the yellow fever mosquito Ae. aegypti, a partial knockdown of the ortholog of GSTe2 led to increased mortality to pyrethroids (deltamethrin), indicating that GSTe2 is also associated with deltamethrin resistance in Ae. aegypti [44]. Moreover, the over-expression of GSTe2 was observed with the high frequency of the 119F resistant allele in Kpome (91% of 119F/F homozygote resistant genotype) [10]. The near fixation of 119F-GSTe2 resistant allele in Kpome, which enlarges the substrate binding sites to increase DDT metabolism [27], coupled with the high overexpression of GSTe2 highlights the key role played by this gene in DDT resistance as previously observed in Pahou. Orthologs of GSTe2 have also been shown to be associated with DDT resistance in other mosquito species such as An. gambiae and Ae. aegypti [33, 44, 45]. Overall, GSTe2 can confer resistance to DDT and permethrin, and this cross-resistance to pyrethroids is of significant concern for malaria control as GSTe2 could protect mosquitoes against the major insecticides used to impregnate LLINs in public health. However, the overexpression of several genes from other gene families in this study highlights the complexity of resistance mechanisms suggesting the involvement of other genes than just GSTe2. In addition, the transcription profile of the duplicated genes CYP6P9a and CYP6P9b (which can metabolise both types I and II pyrethroids as shown by Riveron et al. [27]) in pyrethroid resistance in Kpome is different to that observed for An. funestus population in southern Africa. Indeed these P450 genes were highly upregulated in pyrethroid-resistant laboratory and field population from southern Africa [9, 15, 16, 31, 46, 47] while a lower level of overexpression was recorded in Kpome. The level of expression of CYP6P9a and CYP6P9b recorded in Kpome An. funestus population is similar to what was observed previously in Pahou population (coastal Benin). This result suggests that the two duplicated P450 genes are not strongly associated with permethrin resistance in Benin as observed in Mozambique (southern Africa) [15, 35]. Therefore, the molecular basis of the pyrethroid resistance in Benin is most likely different to that in southern Africa pointing to independent selection events of the pyrethroid resistance across Africa probably under different local selective forces. These different selective pressures could be increased of ITNs coverage across Benin, agricultural use of pesticides [48, 49] and spilt petroleum products [50]. In most African urban areas, insecticides are used for domestic purposes, including the control of mosquitoes in the form of mosquito coils, fumigation bombs or sprays [51]. These insecticides are used in an uncontrolled and heterogeneous manner in term of coverage and doses of insecticides in each household. Such practices may represent an additional selective pressure favouring pyrethroid resistance. The presence of agrochemicals, or industrial pollutants and plant compounds in mosquito breeding sites could also affect insecticide tolerance by modulating mosquito detoxification systems [52]. Concerning DDT resistance, in addition to the GSTe2 gene, two other detoxification genes of GST family were also upregulated: GSTd1-5 and GSTd3. GSTd3 was shown to be upregulated in DDT-resistant An. arabiensis from an urban site in Burkina Faso [53]. GSTd1-5 have been previously implicated in coding for enzymes that directly metabolise DDT or have at least been previously associated with the DDT-resistant phenotype [54, 55]. Further validation of the role of these genes in DDT resistance is required. The two duplicated cytochrome P450 genes, CYP6P9a (FC = 3.7) and CYP6P9b (FC = 3.9), which confer pyrethroid resistance in southern African populations of An. funestus [24] were also upregulated in the Benin population. However, because their encoded proteins are unable to metabolise DDT [24], these genes are not likely involved in DDT resistance observed in Kpome mosquitoes. Nevertheless, further investigation is required to validate this hypothesis. Monitoring the insecticide resistance mechanisms that occur within a population should be an essential component to all insecticide-based vector control programs and improving resistance management involves a better understanding of resistance mechanisms. These data contribute to the growing body of knowledge focussed on pyrethroid and DDT resistance in Benin. This situation emphasises the need for natural resistance and vector monitoring so that adjustments to control programmes can be made timeously and accurately.

Analysis of polymorphisms of the VGSC gene supports a limited role of knockdown resistance

Mechanisms of DDT and permethrin resistance are likely not associated with target site resistance as no kdr mutation was detected in analysed mosquitoes from Kpome. L1014F or L1014S change commonly associated with pyrethroid/DDT resistance in An. gambiae was not detected in this population of An. funestus, as it was also the case for all populations of this species, analysed so far [9, 19,20,21,22]. This suggests that the VGSC is probably evolving neutrally in DDT and permethrin resistance in An. funestus population.

Furthermore, the neutrality tests with Tajima D and Fu and Li D* statistics revealed no signature of directional selection on the sodium channel gene suggesting the limited role of knockdown resistance in both DDT and pyrethroid resistance in An. funestus in Kpome. The neighbour-joining tree revealed that Kpome mosquitoes cluster with Pahou mosquitoes while they are different to Cameroon mosquitoes. This reveals a reduced gene flow between these populations probably through isolation by distance which can also affect the spread of insecticide resistance genes in this species as previously shown for An. funestus populations across the continent [56].

Conclusions

Metabolic resistance is likely driving resistance to both pyrethroids (permethrin) and DDT in the major malaria vector An. funestus in Benin. The glutathione s-transferase gene, GSTe2 is playing a key role in DDT resistance and most likely is responsible for the observed cross-resistance to pyrethroids in An. funestus populations from Kpome and such cross-resistance should be taken into account for the implementation of future insecticide resistance management strategies. Moreover, this study provides knowledge on the resistance profile and underlying resistance mechanisms to the available insecticides in An. funestus, a less studied malaria vector in Benin, in order to develop better insecticide resistance diagnostics. Further investigation should be performed on the expression level of target genes to ascertain the role of metabolic mechanisms in DDT and permethrin resistance in this An. funestus population. Resistance mechanisms detected in this studied population appear to be different from those identified in other African regions showing the need to characterise mosquito populations at country-level for more appropriate and tailored control interventions.

Abbreviations

DDT:

Dichlorodiphenyltrichloroethane

LSTM:

Liverpool School of Tropical Medicine

qRT-PCR:

Quantitative reverse transcriptase polymerase chain reaction

VGSC:

Voltage-gated sodium channel

WHO:

World Health Organisation

References

  1. WHO. World Malaria Report 2016. Geneva: World Health Organisation; 2016. http://www.who.int/malaria/publications/world-malaria-report-2016.

    Google Scholar 

  2. WHO. World Malaria Report 2011. Geneva: World Health Organisation; 2011. http://www.who.int/malaria/world_malaria_report_2011/en/.

    Google Scholar 

  3. Fontenille D, Simard F. Unravelling complexities in human malaria transmission dynamics in Africa through a comprehensive knowledge of vector populations. Comp Immunol Microbiol Infect Dis. 2004;27:357–75.

    Article  Google Scholar 

  4. WHO. World Malaria Report 2014. Geneva: World Health Organisation; 2014. http://www.who.int/malaria/publications/world_malaria_report_2014/report/en/.

    Google Scholar 

  5. 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–96.

    Article  CAS  Google Scholar 

  6. Hamon J, Abonnenc E, Noel E. Contribution à l’étude des Culicidés de l’Ouest du Sénégal. Ann Parasitol Hum Comp. 1955;30:278–308.

    Article  CAS  Google Scholar 

  7. Akogbéto M, Yakoubou S. Resistance of malaria vectors to pyrethrins used for impregnating mosquito nets in Benin, West Africa. Bull Soc Pathol Exot. 1999;92:123–30.

    PubMed  Google Scholar 

  8. Corbel V, N’Guessan R, Brengues C, Chandre F, Djogbenou L, Martin T, et al. Multiple insecticide resistance mechanisms in Anopheles gambiae and Culex quinquefasciatus from Benin, West Africa. Acta Trop. 2007;101:207–16.

    Article  CAS  Google Scholar 

  9. Djouaka R, Irving H, Tukur Z, Wondji CS. Exploring mechanisms of multiple insecticide resistance in a population of the malaria vector Anopheles funestus in Benin. PLoS One. 2011;6:e27760.

    Article  CAS  Google Scholar 

  10. Djouaka R, Riveron JM, Yessoufou A, Tchigossou G, Akoton R, Irving H, et al. Multiple insecticide resistance in an infected population of the malaria vector Anopheles funestus in Benin. Parasit Vectors. 2016;9:453.

    Article  Google Scholar 

  11. Akogbéto MC, Aïkpon RY, Azondékon R, Padonou GG, Ossè RA, Agossa FR, et al. Six years of experience in entomological surveillance of indoor residual spraying against malaria transmission in Benin: lessons learned, challenges and outlooks. Malar J. 2015;14:242.

    Article  Google Scholar 

  12. WHO. Global plan for insecticide resistance management in malaria vectors. Geneva: World Health Organization Press; 2012.

    Google Scholar 

  13. Djouaka RJ, Atoyebi SM, Tchigossou G, Riveron JM, Irving H, Akoton R, et al. Evidence of a multiple insecticide resistance in the malaria vector Anopheles funestus in south-west Nigeria. Malar J. 2016;15:565.

    Article  Google Scholar 

  14. Okoye PN, Brooke BD, Koekemoer LL, Hunt RH, Coetzee M. Characterisation of DDT, pyrethroid and carbamate resistance in Anopheles funestus from Obuasi, Ghana. Trans R Soc Trop Med Hyg. 2008;102:591–8.

    Article  CAS  Google Scholar 

  15. Cuamba N, Morgan JC, Irving H, Steven A, Wondji CS. High level of pyrethroid resistance in an Anopheles funestus population of the Chokwe District in Mozambique. PLoS One. 2010;5:e11010.

    Article  Google Scholar 

  16. Morgan JC, Irving H, Okedi LM, Steven A, Wondji CS. Pyrethroid resistance in an Anopheles funestus population from Uganda. PLoS One. 2010;5:e11872.

    Article  Google Scholar 

  17. Wondji CS, Dabire RK, Tukur Z, Irving H, Djouaka R, Morgan JC. Identification and distribution of a GABA receptor mutation conferring dieldrin resistance in the malaria vector Anopheles funestus in Africa. Insect Biochem Mol Biol. 2011;41:484–91.

    Article  CAS  Google Scholar 

  18. Kawada H, Dida GO, Ohashi K, Komagata O, Kasai S, Tomita T, et al. Multimodal pyrethroid resistance in malaria vectors, Anopheles gambiae s.s., Anopheles arabiensis, and Anopheles funestus s.s. in western Kenya. PLoS One. 2011;6:48.

    Google Scholar 

  19. Samb B, Konate L, Irving H, Riveron JM, Dia I, Faye O, et al. Investigating molecular basis of lambda-cyhalothrin resistance in an Anopheles funestus population from Senegal. Parasit Vectors. 2016;9:449.

    Article  Google Scholar 

  20. Irving H, Charles WS. Investigating knockdown resistance (kdr) mechanism against pyrethroids/DDT in the malaria vector Anopheles funestus across Africa. BMC Genet. 2017;18:76.

    Article  Google Scholar 

  21. Menze BD, Riveron JM, Ibrahim SS, Irving H, Antonio-Nkondjio C, Awono-ambene PH, et al. Multiple insecticide resistance in the malaria vector Anopheles funestus from northern Cameroon is mediated by metabolic resistance alongside potential target site insensitivity mutations. PLoS One. 2016;11:e0163261.

    Article  Google Scholar 

  22. Mulamba C, Riveron JM, Ibrahim SS, Irving H, Barnes KG, Mukwaya LG, et al. Widespread pyrethroid and DDT resistance in the major malaria vector Anopheles funestus in east Africa is driven by metabolic resistance mechanisms. PLoS One. 2014;9:e110058.

    Article  Google Scholar 

  23. Brooke BD, Hunt R, Koekemoer L, Temu E, Taylor M, Small G, et al. Bioassay and biochemical analyses of insecticide resistance in southern African Anopheles funestus (Diptera: Culicidae). Bull Entomol Res. 2001;91:265–72.

    Article  CAS  Google Scholar 

  24. Riveron JM, Ndula M, Barnes KG, Ibrahim SS, Paine MJI, Wondji CS. Directionally selected cytochrome P450 alleles are driving the spread of pyrethroid resistance in the major malaria vector Anopheles funestus. Proc Natl Acad Sci USA. 2013;110:252–7.

    Article  CAS  Google Scholar 

  25. Hargreaves K, Koekemoer LL, Brooke BD, Hunt RH, Mthembu J, Coetzee M. Anopheles funestus resistant to pyrethroid insecticides in South Africa. Med Vet Entomol. 2000;14:181–9.

    Article  CAS  Google Scholar 

  26. Mulamba C, Irving H, Riveron JM, Mukwaya LG, Birungi J, Wondji CS. Contrasting Plasmodium infection rates and insecticide susceptibility profiles between the sympatric sibling species Anopheles parensis and Anopheles funestus s.s.: a potential challenge for malaria vector control in Uganda. Parasit Vectors. 2014;7:71.

    Article  Google Scholar 

  27. Riveron JM, Yunta C, Ibrahim SS, Djouaka R, Irving H, Menze BD, et al. A single mutation in the GSTe2 gene allows tracking of metabolically based insecticide resistance in a major malaria vector Anopheles funestus. Genome Biol. 2014;15:1–20.

    Article  Google Scholar 

  28. Koekemoer LL, Kamau L, Hunt RH, Coetzee M. A cocktail polymerase chain reaction assay to identify members of the Anopheles funestus (Diptera: Culicidae) group. Am J Trop Med Hyg. 2002;66:804–11.

    Article  CAS  Google Scholar 

  29. Crawford J, Guelbeogo W, Sanou A, Traore A, Vernick K, Sagnon N, Lazzaro S. De novo transcriptome sequencing in Anopheles funestus using Illumina RNA-seq technology. PLoS One. 2010;5:e14202.

    Article  CAS  Google Scholar 

  30. Gregory R, Darby AC, Irving H, Coulibaly MB, Hughes M, Koekemoer LL, et al. A de novo expression profiling of Anopheles funestus, malaria vector in Africa, using 454 pyrosequencing. PLoS One. 2011;6:e17418.

    Article  CAS  Google Scholar 

  31. Wondji CS, Irving H, Morgan J, Lobo NF, Collins FH, Hunt RH, et al. Two duplicated P450 genes are associated with pyrethroid resistance in Anopheles funestus, a major malaria vector. Genome Res. 2009;19:452–9.

    Article  CAS  Google Scholar 

  32. Irving H, Riveron J, Ibrahim S, Lobo N, Wondji C. Positional cloning of rp2 QTL associates the P450 genes CYP6Z1, CYP6Z3 and CYP6M7 with pyrethroid resistance in the malaria vector Anopheles funestus. Heredity. 2012;109:383–92.

    Article  CAS  Google Scholar 

  33. David J-P, Strode C, Vontas J, Nikou D, Vaughan A, Pignatelli PM, et al. The Anopheles gambiae detoxification chip: a highly specific microarray to study metabolic-based insecticide resistance in malaria vectors. Proc Natl Acad Sci USA. 2005;102:4080–4.

    Article  CAS  Google Scholar 

  34. Hunt RH, Brooke BD, Pillay C, Koekemoer LL, Coetzee M. Laboratory selection for and characteristics of pyrethroid resistance in the malaria vector Anopheles funestus. Med Vet Entomol. 2005;19(3):271–5.

    Article  CAS  Google Scholar 

  35. Riveron JM, Ibrahim SS, Chanda E, Mzilahowa T, Cuamba N, Irving H, et al. The highly polymorphic CYP6M7 cytochrome P450 gene partners with the directionally selected CYP6P9a and CYP6P9b genes to expand the pyrethroid resistance front in the malaria vector Anopheles funestus in Africa. BMC Genomics. 2014;15:817.

    Article  Google Scholar 

  36. Schmittgen D, Livak J. Analyzing real-time PCR data by the comparative CT method. Nat Protoc. 2008;3:1101–8.

    Article  CAS  Google Scholar 

  37. Thompson JD, Higgins DG, Gibson TJ. CLUSTAL W: Improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res. 1994;22:4673–80.

    Article  CAS  Google Scholar 

  38. Librado P, Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009;25:1451–2.

    Article  CAS  Google Scholar 

  39. Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. MEGA5: Molecular Evolutionary Genetics Analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol. 2011;28:2731–9.

    Article  CAS  Google Scholar 

  40. Duangkaew P, Pethuan S, Kaewpa D, Boonsuepsakul S, Sarapusit S, Rongnoparut P. Characterization of mosquito CYP6P7 and CYP6AA3: differences in substrate preference and kinetic properties. Arch Insect Biochem Physiol. 2011;76:236–48.

    Article  CAS  Google Scholar 

  41. Hemingway J, Hawkes NJ, McCarroll L, Ranson H. The molecular basis of insecticide resistance in mosquitoes. Insect Biochem Mol Biol. 2004;34:653–65.

    Article  CAS  Google Scholar 

  42. Kostaropoulos I, Papadopoulos AI, Metaxakis A, Boukouvala E, Papadopoulou Mourkidou E. Glutathione S-transferase in the defence against pyrethroids in insects. Insect Biochem Mol Biol. 2001;31:313–9.

    Article  CAS  Google Scholar 

  43. Vontas JG, Small GJ, Hemingway J. Glutathione S-transferases as antioxidant defence agents confer pyrethroid resistance in Nilaparvata lugens. Biochem J. 2001;357:65–72.

    Article  CAS  Google Scholar 

  44. Lumjuan N, Rajatileka S, Changsom D, Wicheer J, Leelapat P, Prapanthadara LA, et al. The role of the Aedes aegypti Epsilon glutathione transferases in conferring resistance to DDT and pyrethroid insecticides. Insect Biochem Mol Biol. 2011;41:203–9.

    Article  CAS  Google Scholar 

  45. Strode C, Wondji CS, David JP, Hawkes NJ, Lumjuan N, Nelson DR, et al. Genomic analysis of detoxification genes in the mosquito Aedes aegypti. Insect Biochem Mol Biol. 2008;38:113–23.

    Article  CAS  Google Scholar 

  46. Wondji CS, Morgan J, Coetzee M, Hunt RH, Steen K, Iv WCB, et al. Mapping a Quantitative Trait Locus ( QTL ) conferring pyrethroid resistance in the African malaria vector Anopheles funestus. BMC Genomics. 2007;8:34.

    Article  Google Scholar 

  47. Amenya D, Naguran R, Lo T, Ranson H, Spillings B, Wood O, et al. Over expression of a cytochrome P450 (CYP6P9) in a major African malaria vector, Anopheles funestus, resistant to pyrethroids. Insect Mol Biol. 2008;17:19–25.

    Article  CAS  Google Scholar 

  48. Akogbeto MC, Djouaka R, Noukpo H. Use of agricultural insecticides in Benin. Bull Soc Pathol Exot. 2005;98:400–5.

    CAS  PubMed  Google Scholar 

  49. Yadouleton A, Martin T, Padonou G, Chandre F, Asidi A, Djogbenou L, et al. Cotton pest management practices and the selection of pyrethroid resistance in Anopheles gambiae population in northern Benin. Parasit Vectors. 2011;4:60.

    Article  Google Scholar 

  50. Djouaka RF, Bakare AA, Bankole HS, Doannio JM, Coulibaly ON, Kossou H, et al. Does the spillage of petroleum products in Anopheles breeding sites have an impact on the pyrethroid resistance? Malar J. 2007;6:159.

    Article  Google Scholar 

  51. Elissa N, Mouchet J, Riviere F, Meunier JY, Yao K. Resistance of Anopheles gambiae s.s. to pyrethroids in Côte d’Ivoire. Ann Soc Belg Med Trop. 1993;73:291–4.

    CAS  PubMed  Google Scholar 

  52. Nkya TE, Akhouayri I, Kisinza W, David J-P. Impact of environment on mosquito response to pyrethroid insecticides: facts, evidences and prospects. Insect Biochem Mol Biol. 2013;43:407–16.

    Article  CAS  Google Scholar 

  53. Jones CM, Toé HK, Sanou A, Namountougou M, Hughes A, Diabaté A, D R, Simard F, et al. Additional selection for insecticide resistance in urban malaria vectors: DDT resistance in Anopheles arabiensis from Bobo-Dioulasso, Burkina Faso. PLoS One. 2012;7:e45995.

    Article  CAS  Google Scholar 

  54. Tang A, Tu C. Biochemical characterization of Drosophila glutathione S-transferases D1 and D21. J Biol Chem. 1994;269:27876–84.

    CAS  PubMed  Google Scholar 

  55. Brandt A, Scharf M, Pedra J, Holmes G, Dean A, Kreitman M, et al. Differential expression and induction of two Drosophila cytochrome P450 genes near the Rst(2)DDT locus. Insect Mol Biol. 2002;11:337–41.

    Article  CAS  Google Scholar 

  56. Barnes KG, Weedall GD, Ndula M, Irving H, Mzihalowa T, Hemingway J, et al. Genomic footprints of selective sweeps from metabolic resistance to pyrethroids in African malaria vectors are driven by scale up of insecticide-based vector control. PLoS Genet. 2017;13:e1006539.

    Article  Google Scholar 

Download references

Acknowledgements

We appreciate Kpome community for their cooperation during fieldwork. We thank Innocent Djegbe, Eric Tossou, Gareth Weedall, Claude Gande and Murielle Soglo for their technical assistance and relevant advice.

Funding

The Wellcome Trust supports this work grants Reference 099864/Z/12/Z awarded to RD and a Wellcome Trust Senior Research Fellowship in Biomedical Sciences to CSW (101893/Z/13/Z).

Availability of data and materials

All data generated or analysed during this study are included in this published article and its Additional file.

Author information

Authors and Affiliations

Authors

Contributions

RD and CSW designed the study. GT and RA carried out mosquito collection and GT, RA, AM and JR reared mosquitoes and performed WHO bioassays. GT and HI performed microarray, qRT-PCR analyses and sequencing of resistance genes. GT, JMR and CSW analyzed data. AY gave advise on the study design and contributed to the implementation of the study. GT, RD, AY, JMR and CSW wrote the manuscript. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Genevieve Tchigossou.

Ethics declarations

Ethics approval and consent to participate

Verbal consent was obtained from household heads before mosquito collections.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional file

Additional file 1:

Table S1. The most upregulated genes in Rperm-S, C-S and Rperm-C. Table S2. The most upregulated genes in RDDT-S, Rperm-S and C-S comparisons. (DOCX 56 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Tchigossou, G., Djouaka, R., Akoton, R. et al. Molecular basis of permethrin and DDT resistance in an Anopheles funestus population from Benin. Parasites Vectors 11, 602 (2018). https://doi.org/10.1186/s13071-018-3115-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13071-018-3115-y

Keywords