Population genetics of the malaria vector Anopheles aconitus in China and Southeast Asia

Anopheles aconitus is a well-known vector of malaria and is broadly distributed in the Oriental Region, yet there is no information on its population genetic characteristics. In this study, the genetic differentiation among populations was examined using 140 mtDNA COII sequences from 21 sites throughout southern China, Myanmar, Vietnam, Thailand, Laos and Sri Lanka. The population in Sri Lanka has characteristic rDNA D3 and ITS2, mtDNA COII and ND5 haplotypes, and may be considered a distinct subspecies. Clear genetic structure was observed with highly significant genetic variation present among population groups in Southeast Asia. The greatest genetic diversity exists in Yunnan and Myanmar population groups. All population groups are significantly different from one another in pairwise Fst values, except northern Thailand with central Thailand. Mismatch distributions and extremely significant F s values suggest that the populations passed through a recent demographic expansion. These patterns are discussed in relation to the likely biogeographic history of the region and compared to other Anopheles species.

in China, and south from southern China through Southeast Asia into Indonesia (Harrison, 1980;Chen et al., 2002). It has been incriminated as a vector of malarial parasites in Thailand (Gould et al., 1965) and is considered to be the primary vector in Malaysia and Indonesia (Damar et al., 1981;Rahman et al. 1993; Barcus et al. 2002). Besides malarial parasites, filariae of Wuchereria bancrofti have also been found in this species in Flores, Indonesia (Atmosoedjono and Dennis, 1977). Adults of An. aconitus, An. varuna, An. minimus and An. harrisoni are morphologically very similar to one another, with overlapping ranges of character variation, often resulting in their misidentification (Harrison, 1980;Harbach et al., 2007). Various molecular methods have been developed, as summarized by Van Bortel et al. (2000), to facilitate identification of these and other species of the Myzomyia Series of Anopheles subgenus Cellia, including single-strand conformation polymorphisms (Sharpe et al., 1999), restriction patterns of the ITS2 fragment of rDNA (Van Bortel et al., 2000) and random amplified polymorphic DNA markers (Kengne et al., 2001). Subsequently, Phuc et al. (2003) and Carros et al. (2004) each developed a multiplex PCR assay separately to identify the malaria vector An. minimus and four related species from Southeast Asia, and the major malaria vectors within the African An. funestus and the Oriental An. minimus groups, both using ITS2 locus. Morgan et al. (2011) investigated the impact of pleistocene environmental change in shaping genetic diversity across the Indo-Burma region, by combining COII data of An. aconitus and other 8 Anopheles species. However, they did not deliberately analyze the COII variation and population genetic characteristics of An. aconitus, nor did they sequence other DNA markers to investigate the diversity within the species. Phylogenetic analysis based on nuclear and mitochondrial DNA sequences shows that An. aconitus, An. filipinae, An. mangyanus, An. pampanai and An. varuna form a monophyletic group distinct from other members of the Funestus Group (previously Minimus Group, see Harbach, 2004) of the Myzomyia Series, and the name Aconitus Subgroup was applied to these five closely related species (Chen et al., 2003, Garros et al., 2005. Currently, chemical measures remain central to the control of malarial parasites and their vectors despite the problems of parasite resistance to antimalarial drugs (Meshnick, 1997) and the spread of insecticide-resistance in vector populations (Collins and Paskewitz, 1995;Ranson et al. 2011). The insecticide resistance in Anopheles species has been widely spread in China and Mekong Region (Cui et al., 2006;Van Bortel, et al. 2008). The Mekong Region is still one of the most threatening foci of malaria (Cui et al. 2011), and the absence of knockdown resistance suggests metabolic resistance in the main malaria vectors of the region (Verhaeghen et al. 2009). Knowledge of vector population genetics may allow a more rational management of insecticides to defer the development of resistance and could help explain the role of vectors in relation to geographical variation in malaria transmission. However, the population genetic characteristics of An. aconitus are virtually unknown despite its very wide distribution and epidemiological importance. More importantly, different subspecies might be substantially different in their biology, ecology, behaviour and malaria vector status, and hence malaria control. Despite the wide distribution of this species, whether there is any subdivision at the subspecies level has never been investigated.
Because of its rapid evolution and maternal, non-recombining mode of inheritance, mitochondrial DNA (mtDNA) has been widely used to investigate population differentiation and the evolutionary history of populations and subspecies (Avise et al., 1988;Avise, 1989;Cronin, 1992;Arctander et al., 1996). In particular, the mitochondrial cytochrome oxidase II gene (COII) has been successfully used to study population genetic structure and population history for a wide range of insect species (Krafsur et al., 1999;Sperling et al., 1999;Clements et al., 2000;Franck et al., 2001;Rua et al., 2001;Chen et al., 2004;Chen et al., 2011). In the present study, we sampled populations of An. aconitus in China, Myanmar, Vietnam, Thailand, Laos and Sri Lanka, examined the pattern of variation within and among populations using a large sample of COII sequences, and checked whether there is any subdivision of the morphologically defined species using two mtDNA genes (COII and ND5) and three rDNA sequences (D2, D3 and ITS2 regions) from a smaller set of samples.

Specimens
Mosquitoes were widely collected throughout southern China, Myanmar, Vietnam, Thailand, Laos andSri Lanka in 1996-2003. Most localities visited were areas with relatively high malaria incidence. Adults resting near or landing on humans and tethered bovines were captured with aspirators between 2000 and 2300 h. Mosquitoes were killed by exposing them to ethyl acetate vapour, placed individually in Beem capsules (Agar Scientific Ltd, Stansted, UK) and stored in plastic bags with silica gel. Larvae and pupae were collected in rice fields and irrigation ditches at some sites and reared to adults. The keys of Harrison (1980) were used to identify adults and larvae to An. aconitus and related species.

Molecular identification and DNA sequencing
DNA was extracted from mosquitoes using a phenol-chloroform method (Chen et al., 2002). The single-strand conformation polymorphism (SSCP) method based on the 28S D3 region, developed for the identification of An. minimus, An. harrisoni, An. aconitus, An. varuna, An. pampanai, An. fluviatilis, An. jeyporiensis and An. culicifacies (Sharpe et al., 1999;Chen et al., 2002), was used to exclude samples of species indistinguishable from An. aconitus by morphology. The mtDNA COII gene of those individuals that were identified as An. aconitus by SSCP was sequenced for the study of population structure. A smaller set of these individuals, selected on the basis that they have different COII haplotypes and come from different geographical regions, was also sequenced for the rDNA D2, D3 and ITS2 regions and the mtDNA ND5 gene to confirm the SSCP identification using the D3 locus sequence, and to check the monophyly and unity of the species. Anopheles filipinae and An. pampanai were used as outgroup taxa to analyse the variation of these loci within An. aconitus.
The primers used for PCR amplification and sequencing were taken from Sharpe et al. (2000) for COII, D3 and ITS2 and Krzywinski et al. (2001) for ND5 and D2. Amplifications were performed in 50 μl volume overlaid with two drops of mineral oil on a HYBAID OmniGene cycler (Thermo Hybaid, Ashford, England). Each PCR included 1/100 DNA of a whole mosquito, 5 μl 10× ReddyMix ™ buffer (ABgene, Epsom, England), 200 μM dNTPs, 2 mM MgCl 2 , 600 nM of each primer and 1.3 units of Thermoprime Plus DNA Polymerase (ABgene). Reactions started with denaturation at 95ºC for 5 min, followed by 35 cycles, each cycle consisting of denaturation for 40 s at 95ºC, annealing for 40 s at 55ºC and extension for 1 min at 72ºC, with a final extension at 72ºC for 6 min. PCR products were electrophoresed through ethidium bromide-stained 1% agarose gels in 1× TBE and visualised under UV light to check for successful amplification. The PCR products were then purified using a spin column (Promega Wizard PCR Preps, Madison, Wisconsin, U.S.A.) and sequenced in both directions in an ABI 377 automated sequencer (PE Applied Biosystems, Warrington, England).

Analysis of mtDNA variation
The DNA sequences were edited manually and aligned using Clustal X (Thompson et al., 1997). The COII and ND5 sequences were translated into amino acids using the Drosophila mtDNA genetic code in Transeq (European Bioinformatics Institute: http://www.ebi.ac.uk/). Genetic diversity estimates of the COII locus, θ π , based on the average pairwise number of differences between sequences (Tajima, 1983) and θ s , based on the number of segregating nucleotide sites per sequence (Watterson, 1975), were also estimated and then converted to values per base, at the within-population level and for geographically defined groups of populations ( Fig. 1) using Arlequin, ver 3.5 (Excoffier et al., 2010). Genealogical relationships among haplotypes were constructed using TCS (Clement et al., 2010) with the method described by Templeton et al. (1992). Analysis of population genetic structure was carried out using analysis of molecular variance (AMOVA; Excoffier et al., 2010) also in Arlequin, ver 3.5. AMOVA estimates the proportion of total genetic variation attributable to different hierarchical levels based on the geographical distribution of haplotypes, taking into account the number of molecular differences through squared distances between haplotypes. Correlations among haplotypes are measured in terms of Φ-statistics in a manner similar to F-statistics. The significance of the Φ-statistics was tested by permutation on the original inter-individual squared distance matrix.
Pairwise Fst values used as genetic distances between recently diverged population groups were estimated with the methods of Slatkin (1995) and tested for significance by permutation. Genetic isolation by geographic distance was assessed by a Mantel test. The frequency distributions of the numbers of segregating sites in all possible pairwise comparisons, known as mismatch distributions, were calculated using Arlequin. We also calculated Tajima's D (Tajima, 1989) and performed Fu's Fs test (Fu, 1997) using Arlequin for test of neutrality. LAMARC program (Kuhner 2006) was used to compute effective population size and population growth rate by using the likelihood method. The test involves comparisons of different estimates of the parameter θ, which is defined as 2Nμ for a haploid locus, where N is the effective population size and μ is the mutation rate per sequence per generation (Fu 1997).
Maximum likelihood (ML) phylogenetic analyses were conducted using the version of PAUP* (version 4.0b8) (Swofford, 2001). Each was performed using the heuristic search option employing step-wise addition with 100 random taxon addition sequence replicates and five trees held at each step. The model of DNA substitution for ML was determined using the program Modeltest (version 3.06) (Posada & Crandall, 1998). Modeltest uses hierarchical likelihood ratio tests to determine which of the General Time Reversible (GTR) family of substitution models (sixty-four in all) best fits the data. Node support for MP analysis was assessed using 1,000 bootstrap pseudo-replicates.

COII sequence variation
Anopheles aconitus is difficult to collect because it resides in hilly or mountainous areas and earlier distribution reports are an unreliable guide due to the confusion in morphological identification. Thirty-four sites were visited at least once but only 79 individuals of An. aconitus from 17 sites were successfully collected, and these were sequenced for the entire COII gene with a length of 685 bp. Further 61 COII sequences with a length of 639 bp for An. aconitus were taken from NCBI, reported earlier in Morgan et al. (2011). Totally, 140 sequences with a length of 639 bp from 21 sites were used for the population genetic analysis (Table 1; Fig. 1).
The set of sequences had 46 haplotypes coding for 212 amino acids, and were very AT-rich (75.46% A and T bases). Forty-eight nucleotide substitutions were identified at 47 of the 639 sites, of which 44 were transitions (91.67%) and 4 were transversions (8.33%) ( Table  2). Seven substitutions resulted in amino acid changes (direction of change based on outgroup comparison; Chen et al., 2003): those at site 16 (A→T) for haplotype 24, at site 44 (S→F) for haplotype 44, at site 395 (N→S) for haplotypes 9 and 27, at site 430 (V→I) for haplotypes 6, 11, 15, 16 and 37, at site 463 (V→I) for haplotype 33, at site 496 (V→I) for haplotype 35, and at site 576 (L→F) for haplotypes 46. Twelve amino acid sequences resulted from these amino acid changes do not correspond completely to geographical groups. Sequences have been deposited in GenBank/EMBL (accession numbers in Table 2; AJ512745-AJ512746 for An. filipinae and An. pampanai).
Due to the relatively small sample size for each population, we pooled samples from the 21 collecting sites into 9 population groups based on geographical distribution, but the group N Vietnam with only 4 samples and relatively far geographical separation was not analyzed as a group. Table 1 gives the haplotype distribution in these population groups. Haplotype diversity was relatively low (compared with other mosquito species, see Discussion), with 46 haplotypes found among 140 individuals. This is in accord with the generally low levels of nucleotide diversity estimated by θ π and θ s . The population groups from Yunnan (θ π = 5.80±3.58×10 −3 and θ s = 6.96±3.24×10 −3 per base; Table 1) and from Myanmar (θ π = 4.13±2.54×10 −3 and θ s = 8.70±3.30×10 −3 ) are most diverse followed by the northern and central Thailand groups. The overall estimates of molecular diversity are θ π = 4.82±2.79×10 −3 and θ s = 12.89±3.48×10 −3 , excluding the distinct Sri Lanka sample (see below).
Modeltest selected the general time reversible model GTR+I model as the best fit for ML phylogenetic analysis for each of these 5 data. Maximum-likelihood trees of sequences from each locus were constructed under the GTR+I model using An. filipinae and An. pampanai as outgroup taxa and the Ln likelihood are all quite large (Fig 2). Samples from Sri Lanka form a distinct clade within the species based on COII and ND5 sequence, and have ITS2 and D3 sequences that are not shared with Southeast Asian specimens. Other haplotypes do not show clear geographical differentiation. Figure 3 shows the genealogical network constructed from 46 haplotypes of the COII gene. The number of mutation differences (unit branches) on the network represents the maximum number of mutational connections between pairs of haplotypes justified by the 95% 'parsimony' criterion (Posada & Crandall, 2001). These haplotypes are linked via some missing intermediates that could represent sampling gaps or extinct intermediates. However, because the number of samples greatly exceeds the number of haplotypes, reflecting the low level of genetic variation, these haplotypes should adequately represent COII sequence variation of the species throughout the sampled distribution. Despite the low genetic diversity, a high level of homoplasy is apparent in the network, with some mutational events resulting in closed loops, as observed in other mosquito mtDNA studies (Walton et al., 2000). The network is centred on haplotypes 1, which has the highest outgroup weight (0.09), implying that it may be ancestral. Haplotype 1 is shared by 44 individuals and thus has the highest frequency in the sample, followed by haplotype 2 (20 individuals) and haplotype 3 (10 individuals). Importantly, all three haplotypes from Sri Lanka form a separate branch, well separated from other parts (but haplotype 37) of the network. This is consistent with the distinct D3, ITS2 and ND5 sequences from the island. In addition, the haplotypes of Hainan and southern Thailand were also on separate branches, but these two branches are closely linked to other haplotypes and the sequences of D3, ITS2 and ND5 from these two places do not represent separate lineages.

Genetic structure and genetic distances
In the hierarchical AMOVA, pairwise comparisons were made between 7 population groups (excluding Sri Lanka) and 20 populations (sites). Geographic structure was strong at the large geographical scale of population groups with P<0.001 (Table 3). This is consistent with the strong geographical grouping of haplotypes seen in the network (Fig. 3). Table 4 shows that percent pairwise sequence divergences within groups range from 0.069-0.580%. The each group is significantly different from any other by pairwise Fst value, but the N Thailand from C Thailand. Sri Lanka is the most divergent with an Fst value greater than 0.6 from any of the other groups. It is difficult to define the geographical distance between populations separated by sea. Therefore, a Mantel test for isolation by distance was carried out only for the 8 sites within mainland Southeast Asia with at least 7 samples, ML, MY, MD, CM, LO, VT, PK and NS. The correlation coefficient for the percentage of pairwise sequence divergence between populations with geographic distance was 0.452, but this is not significant based on 1000 permutations (P=0.061).

Pattern of segregating sites and tests of neutrality
The simulated mismatch distribution among haplotypes is closely roughly unimodal (Fig. 4), which is expected in population that have passed through a recent demographic expansion (Rogers and Harpending 1992). The unimodal distribution can also be generated by other processes, notably recombination and high variance in substitution rates. However, the recombination is not expected for these mtDNA sequences, and the possibility caused by high sequence variation is rather low, as the substitutions happen in 48 of 639 sites for the samples investigated. Moreover, the observed mismatch distribution for demographic equilibrium is reliably ragged and often multimodal as it reflects the highly stochastic shape of gene tree (Rogers and Harpending 1992;Harpending et al. 1998). Therefore, we suggest that the unimodal distribution pattern reflects a growing population expansion, and an underlying star-burst-like genealogy in which all of the coalescent events occurred in a narrow time window. No simple statistical test can be made for a fit of the distribution to a Poisson distribution, as the pairwise comparisons are not independent, and the alleles share a common history. Mean distance between haplotypes of the observed mismatch distribution is 3.077 and the variance is 3.934, without much difference with the mean. Despite this, there is a rightward deviation in which the observed distribution differs from a Poisson distribution. There is a deficit of observed values around the mode, a lowness of observed value (compared to the expected) at the lower value and an excess at the upper value (Fig.  4). Neutrality tests showed that Tajima's D and Fu' Fs values are both negative for pooled samples including Sri Lanka (−1.92 and −26.13, respectively) and excluding Sri Lanka (−1.84 and −26.36), which indicates an excess of low-frequency mutations. The extremely significant F s values (P<0.01) in the study most likely indicate a history of population growth or hitchhiking rather than background selection (Fu 1997). The extremely significant Tajima's D values (P<0.01) mean that deviation is more likely to come from selection effects than population expansion, bottleneck or heterogeneity of mutation rates (Tajima 1996). The F s statistic is very sensitive to geographic expansion of populations, which may explain why the F s values are largely negative and extremely significant (P<0.01) also for all population groups investigated.

Effective population size and exponential growth rate
Since evidence earlier in the study indicated population growth, we used the LAMARC program to make simultaneous, phylogenetic estimates of the present-day value of θ, θ 0 , and the population growth rate, g, assuming an exponential modal of growth using a maximum likelihood approach (Kuhner 2006). The result showed that the effective population size, N 0 (4.95 × 10 6 female), is quite large, and indicates a much wider distribution. This agrees with the much wider distribution reported in the literature than represented by the populations studied here. The population growth rate, g (1488.83 ± 75.23), indicates rapid, exponential population growth.

Discussion
In this study, we obtained 140 mtDNA COII sequences of the mosquito species An. aconitus from 21 sites throughout southern China, Myanmar, Vietnam, Thailand, Laos and Sri Lanka (Fig.1). The species has been recorded from 5 provinces in southern China, Yunnan, Guizhou, Guangxi, Hainan and Zhejiang (Lu, 1997), and possibly in Guangdong Province (Baba, 1950). This study confirms its presence in southern Yunnan and Hainan Provinces but no specimens were collected from Guizhou, Guangxi and Guangdong despite visiting many sites there. The species is difficult to distinguish from the An. minimus complex on the basis of morphology. Therefore, it is likely that its reported occurrence in these more northern provinces is due to mistaken identification.
The samples investigated share a single D2 sequence. Despite some intraspecific variation, the sequences from each of the other loci (D3, ITS2, ND5 and COII) are monophyletic with respect to other members of the Myzomyia Series of Anopheles subgenus Cellia (Chen et al., 2003, and unpublished). This supports the current concept of An. aconitus in the area investigated. However, the samples from Sri Lanka have differentiated D3, ITS2, ND5 and COII sequences that each forms a separate clade within An. aconitus (Fig. 2). In the An. minimus complex, the D3 differences between An. minimus, An. harrisoni and An. yaeyamaensis, An. leesoni and An. fluviatilis range between 3 and 5 bp, whereas the difference between Sri Lankan and Southeast Asian An. aconitus reaches 2 bp. Moreover, Sri Lanka population group is the most divergent with an Fst value greater than 0.6 from any of the other groups. These observations suggest that the genetically distinct and geographically separated Sri Lankan populations could be recognised as a subspecies of An. aconitus, but there is a need in sample collecting to cover the distribution area of An. aconitus in further work. Significant isolation by distance was not observed in Southeast Asia based on 1000 permutations (P=0.061), so this pattern may be confirmed by wider sampling. The populations on the island of Hainan and in southern Thailand also show geographically distinct groups of haplotypes for the COII gene (Fig. 3), but the differences are smaller and their sequences for D3, ITS2 and ND5 are not distinct.
The apparently ancestral haplotypes 1 (Fig. 3) has high frequency and occur in several populations. However, they are concentrated in northern and central Thailand, Myanmar and Yunnan, areas that also show higher diversity than any other area. This pattern suggests that it might have been a refuge for the species during the Pleistocene glaciations from which other areas have been colonised.
Interestingly, the Yunnan region has extraordinarily high haplotype diversity with θ π 5.80±3.58 and θ s 6.95±3.24. The samples of the group were collected from Xishuangbanna. This area is in a transitional position from the southern Himalayas to eastern Asia and from tropical Southeast Asia to subtropical China as well as at the junction of the Indian and Burmese plates of Gondwanaland and the Eurasian plate of Laurasia (Audley-Charles, 1987). The Hengduan mountains to the north of the region act as a huge barrier keeping out the cold air from the north in winter, and in valleys and on lower hills below 900 m altitude there are pockets of tropical rain forest, as intermittent tracts or patches in local habitats with a tropical moist climate prevailing due to particular topography (Zhu, 1997). Due to poor access and remoteness from developed and highly populated regions, this area is seldom disturbed by human activity. All these characteristics make the area one of the richest for biodiversity in China, even higher than tropical rain forest of Southeast Asia (Jin and Ou, 1997). For example, 152 and 434 plant species were found to be condensed in 0.0025 and 0.04 km 2 , respectively (Zhu et al., 2004;Jin, 1997), and the arbuscular mycorrhizal fungi spore density can reach 2550 per 100 g soil (Zhao et al., 2003). Previous studies on the biodiversity of these unique tropical rain forests have mainly focused on ecosystem and species biodiversity; however, the present research provides a case of high biodiversity within species. The high haplotype diversity in Yunnan, with no common haplotype, suggests long-term persistence of a locally structured population. Similarly, Sagaing region in Myanmar has also extraordinarily high haplotype diversity with θ π 5.26±3.24 and θ s 8.07±3.55, which was identified as putative refugial region as Yunnan (Morgan et al, 2011).
Hierarchical analysis of molecular variance shows significant variation among population groups, which makes up 36.09% of the total. In this respect, An. aconitus also differs from other Southeast Asian Anopheles species, which show little or no population structure (An. dirus and An. baimaii in Bangladesh, Myanmar and northern, central and southern Thailand, Walton et al., 2000;An. maculatus in northern, central and southern Thailand, Rongnoparut et al., 1999;An. jeyporiensis, Chen et al., 2004).
Although the observed mismatched distribution among haplotypes deviates little from a Poisson fit (Fig. 4), the closely unimodal distributions are consistent with some degree of population expansion, as simulations have shown that stable populations almost never produce this type of profile (Rogers and Harpending 1992;Harpending et al. 1998). The population expansion can be visualised in the star-like haplotype arrays (Fig 3). The extremely significant F s value (P<0.01) also support the inference of population expansion (Table 1). Moreover, the estimate of exponential growth rate, g (1488.83 ± 75.23), in LAMARC indicates a large population expansion as well (Table 5). This result is consistent with that of Morgan et al (2011), which suggests similarly timed population expansions in An. aconitus, An. jeyporiensis, An. baimaii and An. maculatus.
Given an exponential population growth, the coalescence time can be expressed in terms of the number of pairwise differences between samples if a constant mutation rate is assumed (Slatkin and Hudson 1991). The pairwise difference between samples per site are estimated to be 4.82×10 −3 in ARLEQUIN, assuming 10 generations per year and a mutation rate of 10 −8 per site per generation (Powell et al. 1986). From this, it can be inferred that population expansion started approximately 482,000 years ago. The present effective population size, N 0 , is estimated to be 2.47 × 10 7 females ( Table 5).
The immature stages of this species mainly live in rice fields, and therefore Harrison (1980) suggests that the widespread distribution and general abundance of this species were probably due primarily to the spread of the human-rice monoculture system in the Orient. However, the earliest and most convincing archaeological evidence, discovered at Non Nok Tha in northern Thailand, can only trace back to 4000 B.C. for domestication of rice in Southeast Asia (http://www.riceweb.org/history.htm). Originally the larvae and pupae of the species probably occurred primarily in grassy marshes with slow clear running water and along open streams and rivers with grassy margins like An. minimus. In recent history, man has greatly expanded these habitats by cutting forests and exposing streams, by digging ditches for irrigation and by creating artificial grassy marshes in the form of rice fields.
In conclusion, the genetic differentiation analysis of An. aconitus suggested clear genetic structure with highly significant genetic variation present among population groups in Southeast Asia, the closely unimodal mismatch distributions and extremely significant F s values suggest a recent demographic expansion, and the population in Sri Lanka may be considered a distinct subspecies.
The population in Sri Lanka of An. aconitus may be considered a distinct subspecies.
Clear genetic structure was observed with highly significant genetic variation present among population groups in China and Southeast Asia.
Mismatch distributions and extremely significant F s values suggest a recent demographic expansion.  Maximum-likelihood trees of sequences of mtDNA COII and ND5 genes and rDNA D3, D2 and ITS2 loci (An. filipinae and An. pampanai as outgroup taxa), showing their variation within An. aconitus and the subdivision of the Sri Lankan population to subspecies. The sample names indicate their COII haplotypes and localities. Bootstrap percentages of 1000 replicates calculated with the maximum parsimony method are shown above the branches where they exceed 50%. Branch lengths are proportional to the number of character differences. Genealogical relationships among 46 COII haplotypes of An. aconitus estimated by TCS (Clement et al., 2000). The size of a circle corresponds to the COII haplotype frequency, and a unit branch represents one mutation. Small ovals indicate COII haplotypes that were not observed. Observed mismatch distributions among pairwise differences among haplotypes for all populations except Sri Lanka, which is closely unimodal and fits to a Poisson distribution (χ 2 : s 2 = 1.32 and P < 0.001). Table 1 Summary of populations, haplotypes, nucleotide diversity and neutrality tests of An. aconitus. S is the number of segregating sites, and θ π and θ s are the estimates of nucleotide diversity per base.   Genetic distances within and between population groups.  Table 5 FLUCTUATE estimates of effective population size and exponential growth rate for An. aconitus. θ 0 ± SD N 0 g ± SD 0.99 ± 0.145 2.47 × 10 7 1488.83 ± 75.23

Group Hainan Yunnan Myanmar S Vietnam N Thailand C Thailand S Thailand Sri Lanka
Point estimates of θ 0 and g are shown with approximate 95% confidence intervals. The exponential growth rate, g, is in units of 1/μ per generation (θ = 2Nμ in haploid organism). N 0 , the current effective population size of females, is estimated from θ 0 assuming μ of 1 × 10 −8 per site per generation (Powell et al. 1986).