Mapping the Genetic Relatedness of Outdoor-Biting Anopheles Mosquitoes in Zambia

1. Introduction

As African nations strive towards malaria control and elimination, insecticide re-sistance and residual malaria transmission challenge current vector-targeted malaria in-terventions [1,2,3,4]. While the primary vectors Anopheles funestus, An. coluzzii and An. gambiae are the focus of these mitigation strategies [5,6], selection pressure from indoor-focused measures such as indoor residual spraying (IRS) and insecticide treated nets (ITNs) have led to changes in mosquito behavior and populations [7,8,9]. One key change is the recog-nition that secondary Anopheles vectors such as An. rufipes [10,11], An. pharoensis [12], An. squamosus [13] and An. coustani [13,14,15], contribute to sustaining residual transmission in sub-Saharan Africa. These mosquito species are considered largely exophagic and exoph-ilic, behaviors that have allowed them to evade indoor vector control interventions [7,8,16,17].

Despite frequent collection alongside the primary malaria vectors, the genomics, ecology, biology, and behavior of these long-overlooked anopheline mosquito species are poorly understood [9]. Furthermore, the constraints of overlapping and cryptic morpho-logical features [18,19,20,21], together with the absence of available molecular data in genomic databases [17,20,22,23], have made robust morphological and molecular identification of these less studied anopheline species extremely challenging. Studies have also reported evidence of functional heterogeneity in anopheline genomes which influence their behav-ioral plasticity, a crucial characteristic for defining vectorial capacity and adaptability [24,25]. Therefore, the accurate identification and bionomic characterization of under-studied Anopheles species is now critical given their key role as local vectors in driving re-sidual malaria transmission in Zambia, Madagascar, southern Mozambique, Ethiopia and Kenya [13,15,26,27,28,29].

Although there is an extensive list of Anopheles sequences generated using molecular barcodes based on the cytochrome oxidase I (COI) and internal transcribed region 2 (ITS2) genes [12,17,18,20,26,27,28,30,31,32,33,34,35], there still remain limitations in available sequence for these understudied Anopheles species to produce robust differentiation between members of closely related taxa [17,18,31,36]. This includes cryptic species which may be incrimi-nated in residual malaria transmission but have been allotted placeholder names such as An. species 11 [17,18,20], An. species 15 [18] and An. species unknown group 1 [20], par-ticularly in the absence of comprehensive morphological identification to compliment the generated molecular barcode sequences. Furthermore, the use of the single COI gene to validate identification for less studied Anopheles has produced matches with low similari-ties (less than 80%) and weakly supported phylogenies. For instance, this has led to in-conclusive identities for members of the An. coustani group in earlier studies from Zambia [17,31] and Mozambique [26]. Recently, mitochondrial genomes were used to provide conclusive identities and differentiate the cryptic taxa of the An. coustani group into phy-logenetically well-supported taxonomic clades [37].

The acquisition of genomic datasets has become more accessible due to the expan-sion of sequencing and computational technologies, including mitochondrial genomes (mitogenomes) which have shown to be useful in the identification and resolution of phy-logenies for several mosquito species and species groups [38,39,40]. These circular, double stranded DNA molecules encode 37 genes, including 13 protein coding genes (PCGs), 22 transfer RNA (tRNA), 2 ribosomal RNA (rRNA) and a non-coding control region [40,41]. In addition to the 13 PCGs, key characteristics such as low incidence of recombination, high copy number, and maternal inheritance make the mitogenome a more effective tax-onomic tool compared to single barcodes [41,42,43,44]. Expanding mitochondrial genome re-sources to include less-studied mosquito species is essential for accurate species delinea-tion and gaining insights into mosquito ecology and systematics for public health inter-ventions. In this study, we aimed to generate high quality reference sequence data for rep-resentative understudied mosquito species and demonstrate the strength of mitogenomes compared to prior studies that were limited to the COI gene in attempts to resolve phy-logenies.

2. Materials and Methods

2.1. Mosquito collections

Specimens were collected during routine entomological surveillance in Nchelenge, Zambia in 2023-2024 as part of the Southern and Central Africa International Centers of Excellence for Malaria Research (ICEMR) investigations. Miniature CDC Light Traps (John W. Hock Co., Gainesville, Florida), were positioned both indoors and outdoors where people gather in the evening and near to animal pens.

2.2. DNA extraction and sequencing

DNA extractions using a modified extraction method [45] were performed on single mosquito specimens morphologically identified as An. rufipes, An. maculipalpis, An. preto-riensis, An. squamosus and An. pharoensis [21]. The extracted DNA and previously extracted specimens identified using the COI gene as species 11, species 15, unknown group 1, un-known group 2 and unknown group 3 from a previous study [20] were quantified using the Qubit dsDNA assay kit (Thermo Fisher Scientific, Waltham, MA) and shipped to Se-qCenter (Pittsburgh, USA) for library construction and Illumina sequencing. The libraries were 150 bp paired end sequenced to a count of 13.3 million reads per sample.

2.3. Mitochondrial genome assembly and annotation

The mitochondrial genomes were assembled using NOVOPlasty [46] (RRID:SCR_017335) version 4.3.5 with k-mer set at 39 base pairs and An. squamosus (OP_77691) as the seed sequence. Using the MITOchondrial genome annotation (MITOS) [47] galaxy tool, generated contigs were automatically annotated using the invertebrate genetic code under default settings. The start and stop codon positions of the annotated contigs were manually adjusted in Geneious Prime (RRID:SCR_010519) version 2025.1.2 (Biomatters, Auckland, New Zealand) using reference anopheline mitochondrial genomes as a guideline. The generated contigs with corresponding annotations were submitted to the GenBank database for the assignment of accession numbers

2.4. Phylogenetic analysis and tree construction

Using jModelTest (v2.1.10) [48], the best fit base pair substitution model based on the Akaike information criterion (AIC) and the Bayesian information criterion (BIC) was de-termined under default settings using an aligned sequence matrix. This alignment was generated using the MAFFT alignment tool implemented in the Geneious Prime (RRID:SCR_010519) version 2025.1.2 (Biomatters, Auckland, New Zealand from the 13 concatenated PCGs of mitogenomes generated in this study, available mitochondrial ge-nomes of understudied African anopheline species; An. marshallii (NC_064607), An. mou-cheti (NC_064608), An. gibbinsi (OR_569715), An. nili (NC_064610), An. coustani group (PQ_585798, PQ_587039, PQ_587041, PQ_587036, PP_375116), and reference mitoge-nome sequences for the well-studied species An. gambiae (MG_930894) and the An. funestus group (MG_742172, MG_742194, MT_917162, MT_917137, MT_91714, MT_917157, MT_917163). Using Bayesian Evolutionary Analysis by Sampling Trees (BEAST) 2 soft-ware [49], inference analysis was performed on the aligned sequence matrix using tree independent runs and a 20% burn-in rate for tree building purposes under default set-tings. Bayesian analysis was also performed on an alignment generated from COI se-quences available from GenBank complementary to representative species in the mito-chondrial genome tree. Trees were visualized and annotated using FigTree v.1.4.4 (http://tree.bio.ed.ac.uk/software/figtree/).

2.5. Dating time estimation

Divergence time estimations were calculated using the previously mentioned pa-rameters used for Bayesian inference in BEAST2. The frequently referenced Ae-des-Anopheles divergence time of approximately 154.7 million years ago (MYA) [50] was used as the calibration point set to normal distribution.

3. Results

The 10 mitochondrial genomes produced in this study were consistent with other Af-rican anopheline mitogenomes represented in NCBI’s GenBank database as comprised of 13 PCGs, 22 tRNAs and 2 rRNAs, with lengths ranging from 15,534 bp (Unknown group 2) to 15,346 bp (An. pharoensis) and a mean AT content of 77.1 % (Table 1).

Bayesian inference for the mitochondrial genomes resulted in a phylogenetic tree which separated specimen sequences (Figure 1) compared to COI tree which resulted in 4 weakly supported main clades (Figure S1). The most recent common ancestor (MRCA) for An. funestus and An. gambiae with the outdoor biting Anopheles included in this study was 66.45 and 69.65 MYA respectively (Figure 2).

4. Discussion

Bayesian phylogenetic analysis based on the concatenated PCGs from the mitoge-nomes was able to clearly delineate less studied African anopheline taxa compared to phylogenetic analyses using exclusively the COI gene as in earlier studies [17,20]. Several of these species are often misidentified during morphological examination or matched to an unassigned species following molecular barcoding [20,26,33]. This may be attributed at least in part to the unknown diversity of outdoor biting anophelines and relative inexpe-rience in morphological identification of these numerous taxa, some of which may yet to be fully described due to the bias placed on understanding highly anthropophilic, endophagic and endophilic species [6,9,17,28,31,31]. The uncertainty of vector richness of the exophagic anophelines together with their perceived adaptability and undetermined foraging times, habitats, and opportunistic feeding patterns have led to ambiguous spe-cies assignments [9,14,28,33,51]. This combination of characteristics emphasizes the need for the expansion of outdoor surveillance and investigation of the anophelines found in this space, as well as innovative strategies to overcome the shortcomings of in-door-focused interventions.

A growing number of understudied anophelines have been confirmed or implicated in human malaria transmission. Despite historic consideration as a secondary vector [12,17,27,52], An. coustani is now regarded as a primary local vector in regions of Mada-gascar [15,29,53]. Additionally, previous studies have identified An. pretoriensis [34], An. pharoensis [52,54], and cryptic anopheline species [17,18,20] demonstrating opportunistic feeding patterns on humans, some infected with sporozoites of human malaria parasite species. Anopheles squamosus is another species strongly implicated as a malaria vector, with a wide geographical range, and has demonstrated variable foraging behavior to-wards human blood meals [13,55,56]. Related to this are a number of ‘molecular taxa’. Examples include An. species 11 and An. species 15 which are often morphologically keyed as An. squamosus but are differentiated by the COI barcode and even more strongly by the mitogenome sequence (Fig. 1)[17,20]. Others include An. unknown groups 1 – 3 for which morphology and molecular barcoding was inconclusive [20]. Here the mitogenome data provide the most comprehensive insight into the identification of these ‘molecular taxa’, but as with prior studies, without a more extensive sequence database of recognized species, these specimens remain taxonomically unresolved. The fact that many of these exophagic taxa cluster together in the phylogenetic analysis may be an artifact of their shared ancestry and that these share behavioral adaptations may have been reinforced over millennia.

It is clear that full mitogenomes offer much more discriminatory power for a phylo-genetic approach to inquire about shared biological traits and possibly ascertain whether behaviors such as biting preference are due to recent adaptations or reflect the existence of genetically distinct lineages which may have been overlooked when restricted to morpho-logical identification. Dating time estimations from well-recognized malaria vectors fur-ther corroborate the presence of these outdoor biting Anopheles as cryptic lineages with distinct ecological niches, suggestive of understudied species that may maintain trans-mission outdoors, perhaps under certain conditions such as relative absence of non-human hosts, or human behavior that promotes high opportunistic human biting rates. Furthermore, the accurate taxonomic placement of these mosquitoes highlights the relationships between known vectors and putative vector species which may provide fur-ther insights into understanding the differences in biting, foraging, and vectorial capacity of these less studied species. Linking morphological reference specimens to genomic data is key for the accurate identification given the status of unassigned Anopheles species with sporozoites collected in the field [18,19,26,52].

5. Conclusions

Although reference sequences are available for many commonly encountered outdoor biting anopheline species, there remains a paucity of data to accurately identify and taxo-nomically place these species in the wider Anopheles genus. This study contributes valu-able genetic datasets representing exophagic species collected in Zambia and present across the African continent. The generation of mitochondrial genomes for cryptic unas-signed species that are commonly collected has given priority to the use integrative tax-onomy in future research. The linking of molecular data with morphological and type specimens can further strengthen the credibility of species delimitation for the assigned zoological nomenclature of these cryptic taxa. The analyses from this study identified the phylogenetic relationships between the primary malaria vectors and understudied species implicated in malaria transmission, assisting to close the genetic gap of what we know about these anophelines of public health importance.