A fine scale eco-epidemiological study on endemic visceral leishmaniasis in north ethiopian villages

Graphical abstract Hypothetical illustration of Anthroponotic Transmission Cycle of L. donovani in VL Endemic Villages in Northern Ethiopia. Transmission of L. donovani occurs anthroponotically in or close to households near to vertisols. Sand flies breeding in vertisols, are attracted into villages to feed blood on domestic animals, plant communities my serve as resting places for sand fly vectors.


Introduction
Visceral leishmaniasis (VL), also known as Kala-Azar, is a vectorborne parasitic disease caused by protozoan parasites of the Leishmania donovani (Kinetoplastida: Trypanosomatidae) species complex and transmitted by blood-sucking phlebotomine sand flies (Diptera, Psychodidae). It is a systemic disease of the bone marrow, spleen and liver that primarily affects cells of the reticuloendothelial system (Alvar et al., 2012). VL is fatal in over 90% of the cases but timely medical intervention can boost survival rates by 80% (Guerin et al., 2002). Although molecular and serological methods can be indicative of disease, definitive diagnosis is complex and requires invasive medical procedures to obtain lymph node, spleen, or bone marrow aspirates for microscopical and/or molecular detection of Leishmania parasites.
VL is endemic in 98 countries with some 200 million people at risk, 200-500 thousand new cases and 20-40 thousand fatalities annually (Alvar et al., 2012). Over 90% of the cases occur in six countries: India, Bangladesh, Nepal, Sudan, Ethiopia and Brazil (Chappuis et al., 2007).
The species causing VL in the Indian subcontinent and East Africa is L. donovani donovani while L.d. infantum is the causative agent of VL in South America (Chappuis et al., 2007). In East Africa the worst affected countries are Sudan with an estimated 15-20 thousand cases and Ethiopia with 4,500-6,000 cases per year (Herrero et al., 2009) (Alvar et al., 2012). The worst affected regions in Ethiopia are the Humera-Metema lowlands bordering Sudan, which account for approximately 60% of the country's reported VL cases (Leta et al., 2014). Emergence of VL in this region coincided with large-scale deforestation carried out to create vast, intensively cultivated farms producing mostly sesame, cotton and sorghum. The massive scale of these industrial agriculture enterprises necessitates a twice-yearly influx [and efflux] of some 500,000 seasonal workers among which, there is a high prevalence of HIV/VL coinfections (Alvar et al., 2012;Berhe et al., 2001;ter Horst et al., 2008).
Nineteen Phlebotomus species, representing six subgenera, have been documented in Ethiopia (Lewis, 1987). Ph. martini and Ph. celiae transmit VL in the more humid southern foci near the Kenyan border (Gebre-Michael and Lane, 1996) while Ph. orientalis is the predominant vector of VL in Northern Ethiopia and Sudan (Elnaiem, 2011;Gebre-Michael et al., 2004;Gebresilassie et al., 2015c).
In 2006, a new focus of VL emerged in villages close to the town of Sheraro (District of Tahtay Adiyabo, northern Ethiopia, Fig. 1), probably as a result of infected seasonal laborers returning from the VL endemic areas of Humera and Metema (Herrero et al., 2009). Between 2006 and 2011, 209 VL cases and 3 PKDL were diagnosed and treated.
We initiated a multidisciplinary investigation of the Sheraro focus in 2009 performing both epidemiological and entomological studies (Abbasi et al., 2013;Gebresilassie et al., 2015c). The current study comprises part of this project. We hypothesized that transmission L. donovani depended on the population densities of Ph. orientalis vectors and the L. donovani infection rates among these flies. We further reasoned that, these parameters should be influenced by (1) The proximity of sand fly breeding and resting habitats (2) The presence of ample blood sources for sand fly females including infectious hosts. To address these questions we conducted a multilevel fine-scale eco-epidemiological investigation of entomological, epidemiological and ecological factors associated with the transmission of L. donovani to humans. Our analyses were founded in the disciplines of landscapeecology and landscape-epidemiology that provide a conceptual framework for the study of vector-borne disease epidemiology (Emmanuel et al., 2011;Kitron, 1998;Reisen, 2010).

Study area
Tahtay Adiyabo is one of 36 districts (woredas) in the Tigray Region of northern Ethiopia. The district is bordered to the south and southwest by Asigede Tsimbela and the Tekezé River, which separates Tahtay Adiyabo from Wolqayt and Setit Humera. On the north by Eritrea, and on the east by La'ilay Adiyabo. The administrative center of Tahtay Adiyabo is the town of Sheraro (14°23′41"N/37°46′15"E, 1246 m above sea level [masl]). The population of Tahtay Adiyabo is 90,144, with only 6377 (7.07%) living in towns. With an area of 3,841.51 Km 2 , Tahtay Adiyabo has a population density of 23.47 people/Km 2 , (2007 national census, Central Statistical Agency of Ethiopia) (Fig. 1).

Selection of the study villages
In 2010 we conducted a census in the villages of Tahtay Adiyabo comprising data on 11,000 participants (personal data, socio-economic parameters, health status and clinical history of VL). Based on this census, 18 villages were selected for the cohort studies (Abbasi et al.,Fig. 1. The Tahtay Adiabo District comprises the northern most tip of Ethiopia adjacent to the border with Eritrea (red frame in Inset). The district was selected due to high incidence of VL in rural villages close to Sheraro. Pie charts represent the total population at baseline (2011) and the fractional case incidence (2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) in the three villages. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) 2013). The three villages included in the current study were selected based on data from the above census and cohort studies, including the percentage of VL cases; two of the villages (Salamo and Erdwayane) had high incidence of VL and one (Dedevit) had very low incidence (Fig. 2). Baseline parameters such as sand fly species richness and abundance, environmental and geographical variables, accessibility to motorized transportation as well as other logistic hurdles, were also taken in consideration. The three villages cover an area of 7.80 km 2 combined and comprise 367 households, with a total population of 1540 inhabitants. Average household size is 4-5 persons and the average age of the houses is 7.5 years for Erdwayane, 12.3 years for Salamo, and 16.8 years for Dedevit. The population density in Erdwayane was 414 persons/km 2 , in Salamo 99/km 2 , and in Dedevit it was 272 persons/km 2 .
The villages are built on sand-stone hills and outcrops surrounded by vertisol fields, which are used mainly for agricultural cultivation and pasture. The landscape is distinct in the three villages. Salamo and Erdwayane are more rural with open spaces and traditional houses whereas Dedevit is more urbanized with asphalt roads and houses connected to the electric grid. The hills upon which the villages are built are largely deforested due to woodcutting for construction and cooking. The remaining vegetation comprises scattered copses of Acacia sp., Balanites sp., and Zyziphus sp. trees/bushes. The introduced Neem tree (Azadirachta indica) is also abundant around homesteads.

Ethical considerations
Informed consent was obtained from all the adult volunteers who participated in the study. Consent for inclusion of young children, was obtained from parents or guardians. The study protocols were reviewed and approved by the ethical review committee at the Medical Faculty, Addis Ababa University and the National Research Ethics Review Committee (NRERC) at the Ethiopian Ministry of Science and Technology. Each volunteer was treated in accordance with international bio-ethical rules and laws for human sampling. Blood was sampled by certified trained Health officers or nurses. Volunteers requiring medical assistance were treated by project MDs or transported to regional hospitals.

Research design
A cohort study was conducted in the three villages. Volunteers were recruited based on the historical roster of all those previously screened. Former VL Cases were defined as clinically confirmed cases. During the study itself, persons presenting with fever lasting 2 weeks or longer, weight loss and/or splenomegaly and diagnosed either parasitologically or serologically (rk39 dipstick) were considered current VL cases.
To record the data we used an interactive tablet-based poll application QuickTapSurvey (http://www.quicktapsurvey.com Toronto, Canada, 2013). Cases and controls responded to structured interviews summarizing demographic and socio economic characteristics such as age, gender, occupation, education, sleeping habits, bed net use, family size, period living in the current location/house, house construction (materials and condition), ownership of domestic animals (species and number) and farming activities. Respondents were interviewed by trained health professionals in their indigenous languages (Amarignia, Tigrigna and/or Kunamigna).

Blood sampling
Asymptomatic infections were diagnosed by PCR that detected Leishmania DNA in blood as the best indication of infection (Singh et al., 2014). Blood samples were obtained by (i) Finger pricks using disposable blood lancets. (ii) Venous blood was collected using the BD Vacutainer collection system with serum-separating tubes (Gold cap BD Vacutainer ® tubes SST ™ with BD Hemogard™ closure). To identify parasite DNA, four drops of finger prick or venous blood were blotted on Whatman No1 filter papers and kept dry. For serology, tubes were centrifuged at 13,000 rpm for 15 min, then, sera were pulled, aliquoted and stored at −80°C.

DNA extraction and PCR amplification
DNA was extracted with phenol-chloroform from two paper punch disks (r = 3 mm), of Whatman No1 filter paper calculated to have been saturated with approximately 20 μl of blood each, and tested for L. donovani infection via quantitative kDNA real-time PCR (qRT-kDNA/ PCR) as previously described (Abbasi et al., 2013). Samples that tested positive were assessed via conventional PCR to amplify and sequence a Leishmania species-specific segment of the internal transcribed spacer 1 (ITS1) (El Tai et al., 2000).

Enzyme-Linked-Immune sorbent assay (ELISA) for anti-Ph. orientalis saliva antibodies
To estimate the exposure of volunteers to Ph. orientalis bites, antisaliva IgG antibodies in volunteers' sera were measured via ELISA. Microtiter plate wells (Nunc-Immuno™ MicroWell™ 96, Thermo Scientific, Denmark) were coated with Ph. orientalis salivary gland homogenate (SGH, corresponding to 0.2 gland/well, prepared as previously described (Vlkova et al., 2011) and diluted in coating buffer (Carbonate-Bicarbonate, pH 9) and incubated overnight at 4°C. The wells were then washed with 100 μl of PBS-T and incubated with 100 μl of blocking/diluting solution (2% [w/v] low fat dry milk in PBS-T) for 60 min at 37°C to block free binding sites. Thereafter, sera were diluted 1:50 and 1:100 in diluting solution/coating buffer and incubated in duplicate for 90 min at 37°C. Secondary antibody, α-human IgG HRP (Sigma) diluted 1:1000 was added in a total of 100 μl and incubated for 45 min at 37°C. Lastly the unbound antibodies were washed six times using PBS-T three minutes per wash. 100 μl of orthophenylenediamine and H 2 O 2 in McIlwain phosphate-citrate buffer (pH 5.5) were used as a substrate solution. After 5 min incubation in the dark the reaction was stopped by adding of 100 μl of 10% H 2 SO 4 . Absorbance was measured with a Multiskan RC ELISA reader (Labsystems) at 492 nm wavelength.
2.9. Entomology 2.9.1. Trapping methods and sand fly surveillance Sand fly traps were deployed in trapping units (TUs), comprising one CDC light trap (6 V, 150 mA, John W. Hock, Model 512, Gainesville, FL) and two large Sticky Traps (ST) made of white polypropylene boards, measuring 60 × 60 cm, placed horizontally on a square metal frame supporting them approximately 15 cm above ground, with the upper face of the board smeared with sesame oil (Moncaz et al., 2013). The three traps comprising the TU were positioned in a triangular arrangement, 2-3 m apart and deployed either inside the domestic compound or in peridomestic areas. Sand flies were trapped during dry season (March to April 2013). 13 TU were deployed in the village of Salamo, 14 in Erdwayane and 12 in Dedevit (Table 2). TUs operated between 6:00 PM to 6:00 AM the following morning. Sand flies collected by the three traps in each TU were pooled for statistical analyses.
Sand flies were preserved in 100% EtOH and transported to the laboratory for taxonomical identification. All specimens belonging to the genus Phlebotomus and about 5-10% of the Sergentomyia spp. were selected for species identification under a stereo microscope. To identify species, the heads of males and females were removed to expose the internal structure of the pharynx and cibarium. In females, the end of the abdomen was dissected as well, to expose the spermatheca. Flies were mounted in Hoyer's medium on microscope slides. Identification of species was based on cibarial, pharyngeal armature, antennal features and spermathecae of females and external genitalia of males, using different keys (Artemiev, 1980;Kirk and Lewis, 1946;Lewis, 1982;Perfil'ev et al., 1968;Quate, 1964).

Nocturnal activity rhythms
A Collection Bottle Rotator (CBR -John W. Hock. Model 1512) was used to separate the nightly catch into 8 different bottles corresponding to a predetermined period using a programmable timer. The timer was set to rotate every 1.5 h, from 18:00 in the evening until 7:30 the following morning. CBRs were deployed in two different ecotopes (A-Farm field and B -Peridomestic), for 6 consecutive nights, during April 2013.

DNA extraction from sand flies
Prior to dissection, sand flies were placed in a tea strainer, rinsed under tap water, briefly immersed in dilute dish-washing detergent, rinsed again, immersed in 1% sodium hypochlorite (NaClO) solution for 1 min and rinsed in double distilled water to remove all external contamination. DNA was extracted from individual sand flies using phenolchloroform as previously described (Abbasi et al., 2009)

Leishmania infection in sand flies
Leishmania infection rates were determined by qRT-kDNA/PCR as described above for blood samples (Abbasi et al., 2013). Sand flies were pooled in pools of 6 specimens per reaction. Positive samples were reamplified by ITS1 PCR and amplicons were sent for DNA sequencing to identify the Leishmania species (El Tai et al., 2000). Validation of parasite species (i.e. L. donovani) was based on ITS1 sharing > 99% sequence homology with existing L. donovani entries in GenBank.

Determination of sand fly blood sources by Reverse Line Blotting (RLB)
Blood hosts were identified as previously described (Abbasi et al., 2009). Briefly, DNA extracted from abdomens of blood fed Ph. orientalis females (and a small number of other Phlebotomus spp.) was PCR amplified using biotinylated vertebrate species-specific cytochrome b primers. The PCR products were blotted on species-specific oligonucleotides. Specific hybridization was visualized using streptavidin horseradish peroxidase.

DNA sequence analysis
PCR amplicons (from human blood samples and sand flies) were sequenced by The Center for Genomic Technologies at the Hebrew University of Jerusalem using automated DNA Sequencing, based on BigDye Terminator cycle sequencing chemistry from Applied Biosystems (ABI), ABI PRISM 3730xl DNA Analyzer and the ABI's Data collection and Sequence Analysis software. The derived sequences were compared against the GenBank database using NCBI BLAST (https:// blast.ncbi.nlm.nih.gov/Blast.cgi) and default search parameters (Madden, 2013).

Meteorological data
Data loggers with sensors (HOBO Mini-Weather Station, Onset Computer Corporation, Pocasset, MA, and U.S.A) were used to record temperatures (°C) and relative humidity (RH%) during experiments.

Analysis of entomological data
Abundance of sand flies including Ph. orientalis was calculated based on the number of flies captured per trapping unit (see above item 9.1) per night. Data was tested for its distribution and variance by 1-Sample Kolmogorov -Smirnov Z test (K-S normality test) and Shapiro-Wilk test. The nonparametric tests Kruskal-Wallis (K-W) and Mann-Whitney-U (M-W) were applied to analyze non-transformed data. K-W test was employed to compare the median number of sand flies in general or Ph. orientalis specifically collected in the different villages. When K-W was significant Dunn's test for pairwise multiple comparisons was followed. P value was adjusted with Bonferroni correction.
Nocturnal activity, One-way ANOVA was used to compare the hourly activity patterns of male and female Ph. orientalis during the night. Tukey's Studentized test post hoc analysis was utilized for mean separation where ANOVA was significant. For the same analysis, Spearman's rank-correlation analysis (P < 0.05) was also used to compare the effects of average nighttime temperature and humidity on the number of flies captured per hour.
For each categorical variable, the distribution for frequencies across villages was analyzed using a chi-square or Fisher's exact test depending on whether counts were more or less than five, respectively. Pairwise analyses using adjustment by the Benjamini-Hochberg method (i.e., false discovery rate) were conducted for those variables with significant results in the universal test. Continuous variables were described in terms of means and standard deviations.
For the logistic regression models, PCR positivity was defined as ≥1 parasites per ml blood at baseline measurement and clinical VL was defined as developing disease at any time during the study. Unadjusted analyses evaluating the associations between outcome variables and individual-households, were conducted for each village separately using logistic regression.
Multilevel analysis was conducted for data from the entire study sample in aggregate, using a generalized linear mixed model fit by the Laplace approximation.

Spatial analysis
Integration of Geographic information systems (GIS) and Remote Sensing (RS) was used to allow the digitalization and projection of epidemiological and entomological data and permit mapping and evaluation of environmental factors related to the spatial and temporal distribution of the VL cases as related to vector distribution. Spatial data included: the location of houses, land use and land cover as well as the position of the TUs. Data were recorded as waypoints using a GPS device (Garmin Dakota 20), with an accuracy-range of 4-5 m. GPS data were downloaded, projected, codified and rectified into ExpertGPS Map software 5.20 for GARMIN devices (TopoGrafix).
Local distribution maps that depict "hot-spots", areas characterized by high abundance of sand flies and people seropositive for sand fly saliva, were generated using the Inverse Distance Weighting (IDW) interpolation method. Numbers of Ph. orientalis per TU per night and the ODs (optical densities) values of anti-Ph. orientalis saliva obtained by ELISAs were represented on these maps.
Linear regression was applied to evaluate association between the estimated number of Ph. orientalis per TU (as the dependent variable) and the estimated ODs anti-Ph. orientalis saliva ELISAs (as the independent variable). This analysis was done for populated areas only. We represented peridomestic areas as a circular buffer zone around each household (r = 100 m). Furthermore, masks were issued from the results of the produced surface maps to obtain two continuous data sets. The association between the distribution of sand flies and the exposure of people to their saliva and the distance from vertisols was calculated as the shortest Euclidean distance from data collection point (TU or house of exposed volunteer) to the vertisols area. A linear regression analysis was performed to test for association between the number of Ph. orientalis/per TU (as the dependent variable) and the distance from vertisols (as the independent variable). In case of non-linear relationship between the variables, the variable distance was Ln transformed to enable linearization of the variables in the model. For GIS spatial analysis vertisol layers were demarcated by manually generated polygons.

Remote sensing
To study the areas covered by vertisols and their dynamic seasonal changes, we used remotely sensed satellite images of the Sheraro area (NW corner: 15°26′40.20"N, 36°44′32.35"E; SW corner: 38°28′28.42"E, 13°27′34.52"N). Satellite images spanning September 29, 2013 (rainy season, Fig. 7A), January 16, 2014 (early dry season, Fig. 7B), and April 19, 2014 (late dry season, Fig. 7.C) were acquired from Landsat 8 −OLI dataset, Earth Explorer (USGS, http://earthexplorer.usgs.gov/). Supervised classification method was used to cluster pixels in the digital scene dataset into desired spectral classes. Regions of interest (ROI) were generated as defined training classes in the output images corresponding to vertisol areas. Afterwards, maximum likelihood classification was performed to assign to each pixel in the subset image data the trained class of which it had the highest probability of being a member, in this case of being vertisols (Exelis Visual Information Solutions, 2014; Sithiprasasna et al., 2005). The analyses employed data processing based on the following bands: Band 2 (Blue: 0.452-0.512 μm), Band 3 (Green: 0.533-0.590 μm), Band 4 (Red: 0.636-0.673 μm) with 30 m spectral resolution.
Among those individuals that retrospectively reported clinical VL (contracted before 2011) for whom PCR data were available, 6/49 (12.2%) remained PCR positive during the initial sampling in 2011 (3 patients with 1-10 and 3 with > 10 parasites per ml blood). Similarly, 4/48 (8.3%) individuals who developed clinical disease during the study period remained PCR positive 2-3 years later including the final sampling in 2015. Remarkably three VL cases reported before 2011, presented low parasitemia (≥1 parasites per ml blood), and remained PCR positive until the final sampling in 2015. The overall ratio of symptomatic to asymptomatic infected individuals (at any parasitemia and at baseline) across the three villages was calculated as 1:9.
A limited number of volunteers who participated in the past cohort studies was sampled again in 2015 in two of the three villages, Salamo and Dedevit. Blood samples were obtained from 126 volunteers via finger-pricks for PCR. Serum was separated from venous blood for serology. Of the total number of individuals sampled, 44 were past VL cases recorded as of 2007 (VL treated cases), and 82 individuals never had VL (healthy negative) (Table 1 Supp. data).
Anti-L. donovani antibodies were detected in 85 (67.4%) individuals, 40 (31.7%) of whom were past VL cases and 45 (35.7%) healthy individuals without a history of VL. L. donovani DNA was detected in 35(27.7%) volunteers, 20 (15.8%) of which were past VL cases, and 10 (7.9%) healthy negative individuals. Thirty of the volunteers were concurrently positive by ELISA (anti-L. donovani Abs) and PCR (Leishmania DNA). Of these 20 were past VL cases (11 of them had > 10 parasite per ml) and 10 were healthy (of those 5 volunteers with parasitemia of > 10 parasite per ml and 2 presented parasitemia of > 100 parasite per ml). Some 20 (15.8%) past VL cases and 35 (27.7%) healthy volunteers were seropositive by ELISA but negative by PCR (Table 1 Supp. data). One past VL case (0.7%) was seronegative but PCR positive, and 3 (2.3%) were negative for both parameters. Likewise, 4 (3.1%) healthy individuals were positive for k-DNA but seronegative and 33 (26.1%) were negative for both parameters (Table 1 Supp. data). Past VL cases from Dedevit had relatively high titters of anti-L. donovani antibodies (ELISA down to 1:1600 serum dilution). We assume this was the case because all VL cases from Dedevit were relatively recent (since beginning of 2010). In Salamo, one individual that was treated during 2008, still presented with high titters of anti-L. donovani Abs during of the healthy seropositive individuals remained ELISA positive at dilutions downto1:1600.  Most of the asymptomatically infected (qRT-kDNA/PCR positive) population scored 1-10 parasites per ml blood (17.5%), and 6.8% had > 10. A significantly higher number of PCR positive volunteers was found in the village of Salamo (χ 2 test = P < 0.001).

Sand fly bionomics
Sand fly collections were conducted from early February to late May in 2013/14. Seven Phlebotomus spp. and 10 Sergentomyia spp. were identified (Table 2 Supplementary data). The overall population densities of all species of sand flies did not differ significantly between the three villages (Kruskal Wallis test χ2 (2) = 3.027, P = .220). The median TU catch was 158 (111.5, 313.7)/TU/night. The vast majority of sand flies (87.3%) were Sergentomyia spp. All 2483 Phlebotomus spp. were identified showing that Ph. orientalis was by far the predominant species in the genus (93.9%, Table 2 and Table 2 Supplementary data). The male/female ratio of Ph. orientalis was consistently and significantly male biased (males/females (1444/887 = 1.62; Mann Whitney U test, U = 4,192.5, P = .465; Table 2).
Significantly different numbers of Ph. orientalis were recorded in the three villages (Kruskal Wallis test χ2 (2) = 6.417, P = .040) with high numbers in Erdwayane, intermediate in Salamo and the lowest numbers in Dedevit (Table 2)

L. donovani infection rates in Phlebotomus orientalis
Ph. orientalis females (978) were screened for Leishmania infection by qRT-kDNA/PCR. Parasite species was confirmed as L. donovani by sequencing of the ITS1 amplicons in 57 samples. No other Leishmania spp. were identified in any of the kDNA positive Ph. orientalis. The highest rate of L. donovani infections was found in Ph. orientalis from Salamo (23%), followed by Dedevit and Erdwayane with 3 and 2.7% respectively (Table 2).

Host feeding preferences of Phlebotomus orientalis
The source of blood meals ingested by 150 female Ph. orientalis (50 samples from each village) was identified using cyt b PCR followed by vertebrate cyt b specific reverse line blotting (RLB). A high proportion Ph. orientalis was found to have fed only on bovine (79, 52%) or human (22, 14.7%) blood. Forty-nine (32.7%) of the blood-fed flies had mixed blood meals. Of these, the majority were bovine/human (34, 22.7%). Other blood sources were also detected, however always as mixed blood meals (e.g. camelid/bovine avian/human, canine/bovine.

Household characteristics
Across the three study villages, individuals from 324/367 households were included in the analysis. Most of the compounds had a shelter outside the house (Dase) and/or an animal shed (Dambe) (n = 241, 75.3%; n = 201, 62.8%, respectively), with no significant variation across the three villages. Households in Dedevit were less likely to have any animal than those in Erdwayane and Salamo (adjusted p-values: 0.028 for both locations), and more likely to have tin roofs (adjusted p-values: < 0.001 for both locations). Fewer households in Salamo had stone fences compared with households in Erdwayane and Dedevit (adjusted p-values: 0.002 and < 0.001, respectively) ( Table 3).

Unadjusted multi-level associations with asymptomatic Leishmania infection
Bivariate analyses did not reveal a strong association between the odds of infection with Leishmania and any, household, or zone-level factor across all three villages (Table 3, supplementary data). Residing in a household with an external shelter (dase) was associated with a lower L. donovani infection rates in Dedevit (OR: 0.48, 95% CI: 0.25, 0.89), while having any animal led to an 86% increase in odds of infection among residents of Erdwayane (OR 1.86, 95% CI: 1.00, 3.44). In Salamo, the odds of L.donovani infections were significantly reduced among individuals in households with animal sheds (Dambe) outside the compound, as compared to those in households without animal sheds (OR 0.51, 95% CI: 0.26, 0.98) (

Unadjusted multi-level associations with clinical disease
In Erdwayane and Salamo, the odds of developing clinical VL was significantly higher among males than among females (Erdwayane: OR: 2.04, 95% CI: 1.04, 4.02; Salamo: OR: 3.04, 95% CI: 1.57, 5.87) ( Table 4, supplementary data). In Salamo (but not in Erdwayane or Dedevit) individuals residing in households with a stone fence were twice as likely to develop clinical disease as individuals residing in households without a stone fence (OR: 2.00, 95% CI: 1.09, 3.68).

Adjusted multi-level associations
After adjustment for household-level (dase, stone fence, domestic animals) and eco-zone-level factors (vertisol and Acacia trees), significant differences between males and females were calculated. Males had 34% reduced odds of being infected asymptomatically ( To explore the relationship between compositional variables (individualor household level) and contextual variables (eco-zones-levels) versus number of Ph. orientalis, we included houses for which we had data on sand flies that were inhabited by individuals for whom clinical diagnostic data was available (32/39). Results showed that nonsignificant association (Table 5, Supplementary data). Nevertheless, the presence of urban elements and the using of better building quality materials, like corrugated iron roof and walls made of bricks and cement, tend to attract less sand flies, yet, no significant difference was found (K-W test, P = .843; Table 5, Supplementary data).

Ph. orientalis population densities
To assess the densities of Ph. orientalis in different areas, trap yields were interpolated using Inverse Distance Weighting (IDW) and a surface interpolated map was constructed for each village. Areas with predicted high densities of Ph. orientalis were readily discernable ( Fig. 5A-C, dark blue shading). The sand fly hotspots were almost exclusively congruent with vertisol fields located around the villages, or as restricted zones inside the village perimeters (Fig. 5 hatched areas) with scattered trees and shrubs (Balanites sp., Acacia sp. and/or Zizyphus sp.).

Anti-sand fly saliva antibodies in volunteers' blood
We also interpolated the anti-Ph. orientalis serological data (OD levels determined by ELISA) in serum samples obtained from volunteers in Salamo and Dedevit, the two villages with low numbers of Ph. orientalis. Using the IDW method we produced surface interpolated maps representing areas of intense exposure to sand fly bites (Fig. 5A.1 &  B.1). The hotspots of anti-saliva sero-positivity largely overlapped the  areas with dense Ph. orientalis populations close to vertisols. Linear regression predicted the relationship between exposure to sand fly bites (ODs, pink/red circles) and the abundance of Ph. orientalis sand flies (blue shading). A positive significant regression equation was found (Salamo: F (1,111) = 34.535, p < 0.001, with an R 2 of 0.237; Dedevit: F (1,116) = 41.678, p < 0.001, with an R 2 of 0.264)

Distribution of VL cases and asymptomatic infections
We produced village maps showing the spatial distribution of households, VL cases (from 2007 to 2015) and asymptomatic PCR positive, individuals. These data were plotted on the maps representing the density of Ph. orientalis. Geographical; location and number of L. donovani infected Ph. orientalis were also marked on these maps (Fig. 6).
The maps show that VL cases and asymptomatic individuals appear not to be uniformly distributed, and houses with two or more VL cases and/or asymptomatic infected individuals, were more likely to be found relatively close to other households with infected individuals. Moreover, houses of patients and asymptomatic individuals were frequently located in areas with high abundance of Ph. orientalis sand flies. Sand flies infected with L. donovani were most frequently found near houses with infected individuals located a short distance from vertisols (Fig. 6).
In view of these observations it was somewhat surprising that, the prevalence of VL and asymptomatic infections was not significantly associated with the distance form vertisols, the Ph. orientalis population densities or the exposure of persons to their bites (data not shown).

Remote sensing analysis
LANDSAT 8 satellite pictures clearly show the vast expanses of vertisols within which our study sites were located (Fig. 7). Cracked vertisols become exposed during the dry season − black areas in truecolor maps (Fig. 7A, B, and C). In the rainy season, wet soil swell, cracks seal-up and soil becomes partially concealed by annual vegetation (Fig. 7A1). The vegetation cover is easily discerned as red areas in infrared maps (Figs. 7A1, B1, and C1). To emphasize the phenomenon, vertisols were classified in yellow showing minimal areas in the rainy season (Fig. 7A2). Vast areas of vertisols are exposed as the dry season progresses (Figs. 7B2 and C2) because harvesting of agricultural crops and the drying up of seasonal vegetation, promote the soil's visibility in satellite photography.

Discussion
Numerous epidemiological studies have been conducted aiming to identify risk factors that influence the transmission dynamics of VL. However, despite the significant amount of accumulated data, our understanding of the effects of anthropogenic and environmental factors   2) in different seasons. Photo A was acquired during the rainy season when the vegetation cover was dense, largely concealing the vertisols surface. Near infrared bands allows to visualize the vegetation in Red color, which accentuates the dense vegetation cover. Photos B were taken during the early dry season when seasonal vegetation had dried up but was still largely concealing the vertisol surface. Photos C were obtained during the late dry season when the vegetation coverage was relatively sparse revealing the dark vertisols. The three villages (1-Erdwayane, 2-Salamo, and 3-Dedevit) are clearly located within vertisols areas while the town of Sheraro (4) is located in sand stone area (prominently visible in 7A). In the classified images (pictures: A.2, B.2, C.2), vertisols were highlighted in yellow. In the early dry season (B.2), vertisols areas were still partially concealed by vegetation (natural and agricultural). In the late dry season (C.2), cracking vertisol areas become exposed, expanding as the vegetation cover receded posty harvest. Vertisols exposed to the sun, dry up and crack due to accelerated evaporation and subsequent shrinkage. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) well.
In 2011, we screened 4,900 volunteers from 11 villages for VL, asymptomatic infections and seropositivity to L. donovani. We found 14.3% (680/4,757) of the volunteers positive for Leishmania k-DNA and varying incidence of VL during the five preceding years (Table 1, Fig. 2) (Abbasi et al., 2013). Then we selected three of these villages, for further studies. Erdwayane and Salamo, with a relatively high number of VL cases and Dedevit with much fewer (Table 1, Fig. 2). Interestingly, the proportion of asymptomatic infections (qRT-kDNA/PCR) was markedly higher in Salamo compared with Erdwayane and Dedevit (Table 1).
Many infectious diseases affect males more than females. Such biases could potentially arise from physiological differences (e.g. sex hormones interacting with immune effectors) or from gender-related differences in exposure (Guerra-Silveira and Abad-Franch, 2013). Evidence in support of the former hypothesis was documented in hamsters infected with Leishmania (Travi et al., 2002). In humans, male-biased incidence of VL was described in studies in southern Ethiopia and Sudan (Ali and Ashford, 1994;Nackers et al., 2015). In our study too, a higher proportion of males developed VL (males: 69.2%, females: 30.7%), although more females were asymptomatically infected (10.9% of the males and: 13.3% of the females tested). Thus, our study design enabled us to determine that gender-related physiological traits and not behavioral differences between males and females were probably responsible for gender-specific differences in the incidence of VL.
On the other hand, behavioral parameters may have played a role in determining the distinct rates of infection and disease in the three villages (Table 1, Fig. 2). In Salamo and Erdwayane where VL incidence was high, most men are farmers working in the sorghum or sesame fields surrounding their villages and spending the nights in the village. In contrast, most men from Dedevit serve as soldiers in military bases located in the area, where they spend many nights on active duty outside the endemic village. These differences in behavioral profiles may explain the lowered incidence of VL in Dedevit.
The three villages are built on sand-stone outcrops surrounded by vertisols, long suspected as breeding and resting sites of sand flies in Sudan (Elnaiem et al., 1998;Elnaiem et al., 1999) and recently shown as such in Tahtey Adiabo (Moncaz et al., 2014). Despite the closeness of all three villages to vertisols, the densities of Ph. orientalis differed significantly with high densities in Salamo contrasting with significantly lower densities in Dedevit (Table 2). We attribute the lower densities in Dedevit to its more urban character with fewer domestic animals contrasting with the more rural landscape of Salamo and Erdwayane. Dedevit has asphalt roads and houses built closer together of commercial materials. We found a trend showing that lower numbers of Ph. orientalis were associated with houses built with commercial building materials. Nevertheless, we could not find a positive association between compositional variables and the density of Ph. orientalis, most probably because of the small number of traps used (Table 4).
Multi-level associations between compositional (individual orhousehold level) and contextual variables or factors (eco-zone level) that may potentially influence the transmission of VL did not reveal any significant relationships (Tables 2 & 3, Supplementary data). The highest rate of L. donovani infections was found in Ph. orientalis captured in [rural] Salamo and Erdwayane, while flies from [urban] Dedevit were seldom infected (Table 2). These observations may partially explain the relatively high rates of asymptomatic infections and VL among the residents of Salamo and Erdwayane compared with Dedevit (Figs. 2, 3).
We conducted restricted sampling of individuals who participated in the previous cohorts in the village of Salamo and Dedevit, where we found that, almost all of VL treated cases reported before and during the beginning of main cohort study (2011) remained seropisive by ELISA prospectively until 2015, and half of them were positive for Leishmania by kDNA-qRT/PCR.
In line with analysis obtained from the general cohort, performed in the 18 villages in the area, where individuals who were direct agglutination test (DAT) positive at baseline and re-sampled after a year, 29.1% (25/86) remained seropositive (Skrip unpublished). We found that a large number of healthy volunteers were also found to be seropositive by ELISA, demonstrating that there was a high proportion of the population, which had been exposed to L. donovani infections. Moreover, higher titers of anti-L. donovani antibodies were reported in the village of Salamo whereas in Dedevit titers were relatively low. We assume this was the case because all VL cases from Dedevit were relatively recent (since early 2010). A high proportion of the Ph. orientalis females collected in inhabited compounds was found to have fed on bovine (79, 52%) or human (22, 14.7%) blood. The preponderance of bovine blood meals is undoubtedly attributable to the abundance of cows in peridomestic areas compared with other hosts (Gebresilassie et al., 2015a;Gebresilassie et al., 2015c). The abundance of cows close to houses may function as a zooprophylactic barrier, potentially reducing human-vector contact, or it may aggravate the risk of VL infection by attracting flies to houses and supplying them with ample protein-based nutrition required for egg production. Studies in Nepal showed that ownership or proximity of livestock was associated with significant protection from VL, whereas in India risk of VL increased for those living in close proximity to cattle (Bern et al., 2005;Singh et al., 2010). Our findings indicate that having any animal (cow, goat, chicken, camel, donkey, sheep, dog, and cat) increases the probability of being PCR positive for L. donovani in the village of Erdwayane (Table 2A, Supplementary data).
Ph. orientalis is the proven vector of L. donovani in the vast endemic foci of Sudan (Ashford et al., 1992;Elnaiem et al., 1998). The predominance of Ph. orientalis was reported in north-west Ethiopia too (Gebre-Michael et al., 2010;Gebresilassie et al., 2015c;Hailu et al., 1995). Our entomological studies demonstrated that Ph. orientalis was the most prevalent Phlebotomus spp., accounting for more than 93% overall and a relatively large proportion of the females were infected with L. donovani (Table 2) implicating this species as vector of L. donovani.
Previous studies in Northern Ethiopia have identified open vertisol fields as putative breeding grounds for sand flies, where larvae live in the deep vertisol cracks feeding on composting organic material (Moncaz et al., 2014). Here we identify cultivated vertisols as important sources of sand flies foraging in villages during the night (Fig. 5). GISbased analysis allowed us to estimate the dispersion of Ph. orientalis from vertisols to households. Linear regression demonstrated a strong negative association between sand fly abundance and the distance from vertisols. Interpolated surface maps produced by converting pointbased data (number of Ph. orientalis per TU) to surface data, allowed us to estimate the abundance of Ph. orientalis in different areas within each village (Fig. 5). In support of this assertion, we also describe a marked periodicity in the activity patterns of Ph. orientalis, with a peak between 21:00-22:30 in open field and a delayed peak of activity 22:30-00:00, in peridomestic habitat. These finding were also reported in previous studies conducted in Northern Ethiopia (Gebresilassie et al., 2015b;Lemma et al., 2014). We interpret this delay to imply that sand flies foraging in the villages had arrived from the adjacent vertisol fields (Fig. 4).
Furthermore, in additional GIS-based analyses, we measured the levels of sero-positivity (optical density [OD] by ELISA) to-Ph. orientalis saliva as an indicator of exposure to biting sand flies. Careful examination of the resultant surface maps indicates higher levels of α-Ph. orientalis saliva antibodies in volunteers living close to vertisols where sand fly populations are more dense (Fig. 5).
Our GIS-based maps, indicate that VL cases were frequently clustered, with more than one case per house or groups of houses (Fig. 6). In addition, houses of VL patients as well as asymptomatic infected persons were frequently located relatively close to vertisols where sand fly populations were denser (Fig. 6). Hence, proximity to vertisols comprises a risk factor for contracting VL. Taken together, these findings support the hypothesis that, transmission of L. donovani probably occurs anthroponotically (humans serve as parasite reservoir) in or close to households rather than in the fields. In conclusion, sand flies breeding in vertisols, are attracted into villages to feed blood on domestic animals and humans. They are infected with L. donovani as they feed on infected humans (VL patients and asymptomatic cases) and transmit the parasites to naïve individuals living in proximity during subsequent blood meals. Thus, living close to vertisols comprises the single most important risk factor for contracting infections with L. donovani, in foci where the parasite is transmitted by Ph. orientalis.