- Open Access
Ecological niche model of Phlebotomus alexandri and P. papatasi (Diptera: Psychodidae) in the Middle East
International Journal of Health Geographicsvolume 9, Article number: 2 (2010)
The purpose of this study is to create distribution models of two sand fly species, Phlebotomus papatasi (Scopoli) and P. alexandri (Sinton), across the Middle East. Phlebotomus alexandri is a vector of visceral leishmaniasis, while P. papatasi is a vector of cutaneous leishmaniasis and sand fly fever. Collection records were obtained from literature reports from 1950 through 2007 and unpublished field collection records. Environmental layers considered in the model were elevation, precipitation, land cover, and WorldClim bioclimatic variables. Models were evaluated using the threshold-independent area under the curve (AUC) receiver operating characteristic analysis and the threshold-dependent minimum training presence.
For both species, land cover was the most influential environmental layer in model development. The bioclimatic and elevation variables all contributed to model development; however, none influenced the model as strongly as land cover.
While not perfect representations of the absolute distribution of P. papatasi and P. alexandri, these models indicate areas with a higher probability of presence of these species. This information could be used to help guide future research efforts into the ecology of these species and epidemiology of the pathogens that they transmit.
Phlebotomus (Phlebotomus) papatasi (Scopoli) and P. (Paraphlebotomus) alexandri (Sinton) are widely distributed across parts of Europe, Africa, and Asia. Phlebotomus papatasi is a vector of sand fly fever virus and Leishmania major, which causes cutaneous leishmaniasis [1–3]. Phlebotomus alexandri is a vector of L. donovani and is a suspected vector of L. infantum, both of which cause visceral leishmaniasis [4–6]. Though these species are important disease vectors, little is known about the ecology and distribution of each.
Phlebotomus papatasi ranges from Morocco and Spain, across the Mediterranean Basin to India and south to parts of the Sudan and Ethiopia . Phlebotomus papatasi is most abundant in areas with a mean minimum temperature of 16°C and mean maximum temperature of 44°C from May to October . It can be found at elevations ranging from near sea level to over 1100 m .
Phlebotomus alexandri ranges from Spain and Morocco east to the mountains in northwestern China and as far south as southern Ethiopia . This species has been recorded at elevations ranging from sea level to 1500 m above sea level [11, 12]. In Djibouti, this species is found on the coastal plain, inland plateau, and highland valleys .
Characterizing the distribution and ecology of these vector species would be valuable in better understanding the epidemiology of sand fly fever and leishmaniasis. Cross et al. (1996) developed a model of P. papatasi distribution in Southwest Asia based on weather and the normalized difference vegetation index (NDVI); however, the model was not validated. Since the study was completed, powerful presence-only modeling techniques and software have been developed.
These newer modeling methods include ecological niche modeling (ENM). ENM uses presence data in conjunction with environmental data to develop models of habitat range for a given organism . It is often used to examine the distribution of species that have not had intense, methodological sampling. This technique has been used in modeling distribution of diseases, such as dengue, and vectors, such as Anopheles gambiae[15–19]. In addition, niche modeling has been used to examine the distribution and potential distribution of Lutzomyia spp. vectors of leishmaniasis in South America [15, 16, 20].
In this study, we use ENM to develop distribution models for P. papatasi and P. alexandri in the Middle East. Using these models, we attempt to identify environmental factors which influence the distribution of these species.
The MaxEnt model for P. papatasi is shown in Figure 1. The AUC for the training points was 0.944 and for test points was 0.884, with a standard deviation of 0.042. The minimum training presence for a training point was 0.197; therefore, this was set as the threshold for binomial conversion. The fractional predicted area (the area coded as 1 = present) is 0.346 and the omission rate for test points was 0.091. At this threshold, the test points were classified significantly better by the model than by random selection (p < 0.0001).
Jackknife tests of variable importance show that land cover was the most influential variable in model development (Figure 2). The training gain when land cover was the only variable used in model development was high, indicating that it contributes strongly to the model. When land cover was removed from the model, training gain dropped. This indicates that the land cover variable contains unique information that is required for model creation. The land cover types and probabilities associated with the training points are given in Table 1. Points classified as urban, field/woody savanna, and woody savannah coverages have high probabilities of presence. However, the sample size is small for both field/woody savanna and woody savanna. Points classified as bare desert have low probabilities of presence; however the sample size is small. All other classes have either very wide ranges or sample size of one.
The remaining variables contributed much less to model development (Figure 2). Isothermality (Bio3), maximum temperature in the warmest month (Bio5), minimum temperature in the coldest month (Bio6), mean temperature in the coldest quarter (Bio11), precipitation seasonality (Bio15), and precipitation in the coldest quarter (Bio19) all had modest gains when run with only the variable in question. However, the model gain was minimally decreased by the exclusion of the variables from the analysis. This is a reflection of the correlated nature of the bioclimatic variables and an indication that none of these variables in isolation are overwhelmingly contributing to this model.
The MaxEnt model for P. alexandri is shown in Figure 3. The training AUC was 0.942 and the test AUC was 0.844, with a standard deviation of 0.044. The minimum training presence among training points was 0.164. At this threshold, the fractional predicted area was 0.368 and the omission rate for test points was 0.200. The model classifies the test points correctly significantly more than a random model (p < 0.0001).
As in the model for P. papatasi, land cover was the most influential variable in modeling P. alexandri (Figure 4). Jackknife tests show high training gain when land cover is considered alone and a large drop training gain when land cover is omitted from the model. The land cover types and probabilities associated with the training points are given in Table 1. As with P. papatasi, points classified as urban, field/woody savanna, and woody savanna coverages have high probabilities of presence. However, the sample sizes for each of these habitats are small. All other classes have either very wide ranges or small sample size.
The 19 bioclimatic variables and elevation all have very modest training gains when considered in isolation (Figure 4), indicating that none of them strongly contribute to model development on their own. Elevation, minimum temperature in the coldest month (Bio6), mean temperature in the wettest quarter (Bio8), and precipitation in the driest month all show modest decreases in training gain when removed from the model. This indicates that they may contain unique information required for the model.
In this study, land cover contributed strongly to the development of the models for both species. Table 1 shows the different land cover classes associated with the species presence points that were used to develop the models. Though it is tempting draw direct conclusions from these data, one must view these results cautiously. The distribution of model development points in the different land cover classes simply reflects the locations of the points that were randomly selected as training points after the larger data set was entered into the modeling program. While these models may be used to drive hypotheses for further research into the ecology of P. papatasi and P. alexandri, they should not be taken as conclusive evidence of relationships between a certain land cover classification and the relative probability of species presence.
In the models presented here, the urban land cover class was associated with a high probability of presence for both species. This may be partly due to sampling bias, as collections of phlebotomine sand flies tend to be associated with research related to human leishmaniasis. However, in a survey of sand flies in Turkey, P. alexandri was present in the more urban areas of the province, but not in the rural areas . Phlebotomus alexandri is anthropophilic. It is possible that this feeding relationship would drive the species to be more probable in urban environments simply because humans are more readily available. More ecological studies should be conducted to determine if there is a relationship between these species and urban areas or if this is a reflection of sampling bias toward collection of sand flies in areas where there is human disease.
Other non-urban land cover types are also important for these species. In the present study, the logistic probability of presence of both species for the points that fall in the woody savannah and field and woody savanna classes are high (mean 0.9278 and 0.7955, respectively), though the sample sizes are small. The range of probability of presence for both species in relation to semi-desert shrub is extremely wide for both species (Table 1), making it difficult to infer anything about this land cover class from these models. However, research in Morocco has shown that both species are associated with desert, scrub vegetation . For P. papatasi, bare desert appears to be related to a low probability of presence in this study, though only five training points were located in this cover class. In Israel, P. papatasi is more abundant in areas with more humid soils, capable of supporting desert vegetation, than in areas with low soil moisture and less vegetation [23, 24]. Sand flies require a sugar meal, taken from plant material [25, 26]; therefore, they would be expected to be less abundant in areas with little or no vegetation, such as a barren desert. Further field-based ecological research is necessary to better determine what relationships exist between different land cover types and sand flies.
Cross et al. (1996) used normalized difference vegetation index (NDVI) and weather data to develop a predictive model of P. papatasi distribution in southwest Asia; however the models were not validated. We considered including 1-kilometer resolution NDVI data in this model; however, the addition of NDVI to the model development did not affect the outcome. When land cover was excluded and NDVI and the bioclimatic variables were used for model building, the resulting product did not perform as well as the models presented here. While NDVI may help model development in some cases, 1-km resolution NDVI data taken from Advanced Very High Resolution Radiometer (AVHRR) sensors was not as useful as other environmental and climatic variables in this case.
Temperature and precipitation are important for model development, but were not the leading factors in the development of the models presented here. In Morocco, sand flies are most active in the hot, dry season, with P. papatasi most abundant when ambient temperature is 32-36°C range . In Pondicherry, India, P. papatasi reaches its peak abundance at the end of monsoon season . Sand flies in Oman, particularly P. alexandri, are more abundant during periods of low humidity and high temperature (Roberts 1994). In the present study, the bioclimatic variables all contributed toward model development. However, none of these variables were particularly valuable in isolation for either model creation or validation. This is most likely a reflection of the correlated nature of the temperature and precipitation variables.
These models are estimates of the distribution of P. papatasi and P. alexandri, based on the environmental layers chosen in the study and the location of collection records. As such they are not a definitive guide as to whether or not a species will be present in a given area. However, they can be used to estimate to the probability that P. papatasi and/or P. alexandri are present in an area. Since animal reservoir, Leishmania parasite, and/or leishmaniasis disease data are not included in these models, they do not reflect the distribution of disease. Future models incorporating records of Leishmania parasites, disease, and reservoir populations could provide a better understanding of the distribution and risk of cutaneous and visceral leishmaniasis in this area. Finally, further exploration of the relationship between different remote sensing products, such and land cover classifications and vegetations indices, and sand fly populations in the Middle East may help to further refine these models.
Materials and methods
The study area includes part or all of the following countries: Afghanistan, Armenia, Azerbaijan, Cyprus, Egypt, India, Iran, Iraq, Israel, Jordan, Kuwait, Kyrgyzstan, Lebanon, Oman, Pakistan, Qatar, Saudi Arabia, Syria, Tajikistan, Turkey, Turkmenistan, the United Arab Emirates, Uzbekistan, and Yemen (Figure 1). The coordinates delineating the corners of the study area are: northwest corner, N 42.0819, E 25.4443; southwest corner, N 11.3508, E 25.443; northeast corner, N 42.0819, E 75.0586; southeast corner, N 11.3508, E 75.0586.
Presence data for the species were taken from records in the scientific literature dating from 1950 through 2007, collections performed by U.S. military entomologists in Iraq and Afghanistan between 2003 and 2006, and in Turkey in 2006. All coordinates were converted to the decimal degrees format.
The literature search was conducted using PubMed, searching on the terms "Phlebotomus papatasi," "Phlebotomus alexandri," "sand fly," and "sandfly." The search yielded 427 records for "Phlebotomus papatasi," 31 records for "Phlebotomus alexandri," 3026 records for "sand fly," and 3155 records for "sand fly." After exclusion of articles that did not include collection records of either species of interest and articles with records outside the region of interest, there were 31 publications for P. alexandri and 173 for P. papatasi. Records for which there was uncertainty about the location of the sampling or no specific location given were then excluded. This process contributed 98 P. papatasi presence points from 22 articles, by 20 first authors and 79 P. alexandri presence points from 25 publications, by 19 different first authors to the modeling process. For both species, the articles were published between 1968 and 2007.
Once the unpublished Afghanistan, Iraq, and Turkey collection data were added to the literature records, 115 points for P. papatasi and 98 points for P. alexandri were entered into the MaxEnt program for model development. The program was then set to exclude duplicate presence records within the same pixel. For P. papatasi, 25 points were excluded, leaving 90 presence records for model development and validation. Of these, 68 were randomly selected for model development, with the remaining 22 used to test the model. For P. alexandri, 18 records were excluded, leaving 80 presence records for model development and validation. Sixty of these records were randomly selected for model development, with the remaining 20 used to test the model.
Climate and elevation layers were obtained from the WorldClim database, version 1.4 . This database provides climate layers at a spatial resolution of 1 km2 and is derived from weather station data from 1950-2000 . For the purposes of this study, the WorldClim bioclimatic variables were used (Table 2).
Land cover data were obtained from the U.S. Geological Survey's (USGS) Earth Resources Observation and Science (EROS) Data Center . This is a global land cover classification that is broken into 96 land cover classes at 1-km resolution. Of the ninety-six classes, only 60 occur in the study area. The 68 training points for P. papatasi fall into just 10 of these classes (Table 1). The 60 training points for P. alexandri are associated with 10 classes (Table 1).
Model Building and Evaluation
The niche modeling application MaxEnt was used in this analysis . The MaxEnt program develops models of species distribution, subject to environmental variables entered into the model building process, using the principles of the maximum entropy distribution [33, 34]. Models were developed using MaxEnt version 3.2.1.
Seventy-five percent of the data points for each species were randomly selected as training points, used in model building. The remaining 25% of the records were test points, used in model validation. Duplicate presence records were removed by the MaxEnt program prior to model development. The MaxEnt model output was set to logistic, which returns an estimated probability of presence for a given location between the values of 0 (no probability of species presence) and 1 (species is certain to be present). All other parameters were set to the default settings.
The model was evaluated using both threshold-dependent and threshold-independent methods. The area under the curve (AUC) of the receiver operating characteristic (ROC) analysis is a threshold-independent method of evaluating model quality. This technique computes the total area under the curve created by plotting sensitivity against the fractional predicted area for the species [33–35]. The threshold-dependent measure used here is the minimum training presence in which the probabilities are converted to binomial values with 0 being absent and 1 being present [33–35]. Using this method, all pixels with a probability of presence equal to or greater than that of the training point with the lowest probability of presence are classified as present and all pixels with a lower probability of presence are classified as absent. A one-tailed binomial test is then performed with the null hypothesis being that the model does not predict the test points better than random .
In order to determine which variables contribute most to the model development, the MaxEnt program was set to calculate jackknife tests of variable importance. The jackknife procedure produces three different types of models: (1) models created with one variable at a time excluded and all other variables included, (2) models created with only one variable included, and (3) a model created with all variables [33–35]. Variables that are most important to model development are those that decrease the training gain when removed from the model and show gain when the model is developed with only one variable.
Material has been reviewed by the Walter Reed Army Institute of Research. There is no objection to its publication. The opinions or assertions contained herein are the private views of the author, and are not to be construed as official, or as reflecting true views of the Department of the Army or the Department of Defense.
Ashford RW: Leishmaniasis. The Encyclopedia of Arthropod-transmitted Infections. Edited by: Service MW. 2001, New York, NY: CABI Publishing, 269-279.
Ashford RW: Phlebotomus Fevers. The Enclopedia of Arthorpod-transmitted Infections. Edited by: Service MW. 2001, New York, NY: CABI Publishing, 397-401.
Yaghoobi-Ershadi MR, Akhavan AA, Zahraei-Ramazani AR, Jalali-Zand AR, Piazak N: Bionomics of Phlebotomus papatasi (Diptera: Psychodidae) in an endemic focus of zoonotic cutaneous leishmaniasis in central Iran. Journal of Vector Ecology. 2005, 30 (1): 115-118.
Guan L-R, Xu Y-X, Li B-S, Dong J: The role of Phlebotomus alexandri Sinton, 1928 in the transmission of kala-azar. Bulletin of the World Health Organization. 1986, 64 (1): 107-112.
Sukkar F, Yawanis F, Al-Doori NA, Al-Mahdawi SK: The possible vector of infantile visceral leishmaniasis in Iraq. Bulletin of Endemic Diseases. 1985, 26 (1/4): 27-36.
Azizi K, Rassi Y, Javadian E, Motazedian MH, Rafizadeh S, Yaghoobi-Ershadi MR, Mohebali M: Phlebotomus (Paraphlebotomus) alexandri: a probable vector of Leishmania infantum in Iran. Annals of Tropical Medicine and Parasitology. 2006, 100 (1): 63-68. 10.1179/136485906X78454.
Lewis DJ: A taxonomic review of the genus Phlebotomus (Diptera: Psychodidae). Bulletin of the British Museum Natural History (Entomology Series). 1982, 45 (2): 121-209.
Cross ER, Newcomb WW, Tucker CJ: Use of weather data and remote sensing to predict the geographic and seasonal distribution of Phlebotomus papatasi in southwest Asia. American Journal of Tropical Medicine and Hygiene. 1996, 54 (5): 530-536.
Belen A, Alten B, Aytekin AM: Altitudinal variation in morphometric and molecular characterics of Phlebotomus papatasi populations. Medical and Veterinary Entomology. 2004, 18: 343-350. 10.1111/j.0269-283X.2004.00514.x.
Depaquit J: Revision du sous-genre Paraphlebotomus (Phlebotomus - Phlebotominae - Psychodidae - Diptera). Approaches morphologique et moleculaire. Doctoral Dissertation. 1997, Riems, France: L'University de Reims
Kamal HA, Doha SA, El-Hosary SS, Shehata MG, El Sawaf BM: Human zoonotic cutaneous leishmaniasis and associated sand flies (Diptera: Psychodidae) in Sheikh Atiya village, southern Sinai, Egypt. Journal of the Egyptian Society of Parasitology. 2003, 33 (3): 795-803.
Maroli M, Krasnonos L, Gafurov I: Epidemiological and entomological survey in a focus of visceral leishmaniasis in Pap district (Fergana Valley) of Namangan region, Uzbekistan. Acta Tropica. 2001, 80 (3): 223-228. 10.1016/S0001-706X(01)00179-6.
Fryauff D, Cope S, Presley S, Hanafi H, Bailly C, Said-Salah E, Abro M, Dabale D: Sand flies of the Republic of Djibouti: ecological distribution, seasonal population trends, and identification of species. Journal of Vector Ecology. 1995, 20 (2): 168-188.
Peterson AT: Uses and requirments of ecological niche models and related distributional models. Biodiversity Informatics. 2006, 3: 59-72.
Peterson AT, Shaw J: Lutzomyia vectors for cutaneous leishmaniasis in Southern Brazil: ecological niche models, predicted geographic distributions, and climate change effects. International Journal for Parasitology. 2003, 33: 919-931. 10.1016/S0020-7519(03)00094-8.
Peterson AT, Pereira RS, Neves V: Using epidemiolical survey data to infer geographic distributions of leishmaniasis vector species. Revista da Sociedade Brasileira de Medicina Tropical. 2004, 37: 10-14. 10.1590/S0037-86822004000100003.
Peterson AT, Martinez-Campos C, Nakazawa Y, Martinez-Meyer E: Time-specific ecological niche modeling predicts spatial dynamics of vector insects and human dengue cases. Transactions of the Royal Society of Tropical Medicine and Hygiene. 2005, 99: 647-655. 10.1016/j.trstmh.2005.02.004.
Levine RS, Peterson AT, Benedict MQ: Geographic and ecological distributions of the Anopheles gambiae complex predicted using a genetic algorithm. American Journal of Tropical Medicine and Hygiene. 2004, 70 (2): 105-109.
Costa J, Peterson AT, Beard CB: Ecological niche modeling and differentiation of populations of Triatoma brasiliensis Neiva, the most important Chagas' disease vector in northeastern Brazil (Hemiptera, Reduviidae, Triatominae). American Journal of Tropical Medicine and Hygiene. 1911, 67 (5): 516-520.
da Costa SM, Cechinel M, Bandiera V, Zannuncio JC, Lainson R, Rangel EF: Lutzomyia (Nyssomyia) whitmani s.l. (Diptera: Psychodidae: Phlebotominae): geographical distribution and the epidemiology of American cutaneous leishmaniasis in Brazil - mini review. Memorias do Instituto Oswaldo Cruz. 2007, 102 (2): 149-153.
Toprak S, Ozer N: Distribution of sand fly (Diptera: Psychodidae) species and efficieny of capturing methods in Sanliurfa Province, Turkey. Journal of Medical Entomology. 2007, 44 (1): 23-28. 10.1603/0022-2585-44.1.23.
Rispail P, Dereure J, Jarry D: Risk zones of human leishmaniasis in the western mediterranean basin. Correlations between vector sand flies, bioclimatology and phytosociology. Memorias do Instituto Oswaldo Cruz. 2002, 97 (4): 477-483.
Wasserberg G, Abramsky Z, Anders G, El-Fari M, Schoenian G, Schnur L, Kotler BP, Kabalo I, Warburg A: The ecology of cutaneous leishmaniasis in Nizzana, Israel: infection patterns in the reservoir host, and epidemiological implications. International Journal for Parasitology. 2002, 32: 133-143. 10.1016/S0020-7519(01)00326-5.
Wasserberg G, Yarom I, Warburg A: Seasonal abundance patterns of the sandfly Phlebotomus papatasi in climatically distinct foci of cutaneous leishmaniasis in Israeli deserts. Medical and Veterinary Entomology. 2003, 17 (4): 452-456. 10.1111/j.1365-2915.2003.00461.x.
El Said S, Beier JC, El Sawaf BM, Doha S, El Kordy E: Sand flies (Diptera: Psychodidae) associated with visceral leishmaniasis in El Agamy, Alexandria Governorate, Egypt II. Field behavior. Journal of Medical Entomology. 1986, 23 (6): 609-615.
Schlein Y, Jacobson RL: Sugar meals and longevity of the sandfly Phlebotomus papatasi in an arid focus of Leishmania major in the Jordan Valley. Medical and Veterinary Entomology. 1999, 13 (1): 65-71. 10.1046/j.1365-2915.1999.00138.x.
Boussaa S, Guernaoui S, Pesson B, Boumezzough A: Seasonal fluctuations of phlebotomine sand fly populations (Diptera: Psychodidae) in the urban area of Marrakech, Morocco. Acta Tropica. 2005, 95 (2): 86-91. 10.1016/j.actatropica.2005.05.002.
Srinivasan R, Panicker KN, Dhanda V: Population dynamics of Phlebotomus papatasi (Diptera:Phlebotomidae) in Pondicherry, India. Acta Tropica. 1993, 54 (2): 125-130. 10.1016/0001-706X(93)90058-J.
WorldClim - Global Climate Data.
Hijmans RJ, Cameron SE, Parra JL, Jones PG, Jarvis A: Very high resolution interpolated climate surfaces for global land areas. International Journal of Climatology. 2005, 25: 1965-1978. 10.1002/joc.1276.
USGS Earth Resources Observation (EROS) and Science Center.
Phillips SJ, Dudik M, Schapire RE: A maximum entropy approach to species distribution modeling. Proceedings of the 21st International Conference on Machine Learning. 2004, ACM Press, New York, 472-486.
Phillips SJ, Anderson RP, Schapire RE: Maximum entropy modeling of species geogaphic distributions. Ecological Modelling. 2006, 190: 231-259. 10.1016/j.ecolmodel.2005.03.026.
Phillips SJ, Dudik M: Modeling species distributions with Maxent: new extensions and a comprehensive evaluation. Ecography. 2008, 31: 161-175. 10.1111/j.0906-7590.2008.5203.x.
Pearson RG, Raxworthy CJ, Nakamura M, Peterson AT: Predicting species distributions from small numbers of occurence records: a test case using cryptic geckos in Madagascar. Journal of Biogeography. 2007, 34: 102-117. 10.1111/j.1365-2699.2006.01594.x.
We would like to thank all of the U.S. military entomologists deployed to Iraq and Afghanistan for their tireless dedication to their mission. We would like to thank Dr. Ed Rowton and Dr. Katherine Swanson for collections performed in Akbuk, Turkey in 2006. Finally, we would like to thank COL Russell Coleman, Dr. Stephen Davies, COL Richard Johnson and Dr. Michael Sardelis for their support over the course of this project.
The authors declare that they have no competing interests.
MGC conceived and executed this study and drafted the manuscript. PMM obtained the remote sensing imagery used in the study and assisted in model development. JPG helped develop the concept of the study and provided resources for its execution. PMM and JPG assisted in the development of the final manuscript.