Population genetic structure analysis of Anopheles kleini in the Republic of Korea based on the mitochondrial COI gene

Abstract

Background

Anopheles kleini is a competent vector mainly observed in the northern, malaria-risk areas of the Republic of Korea (ROK). In this study, the population genetic structure of An. kleini was analysed for the first time in the ROK using the mitochondrial cytochrome c oxidase subunit I (COI) marker.

Methods

The genetic structure of 249 An. kleini was analysed from three statistically analysable regions, each including more than five mosquitoes.

Results

Network analysis identified 140 haplotypes organized into three clusters. Cluster II was related to An. kleini from eastern Russia and northwestern China. The pairwise genetic distance (F_ST) values among the populations showed regional genetic differences between Gangwon-do and Gyeonggi-do. Analysis of molecular variance (AMOVA) indicated that individual mosquitoes within the population had a significant influence on the total variation. The neutrality test, using three methods (Fu’s Fs, Fu, and Li’s D, and Fu and Li’s F), indicated that all values were negative, suggesting that An. kleini is an expanding population. Anopheles kleini in Yanggu has a significant difference in genetic distance from other regions.

Conclusion

This study provides molecular epidemiologically information for understanding the spatial population structure of An. kleini and is helpful for malaria control in the ROK.

Background

Malaria affects global public health, along with diseases such as acquired immune deficiency syndrome (AIDS), tuberculosis, and amoebiasis. Malaria is caused by an infection by Plasmodium parasites, transmitted to humans via a bite from an infected female Anopheles mosquito and requires interactions between the host, vector, and parasite [1]. Five species of malaria parasites infect humans: Plasmodium falciparum, Plasmodium vivax, Plasmodium malariae, Plasmodium ovale, and Plasmodium knowlesi [2]. Plasmodium falciparum and P. vivax occur in Southeast Asia, while P. knowlesi outbreaks occur in parts of Malaysia and Indonesia [3]. In the Republic of Korea (ROK), P. vivax occurs principally in some northern regions (Gangwon-do, Gyeonggi-do, and Incheon) facing the border of the Democratic People’s Republic of Korea (DPRK) [4]

Approximately 70 species of the genus Anopheles are capable of transmitting malaria worldwide, of which 41 species are known to be the dominant vector species [5]. Eight species of Anopheles mosquitoes (Anopheles belenrae, Anopheles kleini, Anopheles koreicus, Anopheles lesteri, Anopheles lindesayi japonicus, Anopheles pullus, Anopheles sinensis, and Anopheles sineroides) inhabit the ROK, and six species inhabit the DPRK (Anopheles anthropophagus, An. belenrae, An. kleini, An. lesteri, An. sinensis, and Anopheles yatsushiroensis) [6, 7]. In particular, An. kleini, which is found in the both countries, mainly occurs in malaria-risk areas near the demilitarized zone (DMZ) in the ROK [6, 8]. Artificial infection experiments confirmed sporozoite infection in An. kleini [9, 10]. In field-collected Anopheles mosquitoes, An. kleini exhibited a high rate of P. vivax sporozoite infection; therefore, caution is required as a competent vector [11].

Mitochondrial DNA (mtDNA) markers are used to study the population structure and genetic diversity of various species, including mosquitoes. In Thailand, Anopheles minimus was genetically divided into two lineages using cytochrome c oxidase subunit I (COI) and cytochrome c oxidase subunit II (COII) [12]. In Kenya, the genetic structure and diversity of Anopheles funestus were evaluated using a COII marker [13]. Besides the Anopheles mosquitoes, the introduction of Aedes albopictus from abroad was confirmed using a COI marker in Japan [14]. Furthermore, there is a study about the genetic structures of Cyprinidae (fish), Amblyomma tholloni (tick), and Glossina morsitans morsitans (tsetse fly) using COI markers in India, Kenya, Zambia, and Malawi respectively [15,16,17].

Despite various studies using mtDNA markers, no study has been conducted on the genetic structure of An. kleini populations. In a study confirming the genetic structure of An. sinensis population in the ROK, An. kleini was also collected. However, it was excluded from the analysis because it was collected from only one area [18]. This study provides a molecular biological understanding by analysing the population genotypes of An. kleini, which is presumed to be a competent vector of malaria in the ROK.

Methods

Mosquito sampling and species identification

Mosquitoes were collected from approximately 17 civilian houses in malaria-risk areas and seven livestock sheds in non-risk areas from May to October 2022. A total of 266 An. kleini were collected only in malaria-risk areas, and 249 An. kleini were used in population genetic analysis (Table 1). The mosquitoes were collected using a black light trap (Shin-Young Comm. System, Namyangju, ROK) in malaria-risk areas and a light-emitting diode (LED) trap (Biotrap, Gunpo, ROK) in non-risk areas. Mosquito specimens were initially identified using a microscope and morphological characteristics. Species identification was confirmed using multiplex polymerase chain reaction (PCR) assays targeting the internal transcribed spacer 2 region [19, 20]. Total DNA was extracted using DNAzol (Thermo Fisher Scientific, Massachusetts, U.S.A) or the automatic nucleic acid extraction equipment QIAamp 96 DNA QIAcube HT Kit (QIAGEN, Hilden, Germany) following the manufacturer’s protocol.

Table 1 Geographical information on the sampling sites and sizes (n)

COI gene amplification

Total DNA extracted from individual mosquitoes was used as a template to synthesize the 710 bp of COI gene using forward (LCO1490; 5′-GGTCAACAAATCATAAAGATATTGG-3′) and reverse (HCO2198; 5′-TAAACTTCAGGGTGACCAAAAAATCA-3′) primers [21]. DNA was amplified using AccuPower PCR PreMix (Bioneer, Daejeon, Korea). The PCR conditions were as follows: initial denaturation at 95 ℃ for 5 min, 40 cycles of denaturation at 95 ℃ for 30 s, annealing at 53 ℃ for 30 s, elongation at 72 ℃ for 30 s, and final extension at 72 ℃ for 10 min. The PCR products were confirmed using a QIAxcel capillary electrophoresis system (Qiagen, Hilden, Germany), and Sanger sequencing was performed to analyse the DNA sequences.

Data analysis

The amplified sequence of the COI gene was trimmed using DNASTAR Lasergene SeqMan Pro and aligned using MEGA 11 software. The sequence characteristics, such as the number of haplotypes (H), number of segregating sites (S), average number of nucleotide difference (K), average number of mutations per sequence (θ), haplotype diversity (Hd), and nucleotide diversity (Pi), were calculated using DnaSP (version 5.10.01) [22]. Median-joining network of all An. kleini haplotypes were drawn using Network 10.2 software [23]. The pairwise genetic distance (F_ST) between populations was confirmed using the DnaSP program, and gene flow (Nm) values were estimated from the pairwise F_ST [Nm = (1—F_ST)/4F_ST] [24]. Analysis of molecular variance (AMOVA) was performed using Arlequin 3.5.2.2 [25]. The neutrality test was evaluated with three methods (Fu’s Fs, Fu and Li’s D, and Fu and Li’s F) using the DnaSP program [22].

Results

Population sampling

Mosquitoes were collected from 24 locations in the ROK; however, 266 An. kleini were identified at only 10 sites, which are malaria-risk areas (Table 1). All An. kleini samples were confirmed using multiplex PCR assays. Only 249 samples from three regions, Baekyeon-ri (BY), Josan-ri (JS), and Yanggu (YG), were used for population genetic analysis since the sample sizes from these regions had more than five individuals (Fig. 1).

Fig. 1

Anopheles kleini sampling sites in the Republic of Korea (ROK) in 2022. The site names are abbreviated, and the full names and locations are listed in Table 1. The sites where An. kleini were collected and marked with blue circles

Sequence characteristics

The mtDNA COI gene was successfully amplified from individual An. kleini. Compared with the sequences from a previous study [26,27,28], 594 bp out of a total of 710 bp were analysed. There were 499 conserved and 95 variable sites.

The number of segregating sites (S), JS was the highest at 85, followed by YG at 48 and BY at 34. The average number of nucleotide differences (K) was 10.41218 ± 1.39397 and the average number of mutations per sequence (θ) was 13.39431 ± 3.05977. Haplotype diversity (Hd) was high in all three analysis regions, from 0.92105 (YG) to 0.97895 (BY), whereas nucleotide diversity (Pi) was low, from 0.01543 (JS) to 0.02006 (YG) (Table 2).

Table 2 Sequence characteristics of Anopheles kleini in the three malaria-risk areas

Haplotype network analysis

A total of 140 haplotypes were identified in 249 An. kleini individuals. The haplotype percentages for each collection area were BY 85.00% (17/20), YG 75.00% (15/20), and JS 56.94% (119/209) (Table 2). There were three dominant haplotypes with more than ten individuals: H_12 (n = 28), H_2 (n = 25), and H_1 (n = 15). Haplotypes with more than two individuals were confirmed as 25 (17.86%), and the rest were found from only one individual. Region-specific haplotypes from one individual were detected in all three analysis regions but were especially abundant in JS.

Three clusters were identified in the COI haplotype network and some haplotypes (median vectors) were not detected (Fig. 2). Cluster I contained 75 haplotypes, including the dominant haplotypes (H_12, H_2, and H_1). Furthermore, 59 haplotypes belonged to Cluster II and six haplotypes belonged to Cluster III. When all haplotypes were compared with the existing foreign reference sequence, H_101 in Cluster II was identical to the Chinese sequence (OP311323), and H_105 was identical to the Russian sequence (KC855655).

Fig. 2

The haplotype network of the cytochrome c oxidase subunit I (COI) gene was constructed using Network 10.2. Each circle represents a haplotype and its size is proportional to the number of individuals included in the haplotype. Pink, green, and blue represent individuals from Baekyeon-ri (BY), Josan-ri (JS), and Yanggu (YG), respectively. White dots (median vectors) represent undetected hypothetical haplotypes. Black circles indicate haplotypes with the same sequence as the reference (H_101: OP311323, Heilongjiang, China; H_105: KC855655, Khabarovsk, Russia)

Genetic structure among population

Pairwise genetic distance (F_ST) was used to evaluate the genetic distance between populations in the analysis region. The F_ST values between YG and BY were 0.08899 and 0.74900 for JS.

In other words, two groups, Gangwon-do (YG) and Paju, Gyeonggi-do (BY and JS), had high genetic divergence and were regionally differentiated [29]. However, the genetic distinction of populations was difficult to analyse because the F_ST value between BY and JS was negative, and gene flow (Nm) was not available (Table 3).

Table 3 Pairwise genetic distance (F_ST) and gene flow (Nm) values among populations

AMOVA was analysed by dividing into two groups as Gangwon-do (YG) and Gyeonggi-do Paju (BY, JS). As a result, “within populations” was noticeably higher than “among groups” and “among populations within groups” (Table 4), which implies that individuals within a population have a substantial influence on the total variation than “among groups” and “among populations within groups.”

Table 4 Analysis of molecular variance (AMOVA) for geographic variation of Anopheles kleini in three malaria-risk areas

Neutrality test

The neutrality test was calculated using DnaSP with three methods: Fu’s Fs, Fu and Li’s D, and Fu and Li’s F. All data showed negative values, suggesting that the population was expanding and had many low-frequency mutations (Table 5). In particular, strong negative values (P < 0.05) of Fu and Li’s D and Fu and Li’s F were identified in the JS population, and the neutrality deviation was more sensitive than that in other regions.

Table 5 Neutrality test of Anopheles kleini in the three malaria-risk areas

Discussion

Population genetics research is useful for understanding and predicting the epidemiology of vector-borne diseases and provides information on the spatial limits of natural populations and the characteristics of gene flow between populations [30]. The COI is a highly representative marker used to study the genetic diversity of insects, including Anopheles mosquitoes [31].

Anopheles kleini belongs to the Hyrcanus group and mainly inhabits continental monsoon climates, such as far-east Russia and northwest China [26, 28]. In the ROK, An. kleini is rarely identified in the southern region, but most specimens have been collected near the demilitarized zone (DMZ), especially in the northern region of Gyeonggi-do [8, 32,33,34]. In this study, a population genetic analysis of An. kleini in the ROK is conducted using a mitochondrial COI marker. All An. kleini samples were collected in malaria-risk areas near the DMZ, and regions with less than five individuals were excluded from the analysis to facilitate statistical analysis (Table 1). The biased collection of An. kleini populations may be because of their ecological preference for cool temperatures [35]. Most specimens were collected in Josan-ri (JS), the area closest to the Democratic People’s Republic of Korea (DPRK).

Haplotype diversity (Hd) and nucleotide diversity (Pi) indicate genetic diversity among populations. When analysing Baekyeon-ri (BY), JS, and Yanggu (YG), we demonstrated the characteristics of a migratory group with high Hd and low Pi. This indicated that An. kleini populations recently expanded to a small effective population after experiencing a bottleneck [31]. In addition, high levels of Hd may occur because of the various environments and lifestyles to which An. kleini are exposed during the process of adapting to rapid natural development [36].

Three clusters were identified in the COI haplotype network of the An. kleini population in the ROK of three analysis regions (Fig. 2). Cluster I contained three dominant haplotypes (H_12, H_2, and H_1), and when compared with the existing reference, H_12 was confirmed to have the same sequence as OP150362 in the National Center for Biotechnology Information (NCBI) database. OP150362 is an An. kleini collected in 2021 from the Neutral Nations Supervisory Commission camp (< 10 m from the DMZ) [27]. H_12 had a large population size in the network but was not a central haplotype. Therefore, H_12 considered a migratory group. Cluster II consisted of small haplotypes with no more than three individuals. H_101 and H_105 in Cluster II had sequences identical to OP311323 from China and KC855655 from Russia, respectively [26, 28]. In particular, H_101 and H_105 were linked to various haplotypes in Cluster II; therefore, many individuals were presumed to have immigrated. Because of migratory insect pest populations, such as Nilaparvata lugens, Sogatella furcifera, and Laodelphax striatellus which entered the ROK from China via the jet stream [37], the introduction of mosquitoes is also plausible. Cluster III had only six haplotypes, most of which were identified only in YG (H_78, H_79, H_80, and H_124). This supports the genetic structure of the An. kleini individuals differed regionally, as the F_ST value of YG in Gangwon-do was distinct by more than 0.05 from both BY and JS in Gyeonggi-do (Table 3). These regional differences may have been influenced by environmental factors, such as ecological isolation since the YG is surrounded by mountain ranges.

In total, 140 COI haplotypes were identified in 249 An. kleini and most haplotypes had only one individual. Based on the results of the AMOVA, individuals within the populations had a significant impact on the overall variation of the population (Table 4). However, if there are too many low-frequency mutations, the genetic characteristics of population expansion may persist for a long time, obscuring the genetic structure of the true ecological populations that exist in nature [38].

Additionally, in the neutrality test using the three methods, all values were negative (Table 5), confirming that An. kleini populations are expanding [39]. The high number of haplotypes and population expansion suggested the excellent adaptability of An. kleini, which may improve the vector’s ability to respond to ecological environments and maintain malaria transmission [13].

Plasmodium vivax was detected in An. kleini by targeting the small subunit ribosomal RNA gene [40]; however, all samples were negative in this study. According to the 2023 malaria vector mosquito surveillance by the Korea Disease Control and Prevention Agency, P. vivax was detected in Anopheles mosquitoes collected from JS [41]. Based on internal data, P. vivax was consistently confirmed in Anopheles kleini was collected using a high-altitude insect net installed near the DMZ in 2024. Most belonged to Cluster II, which was related to northern populations in the haplotype network. Since approximately 2,136 patients have identified in the DPRK in 2022 [3], P. vivax-positive An. kleini likely migrated from the north.

To our knowledge, this is the first study to analyse An. kleini population genetics at malaria-risk areas in the ROK. Although this study is limited by the fact that the population is biased toward one region (JS) and the population size is different, the data distinguishes between populations of An. kleini using regional characteristic analysis. Following the continuous surveillance of Anopheles mosquitoes as vectors, the accumulation of data related to the genetic structure can be used to understand the characteristics of An. kleini and malaria eradication in the ROK.

Conclusions

This study analysed the population genetic structure of An. kleini for the first time in the ROK based on the COI gene. Anopheles kleini population was divided into three clusters and showed characteristics of a migratory population. The An. kleini population is expanding and is estimated to be closely related to the northern reference sequences. Such information provides insight into the malaria vector and can be provided as basic information for malaria control.