The field surveys were carried out as a part of the Zambian Ministry of Health (MoH) Lymphatic Filariasis Control Programme (2003–2005) and Programme for Integrated Control of Neglected Tropical Diseases (2009–2011), and followed protocols approved by the MoH for these programmes. The selected survey populations were called for meetings during which they were given detailed information about LF and the background, purpose and implications of the survey. Individuals volunteering to be examined provided oral informed consent under observation of both project staff and village authorities (parents/guardians consented on behalf of children below 15 years). Oral consent is the traditional way for making agreements in the survey areas, where written consent is unfamiliar and would cause suspicion and refusal to participate.

All 72 districts of Zambia existing at the start of the activity in 2003 (some have later been split and/or reorganized) were targeted for LF mapping. Based on previous reports and hospital records indicating possible cases of LF, 14 districts located in eight provinces were first selected. These were Choma and Sinazongwe (Southern Province), Mpongwe (Copperbelt Province), Kalabo, Sesheke and Senanga (Western Province), Mbala and Chinsali (Northern Province), Chama and Lundazi (Eastern Province), Luangwa and Kafue (Lusaka Province), Serenje (Central Province) and Zambezi (North-Western Province). In each of these districts, three chiefdoms were selected to provide 100 volunteers each to be tested for circulating filarial antigen (CFA) during 2003–2005.

In the remaining 58 districts, which were considered less likely to have LF, one site was identified for the mapping exercise and 100 volunteers were targeted for CFA testing at each site during 2009–2011. Selection of the sites was facilitated by local health personnel who led the survey team to areas where the population of people was high enough to allow the required number of people to be tested.

Members of the community were usually called to one central place for the CFA test. A clinic or health centre was found to be convenient for the purpose. Local health personnel were requested to assist in the exercise, and their presence brought confidence and trust, or less suspicion, from the community members. Geographical coordinates (longitude, latitude and elevation) were taken at the survey sites using a hand held GPS receiver (eTrex Summit, Garmin Corporation, Taiwan).

Following WHO guidelines [8], [11], about 100 volunteers above the age of 15 years were tested for CFA from each survey site. At few sites, however, volunteers down to the age of 12 years were allowed to participate due to low numbers coming forward for the test. From each individual, 100 µl finger-prick blood was collected using a heparinized capillary tube. The blood was applied to the specimen pad of a rapid immunochromatographic test card (ICT card, Binax Inc., USA). The result was read as positive or negative ten minutes after the card was closed and was recorded on a survey form together with the name, sex and age of the volunteer. The data were entered in Excel, and later transferred to SPSS for exploratory analysis.

Proxy environmental variables that may potentially influence the distribution of the filarial parasite-host-mosquito system and hence LF transmission [9] were extracted from freely accessible Remote Sensing (RS) sources at spatial and temporal resolutions shown in Table 1. Daytime land surface temperature (LST day), night time land surface temperature (LST night) and the Normalized Difference Vegetation Index (NDVI) were averaged over the period 2001–2010 representing the climatic period of the LF survey, according to Zambia's three distinct climatic seasons: i) cold/dry season (May–August), ii) hot/dry season (September–November) and iii) hot/rainy season (December–April). Land cover data contained 23 different land cover classes for Zambia, which were re-classified into 7 categories; water bodies, wetlands, forests, urban areas, shrublands, grasslands and croplands and re-sampled to 1 km resolution. Rainfall estimates averaged over the climatic normal period 1950–2000 were obtained from the Worldclim project [12]. As a proxy for changes in the environment due to anthropometric activities, the human influence index (HII) [13] was used. A selection of maps of environmental predictors can be seen in Figure S1 in the supplementary material.

The MODIS Reprojection Tool (USGS) was used to convert the RS data to geo-referenced maps. Further processing of the environmental data and distance calculation to the nearest water bodies was carried out in ArcMap v. 10.0 (ESRI). Additional data processing was performed in Revolution R Enterprise version 4.0 (Revolution Analytics; Palo Alto, USA) and Stata/SE 10 (StataCorp LP; College Station, USA). To elucidate potential co-linearity among the environmental variables, a correlation (Pearson's test) matrix was constructed based on 10,000 randomly extracted pixel values for each of the environmental predictors, with variables above a threshold of r>0.75 not allowed to enter the same model.

To explore the ecological niche of the LF parasite-vector-host biocoenose in Zambia, a species distribution modeling approach was deployed. Species distribution models, also referred to as ecological niche models [14], are commonly used to predict the geographic range of a species by extracting associations between point presence data and environmental data layers. The relationships are then used to characterize the environmental requirements of the species, and finally to predict suitable habitats across unsurveyed areas.

Here, species distribution modeling was implemented using the MaxEnt approach [15], commonly used to explore and predict environmental suitability for species, and has been shown to perform well compared to other predictive algorithms in comparative studies [16], [17]. A brief explanation of MaxEnt modeling is given in Methods S1 in the supplementary material. Recently, the MaxEnt approach has also been applied in mapping the Africa continent-wide current and potential future distribution of LF [18] and schistososmiasis [19]. Specifically, MaxEnt, which builds on the principles of maximum entropy, was chosen as it allows a flexible modelling of the often complex non-linear associations of infection presence with environmental variables [16], [18]. This flexibility can help facilitate an improved understanding of the ecological niche of a species, a prerequisite for a more reliable mapping of the potential distribution [20]–[22]. Furthermore, the MaxEnt method does not require absence data for the species being modeled; instead it takes advantage of the background environmental data for the entire study area through the background sampling procedure (see supporting information for more details). An advantage of this is dealing with the risk of including ‘false’ absence records in the model that can arise from limitation of parasite detectability [23] and hence falsely indicate non-suitability of a location. Finally, MaxEnt copes relatively well with correlated variables (which environmental variables often are) through the inbuilt method for regularization (L1-regularization) known to be well-performing [24], making it possible to explore a wider breath of potential environmental dependencies.

Two separate models were explored, based on different prevalence value cut-offs: Model 1 was based on survey sites that had at least 5% prevalence, and model 2 used survey sites with at least 15% CFA prevalence as MaxEnt model input data. This was done to get an indication of the drivers of both the general distribution of endemic LF in Zambia (represented by the distribution of at least 5% CFA prevalence), as well as the distribution of medium to high levels of infection prevalence (at least15% prevalence).

The spatial output of the MaxEnt model consists of a continuous range of relative probabilities indicating, in the case of this study, presence of the host–parasite system at the given prevalence threshold. The default logistic model that gives predicted estimates between 0 and 1 of the probability of infection presence for each pixel in the map was used. It was chosen to fit only linear, quadratic and product relationships, since more complex models can be difficult to specify a priori based on ecological theory [25]. Other parameterizations (maximum number of iterations and convergence threshold) followed recommendations by the model developers [15], [26].

The importance of the environmental variables was evaluated by comparing estimates of the relative contribution of environmental factors to overall model training gain. The gain is a measure closely related to deviance, the goodness of fit measure used in generalized additive and generalized linear models [15]. Furthermore, the explanatory information in each variable when used in isolation and the information lost when omitted from the model was quantified using a jackknife cross-evaluation procedure.

The continuous probability maps were furthermore converted into binary presence/absence maps of the LF host–parasite system, using the threshold indicating maximum training sensitivity plus specificity (i.e., that threshold which maximizes the sum of sensitivity and specificity for the training data). This is one of 11 thresholds calculated by MaxEnt, and is in considered one of the more robust of several standard thresholds for converting continuous probability surface to presence/absence surface [27], [28]. This distribution was then used to define the spatial limits for each of the two categories (≥5% and ≥15%) of infection prevalence.

A validation procedure was implemented by randomly dividing the occurrence data in training and test data sets (based on a 80–20% splitting of the data set). The evaluation focused on predictive performance at sites. Three statistics were applied; 1) the Area under the receiver operating characteristic Curve (AUC), 2) correlation (COR) and 3) sensitivity and specificity, to assess the agreement between the prevalence recorded at sites and the predictions.

AUC ranges from 0 to 1, where an AUC≤0.5 indicates that model performance is equal to or worse than that of a random prediction while an AUC above 0.75 is normally considered useful [17]. COR was calculated as the Pearson correlation coefficient between the full range of prevalence values in the test dataset (including negative sites) and the model logistic prediction [16], [26]. Sensitivity was calculated as the proportion of true positives/negatives (‘presence/absence’ points) falling within the predicted presence/absence area, and specificity as the proportion of true negatives falling within the predicted absence area.

A total of 10193 volunteers from 108 survey sites located in all 72 districts and 9 provinces of Zambia were surveyed for CFA. Among these, 9964 (97.8%) had a valid test card result and comprise the study population of examined individuals analyzed in this study. An overview of the survey sites, and the number, positivity for CFA, age and sex of the study population, is presented in Table 2. A list of geographical coordinates for the study sites is given in Table S1 in the supplementary material.

Most survey sites (83 or 76.9%) had more than 90 examined individuals, whereas 14 sites (13.0%) had less than 70. The highest mean number of examined individuals per site (100.9) was in Copperbelt Province, whereas the lowest (55.9) was in Luapula Province. The age of examined individuals ranged from 12 to 96 years. The mean age for the survey sites ranged from 21.2 to 46.0 years, and the overall mean age was 34.0 years. Many more females than males were examined (6376 vs. 3585), and the great majority of sites had more examined females than males (94 or 87.0%).

CFA positive cases were identified at 84 (77.8%) of the survey sites, where the prevalence ranged from 1.0 to 53.9%. The prevalence was ≥5% at 49 sites and ≥15% at 14 sites. The highest mean CFA prevalences were seen in Western (19.0%) and Lusaka (18.8%) provinces, whereas the lowest were in Copperbelt (3.4%) and North-Western (2.5%) provinces. The overall mean CFA prevalence for all examined sites was 7.4%.

A graphical presentation of the measured CFA prevalence at the different survey sites is shown in Figure 1, which thus gives an overview of the distribution pattern of LF in Zambia. All provinces had sites with a low CFA prevalence below 15%. However, six foci with CFA prevalences above 15% are clearly identified from the figure. Named by the district of location, these are the Kalabo and Senanga focus (both in Western Province), the Luangwa and Kafue focus (both in Lusaka Province), the Serenje focus (Central Province) and the Lundazi focus (Eastern Province). The first four of these foci had sites with particularly high CFA prevalences of >25%, and among these the Kalabo focus had sites with >50% CFA prevalence.

The importance of the environmental determinants of LF distribution in Zambia, as measured by their contribution to overall model training gain, varied substantially between the model based on ≥5% and the model based on ≥15% CFA prevalence data. The relative contribution of the 7 most important (of a total of 16) of the environmental predictor variables is given in Table 3. Between them, these 7 predictors were ranked in the top three of at least one of the two models. Overall the most important predictor was land cover, which in particular for model 1 (CFA≥5%) contributed to a significant part of model gain. In particular croplands and grasslands were associated with high probabilities of presence of infection, whereas forested areas were predicted as the least suitable of the land cover classes. The second most important predictor variable was day land surface temperature (hot/dry season). It was especially important in the CFA≥15% model, where it contributed 22.4% of total training gain, which is in accordance with the jackknife procedure that indicated that it was also the variable with the highest model gain when used in isolation. The human influence index, HII, was also an important predictor in model 1 where it contributed 20.9% of total training gain, whereas it did not play a significant role in model 2, contributing only 1.5% of total training gain.

The least important environmental factors for both models, as judged from the total gain, were rainfall and night time LST. The environmental variable that decreased the gain most when omitted was the distance to surface water bodies, which therefore appeared to have the most information not present in the other variables.

The functional relationship between the most important continuous predictor variables and the predicted probability of presence of either ≥5% or ≥15% CFA is depicted in the response curves in Figure 2. Each curve is made by generating a MaxEnt model using only the corresponding predictor variable, disregarding all other variables.

Maps of the MaxEnt predicted distributions of low (≥5% CFA) and medium-high LF infection prevalence (≥15% CFA) categories are presented in Figure 3a and 3b respectively. The heat map values represent the probabilities of ‘presence’ of each prevalence category, with relative probability values ranging from 0 (green colors) to 1 (red colors). The scale is defined for each map so that red areas correspond to ‘presence areas’ as defined by the threshold indicating maximum training sensitivity plus specificity.

Both maps indicate that LF infection potentially is present across Zambia with a somewhat patchy distribution, but with particularly high probability of presence in the floodplains of Western Province, the western part of North-western Province, the flood plain areas surrounding Zambezi River and its tributaries, the areas along Lake Kariba, the Kafue plains and the low plateau and river floodplains of Luangwa River. The most notable difference between the two maps is the much more confined presence areas predicted for the ≥15% prevalence category in the Northern and Luapula Provinces as compared to the relatively large areas predicted as potential ≥5% prevalence presence in these provinces.

Superimposing the binary presence/absence maps to produce one risk map (Figure 4) furthermore highlighted differences and similarities between the two model predictions. The orange color in this figure represents areas where only model 1 (CFA≥5%) predicts presence and the dark-brown color shows where model 2 (CFA≥15%) predicts presence (nested within model 1 predicted presence areas). Finally, the light yellow color in the map delineates areas where none of the models predict presence, i.e. areas expected to have no or less than 5% infection prevalence.

Measures of model accuracy are presented in Table 3. AUC values ranged from 0.866 to 0.892, indicating that the ‘suitability’ for LF infection was correctly ranked for 87–89% of the evaluated map pixels.

The correlation (COR) between the MaxEnt model predicted suitability and the observed full range of CFA prevalences at all 108 localities ranged from 0.117–0.355, and increased with CFA prevalence cut-off level (Table 3). This indicates that MaxEnt modeled the ‘true’ prevalence pattern of LF infection in Zambia better when using medium to high prevalence localities only (model 2), rather than the more general presence of infection (≥5%) which showed a non-significant correlation to the observed CFA prevalences at survey sites. Based on the presence/absence map, model 2 also had the best predictive positive and negative performance as evaluated by its sensitivity (76.9%) and specificity (64.5%) meaning that 76.9% of the ≥15% prevalence data points were correctly identified within the predicted ‘≥15% prevalence zone’, and that 64.5% of the true negatives were correctly identified within the ‘<15% prevalence zone’.