Genomic signatures of population decline in the malaria mosquito Anopheles gambiae
© O’Loughlin et al. 2016
Received: 29 October 2015
Accepted: 5 March 2016
Published: 24 March 2016
Population genomic features such as nucleotide diversity and linkage disequilibrium are expected to be strongly shaped by changes in population size, and might therefore be useful for monitoring the success of a control campaign. In the Kilifi district of Kenya, there has been a marked decline in the abundance of the malaria vector Anopheles gambiae subsequent to the rollout of insecticide-treated bed nets.
To investigate whether this decline left a detectable population genomic signature, simulations were performed to compare the effect of population crashes on nucleotide diversity, Tajima’s D, and linkage disequilibrium (as measured by the population recombination parameter ρ). Linkage disequilibrium and ρ were estimated for An. gambiae from Kilifi, and compared them to values for Anopheles arabiensis and Anopheles merus at the same location, and for An. gambiae in a location 200 km from Kilifi.
In the first simulations ρ changed more rapidly after a population crash than the other statistics, and therefore is a more sensitive indicator of recent population decline. In the empirical data, linkage disequilibrium extends 100–1000 times further, and ρ is 100–1000 times smaller, for the Kilifi population of An. gambiae than for any of the other populations. There were also significant runs of homozygosity in many of the individual An. gambiae mosquitoes from Kilifi.
These results support the hypothesis that the recent decline in An. gambiae was driven by the rollout of bed nets. Measuring population genomic parameters in a small sample of individuals before, during and after vector or pest control may be a valuable method of tracking the effectiveness of interventions.
Many population genomic parameters depend upon population size. These include genetic diversity (θ or π), linkage disequilibrium, the population recombination parameter ρ, and runs of homozygosity [1–4]. Changes in population size can also lead to transient changes in the allele frequency spectrum and statistics based on it (e.g. Tajima’s D; [5, 6]). Thus population genomic data contain information on the past demographic history and might, therefore, show the effects of efforts to suppress the population, as has been observed in the Plasmodium genome after malaria control [7, 8].
In many parts of Africa there have been concerted efforts to control malaria transmission by controlling the mosquito vector using insecticide treated bed nets (ITNs) and indoor residual spraying (IRS). In some places these efforts have been successful, resulting in substantial reductions in the numbers of Anopheles gambiae (the most important vector species in sub-Saharan Africa) and in malaria transmission [9, 10]. An. gambiae is particularly susceptible to ITNs and IRS because of its propensity to bite and rest indoors. Other vectors may be less susceptible to these control methods, including the sibling species Anopheles arabiensis, which is able to bite earlier and outdoors due to increased resistance to desiccation [10, 11].
One place where there has been a particularly striking reduction in An. gambiae abundance is in the Kilifi district of coastal Kenya: entomological surveys have revealed an overall reduction in density of An. gambiae s.l. and Anopheles funestus, accompanied by a shift in the proportions of different species, with An. arabiensis and Anopheles merus replacing An. gambiae s.s. and An. funestus as the major vectors . The authors attribute this shift in species composition to the widespread distribution and use of ITNs from 2006 onwards, although they do not rule out land-use change and improvements in house construction as contributing factors.
A previous study by O’Loughlin et al.  reported a RADseq analysis of the An. gambiae s.l. species complex from three locations approx. 200 km apart in East Africa: Moshi, Muheza, and Kilifi. Although mosquito control is also being carried out in Muheza and Moshi, at the time of sampling there had been no reported decline in mosquito numbers. In this study it was found that genetic diversity in An. gambiae s.s. was slightly but significantly lower in Kilifi than in Muheza (~5 % lower π and ~15 % lower θW averaged across all chromosomes). Diversity did not differ among the three An. arabiensis populations. The study also found that An. gambiae in Kilifi was the only population with a positive value for Tajima’s D, reflecting a deficit in low frequency polymorphisms, consistent with a recent decline in population size . Modelling of the allele frequency spectra showed that An. arabiensis and An. merus fitted a simple model of modest population expansion, whereas the An. gambiae populations showed a more complex history of past population expansion followed by population decline. In the case of An. gambiae from Kilifi, the present population size was inferred to be smaller than the historical, pre-expansion size.
These results appear to be consistent with a recent population reduction for An. gambiae in Kilifi, perhaps due to control efforts . To investigate this hypothesis more closely, here the analysis is expanded to consider linkage disequilibrium and the population recombination parameter ρ. RADseq data consists of a small fraction of the genome so is not suitable for some linkage-based methods of inferring population history such PSMC and MSMC [15, 16]. However the number of SNPs and their location throughout the genome make it ideal for calculating ρ. In the standard neutral model at equilibrium ρ has an expected value (or is defined as) 4Ner, where Ne is the effective population size and r is the recombination rate per base per generation. ρ is inversely related to the levels of linkage disequilibrium in a sample. It has previously been observed, both empirically and by simulation, that ρ is strongly affected by non-equilibrium demographics and selection [3, 17]. Although it is well established that ρ decreases after population bottlenecks followed by recovery (e.g. [3, 17, 18]), these studies did not explore the effect of very recent population declines without recovery, such as after successful vector control. Therefore, in this study, simulations are used to study the time-scale over which the different population genomic parameters are expected to change in response to reductions in population size.
Simulations: genomic signatures of successful control
To investigate the expected effect of a recent population crash on ρ, θW, π and Tajima’s D, sequences were simulated under different demographic scenarios using Hudson’s ms . A sample size of n = 26 was used (equivalent to 13 diploid individuals), chosen to match the sample size that were analysed with RADseq, and simulated sequences of 50 kb in length. Populations were simulated with an ancestral size of 2 million (the estimated long term Ne for An. gambiae population from ), a mutation rate of 1.1 × 10−9 per generation (estimated from divergence of Drosophila lineages  and assuming ten generations per year), and a ρ of 10 times present θ, which is the neutral expectation of ρ/θ calculated using the recombination rate for chromosome 3L  and is within the range of values seen at selective and demographic neutrality in Drosophila populations [22, 23]. Throughout the simulations parameters for the 3L chromosome arm were used, because 2L and 2R contain polymorphic inversions in An. gambiae and similarly 2R and 3R in An. arabiensis. Population crashes were simulated in which the population size after the crash was 10−2, 10−3, 10−4 or 10−5 of the ancestral population, and occurred 10, 102, 103 or 104 generations in the past. Ms commands are given in Additional file 1.
Empirical data from East Africa
Metrics of data sets containing variant and invariant sites
No. of tag
No. of SNPsa
SNPs per tag
Another potential sign of population crash is extended runs of homozygosity caused by mating between related individuals [4, 29]. To look for runs of homozygosity, heterozygosity at every site was plotted across the genome for each mosquito individually.
Simulations: genomic signatures of successful control
Empirical data from East Africa
Population recombination parameter ρ
The results from the Kilifi population of An. gambiae are markedly different, with ρ low in all arms. Comparison with the Muheza population demonstrates a statistically significant difference (paired t test, t = 2.39, p = 0.04). The reduction in ρ in Kilifi vs Muheza is between 11 % (2L) and 99.9 % (3R), with an average across chromosomes of 82 %. This difference is substantially greater than the reductions seen in diversity measures π and θW. As a result, the ratio ρ/θW for An. gambiae from Kilifi is more than 1000-fold lower then neutral expectations. The comparison among populations and species is perhaps clearest for chromosome arm 3L, which does not have segregating inversions in any of these three mosquito species (hatched bars in Fig. 3).
Estimating the timing and severity of population crash from ρ
Runs of homozygosity
Vector control is an important tool in the fight against malaria and other vector borne diseases, but as yet measuring the entomological impact of these methods has been largely anecdotal or ad-hoc, with even large shifts in abundance or species proportions being difficult to quantify. This is because monitoring changes in mosquito numbers is not an easy activity; entomological surveying is labour-intensive and expensive, and prone to variation caused by seasonal fluctuations, different collection methods and degree of collection effort (e.g. [31, 32]). One alternative for detecting changes in population size may be monitoring the genome.
In the Kilifi district of Kenya, ITNs were first introduced in a randomized trial in 1993, when their effectiveness in reducing malaria incidence in children was proven . Since then, ITN coverage increased gradually until a large-scale distribution program in 2006 resulted in coverage rising to 67 % . Over a similar time period, entomological surveys in the Kilifi region have detected a decline in vector density (from 1990 to 2010), with an accompanying change in An. gambiae s.l. species composition . Anopheles gambiae declined from 79 % of An. gambiae s.l. in 1997–1998 to an undetectable level in 2007–2008. The dominant An. gambiae s.l. species in 2007–2008 was An. arabiensis (93 %), with An. merus contributing 5 %.
This change in mosquito abundance makes it a good model system to test whether the effect of population decline can be detected in the genome. Previously it has been shown that An. gambiae from Kilifi have lower π and θ and higher Tajima’s D compared to those from Muheza, where there has been no reported change in mosquito abundance . Anopheles arabiensis from the same locations did not show any differences. These results were consistent with the observed recent changes in abundance.
Here, simulations are presented showing that analysis of linkage disequilibrium and ρ can be a more sensitive test for population decline, as large signals are seen more quickly. Therefore, these statistics were compared among populations and species, and indeed linkage disequilibrium extends 100–1000 times further in Kilifi than in comparator populations, and ρ is 100–1000 times lower. Population genetic simulations indicate that the observed ρ value implies that the sampled population is no more than a 5 × 10−4 (1/2000th) of the ancestral population. The observed values would also be consistent with an even greater, recent population decline. The unusual runs of homozygosity in Kilifi An. gambiae also suggest that a recent and severe population crash has occurred, which is resulting in signs of inbreeding in some individuals.
One weakness of the simulations is that they do not take into account spatial structure, which is known to effect the extent of linkage disequilibrium, especially when all the samples are from a single sub-population in a large “stepping-stone” array [3, 30]. In general, there is little differentiation among these East African populations of the same species , but since the extent of differentiation is determined by the product of the migration rate and the population size, there may be some interaction between a population crash and population structure affecting linkage disequilibrium and ρ. The inferences would also be improved by analysing a time series of samples from before and after the population crash, but unfortunately pre-control samples were not available for analysis.
The fact that a similar pattern of increased linkage disequilibrium and reduced ρ was not observed in An. arabiensis and An. merus collected from Kilifi at the same time as the An. gambiae samples suggests that these species have not been affected by the same population decline, and is consistent with the entomological observations of a shift in species composition . This difference among species supports the hypothesis that the population decline has largely been due to the use of ITNs, as these are expected to have a larger impact on the highly anthropophilic and indoor-biting An. gambiae, compared with partially zoophilic and outdoor biting species such as An. arabiensis and An. merus .
The results of the simulations and observed data from Kilifi show promise for the prospective monitoring of vector control efforts. Mosquitoes have ~10 generations per year, so measurements of ρ, θ, π and Tajima’s D taken pre-intervention and at one-yearly intervals should detect whether control is succeeding. Whole-genome sequencing is not necessary for measuring ρ. The simulated data suggests it is possible to get reliable estimates of ρ from as few as 300 segregating sites, so a medium-throughput SNP genotyping platform such as RADseq or Golden-Gate assay would be sufficient for monitoring. For Anopheles species it is important to use SNPs that are not on chromosome arms containing segregating inversions. The sample here of just 11–13 mosquitoes per population was sufficient to distinguish clear differences in ρ between populations and chromosome arms.
Observations of genomic diversity and linkage disequilibrium in a small sample of just 13 mosquitoes provide compelling evidence that An. gambiae in Kilifi has undergone a recent population crash. In practical terms, this means that regular monitoring of a small number of genomes could allow rapid detection of whether a control intervention is succeeding. In nature, there may be complicating factors such as seasonal variation and immigration from non-treated areas, but the results presented here from Kilifi suggest that a severe population crash will be detectable despite these factors. Given the practical difficulties of measuring mosquito abundance by direct surveying, genotyping a small number of mosquitoes could be an attractive alternative for assessing the entomological impact of vector control.
Data available in Dryad Digital Repository: doi:10.5061/dryad.hm6tt.
AB and SO devised the study and wrote the manuscript. SO performed molecular work, simulations and data analysis. SM, CM, JM and FM facilitated field collections of mosquitoes. All authors read and approved the final manuscript.
Mosquito collections were undertaken by staff of the NIMR Amani Research Centre (Muheza, Tanzania), Kilimanjaro Christian Medical Centre (Moshi, Tanzania) and KEMRI (Kilifi, Kenya). This manuscript is published with the permission of the director, KEMRI. This work was supported by a grant from the Foundation for the National Institutes of Health through the Vector-Based Transmission of Control: Discovery Research (VCTR) program of the Grand Challenges in Global Health Initiative, and from the European Union’s Seventh Framework Programme (FP7/2007–2013) under grant agreement no 228,421- INFRAVEC.
The authors declare that they have no competing interests.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
- Watterson GA. On the number of segregating sites in genetical models without recombination. Theor Popul Biol. 1975;7:256–76.View ArticlePubMedGoogle Scholar
- Hill WG. Estimation of effective population-size from data on linkage disequilibrium. Genet Res. 1981;38:209–16.View ArticleGoogle Scholar
- Andolfatto P, Przeworski M. A genome-wide departure from the standard neutral model in natural populations of Drosophila. Genetics. 2000;156:257–68.PubMedPubMed CentralGoogle Scholar
- Kirin M, McQuillan R, Franklin CS, Campbell H, McKeigue PM, Wilson JF. Genomic runs of homozygosity record population history and consanguinity. PLoS One. 2010;5:e13996. doi:10.1371/journal.pone.0013996.View ArticlePubMedPubMed CentralGoogle Scholar
- Tajima F. Statistical-method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989;123:585–95.PubMedPubMed CentralGoogle Scholar
- Gutenkunst RN, Hernandez RD, Williamson SH, Bustamante CD. Inferring the joint demographic history of multiple populations from multidimensional SNP frequency data. PLoS Genet. 2009;5:e1000695. doi:10.1371/journal.pgen.1000695.View ArticlePubMedPubMed CentralGoogle Scholar
- Nkhoma SC, Nair S, Al-Saai S, Ashley E, Mcgready R, Phyo AP, et al. Population genetic correlates of declining transmission in a human pathogen. Mol Ecol. 2013;22:273–85. doi:10.1111/mec.12099.View ArticlePubMedPubMed CentralGoogle Scholar
- Daniels RF, Schaffner SF, Wenger EA, Proctor JL, Chang HH, Wong W, et al. Modeling malaria genomics reveals transmission decline and rebound in Senegal. Proc Natl Acad Sci USA. 2015;112:7067–72. doi:10.1073/pnas.1505691112.View ArticlePubMedPubMed CentralGoogle Scholar
- Bayoh MN, Mathias DK, Odiere MR, Mutuku FM, Kamau L, Gimnig JE, et al. Anopheles gambiae: historical population decline associated with regional distribution of insecticide-treated bed nets in western Nyanza Province, Kenya. Malar J. 2010;9:62. doi:10.1186/1475-2875-9-62.View ArticlePubMedPubMed CentralGoogle Scholar
- Russell TL, Govella NJ, Azizi S, Drakeley CJ, Kachur SP, Killeen GF. Increased proportions of outdoor feeding among residual malaria vector populations following increased use of insecticide-treated nets in rural Tanzania. Malar J. 2011;10:80. doi:10.1186/1475-2875-10-80.View ArticlePubMedPubMed CentralGoogle Scholar
- Lindblade KA, Gimnig JE, Kamau L, Hawley WA, Odhiambo F, Olang G, et al. Impact of sustained use of insecticide-treated bednets on malaria vector species distribution and culicine mosquitoes. J Med Entomol. 2006;43:428–32. doi:10.1603/0022-2585.View ArticlePubMedGoogle Scholar
- Mwangangi JM, Mbogo CM, Orindi BO, Muturi EJ, Midega JT, Nzovu J, et al. Shifts in malaria vector species composition and transmission dynamics along the Kenyan coast over the past 20 years. Malar J. 2013;12:13. doi:10.1186/1475-2875-12-13.View ArticlePubMedPubMed CentralGoogle Scholar
- O’Loughlin SM, Magesa S, Mbogo C, Mosha F, Midega J, Lomas S, et al. Genomic analyses of three malaria vectors reveals extensive shared polymorphism but contrasting population histories. Mol Biol Evol. 2014;31:889–902. doi:10.1093/molbev/msu040.View ArticlePubMedPubMed CentralGoogle Scholar
- Tajima F. The effect of change in population-size on DNA polymorphism. Genetics. 1989;123:597–601.PubMedPubMed CentralGoogle Scholar
- Li H, Durbin R. Inference of human population history from individual whole-genome sequences. Nature. 2011;475:493–6. doi:10.1038/nature10231.View ArticlePubMedPubMed CentralGoogle Scholar
- Schiffels S, Durbin R. Inferring human population size and separation history from multiple genome sequences. Nat Genet. 2014;46:919–25. doi:10.1038/ng.3015.View ArticlePubMedPubMed CentralGoogle Scholar
- Haddrill PR, Thornton KR, Charlesworth B, Andolfatto P. Multilocus patterns of nucleotide variability and the demographic and selection history of Drosophila melanogaster populations. Genome Res. 2005;15:790–9. doi:10.1101/gr.3541005.View ArticlePubMedPubMed CentralGoogle Scholar
- Pritchard JK, Przeworski M. Linkage disequilibrium in humans: models and data. Am J Hum Genet. 2001;69:1–14. doi:10.1086/321275.View ArticlePubMedPubMed CentralGoogle Scholar
- Hudson RR. Generating samples under a Wright–Fisher neutral model of genetic variation. Bioinformatics. 2002;18:337–8. doi:10.1093/Bioinformatics/18.2.337.View ArticlePubMedGoogle Scholar
- Tamura K, Subramanian S, Kumar S. Temporal patterns of fruit fly (Drosophila) evolution revealed by mutation clocks. Mol Biol Evol. 2004;21:36–44. doi:10.1093/molbev/msg236.View ArticlePubMedGoogle Scholar
- Pombi M, Stump AD, Della Torre A, Besansky NJ. Variation in recombination rate across the X chromosome of Anopheles gambiae. Am J Trop Med Hyg. 2006;75:901–3.PubMedGoogle Scholar
- Bachtrog D, Andolfatto P. Selection, recombination and demographic history in Drosophila miranda. Genetics. 2006;174:2045–59. doi:10.1534/genetics.106.062760.View ArticlePubMedPubMed CentralGoogle Scholar
- Thornton K, Andolfatto P. Approximate Bayesian inference reveals evidence for a recent, severe bottleneck in a Netherlands population of Drosophila melanogaster. Genetics. 2006;172:1607–19. doi:10.1534/genetics.105.048223.View ArticlePubMedPubMed CentralGoogle Scholar
- Li H, Durbin R. Fast and accurate short read alignment with Burrows–Wheeler transform. Bioinformatics. 2009;25:1754–60. doi:10.1093/bioinformatics/btp324.View ArticlePubMedPubMed CentralGoogle Scholar
- Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25:2078–9. doi:10.1093/bioinformatics/btp352.View ArticlePubMedPubMed CentralGoogle Scholar
- McVean G, Awadalla P, Fearnhead P. A coalescent-based method for detecting and estimating recombination from gene sequences. Genetics. 2002;160:1231–41.PubMedPubMed CentralGoogle Scholar
- Stump AD, Pombi M, Goeddel L, Ribeiro JMC, Wilder JA, Torre AD, et al. Genetic exchange in 2La inversion heterokaryotypes of Anopheles gambiae. Insect Mol Biol. 2007;16:703–9.View ArticlePubMedGoogle Scholar
- Gillespie JH. Population genetics: a concise guide. Baltimore: John Hopkins University Press; 2004.Google Scholar
- McQuillan R, Leutenegger AL, Abdel-Rahman R, Franklin CS, Pericic M, Barac-Lauc L et al. Erratum: runs of homozygosity in European populations (Am J Hum Genet. 2008;83:359-72, 2008). Am J Hum Genet. 2008;83:658. doi:10.1016/j.ajhg.2008.10.009.
- De A, Durrett R. Stepping-stone spatial structure causes slow decay of linkage disequilibrium and shifts the site frequency spectrum. Genetics. 2007;176:969–81. doi:10.1534/genetics.107.071464.View ArticlePubMedPubMed CentralGoogle Scholar
- Govella NJ, Chaki PP, Geissbuhler Y, Kannady K, Okumu F, Charlwood JD, et al. A new tent trap for sampling exophagic and endophagic members of the Anopheles gambiae complex. Malar J. 2009;8:157. doi:10.1186/1475-2875-8-157.View ArticlePubMedPubMed CentralGoogle Scholar
- Chaki PP, Mlacha Y, Msellemu D, Muhili A, Malishee AD, Mtema ZJ, et al. An affordable, quality-assured community-based system for high-resolution entomological surveillance of vector mosquitoes that reflects human malaria infection risk patterns. Malar J. 2012;11:172. doi:10.1186/1475-2875-11-172.View ArticlePubMedPubMed CentralGoogle Scholar
- Nevill CG, Some ES, Mungala VO, Mutemi W, New L, Marsh K, et al. Insecticide-treated bednets reduce mortality and severe morbidity from malaria among children on the Kenyan coast. Trop Med Int Health. 1996;1:139–46.View ArticlePubMedGoogle Scholar
- Okiro EA, Hay SI, Gikandi PW, Sharif SK, Noor AM, Peshu N, et al. The decline in paediatric malaria admissions on the coast of Kenya. Malar J. 2007;6:151. doi:10.1186/1475-2875-6-151.View ArticlePubMedPubMed CentralGoogle Scholar
- Kipyab PC, Khaemba BM, Mwangangi JM, Mbogo CM. The bionomics of Anopheles merus (Diptera: culicidae) along the Kenyan coast. Parasit Vectors. 2013;6:37. doi:10.1186/1756-3305-6-37.View ArticlePubMedPubMed CentralGoogle Scholar