Mapping a Quantitative Trait Locus (QTL) conferring pyrethroid resistance in the African malaria vector Anopheles funestus

Linkage mapping

Data from the eleven families were used to build a combined genetic linkage map by genotyping a total of 349 F₂ individuals at 56 loci (Figure 1). Three linkage groups were resolved instead of four in the previous map [18], with each linkage group corresponding to one of the three chromosomes. The linkage map of the X chromosome with 6 markers was 26 cM in length compared to 44 cM for five markers previously. 25 markers were mapped to chromosome 2 for a genetic distance of 64 cM compared to 158 cM for 16 markers previously. One linkage group was resolved for chromosome 3 instead of two previously with 22 markers giving a genetic distance of 55 cM which is well below the total distance of 189 cM for 18 markers observed in the first map. The average resolution for each chromosome is 4.3, 2.5 and 2.5 cM respectively for chromosome X, 2 and 3. The total length of the current map is 145 cM with an average resolution of 2.7 cM/marker.

QTL mapping

Genotype-phenotype associations

We used χ² goodness-of fit tests to identify loci statistically (P < 0.05) associated with pyrethroid resistance in each family, the null hypothesis being that resistance or susceptibility to pyrethroids is equal in each genotype class. Significant loci are indicated in Table 1.

Significant loci	Chi-square	P-values
Family 1
7P6P4	5.92	0.015
AChE	5	0.03
Family 2
7P6P4	6.02	0.013
Family 3
7P6P4	6.32	0.011
Family 4
7P6P4	7.33	0.019
3P6P3	5.9	0.032
AChE	5.5	0.035
AFND6	5.3	0.039
Family 6
3BU82	5.87	0.03
6BU40	5.69	0.034
Family 7
7P6P4	7	0.02
6BU40	5.8	0.034
Family 10
No significant marker
Family 10X
7P6P4	6.5	0.01
3P6P3	6	0.022
Family 11
7P6P4	9.78	0.004
3P6P3	7.46	0.015
AFND6	6.53	0.028
Family 12
No significant marker
Family 16
7P6P4	7.45	0.01
3P6P3	6.9	0.029

For family 4, we found significant associations between four loci and resistance to pyrethroids. Two of these loci (7P6P4 and 3P6P3) were semi-informative (one of the P₁ parent being heterozygote while the other is homozygote for the locus) for this family and the others (AChE and AFND6) were fully informative. For 7P6P4, 83% of F₂ individuals with an "aa" genotype (genotype of the parental susceptible female) and 30% of F₂ heterozygous "ab" died after 24 h exposure to permethrin. These values are significantly different and also oriented in the anticipated direction. This trend is repeated for marker 3P6P3 where 100% of homozygous "aa" and 44% of heterozygous "ab" F₂ individuals died after exposure to permethrin. For locus AChE, mortality rates of 80% for F₂ "aa" homozygous individuals, 64% for heterozygous "ab" and 23% for homozygous "bb" (genotype of the parental resistant male) were obtained. A similar pattern is seen for locus AFND6 (Figure 2).

For family 11, statistically significant associations were found between genotypes at three loci (7P6P4, 3P6P3 and AFND6) and the resistance trait. These loci showed a similar pattern of association with pyrethroid resistance as in family 4. A high mortality rate was observed in "aa" homozygous F₂ individuals (100, 60 and 100% respectively for 7P6P4, 3P6P3 and AFND6) while mortality was low in "bb" individuals (19 and 39%, respectively for 7P6P4 and AF6).

For family 6, only two markers (3BU82 and 6BU40) were statistically associated with pyrethroid resistance. Again, the mortality rate was higher in F₂ "aa" homozygotes (65 and 67% respectively for 3BU82 and 6BU40) than in F₂ "bb" homozygotes (25 and 35% respectively for 3BU82 and 6BU40).

This analysis indicates some important factors linked to the inheritance of pyrethroid resistance trait. All four loci associated with the resistance trait are located on chromosome 2R (between divisions 9 and 12) suggesting that a locus or loci that affect resistance to pyrethroids is located within this region of the genome. Furthermore, plotting resistance against genotypes (Figures 2 and 3) also suggests that alleles at these loci are additive in their effect on resistance.

We next used interval mapping (IM), composite-interval mapping and multiple-interval mapping to predict the location of this QTL. Three families with a sample size of more than 40 F₂ progeny (families 4, 6 and 11) were independently analyzed for the presence of QTL associated with pyrethroid resistance.

For family 4, a QTL was detected at the end of chromosome 2R with both IM and CIM (Figure 4) between markers 7P6P4 and 3P6P3 already found to be associated with pyrethroid resistance by χ² goodness-of fit tests. The LOD scores associated with this QTL are 3.3 and 4 by IM and CIM respectively. No QTL were detected on chromosome X and chromosome 3. Multiple-interval mapping confirmed the presence of this QTL (Table 2) and estimated that its genetic variance (σ_g²) accounted for 27.9% of the phenotypic variance (σ_p²) for pyrethroid resistance (27.9% resulting from 31.1% for additive effect and -3.2% for dominance effect). Most of the genetic variance (σ_g²) was attributable to additive effects (31.1%, Table 2). Two additional QTLs were detected in this family after refining the MIM model by searching and testing for new QTLs. These were located at 20 cM on chromosome X between markers G17 and FUNQ and at position 75.2 cM on chromosome 3 between markers ND10 and FUB1 (Table 2). These QTL were not detected by any other analysis method, nor were they found in other single family analysis.

σg2(% σp2)	σe2(% σp2)	σp2	Nearest marker	Genetic distance (cM)	LOD	Effect	% σg2
Family 4
14.72 (58.87%)	10.28 (41.13%)	25	7P6P4	0.1 (0.0–1.4)	1.09 0.18	A, -5.589 D, -0.660	31.1 -3.2
			G17	21 (20–22)	0.77 1.09	A, 2.1649 D, 3.6980	14.8 7.8
			ND10	75.2 (75.0–76.2)	0.27 0.95	A, 1.0312 D, 3.2745	1.4 7.1
Family 11
15.86 (63.4%)	9.14 (36.6%)	25	7P6P4	7.3 (7.0–9.8)	2.26 0.02	A, -4.825 D, +2.762	49.8 13.6
Family 6
2.428 (9.71%)	22.57 (90.29%)	25	3BU82	14.1 (14–16)	0.38 0.12	A, -2.4513 D, -0.7513	10.4 -0.8

σ_g², σ_p², σ_e² respectively for genetic, phenotypic (in parentheses) and environmental variance; A for additive; D for dominance; confidence intervals for QTL position are in parentheses below the position estimate.

For family 11, using IM and CIM, we detected a QTL on chromosome 2R between markers 7P6P4 and 3P6P3, at a similar position to that in family 4. This QTL was located at 7.3 cM with 7P6P4 as the nearest marker. The LOD scores for this QTL were 5 and 5.8 with IM and CIM respectively (Figure 5). With MIM, only one significant QTL was detected which, in this family, accounted for 63.4% of phenotypic variance (σ_p²) (Table 2). Additive effects accounted for the majority of the genetic variance. We have named this QTL rp1 for r esistance to p ermethrin 1. No QTL was detected on chromosome 3 and no genetic linkage was constructed for chromosome X in this family due to a lack of informative markers on this chromosome.

In families 6 and 10 (data not shown in Table 2 for family 10 because of the low sample size), the markers 7P6P4 and 3P6P3 were not informative and only low LOD scores were observed at markers at the end of 2R chromosome. However when using multiple-interval mapping, we detected the same QTL as that observed in families 4 and 11. The LOD score for this QTL in family 6 and 10 was 1.9 and 1.6 respectively with composite-interval mapping, which is considerably smaller than that observed for families 4 and 11. This is probably due to the small sample size (family 10) and the lack of informative markers in this region of the genome. In families 4 and 11, as in family 6 (Table 2), additive effects accounted for the majority of the genetic variance.

Physical mapping of rp1 QTL

The genes implicated in insecticide resistance in other species (cytochrome P450s, kdr, AChE) were physically mapped to the An. funestus polytene chromosomes in order to locate possible candidate resistance genes within the boundaries of the QTL. The location of these genes is shown in Figure 6 and the map position of these genes and of their orthologs in An. gambiae is shown in Table 3. There is a high level of synteny between An. funestus and An. gambiae in the location of these genes. With the exception of CYP6M1 and CYP6P1, all of the An. funestus probes hybridized to the expected chromosome arms according to the established synteny between the two species [19]. CYP6M1 in An. funestus is found on chromosome 2R instead of the predicted 2L and CYP6P1 was located on 2L arm and not 2R as expected.

Five CYP6 P450s physically mapped to chromosome arm 2R, within the region of rp1 QTL (Figure 6). These included CYP6P4 and CYP6P3 from which two SNP loci, (7P6P4 and 3P6P3) were identified. In agreement with the physical mapping data, we found that these SNPs were closely linked to rp1 in all families for which these markers were informative.

The kdr locus is located on chromosome 3R, division 36, and is not linked to pyrethroid resistance in the FUMOZ-R strain of An. funestus. Again this agrees with the chromosomal location determined for this gene in A. gambiae (3R in An. funestus is equivalent to 2L in A. gambiae). No amino acid substitutions were found in the amplified 1342 bp fragment of the voltage gated sodium channel gene in the FOMUZ-R strain compared to the field samples of An. funestus and the classical kdr leucine substitution at position 1014, associated with pyrethroid resistance in many insect species, was not found in An. funestus.

The fragment of the AChE gene sequenced encompasses the common ace-1 mutation site associated with insensitive AChE in An. gambiae [20]. We did not detect this mutation in the An. funestus AChE.

Mosquito strains and bioassays

The two strains of An. funestus used in this study were colonized by the Vector Control Reference Unit of the National Institute for Communicable Diseases in South Africa [17]. The FUMOZ strain originated from southern Mozambique in 2000 and was initially a heterozygous permethrin resistant population. After selection of the parental strain with 0.1% lambda cyhalothrin [17], a highly resistant strain called FUMOZ-R was generated. The FANG colony which is completely susceptible to all pyrethroids was colonized from Calueque in southern Angola in 2002.

Mosquito crosses

Reciprocal crosses using virgin FANG and FUMOZ-R females with males from the alternative strain were set-up. Blood fed, mated females were left to oviposit singly and the F1 were crossed to generate F₂ progeny, producing iso-female lines. Bioassays were carried out on 1–3-day-old adults using the standard WHO bioassay procedure [29]. Mosquitoes were exposed to 0.75% permethrin for 1 hour and mortality was recorded 24 h post-exposure. Dead mosquitoes were considered as susceptible to permethrin and separated from the living progeny which were considered as resistant. Parental females (P₁), F₁ progeny as well as F₂ progeny from each family were collected and kept in silica gel Eppendorf tubes for later DNA extraction and genotype determination.

Genotyping of molecular markers

Single Nucleotide Polymorphisms (SNPs) genotyping

Fifteen SNP markers were previously identified in An. funestus [18]. These were scored in the P₁ and F₁ parents in each of the 11 families and informative SNPs were scored either by the heated ligation oligonucleotide assay (HOLA) [31] or by single-base pair extension reaction using terminator dyes and CEQ8000 software or using the pyrosequencing method, as described below. An additional 11 informative SNP markers were identified by direct sequencing of PCR products from cytochrome P450 genes or from genes that had been previously mapped to An. funestus polytene chromosomes [19].

Pyrosequencing reactions

Single-stranded biotinylated PCR products were prepared for the pyrosequencing reaction using a Vacuum Prep Tool (Biotage AB). The biotinylated PCR products were immobilized onto Streptavidin Sepharose high performance beads. For a single sample, 3 μl of Streptavidin Sepharose™ HP (Amersham) were added to 37 μl binding buffer (10 mM Tris-HCl pH 7.6, 2 M NaCl, 1 mM EDTA, 0.1% Tween 20) and mixed with 25 μl PCR product and 15 μl water on a mechanical shaker for 5 min at room temperature in a 96-well plate. The beads containing the immobilised templates were captured onto the filter probes after applying the vacuum and then washed with 70% ethanol for 5 sec, denaturizing solution (0.2 M NaOH) for 5 sec and washing buffer (10 mM Tris-Acetate pH 7.6) for 5 sec. The vacuum was switched off and the beads were released into a PSQ 96 well plate containing 38.4 μl annealing buffer (20 mM Tris-Acetate, 2 mM MgAc2 pH 7.6) and 1.6 μl of the corresponding sequencing primer (10 μM). Annealing was achieved by heating the samples to 80°C for 3 min and then allowing them to cool to room temperature. Pyrosequencing reactions were performed according to the manufacturer's instructions using the PSQ 96 SNP Reagent Kit (Biotage AB) and the genotype was determined using SNP Software (Biotage AB).

Linkage mapping

The JoinMap 2.0 package [32] was used to build a genetic linkage map for each individual family and for the combined genotype data from all families. Genotype data for each marker were tested for conformity to Mendelian ratio with a χ² goodness-of-fit analysis using the JMSLA procedure. Loci were separated into linkage groups using the JMGRP1 and JMSPL procedures with minimum and maximum LOD thresholds of 0.0 and 6.0 respectively and LOD increments of 0.1. The JMREC program estimated pairwise cM distances between all pairs of informative loci in each linkage group. The JMMAP program estimated the maximum likelihood map using the Kosambi distances. DrawMap 1.1 [33] software was then used to plot the genetic map.

QTL mapping

Associations between genotypes at each locus and the resistance phenotype (dead or alive) were assessed using a contingency χ² analysis. The null hypothesis was that susceptibility to pyrethroid is equal in each genotype class. The marginal probabilities were the frequencies of each genotype at a locus and the mortality and survival rates after pyrethroid exposure. For loci with a significant χ² we analyzed the inheritance of the alleles at these loci. The a priori hypothesis was that higher mortality rate would occur among F₂ individuals with one or both alleles inherited from the susceptible parent while lower mortality rate will be observed among individuals with alleles inherited from the resistant parent.

The JoinMap linkage map and the genotype/phenotype data were entered into Windows QTL Cartographer 2.5 [34]. Interval mapping [22], composite-interval mapping [24] and multiple-interval mapping [35] procedures were implemented for each family. Interval mapping uses two observable flanking markers to construct an interval within which to search for QTL. The optimum LOD thresholds were estimated by permutation of trait and marker data 1000 times with a walking speed of 2 cM. Composite-interval mapping (CIM) tests whether an interval between two markers contains a QTL while simultaneously controlling for the effect of proximal QTLs located outside the interval. CIM was performed using the standard model with a control marker number of n _p = 5 and a window size of w _s = 10 cM. We used 1000 permutation analysis to determine the optimum significance threshold of the LOD. Multiple-interval mapping (MIM) analyzes multiple marker intervals simultaneously to fit multiple putative QTLs. An initial MIM model was estimated by forward and backward marker selection with a probability of a partial R² set to 0.01. We then optimized QTL positions, searched for new QTLs, tested for existing QTL before saving the QTL model. The MIM model summary procedure estimated additive and dominance effects, the QTL likelihood ratio, the confidence interval around QTL positions, and partitioned the variance explained by QTL.

Cloning of partial sequences of sodium channel (kdr) and acetylcholinesterase (AChE) genes in An. funestus

We cloned partial fragments of the two major target sites of insecticides used in malaria control, the acetylcholinesterase gene (AChE) which is the target site of carbamate and organophosphate insecticides and the voltage gated sodium channel gene (kdr) which is the target site of pyrethroids and DDT. Primers designed from a conserved region of the S6 region of the voltage-gated sodium channel gene AAK2 (5'-TGC GGA GAA TGG ATY GAA TC-3') and Ask2 (CTT AGC CTT GCT TTT GTC AAA) were used to amplify a 200 bp fragment of this gene from An. funestus genomic DNA. The PCR conditions were 96°C for 1 min followed by 30 cycles of 94°C for 30 s, 60°C for 30 s and 72°C for 30 s with a final 10 min extension at 72°C. The PCR product was ligated into pGEM-T easy vector and sequenced with M13 forward and reverse primers. From this 200 bp product, the Adf2 primer (5'-TGC AAA ATA GAG TCA TTG GTG AA-3') was designed and used with the Dip3 primer (5'-ATC ATC TTC ATC TTT GC-3') [36] to amplify a larger fragment of 1342 bp. This fragment was amplified and sequenced from individuals from the resistant FUMOZ-R strain and from field collected populations from Kenya, Mozambique and Malawi. A similar protocol was used to amplify a 1487 bp partial fragment of the acetylcholinesterase gene (AChE) in An. funestus using primers ANace-F1 (5'-CCG GGG GCG ACT ATG TGG AAC-3') and ANace-R3 (5'-GTT GCT GTT CGG GTT GTC CG-3'). The partial sodium channel and acetylcholinesterase of An. funestus sequences have been deposited in GenBank, accession numbers DQ534436 and DQ534435 respectively for kdr and AChE.

Polytene chromosome in situ hybridization

A subset of 12 partial P450 Cytochrome P450 genes were chosen for physical mapping as increased activity of the P450 enzymes has been demonstrated in the FUMOZ-R resistant strain [17] and because this is an important pyrethroid resistance mechanism in other insects [27, 37]. Half-gravid An. funestus females were collected from Kela village in Chikwawa district in southern Malawi. Ovaries were removed and preserved in Carnoy's fixative solution (3 ethanol: 1 glacial acetic acid by volume) and stored at -20° until use. The identity of specimens used in this study was confirmed by a rDNA diagnostic PCR which distinguishes four members of the An. funestus group [38]. Probes were prepared with the GIBCO BRL in situ hybridization and detection system, using the manufacturer's recommended protocol. Using fluorescence in situ hybridization (FISH), An. funestus clones were mapped on the five arms of the polytene chromosome. After hybridization, chromosomes were washed in 0.2× SSC and counterstained with YOYO-1 and mounted in DABCO antifade solution (Sigma) [19]. Fluorescent signals were detected using Zeiss LSM Pascal scanning confocal microscope.