Mosquito abundance, biting rate and morphological identifications
A total of 3920 Anopheles sp. females, 1167 and 2753 during the dry and rainy seasons respectively, were captured on a total of 60 collection days. Overall 81% (3187/3920) of the samples were collected in the cow odor-baited double net traps (CBNTs) and while this relative abundance was rather consistent between different collection sites for the CBNTs, 67% (490/733) of the Anopheles from the human odor-baited double net traps (HBNTs) were collected in the forest sites (Table S1).
The biting rate (# of females/trap/day) for the HBNTs was consistently higher in the forest sites compared to all other locations during both the rainy and the dry seasons (Table 1). However, for the CBNTs, while the biting rate was the highest in the forest sites during the dry season, the tendency changed during the rainy season with a higher biting rate in the villages and the forests near the villages compared to the forest sites (Table 1).
A total of 3131 females were morphologically identified as 14 different Anopheles species or complexes of morphologically indistinguishable sibling species. Based on these morphological identifications, species thought to be primary vectors comprised only 10.2% of the collected mosquitoes: Anopheles dirus s.l. (8.1%, n = 319), A. minimus s.l. (0.4%, n = 15) and A. maculatus s.l. (1.7%, n = 67). The most abundant species (represented by more than a hundred individuals in our collection) constituted 75.8% of the collected Anopheles mosquitoes and were represented by 6 species complexes: A. barbirostris (21.2%, n = 831), A. philippinensis (14.6%, n = 571), A. hyrcanus (13.6%, n = 535), A. kochi (10.5%, n = 412), A. dirus (8.1%), A. aconitus (7.7%, n = 303).
Molecular determination of mosquito species
A total of 844 females were molecularly characterized for species in the random subset and represent 26 distinct Anopheles species as determined by ITS2 and CO1 (Table S2). The most abundant species (representing ≥ 5% of the samples; n ≥ 42) comprise 77.8% of the molecularly typed individuals and represent 8 species from 6 different species complexes. These most abundant species included A. dirus (13.2%, n = 112) from the Dirus complex. From the Barbirostris complex; A. dissidens (13.2%, n = 112), and A. campestris-wejchoochotei (8.1%, n = 69). From the Hyrcanus Group, A. peditaeniatus (12.8%, n = 108), and A. nitidus (5.7%, n = 48). The Annularis, Funestus, and Kochi Groups were each represented by a single species A. nivipes (9.2%, n = 78), A. aconitus (6.7%, n = 57), and A. kochi (8.7%, n = 74), respectively. The 18 less abundant species, represent by fewer than 42 samples and in many cases just a handful of samples included A. philippinensis (n = 17) and A. annularis (n = 1) from the Annularis Group, A. jamesii (n = 16), A. pseudojamesi (n = 1), and A. splendidus (n = 1) from the Jamesii Group and A. saeungae (n = 29) and A. barbirostris (n = 2) from the Barbirostris Group. From the Hyrcanus Group A. crawfordi (n = 40), A. argyropus (n = 1), An. nigerrimus (n = 28), and A. sinensis (n = 3) were sampled. Anopheles maculatus (n = 22), A. sawadwongporni (n = 4), and A. rampae (n = 2) from the Maculatus Group. Anopheles tessellatus from the Tessellatus Group and A. interruptus from the Asiaticus Group were each sampled once. A. vagus (n = 12) and A. karwari (n = 3) were also present. There were 2 mosquitoes that had 99.9% identical ITS2 and 99.4% identical CO1 sequences but matched no species in the NCBI database. In addition to the random subset, 79 Plasmodium sp. infected samples were molecularly characterized for species which resulted in a total of 29 Anopheles species as determined by ITS2 and CO1.
Day biting rate
Overall 20.2 ± 1.2% of the Anopheles females were captured during the daytime (between 06:00 and 18:00). Indeed, while the majority of Anopheles mosquitoes bite at night, an important proportion was active during the day (Fig. 2). Excluding species with extremely low sample sizes and unidentified samples, day biting behaviour was observed for all the Anopheles species and varied from 13 to 68% (Table S3).
Average number of Anopheles females collected per trap per hour in the different collection sites in the HBNTs and the CBNTs.
The day biting rate in the HBNTs was not different across collection sites (19.6 ± 2.9%; χ2 = 3.6, df = 3, p = 0.3; Fig. 3a) but was higher during the dry season (25.9 ± 4.6%) compared to the rainy season (13.8 ± 3.5%; χ2 = 19.08, df = 1, p < 0.001). There was an overall trend for a higher day biting rate during the dry season compared to the rainy season for all sites except in the forest near the villages where an opposite tendency was observed, however the interaction between collection sites and seasons was not significant (χ2 = 7.5, df = 3, p = 0.057).
Day biting rate (proportion of Anopheles females collected between 6am and 6 pm) of mosquitoes captured in the different collection sites and seasons in (a) HBNTs and (b) CBNTs. The day biting rate was significantly higher during the dry season compared to the rainy season in both kind of traps (HBNTS χ2 = 19.08, df = 1, p < 0.001 and CBNTS χ2 = 115.53, df = 1, p < 0.0001). Data show average proportion ± 95% confidence intervals.
The day biting rate in the CBNTs was higher during the dry season (37.7 ± 3.3%) compared to the rainy season (14.5 ± 1.4%; χ2 = 115.53, df = 1, p < 0.0001). There was no significant difference across the collection sites (20.4 ± 1.4%; χ2 = 0.65, df = 2, p = 0.72; Fig. 3b) nor of the interaction between collection sites and seasons (χ2 = 1.97, df = 2, p = 0.37).
Malaria infection rates
Samples were screened to detect Plasmodium sp. in the head and thorax, thus indicating mosquitoes likely to be infectious due to carriage of sporozoites in the salivary glands (Fig. 4). Over the 3919 mosquitoes tested, 141 positive samples were found (3.6%; 122 P. vivax, 13 P. falciparum, 4 mixed P. falciparum and P. vivax, 1 mixed P. falciparum and P. vivax and P. ovale, 1 mixed P. falciparum and P. inui). In addition, 2 mosquitoes were found infected by P. cynomologi and 3 by P. inui, these non-human malaria parasites were not included in the following analyses. Malaria infection rate in HBNTs (4.77 ± 1.5%) was not significantly higher than in CBNTs (3.33 ± 0.62%; χ2 = 3.19, df = 1, p = 0.07). In the HBNTs, malaria infection rate was significantly higher during the rainy season (8.2 ± 2.8%) compared to the dry season (1.13 ± 1.1%; χ2 = 8.57, df = 1, p = 0.003; Odds Ratio = 5.04 [1.53–16.53]). Plasmodium sp. prevalence was not significantly correlated with collection site (χ2 = 1.89, df = 3, p = 0.60), nor with time of sampling (χ2 = 0.22, df = 1, p = 0.64).
Malaria parasite prevalence in mosquitoes captured across different collection sites and seasons in (a) HBNTs and (b) CBNTs. Data show average proportion ± 95% confidence intervals. Numbers below the bars indicate the total number of mosquitoes. Zero in the graph indicates no infected mosquitoes. Malaria parasite prevalence was significantly higher during the dry season compared to the rainy season in the HBNTS (χ2 = 8.57, df = 1, p = 0.003) whereas the opposite was observed in the CBNTS (χ2 = 5.44, df = 1, p = 0.02). In the CBNTS, malaria parasite prevalence was significantly affected by daytime (χ2 = 8.13, df = 1, p = 0.004) and the interaction between collection sites and seasons (χ2 = 6.94, df = 2, p = 0.03).
In the CBNTs, the prevalence of Plasmodium-positive mosquitoes was significantly higher during the dry season (5.66 ± 1.6%) compared to the rainy season (2.53 ± 0.6%; χ2 = 5.45, df = 1, p = 0.02). Malaria infection rate was not significantly different across collection sites (χ2 = 5.87, df = 2, p = 0.053), although more infectious mosquitoes were collected in the forests (5.55 ± 1.32%) compared to the forest near the villages (2.86 ± 1%) and the villages (1.1 ± 0.7%). Malaria infection rate was higher in mosquitoes collected during daytime (6 ± 1.8%) compared to night time (2.64 ± 0.6%; χ2 = 8.13, df = 1, p = 0.004) regardless of the collection sites (daytime: collection site interaction: χ2 = 3.28, df = 2, p = 0.19). Malaria infection rate was significantly influenced by the interaction between collection sites and seasons (χ2 = 6.89, df = 2, p = 0.03) with opposite patterns in the malaria infection rates observed during the dry and the rainy seasons in the forests and the forest near the villages (χ2 = 6.51, df = 1, p = 0.032, all other post-hoc comparisons were non-significant; Fig. 4).
The infectious mosquitoes collected in the HBNTs belonged to 8 different mosquito species (molecular identification; Table 2). Interestingly, an additional 12 Anopheles species were also found infected by human malaria parasites among the mosquitoes collected in the CBNTs (A. aconitus, A. campestris-wejchoochotei, A. crawfordi, A. jamesii, A. karwari, A. minimus, A. nivipes, A. pallidus, A. rampae, A. saeungae, A. sawadwongporni, A. vagus), bringing to 20 the total number of malaria vectors collected in this study among the 29 species molecularly identified (Table 2).
Mosquito host preferences
We assessed the inherent mosquito host preference by calculating the anthropophily index (AI) for the different species, defined as the number of Anopheles sp. individuals caught in the human-baited trap over the total number of mosquitoes caught in both human- and calf- baited traps (see Methods). Most of the collected mosquitoes were zoophilic (AI < 50) or generalist (AI ~ 50) including A. dirus or A. maculatus mosquitoes (Table S4).
Anopheles dirus—A total of 319 An. dirus females were collected, mainly in the forest sites (n = 303). Human Plasmodium sp. prevalence was 8.46 ± 3% overall and 7.9 ± 3% in the forest sites. The anthropophily index was calculated in the forest sites only and was 48.51 ± 5.63%. Infectious mosquitoes were significantly more anthropophilic (87.5 ± 13.23%) than uninfected mosquitoes (45.16 ± 5.84%; χ2 = 21.56, df = 1, p < 0.0001). Anopheles dirus females were more anthropophilic during the dry season (51.83 ± 7.65%) compared to the rainy season (44.6 ± 8.26%; χ2 = 7.19, df = 1, p = 0.008). There was no difference in the anthropophily index of A. dirus collected during the day or night time (χ2 = 3.39, df = 1, p = 0.066).
Anopheles kochi—A total of 412 A. kochi females were collected with more than half collected in the forest sites (n = 257). Human Plasmodium sp. prevalence was 3.15 ± 1.69% overall and 4.67 ± 2.58% in the forest sites. There was no significant effect of season (χ2 = 1.98, df = 1, p = 0.16) or infection status (χ2 = 0.67, df = 1, p = 0.41) on the anthropophily index of A. kochi. Anopheles kochi females were significantly more anthropophilic during the daytime (25.64 ± 9.69%) compared to night time (8.98 ± 3.07%; χ2 = 8.53, df = 1, p = 0.003).
Out of the random subset 26 species were molecularly identified, among which 19 species were collected in both HBNTs and CBNTs although in variable proportions (Fig. 5). Indeed, A. annularis (n = 2), A. argyporus (n = 1), A. splendidus (n = 1) and A. tessellatus (n = 1) were only observe in the CBNTs, whereas A. barbirostris (n = 2), A. interruptus (n = 1) and A. pseudojamesi (n = 1) were only observed in the HBNTs. Among the 19 species collected in both type of traps, 3 species were enriched in the CBNTs—A. nivipes (n = 78, χ2 = 8.332, df = 1, p = 0.039), A. peditaeniatus (n = 108, χ2 = 17.334, df = 1, p < 0.0001), A. kochi (n = 74, χ2 = 4.944, df = 1, p = 0.026)—whereas 2 species were enriched in the HBNTs—A. campestris-weichoochetei (n = 69, χ2 = 40.678, df = 1, p < 0.0001), A. dirus (n = 110, χ2 = 16.413, df = 1, p < 0.0001).
Proportions of the molecularly identified samples in the CBNTS and HBNTS.
Insecticide resistance mutations
192 mosquitoes were genotyped by amplicon sequencing for known mutations associated with resistance to different classes of insecticides: pyrethroids, carbamates, organochlorides, and organophosphates. We genotyped previously characterized mutations in ace-1, kdr, and rdl genes, G119S45, L1014F45 and A296S52, respectively. Of the 192 individuals, 187 were successfully genotyped for G119S, 189 for L1014F, and 191 for A296S. The few missing genotypes were due to PCR amplification failure. Twenty-seven individuals (14%) carried the G119S mutation, none of which were heterozygotes. Nineteen (10%) individuals carried the L1014F mutation, five of which were heterozygotes. Twenty-nine (15%) individuals carried the A296S mutation, seven of which were heterozygotes. There were four species (of 109 total) that carried resistant mutations. Two A. dirus (of 31 total) carried the rdl A296S mutation, a single A. kochi individual carried the ace-1 G119S mutation and an additional A. kochi individual carried the rdl A296S mutation (18 A. kochi individuals total). Two A. nivipes individuals carried the rdl A296S mutation (21 total). Interestingly, a single species, A. peditaeniatus was significantly overrepresented for all resistant genotypes (ace-1: χ2 = 179, df = 1, p = 8.01 × 10–41, kdr: χ2 = 132.4, df = 1, p = 1.22 × 10–30, rdl: χ2 = 168.1, df = 1, p = 1.92 × 10–38). A. peditaeniatus accounted for 96% of individuals carrying the ace-1 mutation (26/27), 100% of individuals with the kdr mutation (19/19), and 83% of individuals with the rdl mutation (24/29); despite the fact that it comprised only 13% (26/192) of the sample set genotyped. To gain a better understanding of this apparent enrichment of known insecticide resistance mutations in A. peditaeniatus, we returned to the larger sample set (n = 844) and genotyped the additional 82 individuals not included in the representation plates, for a total of 108 genotyped A. peditaeniatus. When resistant allele frequencies were compared between the A. peditaeniatus in the representation plates (n = 26) and all other A. peditaenitus (n = 82), there was no statistically significant difference (ace-1: χ2 = 1.69, df = 1, p = 0.2, kdr: χ2 = 0.01, df = 1, p = 0.91, rdl: χ2 = 0.78, df = 1, p = 0.38), indicating that the representation plates accurately represented the larger dataset.
Of the 108 total A. peditaeniatus, 106 carried the ace-1 G119S mutation (3 heterozygotes; 103 homozygotes), 78 carried the kdr L1014F mutation (22 heterozygotes; 56 homozygotes), and 104 carried the rdl A296S mutation (12 heterozygotes, 92 homozygotes; Fig. 6). The presence of the ace-1 mutation, either homozygote or heterozygous was also a strong predictor of whether that individual was also resistant to other insecticides as 97% of individuals that carried the ace G119S mutation also had at least one other resistance mutation, homozygote or heterozygote. For the L1014F mutation, the distribution of alleles by ecological site was not uniform (Fig. 7). The L1014F kdr resistant allele was significantly more common in the Forests Near Village 2 and 3 (χ2 = 4.0, df = 1, p = 0.04 and χ2 = 5.1, df = 1, p = 0.02 respectively) as compared to all other sites. Conversely, Forest 2 had a significantly greater proportion of the susceptible kdr L1014F allele (χ2 = 4.5, df = 1, p = 0.03). The frequency of the ace-1 G119S mutation was significantly greater than the other mutations (ace-1 vs kdr: χ2 = 76, df = 1, p = 2.76 × 10−18, ace vs rdl: χ2 = 7.8, df = 1, p = 5.3 × 10−3). Additionally, the frequency of the A296S rdl resistant allele was greater than the frequency of the L1014F kdr resistant allele (χ2 = 45.2, df = 1, p = 1.8 × 10−11). All heterozygotes were verified by amplifying and resequencing the necessary fragment and visualizing the resultant sequence trace files.
Insecticide resistance mutations are very frequent in Anopheles peditaeniatus. All A. peditaeniatus (n = 108) in the molecularly typed set of 844 are depicted here sorted by their collection location and within collection site by their insecticide resistance genotypes. Each row of three boxes represents genotype data from a single A. peditaeniatus individual. Out of 108 A. peditaeniatus, 106 carried the G119S mutation (3 heterozygotes, 103 homozygotes), 78 carried the L1014F mutation (22 heterozygotes, 56 homozygotes), and 104 carried the A296S mutation (12 heterozygotes, 92 homozygotes). Forty nine percent (53/108) of individuals have homozygous resistant genotypes for all 3 loci and an additional 23% (25/108) carry at least one resistant allele for each of the 3 insecticide loci typed. Of the remaining 28% of individuals (30/108), 28 carried the homozygous susceptible allele at the L1014F kdr mutation (2 individuals with no available sequence), 17 were homozygous for the resistant allele at both ace-1 and rdl, another 7 individuals were homozygous for ace G119S and heterozygous for rdl A296S. A single individual was homozygous for the rdl resistant mutation and heterozygous for the ace-1 mutation. Five individuals carried a resistance mutation at only one of the 3 loci tested, 3 were ace-1 mutation homozygotes, 1 was heterozygous for rdl resistance and 1 was homozygous for rdl resistance with no sequence available for ace-1 resistance.
Geographic distribution of insecticide resistance mutations in A. peditaeniatus. Allele frequencies for ace-1 G1119S (top), kdr L1014F (middle) and rdl A196S (bottom) are shown as a function of sampling location. For kdr there is a detectable geographic pattern with the Forest site 2 exhibiting significantly (p = 0.03) higher susceptible allele frequency than all other sites combined, while Forest near Village 2 (p = 0.04) and 3 (p = 0.02) exhibit significantly higher resistant allele frequencies than all other sites combined. A similar pattern is observed for rdl, but here the total number of susceptible alleles does not provide the power to detect significant differences. For locations: F = Forest, FNV = Forest Near Village, P = Plantation, and V = Village. Only those 8 sites with at least 5 A. peditaeniatus (10 alleles) are shown here. Resistant (res) and susceptible (sus) alleles are colored according to the legends, with resistant alleles in the darker shade.
We also examined the frequency of the kdr L1014F, ace-1 G119S and rdl A296S insecticide resistance mutations in rarer mosquito species that were detected in the random subset of 844 samples molecularly typed for species. These 68 mosquito samples of rarer species represent 14 species and 2 unknown Anopheline samples (each species represented by 1–17 individuals; median = 2; Supplementary Table S5). In these rare species, rates of resistance alleles were high for rdl A296S with resistance being present in 8 of the 14 species and in the unknown Anopheline. In contrast no mosquitoes carried the L1014F mutation in kdr and only A. sinensis carried the ace-1 G119S mutation.
Landscape indices
A 2 km buffer zone around each mosquito collection site was used to determine the percentage of the different land cover classes present and showed high variability between sites (Fig. 8a). Patches of evergreen forest, forest and plantation forest were observed in all sites. Croplands were observed in all buffer zones, from 0.16% in Forest 3 up to 40.2% in Plantation 1, however the largest human-made land cover class was represented by plantation forests, from 7.5% in Forest 4–89.7% in Forest 1. While more A. dirus were observed in Forests 3 and 4 (n = 122 and n = 92 respectively), it was still present in the more degraded landscapes of Forests 1 and 2 (n = 71 and n = 18). The most conserved forest sites (Forests 3 and 4) were also associated with lower densities and lower diversities of mosquitoes as also observed in Plantations 2 and 3 (Fig. 8b).
(a) Percentage of each land cover category per collection site determined by a 2 km buffer zone and based on Landsat imagery (Landsat 4–8, 30 m resolution, 2017). Land cover classes defined as: cropland (land with herbaceaous and shrubby crops with harvesting and bare soil periods), evergreen broadleaf (land dominated by trees > 60% canopy cover with a tree above 5 m, evergreen broadleaf trees make up > 60% of the total tree cover), forest (land spanning more than 0.5 hectares with trees higher than 5 m and a canopy cover of more than 10%), mixed forests (land with > 60% tree canopy cover, tree height is greater than 5 m, and the forest composition is mixed such that no single forest type makes up > 60% of the total tree cover), plantation forests (land cultivated with perennial crops that reach heights above 5 m and occupy the land for long periods). The table is indicating patch density (PD), edge density (ED), Shannon’s diversity index (SHDI), and Simpson’s diversity index (SIDI) for each buffer zone. (b) Percentage of each mosquito species (molecular identification) sampled in each collection site. Sample sizes are given in parentheses.
Although Forests 3 and 4 had high proportions of natural forests, their high PD and ED values indicate more fragmented landscapes than Forests 1 and 2. Patch richness only varies between 4 and 5 types of land cover classes therefore an increase value of SHDI and SIDI corresponded to an increase evenness between patch types, with the highest patch evenness observed in Forests 2, 3 and 4 as well as Plantations 2 and 3. Villages and forest near the villages exhibited uneven patch distribution associated with fragmented landscapes.
Human behaviour indicators
More than half of the villagers reported going to bed between 19:00 and 20:00 pm and using a bed net the previous night (92%; Table S6). However, only 41% were insecticide-treated bed nets. More than three quarters of the villagers were working in cashew plantations and cassava fields at the time of the study (77%, Table S6). Interestingly, sleeping overnight in the field sites was reported in 3.7–7.81% of the cases, with a vast majority reporting the use of bed nets, whereas 77.9% of people working in the forest sites also slept there with only 25% of them using a bednet (Table S6).
Source: Ecology - nature.com
