Warning: fopen(/home/virtual/epih/journal/upload/ip_log/ip_log_2024-05.txt): failed to open stream: Permission denied in /home/virtual/lib/view_data.php on line 95 Warning: fwrite() expects parameter 1 to be resource, boolean given in /home/virtual/lib/view_data.php on line 96 Spatial modeling of cutaneous leishmaniasis in Iranian army units during 2014-2017 using a hierarchical Bayesian method and the spatial scan statistic

Spatial modeling of cutaneous leishmaniasis in Iranian army units during 2014-2017 using a hierarchical Bayesian method and the spatial scan statistic

Article information

Epidemiol Health. 2018;40.e2018032
Publication date (electronic) : 2018 July 13
doi : https://doi.org/10.4178/epih.e2018032
1Modern Epidemiology Research Center, AJA University of Medical Sciences, Tehran, Iran
2Infectious Diseases Research Center, AJA University of Medical Sciences, Tehran, Iran
3Department of Preventive Medicine, Deputy of Health, AJA University of Medical Sciences, Tehran, Iran
Correspondence: Ali Reza Khoshdel   Modern Epidemiology Research Center, AJA University of Medical Sciences, Fatemi St., Shahid Etemadzadeh St., Tehran 1411718541, Iran  E-mail: alikhoshdel@yahoo.com
Received 2018 June 3; Accepted 2018 July 13.



This study aimed to map the incidence of cutaneous leishmaniasis (CL) in Iranian army units (IAUs) and to identify possible spatial clusters.


This ecological study investigated incident cases of CL between 2014 and 2017. CL data were extracted from the CL registry maintained by the deputy of health of AJA University of Medical Sciences. The standardized incidence ratio (SIR) of CL was computed with a Besag, York, and Mollié model. The purely spatial scan statistic was employed to detect the most likely high-and low-rate clusters and to obtain the observed-to-expected (O/E) ratio for each detected cluster. The statistical significance of the clusters was assessed using the log likelihood ratio (LLR) test and Monte Carlo hypothesis testing.


A total of 1,144 new CL cases occurred in IAUs from 2014 to 2017, with an incidence rate of 260 per 100,000. Isfahan and Khuzestan Provinces were found to have more CL cases than expected in all studied years (SIR>1), while Kermanshah, Kerman, and Fars Provinces were observed to have been high-risk areas in only some years of the study period. The most significant CL cluster was in Kermanshah Province (O/E, 67.88; LLR, 1,200.62; p<0.001), followed by clusters in Isfahan Province (O/E, 6.02; LLR, 513.24; p<0.001) and Khuzestan Province (O/E, 2.35; LLR, 73.71; p<0.001), while low-rate clusters were located in the northeast areas, including Razavi Khorasan, North Khorasan, Semnan, and Golestan Provinces (O/E, 0.03; LLR, 95.11; p<0.001).


This study identified high-risk areas for CL. These findings have public health implications and should be considered when planning control interventions among IAUs.


Cutaneous leishmaniasis (CL), despite being a neglected disease, is one of the most important parasitic diseases in Iran due to several factors, including areas of high endemicity [1-3], the annual incidence of more than 24,000 cases [4], high outbreak rates [5-7], and emergence in new foci [8]. Zoonotic CL, caused by Leishmania major, is distributed mainly in the central [9] and southwest regions [10], while anthroponotic CL, caused by Leishmania tropica, affects almost all urban areas [11]. The spatial inequality of CL in Iran makes it a threat to public health and poses major challenges to control strategies [12]. Studying the spatial and geographic patterns of CL is important because the components of the chain of infection, including the parasite, host vector, and required environmental conditions, are spatially distributed [13-15]. Although several studies have been published on the epidemiological and spatial patterns of leishmaniasis in the general population in Iran, the spatial distribution of CL among special populations, such as the army and military personnel, has not been adequately studied [16,17]. The army is the most vulnerable group, with the highest incidence of this disease, due to its distinct spatial status [18], and immunologically naive troops are sometimes deployed to endemic areas for training or operational activities [19,20]. Knowledge of the spatial patterns of CL in at-risk populations can help direct control strategies for truly needy areas and reduce the burden of the disease.

When studying the spatial distribution of a disease on a defined map, however, there are some methodological challenges that must be addressed. Overlooking the effect of spatial aggregation when attempting to detect spatial clusters can increase the risk of false detection [21]. Kulldorff et al. [22]’s spatial scan statistic was introduced as an efficient method to detect spatial clusters and disease outbreaks. There is degree of population heterogeneity in the various areas of a defined map, and estimating the rate of a disease in low-population areas is accompanied by a degree of uncertainty, a phenomenon known as small area estimation or the sparse data bias; however, the rates in various areas have a degree of spatial autocorrelation, and overlooking this in spatial modeling can also induce serious bias [23-25]. Spatial methods, such as the Besag, York, and Mollié (BYM) spatial model [25,26], have been developed to address these challenges. Therefore, using up-to-date data on CL incidence from 2015 to 2017 in Iranian army units (IAUs), the objectives of the present study were to investigate the spatial distribution of CL using a BYM model and to detect significant spatial clusters of CL by applying Kulldorff et al. [22]’s spatial scan statistic.


Study area and setting

Iran is situated in southwestern Asia and is the second largest country in the Middle East, with an area of 1,648,195 km². The country is located between the Caspian Sea in the north and the Persian Gulf and Gulf of Oman in the south. Iran is bordered by Afghanistan, Armenia, Azerbaijan, Iraq, Pakistan, Turkey, and Turkmenistan. The latitude and longitude of Iran are 32.4279°N and 53.6880°E, respectively. Iran has great diversity in topographical features, with elevation extremes ranging from -28 m at the Caspian Sea to 5,610 m at Mount Damavand. The great variations in altitude result in a correspondingly great diversity of environmental conditions, geographical factors, and climate indicators. In Iran, on some days temperatures can easily reach to lower than -30°C and some days reach to more than 40°C; January is the coldest month with a mean temperature of 5 to 10°C and August is the hottest month with mean temperatures of 20 to 30°C. Rainfall has a seasonal pattern, with great variation in the amount of rainfall (e.g., the mean annual rainfall ranges from 100 to 2,000 mm) [27].

Study design and data source

An ecological study was conducted. The study included confirmed cases of CL registered by the deputy of health of AJA University of Medical Sciences between March 2014 and March 2017. The surveillance system for diseases in army units involves the systematic collection, analysis, and interpretation of health-related outcomes for use in planning, implementing, and evaluating public health and prevention interventions. This surveillance system provides a capacity for the timely referral of confirmed cases to the next higher level and finally to the deputy of health of AJA University of Medical Sciences. Reporting systems vary depending on the type of health-related outcome data and information being reported; for example, confirmed CL cases are reported monthly. The data is registered for each Anno Hegirae Solaris (A.H.S.) or solar Hijri year. A.H.S. is the unit of years used in the official calendar of Iran, with 12 months that roughly range from April of one year to March of the following year. Therefore, each year in this study was defined from April of one year to March of the following year (e.g., April 2014 to March 2015).

The data were extracted from regular surveillance reports on CL cases and transferred to Microsoft Excel (Microsoft, Redmond, WA, USA). Then, the data were entered into the Iran ArcGIS shapefile map layer for spatial analysis and to generate spatial visualizations. A queen spatial weight matrix (first order) was defined using GeoDa [28] for calculating spatial autocorrelation and for hierarchical Bayesian analysis. Weighting based on the queen order indicates whether spatial units share a boundary, so that all spatial units in contact with spatial unit i are considered as its neighbors.

Statistical analysis

Spatial autocorrelation

The global spatial autocorrelation of CL incidence in IAUs at the province level during the study period was checked using the Moran’s index (MI). This index varies between -1.0 and 1.0; values near -1.0 indicate negative spatial autocorrelation (clustering of dissimilar values), while those close to 1.0 indicate positive spatial autocorrelation (clustering of similar values) and 0.0 indicates randomness. GeoDa [28] was used to calculate the MI and to test its statistical significance.

Spatial mapping

The crude incidence rate for area i can be calculated through dividing the number of new cases by the number of people residing in area i. An important goal of spatial epidemiology is to compare the incidence rates of a disease in different areas. However, the crude incidence rate does not allow valid and robust comparative inferences across areas because it is not corrected for the number of people at risk or the sex and age structure of each area, and it faces several other methodological issues in a spatial setting. In order to derive valid inferences, the incidence rates must be standardized, and indicators such as the standardized incidence ratio (SIR) must be estimated [29]. In the present study, the incidence rate was corrected only for the population at risk. Since the number of cases in each area follows a Poisson distribution, the SIR is calculated as the ratio of the observed number of incident cases of a disease (Oi) in the study area i to the number of cases that would be expected (Ei) if the study area i had the standard or total incidence rate. Ei is calculated as follows:

Ei = nii0iini,i = 1,2, ... . I

Where ni is the number of at-risk individuals in area i. Thus, the SIR is calculated as follows:

SIRi = OiEi

Crude SIRs are not corrected for spatial autocorrelation or small area estimation. Ignoring these methodological issues can result in SIRs that are biased due to overdispersion or extra-Poisson variation [24]. To obtain maps that remove the sources of overdispersion, hierarchical Bayesian approaches, such as the BYM model, have been introduced [26]. The BYM model and how the prior and hyperprior are defined for the parameters of interest have been explained in detail elsewhere [30]. Briefly, the smoothed SIRs (posterior distribution of the parameters) in a BYM model are calculated by combining the predefined information about sources of overdispersion (prior distribution) and the observed data. Sources of overdispersion are explained in the BYM model through 2 random effects: (1) the exchangeable (non-spatial) random effect ui , which refers to how the SIR in area i is shrunk toward the global mean of the study area (the prior distribution of ui follows a normal distribution, ui~N0, τu2; and (2) the spatial random effect vi , which refers to how the SIR in area i is shrunk towards the local mean of neighboring areas and is defined by a conditional autoregressive model Nνi¯, τi2. τu and τv are the parameters of precision, and the hyperprior will be defined for them in the next step. We used a gamma (0.01, 0.01) hyperprior for τu2, a non-informative prior. It has been suggested that this choice involves an empirical trade-off between (a) not using a more than sufficiently informative prior when there is no robust knowledge about the spatial pattern of the outcome (e.g., CL among IAUs) and (b) enough assurance for convergence of the Markov-chain Monte Carlo (MCMC) algorithm when a flat prior is used [31]. For the spatial random effect, a gamma (0.5, 0.005) hyperprior, which is a diffuse prior, was used for τv2. It has been proven that with these values, an artificial spatial structure will not be imposed on the estimates of relative risks (RRs)[31].

For this study, the posterior estimates of the SIRs were obtained by simulating from the joint posterior by means of 50,000 MCMC iterations on 2 parallel chains, with the initial 4,999 discarded in the burn-in phase. The MCMC convergence was checked in several ways. First, the Brooks-Gelman-Rubin (BGR) diagnostic was used, as values of BGR near 1 indicate the convergence of the model [32]. Second, the integrated autocorrelation time was used to quantify the effects of sampling error on the results to determine whether the samples in the chain were independent. The integrated autocorrelation time directly quantifies the Monte Carlo (MC) error. As rule of thumb, an MC error lower than 5% of the posterior standard deviation (SD) guarantee convergence. Thinning was used to reduce autocorrelation in the samples in the chains, as some samples were discarded from the simulation [33]. We specified thin= 1 so that the first parameter value would be retained. We used the open-access software OpenBUGS 3.2.3 (http://openbugs.net/w/FrontPage) to conduct hierarchical Bayesian analysis. The posterior effect size estimate, median SIR, and Bayesian 95% credible intervals (2.5th and 97.5th sample percentiles) were reported. The code for the BYM model is presented in Supplementary Material 1.

Spatial cluster analysis

Kulldorff et al. [22]’s purely spatial scan statistic was used to identify spatial clusters with significantly high and low rates. A high-rate cluster was defined as an area where the observed number of CL cases exceeded the expected number, and a low-rate cluster was defined as an area where the observed number was lower than the expected number of cases. The spatial scan statistic imposes thousands or millions of overlapping circular windows to scan the study area. The windows (potential clusters) are centered in each study region. The radius of the window varies continuously in size from zero to a specified upper limit. By default, the maximum spatial window size (MSWS) involved 50% or less of the total population at risk. Often, the set of potential clusters is identified hierarchically, with no geographical overlap. However, the results obtained when using the hierarchical approach to define a set of non-overlapping clusters can be misleading. When the MSWS is set to 50% of the population, the hierarchical no-geographical-overlap approach provides overly large clusters with relatively small RR values, whereas when the upper limits of the MSWS are set with small sizes (e.g., 2%), small clusters with unexpectedly large RRs result [34]. It has been shown that the Gini coefficient can be used to determine the optimal MSWS [34]. To create a set of potential clusters based on the Gini index, first, the MSWS is defined with different upper limits for the size of the cluster, and then the hierarchical no-geographical-overlap approach is used to create a collection of clusters for each different upper limit. After that, the Gini index is calculated for a set of non-overlapping clusters of each upper limit, and the collection that maximizes the Gini index is the optimal upper limit [34].

Since the number of CL cases in each potential cluster has a Poisson distribution, the exponential model-based spatial scan statistic was used as a probability model to perform cluster analysis and testing. For each circular window (identified potential cluster), the likelihood ratio statistic under the Poisson distribution was calculated, as follows:


Where C is the total number of CL cases, c is the observed number of CL cases within the window, E[c] is the crude expected number of CL cases within the window under the null hypothesis, and C-E[c] is the expected number of CL cases outside the window. The cluster with the maximum log likelihood ratio was taken as the most likely cluster. The maximum likelihood method test was used for deviations from the null hypothesis that the number of CL cases inside and outside of the window would be equal (Ho: θinout vs. the alternative hypothesis Ha: θin≠θout). The statistical significance of the detected clusters was tested using MC hypothesis testing with 999 permutations. SaTScan version 9.4.2 (https://www.satscan.org/) developed by Kulldorff [35] was used for spatial cluster analysis. All cartographic manipulations and displays were performed in ArcGIS version 10.3 (Esri, Redlands, CA, USA).


In total, 1,144 cases of CL were identified, with an incidence rate of 260 per 100,000 (86.66 per 100,000 in each year) during the study period. The MI for CL incidence showed a degree of negative spatial autocorrelation, although it was non-significant (MI, -0.12; p= 0.20 for 999 permutations).

The results of the crude SIRs (a frequentist analysis), as well as estimates from the BYM model (a Bayesian analysis), are presented in detail in Supplementary Materials 2-5. As expected, the absence of any CL cases in some provinces resulted in a crude SIR of 0; however, by borrowing information from neighbors, the BYM model provided smoothed results and the SIRs of 0 disappeared. There was little difference between the SIRs estimated in the BYM model and the crude SIRs. The value of MC error (less than 5% of posterior SD) for the estimated smoothed SIRs indicated convergence in the MCMC algorithm.

Figure 1 illustrates the smoothed SIRs from the BYM model for the time period of 2014 to 2017. Figure 1 suggests that the incidence of CL was highest in the provinces of Kermanshah (SIR, 67.60), Isfahan (SIR, 6.00), and Khuzestan (SIR, 2.34). Figure 2 illustrates the results of the SaTScan spatial analysis for the most likely cluster and secondary clusters 1 and 2 for the time period of 2014 to 2017. The results of the SaTScan spatial analysis are presented in detail in Table 1. The Gini index of 0.84 indicates great statistical dispersion among the detected clusters. The maximum size cluster was 10%. A total of 10 non-overlapping SaTScan highrate and low-rate clusters were identified, only 1 of which was not a statistically significant circular window. The SaTScan spatial analysis revealed that the most likely primary spatial high-rate cluster (observed to expected [O/E], 67.88; LRR, 1,200.62; p< 0.001) was in Kermanshah Province, implying that the incidence of CL was 67.88 times greater within this cluster than in the rest of the study area during the study period. The most likely low-rate cluster (O/E, 0.03; LRR, 95.11; p< 0.001) included Razavi Khorasan, North Khorasan, Semnan, and Golestan Provinces, implying that the incidence of CL was 33.33 times lower within this cluster than in the rest of the study area.

Figure 1.

The smoothed standardized incidence ratios (SIRs) of cutaneous leishmaniasis in Iranian army units based on a Besag, York, and Mollié model (2014-2017).

Figure 2.

Clusters with a statistically significantly higher or lower than expected incidence of cutaneous leishmaniasis in Iranian army units (2014-2017).

Characteristics of detected clusters with a higher or lower than expected incidence of CL in Iranian army units (2014-2017)

The spatial patterns of the smoothed SIRs from the BYM model in the 3 distinct time periods of 2014-2015, 2015-2016 and 2016-2017 are shown in Figure 3. The spatial pattern of the most likely high-rate and low-rate clusters from the SaTScan spatial analysis for the 3 time periods of 2014-2015, 2015-2016, and 2016-2017 are shown in Supplementary Material 6. The estimation performance of the BYM model and the SaTScan spatial analysis was similar. In 2014-2015, the highest rates occurred in Kermanshah Province (SIR, 95.93), followed by Isfahan and Khuzestan Province (SIR, 2.00-4.00) (Figure 3A). In 2015-2016, the number of CL cases exceeded the expected values in the 4 provinces of Isfahan (SIR, 10.11), Khuzestan (SIR, 1.98), Kerman (SIR, 1.75), and Fars (SIR, 1.12) (Figure 3B). However, in 2016- 2017, only 2 provinces (Isfahan and Khuzestan) were at an elevated risk for CL (e.g., the number of CL cases in Isfahan Province was 12.27 times higher than was expected) (Figure 3C).

Figure 3.

The smoothed standardized incidence ratios (SIRs) of cutaneous leishmaniasis in Iranian army units based on a Besag, York, and Mollié model (A) 2014-2015, (B) 2015-2016, and (C) 2016-2017.


Spatial analysis at different geographical levels can provide valuable information for developing policies to support the equitable distribution of healthcare and prevention resources. This study aimed to explore and map the spatial distribution of CL in IAUs, with cluster detection of CL incidence at the province level. The incidence of CL was 260 per 100,000 among IAUs during 2014- 2017 (86.66 per 100,000 in each year). A hierarchical Bayesian analysis showed a degree of inequality in the spatial distribution of CL, as the provinces of Isfahan, Khuzestan, Kermanshah, Kerman, and Fars had a higher than expected risk. The spatial scan statistic identified significant clusters that indicated that the most likely clusters of CL were primarily located in Kermanshah, Isfahan, and Khuzestan Provinces. In contrast, the most significant low-rate clusters of CL were identified in the northeastern provinces of Razavi Khorasan, North Khorasan, Golestan, and Semnan.

In the present study, the annual incidence rate of CL in IAUs was found to be 86.66 per 100,000, which is in accordance with previous studies that have shown the incidence of CL in the military population to be higher than its incidence in the general population in Iran [16,36]. For example, the study of Holakouie-Naieni et al. [36], which aimed to conduct spatial modeling of CL in Iran from 1983 to 2013, showed that the annual incidence of CL was 30.90 per 100,000 in the Iranian population. The study of Pakzad et al. [16], which compared the incidence of CL between IAUs and the general population during 2005-2014, showed that the incidence of CL in IAUs (143.68 per 100,000) was higher than the incidence in the general population (25.86 per 100,000). Based on a comparison with the study of Pakzad et al. [16] it can be deduced that a considerable decline in CL incidence had occurred among IAUs by 2014-2017.

Military forces are one of the most vulnerable groups to CL, and special attention to this group is very important for the prevention and control of this disease. The employment conditions of the armed forces are unique, in that soldiers are forced to settle in rural, border, and marginal urban areas and to live in areas where leishmaniasis is endemic; moreover, due to occasional large-scale transportation (e.g., military operations), they transfer the disease to other parts of the country [17,19,37].

Our results showed that the estimation performance of the BYM model and the Poisson model was similar. This implies that there was no underlying extra-Poisson variability in the data and that the introduction of spatial and non-spatial structured random effects in the BYM model fitted to the dataset was of little value for mapping modeling. In other words, modeling a spatial structure when it does not actually exist in the data may lead to a phenomenon known as over-fitting. The main reason for applying the BYM model in the present study was to explain all residual spatial variation, so we focused on providing robust and valid RRs. Our results showed that there was a degree of negative spatial autocorrelation in the distribution of CL among IAUs, although it was non-significant. Moreover, we considered a diffuse prior (e.g., a gamma [0.5, 0.005] hyperprior) for the precision parameter of the spatial random effect, as this choice precludes imposing an artificial spatial structure on the data. Although the differences between the SIRs from the Poisson model and those from the BYM model were negligible, explaining structured and unstructured spatial residuals with the BYM model may lead to much less bias in the SIRs. Moreover, a simulation study suggested that the results of ecological associations from the BYM model may be the least biased, even when true extra-Poisson variability does not exist in the dataset [38].

Our study showed that the provinces located in the western (Kermanshah), central (Isfahan), and south and southwestern areas (Khuzestan) were at a higher risk for CL during the study period. The studies in this field have shown that provinces located in cold regions, such as northern and northwestern Iran, mostly along the Zagros and Alborz mountains, had a lower incidence of CL. Due to the high altitude of those regions, the common vector of CL cannot transmit the pathogen. It has been concluded that the incidence of CL has an inverse correlation with altitude, because cold has a negative effect on the survival of the vector and the parasite of CL [39-41]. Similarly, other studies have shown that dry and desert areas had higher incidence rates due to having appropriate temperatures for carriers and parasites [42].

Due to variability in climactic conditions and the distinct activities of military forces, the risk of outbreaks should always be considered. The SaTScan spatial analysis showed that the most likely cluster of higher than expected incidence of CL was located in Kermanshah Province. Kermanshah Province is located in a mountainous region along the Zagros mountains in Iran, but within this province, there exists great variability in climactic and ecological conditions, as some counties (e.g., Sarpol-e Zahab and Qasr-e-Shirin) are located in arid and tropical regions. Studies have shown that the majority of CL cases in Kermanshah Province could be attributed to the 2 counties of Sarpol-e Zahab and Qasr-e-Shirin, which have warm and dry weather that provides a suitable setting for the survival of the vector and parasite of CL [43,44]. Thus, from the O/E ratio of 67.88 for the Kermanshah cluster, it can be concluded that when environmental conditions permit the chain of infection of CL to be established, outbreaks of CL can occur among army personnel to a remarkable extent.

This study also showed that secondary and tertiary high-risk clusters were located in Isfahan and Khuzestan Provinces, respectively. These findings are consistent with studies conducted in this field [16,17]. The existence of ideal environmental conditions for the vector, reservoir, and parasite of CL in these provinces has made them active foci for CL [45-47]. Additionally, inappropriate and inadequate sanitary facilities increase the chance of growth and proliferation of reservoirs (rodents) and vectors of CL in these 2 provinces. Our study has several limitations that should be considered. First, an obvious weakness of ecological studies, such as our study, is the ecological fallacy, whereby inferences on the aggregate level may not be true on the individual level. Second, the observed results may have been influenced by the modifiable areal unit problem, which is a well-known problem in geography in which results may be the function of the defined area units. Third, the estimated SIRs for CL may have been influenced by climactic variables, because the components of the chain of CL infection are all environmentally sensitive. However, data on climactic variables were not available in our study. Finally, we used a type of administrative and surveillance data, but results derived from such data should be interpreted with caution, because these data may not have been gathered for research purposes.

In conclusion, our study showed that the incidence of CL in IAUs was distributed differently at the province level. Due to the distinct conditions of army personnel, outbreaks of CL in non-endemic areas such as Kermanshah are expected. Additionally, army personnel in regions with clusters of CL, such as Isfahan and Khuzestan Provinces, must be targeted for further prevention and control interventions.


We would like to thank the staff of the deputy of health of AJA University of Medical Sciences for providing the data.


The authors have no conflicts of interest to declare for this study.


Supplementary Material 1: OpenBUGS code for BYM model is available at http://www.e-epih.org/.


Supplementary Material 2: Table S1 is available at http://www.e-epih.org/.


Supplementary Material 3: Table S2 is available at http://www.e-epih.org/.


Supplementary Material 4: Table S3 is available at http://www.e-epih.org/.


Supplementary Material 5: Table S4 is available at http://www.e-epih.org/.


Supplementary Material 6: Figure S1 is available at http://www.e-epih.org/.



1. Bamorovat M, Sharifi I, Aflatoonian MR, Sharifi H, Karamoozian A, Sharifi F, et al. Risk factors for anthroponotic cutaneous leishmaniasis in unresponsive and responsive patients in a major focus, southeast of Iran. PLoS One 2018;13e0192236.
2. Ramezankhani R, Sajjadi N, Nezakati Esmaeilzadeh R, Jozi SA, Shirzadi MR. Climate and environmental factors affecting the incidence of cutaneous leishmaniasis in Isfahan, Iran. Environ Sci Pollut Res Int 2018;25:11516–11526.
3. Saghafipour A, Vatandoost H, Zahraei-Ramazani AR, YaghoobiErshadi MR, Jooshin MK, Rassi Y, et al. Epidemiological study on cutaneous leishmaniasis in an endemic area, of Qom Province, central Iran. J Arthropod Borne Dis 2017;11:403–413.
4. Alvar J, Vélez ID, Bern C, Herrero M, Desjeux P, Cano J, et al. Leishmaniasis worldwide and global estimates of its incidence. PLoS One 2012;7e35671.
5. Aflatoonian MR, Sharifi I, Hakimi Parizi M, Fekri AR, Aflatoonian B, Sharifi M, et al. A prospective cohort study of cutaneous leishmaniasis risk and opium addiction in south eastern Iran. PLoS One 2014;9e89043.
6. Karimi A, Hanafi-Bojd AA, Yaghoobi-Ershadi MR, Akhavan AA, Ghezelbash Z. Spatial and temporal distributions of phlebotomine sand flies (Diptera: Psychodidae), vectors of leishmaniasis, in Iran. Acta Trop 2014;132:131–139.
7. Saghafipour A, Rassi Y, Abai MR, Oshaghi MA, Farzinnia B, Mostafavi R, et al. Outbreak of zoonotic cutaneous leishmaniasis: a report. Arch Hyg Sci 2013;2:48–54.
8. Ghatee MA, Haghdoost AA, Kooreshnia F, Kanannejad Z, Parisaie Z, Karamian M, et al. Role of environmental, climatic risk factors and livestock animals on the occurrence of cutaneous leishmaniasis in newly emerging focus in Iran. J Infect Public Health 2018;11:425–433.
9. Karami M, Doudi M, Setorki M. Assessing epidemiology of cutaneous leishmaniasis in Isfahan, Iran. J Vector Borne Dis 2013;50:30–37.
10. Ostad M, Shirian S, Pishro F, Abbasi T, Ai A, Azimi F. Control of cutaneous leishmaniasis using geographic information systems from 2010 to 2014 in Khuzestan Province, Iran. PLoS One 2016;11e0159546.
11. Tashakori M, Ajdary S, Kariminia A, Mahboudi F, Alimohammadian MH. Characterization of leishmania species and L. major strains in different endemic areas of cutaneous leishmaniasis in Iran. Iran Biomed J 2003;7:43–50.
12. Norouzinezhad F, Ghaffari F, Norouzinejad A, Kaveh F, Gouya MM. Cutaneous leishmaniasis in Iran: results from an epidemiological study in urban and rural provinces. Asian Pac J Trop Biomed 2016;6:614–619.
13. Bounoua L, Kahime K, Houti L, Blakey T, Ebi KL, Zhang P, et al. Linking climate to incidence of zoonotic cutaneous leishmaniasis (L. major) in pre-Saharan North Africa. Int J Environ Res Public Health 2013;10:3172–3191.
14. Purse BV, Masante D, Golding N, Pigott D, Day JC, Ibañez-Bernal S, et al. How will climate change pathways and mitigation options alter incidence of vector-borne diseases? A framework for leishmaniasis in South and Meso-America. PLoS One 2017;12e0183583.
15. Talmoudi K, Bellali H, Ben-Alaya N, Saez M, Malouche D, Chahed MK. Modeling zoonotic cutaneous leishmaniasis incidence in central Tunisia from 2009-2015: forecasting models using climate variables as predictors. PLoS Negl Trop Dis 2017;11e0005844.
16. Pakzad R, Dabbagh-Moghaddam A, Mohebali M, Safiri S, Barati M. Spatio-temporal analysis of cutaneous leishmaniasis using geographic information system among Iranian Army Units and its comparison with the general population of Iran during 2005-2014. J Parasit Dis 2017;41:1114–1122.
17. Salehi-Moghadam A, Barati M, DabbaghMoghadam A, Khoshdel A, Totonchian M, Noorifard M. Temporal changes and mapping leishmaniasis in military units of IRI army. Hormozgan Med J 2014;18:91–98. (Persian).
18. Weina PJ, Neafie RC, Wortmann G, Polhemus M, Aronson NE. Old world leishmaniasis: an emerging infection among deployed US military and civilian workers. Clin Infect Dis 2004;39:1674–1680.
19. Oré M, Sáenz E, Cabrera R, Sanchez JF, De Los Santos MB, Lucas CM, et al. Outbreak of cutaneous leishmaniasis in Peruvian military personnel undertaking training activities in the Amazon Basin, 2010. Am J Trop Med Hyg 2015;93:340–346.
20. Patino LH, Mendez C, Rodriguez O, Romero Y, Velandia D, Alvarado M, et al. Spatial distribution, leishmania species and clinical traits of cutaneous leishmaniasis cases in the Colombian army. PLoS Negl Trop Dis 2017;11e0005876.
21. Ozonoff A, Jeffery C, Manjourides J, White LF, Pagano M. Effect of spatial resolution on cluster detection: a simulation study. Int J Health Geogr 2007;6:52.
22. Kulldorff M, Heffernan R, Hartman J, Assunção R, Mostashari F. A space-time permutation scan statistic for disease outbreak detection. PLoS Med 2005;2e59.
23. Elliott P, Wartenberg D. Spatial epidemiology: current approaches and future challenges. Environ Health Perspect 2004;112:998–1006.
24. Lawson AB, Browne WJ, Rodeiro CL. Disease mapping with WinBUGS and MLwiN New Jersey: John Wiley & Sons; 2003. p. 10–13.
25. Waller LA, Carlin BP. Disease mapping. Chapman Hall CRC Handb Mod Stat Methods 2010;2010:217–243.
26. Besag J, York J, Mollié A. Bayesian image restoration, with two applications in spatial statistics. Ann Inst Stat Math 1991;43:1–59.
27. Mostafavi E, Haghdoost A, Khakifirouz S, Chinikar S. Spatial analysis of Crimean Congo hemorrhagic fever in Iran. Am J Trop Med Hyg 2013;89:1135–1141.
28. Anselin L, Syabri I, Kho Y. GeoDa: an introduction to spatial data analysis. Geogr Anal 2006;38:5–22.
29. Gordis L. Epidemiology 3rd edth ed. Philadelphia: Saunders; 2004. p. 64–65.
30. Ayubi E, Mansournia MA, Motlagh AG, Mosavi-Jarrahi A, Hosseini A, Yazdani K. Exploring neighborhood inequality in female breast cancer incidence in Tehran using Bayesian spatial models and a spatial scan statistic. Epidemiol Health 2017;39e2017021.
31. Kelsall JE, Wakefield JC. Discussion of ‘Bayesian models for spatially correlated disease and exposure data’, by Best et al. Bayesian Stat 1999;6:151.
32. Brooks SP, Gelman A. General methods for monitoring convergence of iterative simulations. J Comput Graph Stat 1998;7:434–455.
33. Hamra G, MacLehose R, Richardson D. Markov chain Monte Carlo: an introduction for epidemiologists. Int J Epidemiol 2013;42:627–634.
34. Han J, Zhu L, Kulldorff M, Hostovich S, Stinchcomb DG, Tatalovich Z, et al. Using Gini coefficient to determining optimal cluster reporting sizes for spatial scan statistics. Int J Health Geogr 2016;15:27.
35. Kulldorff M. SaTScan TM user guide for version 9.4 [cited 2015 Sep 15]. Available from: http://www.satscan.org/.
36. Holakouie-Naieni K, Mostafavi E, Boloorani AD, Mohebali M, Pakzad R. Spatial modeling of cutaneous leishmaniasis in Iran from 1983 to 2013. Acta Trop 2017;166:67–73.
37. Hamidi FR, Tajik AR, Keshavarz A, Rajabi J, Hosseini SJ. Evaluating outbreak of cutaneous leishmaniasis in military personnel referring to the 554 military Zanjan Hospital in 2008. Ann Mil Health Sci Res 2011;9:105–109.
38. Latouche A, Guihenneuc-Jouyaux C, Girard C, Hémon D. Robustness of the BYM model in absence of spatial variation in the residuals. Int J Health Geogr 2007;6:39.
39. Salahi-Moghaddam A, Khoshdel A, Hanafi-Bojd AA, Sedaghat MM. Mapping and review of leishmaniasis, its vectors and main reservoirs in Iran. J Kerman Univ Med Sci 2015;22:83–104. (Persian).
40. Ali-Akbarpour M, Mohammadbeigi A, Tabatabaee SH, Hatam G. Spatial analysis of eco-environmental risk factors of cutaneous leishmaniasis in southern Iran. J Cutan Aesthet Surg 2012;5:30–35.
41. Nadim A, Javadian E, Mohebali M, Momeni A, Ardehali S, Yaghoobi-Ershadi MR, et al. Leishmania and leishmaniasis 3rdth ed. Tehran: Academic Publication Center; 2008. p. 45–46. (Persian).
42. Chelbi I, Kaabi B, Béjaoui M, Derbali M, Zhioua E. Spatial correlation between Phlebotomus papatasi Scopoli (Diptera: Psychodidae) and incidence of zoonotic cutaneous leishmaniasis in Tunisia. J Med Entomol 2009;46:400–402.
43. Goldust K. Prevalence of cutaneous leishmaniasis in the patients referred to the skin section in the center for education and treatment in Kermanshah district 1990-1994 [dissertation]. Kermanshah: Kermanshah University; 1994. (Persian).
44. Hamzavi Y, Sobhi SA, Rezaei M. Epidemiological factors of cutaneous lieshmaniasis in patients reffered to health centers in Kermanshah Province (2001-2006). J Kermanshah Univ Med Sci 2009;13:151–161. (Persian).
45. Yaghoobi-Ershadi M, Javadian E, Tahvildare-Bidruni GH. Leishmania major MON-26 isolated from naturally infected Phlebotomus papatasi (Diptera: Psychodidae) in Isfahan Province, Iran. Acta Trop 1995;59:279–282.
46. Parvizi P, Baghban N, Novin EA, Absavaran A. Detection, identification and molecular typing of Leishmania major in Phlebotomus papatasi from a focus of zoonotic cutaneous leishmaniasis in central of Iran. Exp Parasitol 2010;124:232–237.
47. Talari SA, Talaei R, Shajari G, Vakili Z, Taghaviardakani A. Childhood cutaneous leishmaniasis: report of 117 cases from Iran. Korean J Parasitol 2006;44:355–360.

Article information Continued

Figure 1.

The smoothed standardized incidence ratios (SIRs) of cutaneous leishmaniasis in Iranian army units based on a Besag, York, and Mollié model (2014-2017).

Figure 2.

Clusters with a statistically significantly higher or lower than expected incidence of cutaneous leishmaniasis in Iranian army units (2014-2017).

Figure 3.

The smoothed standardized incidence ratios (SIRs) of cutaneous leishmaniasis in Iranian army units based on a Besag, York, and Mollié model (A) 2014-2015, (B) 2015-2016, and (C) 2016-2017.

Table 1.

Characteristics of detected clusters with a higher or lower than expected incidence of CL in Iranian army units (2014-2017)

Optimal Gini coefficient MSC Clusters detected Involved province(s) Observed cases Expected cases O/E RR1 LLR p-value
Area rate
 High 0.84 0.1 1 Kermanshah 353 5.20 67.88 97.73 1,200.62 <0.001
2 Isfahan 454 75.40 6.02 9.33 513.24 <0.001
3 Khuzestan 232 98.8 2.35 2.69 73.71 <0.001
 Low 1 Razavi Khorasan, North Khorasan, Semnan, Golestan 3 104.00 0.03 0.03 95.11 <0.001
2 Zanjan, Qazvin, Hamadan, Kurdistan, Gilan 1 96.20 0.01 0.01 94.83 <0.001
3 Kerman, South Khorasan, Hormozgan 15 114.40 0.13 0.12 73.57 <0.001
4 East Azerbaijan 13 91 0.14 0.13 55.52 <0.001
5 Bushehr, Fars 28 104 0.27 0.25 41.97 <0.001
6 West Azerbaijan 2 31.20 0.06 0.06 24.08 <0.001
7 Qom, Markazi 1 7.80 0.13 0.13 4.76 0.09

CL, cutaneous leishmaniasis; MSC, maximum size cluster; O/E, observed-to-expected; RR, relative risk; LLR, log likelihood ratio.


RR was calculated as the observed number of cases divided by the expected number of cases within the cluster divided by the corresponding ratio outside the cluster.