Genetic diversity and structure in Leishmania infantum populations from southeastern Europe revealed by microsatellite analysis

Background The dynamic re-emergence of visceral leishmaniasis (VL) in south Europe and the northward shift to Leishmania-free European countries are well-documented. However, the epidemiology of VL due to Leishmania infantum in southeastern (SE) Europe and the Balkans is inadequately examined. Herein, we aim to re-evaluate and compare the population structure of L. infantum in SE and southwestern (SW) Europe. Methods Leishmania strains collected from humans and canines in Turkey, Cyprus, Bulgaria, Greece, Albania and Croatia, were characterized by the K26-PCR assay and multilocus enzyme electrophoresis (MLEE). Genetic diversity was assessed by multilocus microsatellite typing (MLMT) and MLM Types were analyzed by model- and distance- based algorithms to infer the population structure of 128 L. infantum strains. Results L. infantum MON-1 was found predominant in SE Europe, whilst 16.8% of strains were MON-98. Distinct genetic populations revealed clear differentiation between SE and SW European strains. Interestingly, Cypriot canine isolates were genetically isolated and formed a monophyletic group, suggesting the constitution of a clonal MON-1 population circulating among dogs. In contrast, two highly heterogeneous populations enclosed all MON-1 and MON-98 strains from the other SE European countries. Structure sub-clustering, phylogenetic and Splitstree analysis also revealed two distinct Croatian subpopulations. A mosaic of evolutionary effects resulted in consecutive sub-structuring, which indicated substantial differentiation and gene flow among strains of both zymodemes. Conclusions This is the first population genetic study of L. infantum in SE Europe and the Balkans. Our findings demonstrate the differentiation between SE and SW European strains; revealing the partition of Croatian strains between these populations and the genetic isolation of Cypriot strains. This mirrors the geographic position of Croatia located in central Europe and the natural isolation of the island of Cyprus. We have analysed the largest number of MON-98 strains so far. Our results indicate extensive gene flow, recombination and no differentiation between MON-1 and MON-98 zymodemes. No correlation either to host specificity or place and year of strain isolation was identified. Our findings may be associated with intensive host migration and common eco-epidemiological characteristics in these countries and give valuable insight into the dynamics of VL.


Background
Leishmania infantum is the main causative agent of human visceral leishmaniasis (HVL) and canine visceral leishmaniasis (CanL) and also responsible for sporadic cases of cutaneous leishmaniasis (CL) across the Mediterranean basin [1,2]. The concurrent species Leishmania tropica extending rarely to Europe via Greece [1,[3][4][5] and Leishmania major present in North Africa and Middle East, are also accountable for anthroponotic and zoonotic CL, correspondingly [1].
Dogs constitute the main reservoir host of L. infantum MON-1 [9] and probably MON-98 [5], while cats could be a secondary domestic L. infantum MON-1 reservoir [17,18]. Feline leishmaniasis is reported in countries highly endemic for CanL [19] and cats were recently found capable of transmitting MON-1 parasites to a proven L. infantum vector, Phlebotomus perniciosus [17]. The large number of sub-clinically infected dogs is a major veterinary and public health problem in south Europe [1]. Although clinically affected dogs are more prone to infect sand fly vectors, those being sub-clinically infected and seronegative can also transmit Leishmania parasites to sand flies and, therefore, contribute to the parasites' maintenance [20][21][22][23][24].
Further to this, the distribution range of CanL has surpassed the limits of south Europe and is spreading toward northern European countries, presently considered non-endemic. In Germany, there is an estimate of 20,000 infected dogs that were either imported from the Mediterranean basin or infected after travelling to this region [25,26]. Furthermore, the presence of a L. infantum vector (P. perniciosus) in Germany and Switzerland has been recently reported [27]. The northward spread of CanL is also apparent in northern Spain, where new foci of CanL and the presence of its vectors (Phlebotomus ariasi and P. perniciosus) have been reported. In Italy, canine incidence has reached 30.3% in West Liguria [28], and in France seroprevalence fluctuates between 3.2% in Alpes-Maritimes and 26.5% in Corsica [29].
Concurrently, a considerable increase in VL cases due to L. infantum MON-1 is documented in traditionally endemic southeastern (SE) European countries. Since the first reported HVL cases in the islands of Spetses (1879) and Hydra (1881) [30], leishmaniasis is endemic in most islands and coastal regions of Greece and constitutes a significant veterinary and public health concern.
During 1994-2001, HVL has re-emerged in inland regions of Greece, where seroprevalence reached 33.4% in Attica and 12.6% in Epirus [31]. During the same 8-year survey, the highest CanL seroprevalence in clinically affected dogs was found in Attica (48.4%) followed by Epirus (45.4%) and central Macedonia (40%) [31]. Of additional concern are the high seropositivity rates found in clinically healthy dogs across the mainland of Greece. According to recent studies, seroprevalence in healthy dog populations ranged between 22.4-30.1% in Attica, 21-24.4% in Epirus, 12.3%-17% in central Greece and 10.8%-20.5% in central Macedonia [31,32]. Among all Greek islands most dramatic changes have been witnessed in Crete, where CanL increased from 0.3% to 19.8% during 1990-2009 and HVL prevalence from 0 to 1.33 per 100,000 inhabitants [4,5].
In south Cyprus, HVL and CanL do co-exist, but show different epidemiological patterns. In contrast to other endemic Mediterranean countries, in Cyprus two separate concurrent disease patterns involve CanL due to L. infantum MON-1 [11,33] and HVL due to Leishmania donovani MON-37 [34][35][36]. CanL was widespread by 1935 and constituted a major veterinary health problem until 1945 [33]. Since the completion of subsequent eradication campaigns against malaria and echinococcosis, CanL gradually re-emerged in coastal areas and is now widespread and highly endemic. In particular, the overall canine incidence has risen from 1.7% to 14.9%-19.6% [11,33]. Currently, the highest prevalence rates found in all prefectures of south Cyprus are 33.3% (Larnaca), 32.4% (Limassol), 30.8% (Paphos) and 20% (Nicosia), while 24.3% of the PCR positive and 9.2% of the seropositive dogs are subclinically infected [11]. HVL due to L. infantum MON-1 has been absent since 1987 [33] and there are no official reports of infected natives even in areas highly endemic for CanL in south Cyprus. Yet, an increased incidence of human leishmaniasis possibly due to L. infantum has been reported in the northern part of the island [11,37].
In Croatia, HVL and CanL existed across the Dalmatian coasts (Split, Makarska, and Dubrovnik areas) until 1957. The majority of reported HVL cases (90.4%) involved children under 10 years old and 4% of these cases were fatal [40]. The disease occurred only sporadically for 30 years due to mass spraying with insecticides for malaria control, but dramatically re-emerged during the war and postwar periods. During the past two decades, human incidence has reached 0.4/10 5 of the population in southern Dalmatia [40] and is still considered a paediatric disease as 64% of cases involve children [41]. Reports on HVL cases suggest that existing differences regarding the clinical course of infection and/or disease outcome may be associated, amongst other factors, with variations of parasite species [42]. CanL is endemic in the city of Split (14.7%) and across the enzootic Split-Dalmatia County (7.1%-42.8%), while 45.7% of seropositive dogs were found to be sub-clinically infected [43].
In Albania, the incidence and geographical spread of HVL and CanL has gradually increased since World War II. However, the disease sharply re-emerged during 1997-2001 with an average 173 HVL cases/year and 2.8/ 10 4 population morbidity, which exceeded 1/10 3 lethality in active foci across north and central Albania [44]. Endemic zones are no longer restricted to northern regions or coastal and hilly areas. Children comprise 70-87.7% of HVL cases with an annual incidence of 25/10 5 (0-6 year age group) [44][45][46]. Adult HVL cases and HIV/Leishmania co-infections are also recorded [45]. Since 1977, CanL is widespread across the country with 16-17% seroprevalence recorded in Tirana and several other districts [44,46].
In Bulgaria, canine seroprevalence has gradually decreased from 81% in 1941 to 13.8% in 1993 [47,48]. However, during 1988-2004 an HVL outbreak occurred in south Bulgaria and an unexpectedly high morbidity rate (50.7%) was noticed in Stara Zagora. In 2005, two autochthonous HVL cases were reported in the same region, but all stray dogs tested were found seronegative [49]. Conversely, in 2006 two breed dogs, which were not imported and had never travelled outside Bulgaria, were diagnosed with CanL [47]. L. infantum is incriminated for human and canine disease [49], underreporting is considered mild and the estimated annual VL incidence ranged between 8 and 12 cases/year during 2004-2008 [50].
The epidemiology of VL in these SE European and Balkan countries has been either inadequately or never examined so far. Multilocus microsatellite typing (MLMT), a powerful tool for phylogenetic and population genetic studies in Leishmania [51][52][53][54][55], previously showed that L. infantum MON-1 strains from Greece differentiate from those of Western Mediterranean [15] and North Africa [56,57]. In addition, there is no indication of subdivision between Turkish and Greek strains or within Greek strains as well as no differentiation between MON-1 and MON-98 strains from Greece [15]. These findings prompted us to initiate a comprehensive population genetic study in SE Europe. For this, 109 Leishmania strains (25 human and 84 canine isolates) were collected from Turkey, Cyprus, Bulgaria, Greece, Albania and Croatia, and characterized by MLEE and K26 typing. The analysis of L. infantum strains by MLEE has demonstrated a broad enzymatic polymorphism within this species but has limited epidemiological applicability in Europe due to the predominance of MON-1 [6]. The K26-PCR assay, which is specific for the L. donovani complex, allows simultaneous discrimination between L. infantum MON-1 and more than nine L. infantum non MON-1 zymodemes that are known to exist in Europe [14,58,59]. Notably, the 626 bp K26 amplicon corresponds only to zymodeme MON-1, while the 940 bp K26 amplicon may correspond either to zymodeme MON-1 or MON-98. For microsatellite typing, a standard set of highly discriminative markers [52] was employed for all strains from SE and their MLMT profiles were compared to 19 previously analyzed strains of L. infantum MON-1 from southwestern (SW) European countries (France, Spain and Portugal) aiming to re-assess the genetic diversity and population structure of L. infantum in SE Europe.

Methods
Leishmania strains and study area Table 1 summarizes the characteristics of the 128 L. infantum strains included in the study. Fifty-four strains were isolated from 20 human and 34 canine VL cases in diverse geographic regions of the mainland and islands of Greece. Sixteen human and canine isolates were obtained from Turkey, mainly from the Aegean region. Of these, nine MON-1 or MON-98 strains from Greece and two Turkish strains of unknown zymodeme type (strains EP3 and EP16) were previously analyzed in different studies [15,52]. Twenty-three canine isolates were collected from five prefectures of southern Cyprus and sixteen strains were isolated in Bulgaria, Albania and Croatia, mostly from canines. Nineteen previously analyzed L. infantum MON-1 strains from France, Spain and Portugal [15,52] were incorporated for comparison.

Parasite cultures, cloning and DNA extraction
Leishmania promastigotes were cultured at 25°C in RPMI 1640 or M199 cell growth medium supplemented with 10-20% heat-inactivated FBS, 5U/ml penicillin, 5 μg/ ml streptomycin and 10 mM HEPES buffer. The number of in vitro passages was limited to a minimum to avoid contamination and genetic drift due to serial sub-culturing. Clonal populations from the mixed L. infantum MON-1/L. donovani MON-37 culture of the CD44 strain were obtained as described elsewhere [36]. Leishmania DNA was extracted from parasite cultures using a classical phenol/ chloroform protocol [60] or the standard protocol for cultured cells using the NucleoSpin® Tissue kit (Macherey-Nagel). Extracted DNA was respectively dissolved in dH 2 O    or eluted in TE buffer and stored at −20°C, whereas working solutions were kept at 4°C until use.

K26-PCR, MLEE and MLMT analysis
The K26-PCR assay was applied to discriminate MON-1 from other L. infantum zymodemes according to the length polymorphism of the K26 gene, as described previously [14]. MLEE analysis based on starch gel electrophoresis of fifteen enzymatic loci [6] was performed for all but six Turkish strains, at the University of Montpellier (France), to complement K26 typing.
Multilocus microsatellite typing (MLMT) was applied using fourteen microsatellite markers previously found highly discriminative within the L. donovani complex [52,61]. Microsatellite fragments were amplified using established PCR protocols [52]. Amplicons were resolved on high-resolution MetaPhor® gels (3.5-4.5%) using a large gel horizontal electrophoresis system with built-in recirculation (Owl Large Gel System, Thermo Scientific) and screened for length polymorphisms using the AlphaIma-ger® imaging system (Alpha Innotech Corp.) as previously described [36]. Alternatively, PCR assays were performed using fluorescence conjugated primers (WellRed dyes, Proligo) and analyzed by capillary electrophoresis (CEQ™ 8000, Beckman Coultier), as described elsewhere [36]. Control strains were included in both PCR protocols and reading methods to confirm the reproducibility of the results.

Population structure and phylogenetic analysis
The MLMT dataset of 128 strains, including that of 30 previously analyzed strains [15,36,52,55,61], was processed by models based on Bayesian statistics and genetic distances. Model-based analysis was performed by STRUC-TURE v2.3.1 [62] to infer the population structure and identify genetically distinct L. infantum clusters. Based on multilocus genotypes and independent of a mutation model, this algorithm was used to determine patterns of allele frequencies per analyzed locus and estimated the most probable degree of differentiation. Individuals were assigned in discrete genetic populations and exposed their membership to each estimated population by genotypic fractions. The allele frequencies among populations were correlated by admixture modeling for a series of runs using a 'burn-in' period of 2x10 5 iterations and probability estimates based on 2×10 6 of MCMC repeats. For each possible number of populations between K = 1 and K = 10, ten independent simulations were conducted to estimate the delta (Δ)K values [63]. The likelihood of population number was calculated by STRUCTURE HARVESTER [64], implementing the Evanno method [63]. Subsequent STRUCTURE analyses were performed to achieve assignment at different hierarchies.
Phylogenetic analysis was performed based on microsatellite genetic distances, which were calculated by MICRO-SAT software using the proportion of shared alleles distance measure (D AS ). The resulting distance matrix input was processed by PHYLIP v3.6 to construct a Neighbour-joining (NJ) tree with confidence intervals by bootstrapping (1000 replications). A consensus tree with the corresponding bootstrap values and branch lengths was re-constructed using GENEIOUS v5.4 [65] and the derived mid-point rooted NJ tree was edited by FIGTREE (http://tree.bio.ed.ac.uk/software/figtree). Genetic relatedness and recombination were evaluated through Neigh-borNet networks constructed by SplitsTree4 software and through the population membership coefficients estimated by STRUCTURE. The MLMT data were analyzed by GDA (http://hydrodictyon.eeb.uconn.edu/people/plewis/ software.php) with respect to the proportion of polymorphic loci (P), allelic variation (A), inbreeding coefficient (F IS ), expected (H e ) and observed (H o ) heterozygosity. GDA analysis was performed per microsatellite locus and per population/subpopulation defined by STRUCTURE. F-statistics were estimated for all populations/subpopulations to evaluate genetic differentiation and gene flow. For this, the F ST values (pairwise Wright's fixation index) with their corresponding p-values (confidence test) were calculated using MSA [66]. Moreover, the total strain-set was sorted by geographical source to calculate the descriptive statistics per endemic country. Similarly, a separate descriptive analysis was performed for Turkish/Greek MON-1 and MON-98 strains based on their zymodeme type and origin.
The spatial distribution of the observed genetic structure at population and individual level was assessed based on allelic similarity by factorial correspondence analysis (FCA) using GENETIX v4.03 [67].

Ethics statement
The Leishmania strains designated with an asterisk (*, Table 1) and those originating from SW Europe ( R , Table 1), which were used in this microsatellite analysis for comparison, have been the object of previous publications [14,15,36,52,55,61,68].
Strains from Albania were obtained from the National Reference Centre of Leishmania at Montpellier, France. Ethical considerations are described in Petrela et al. [45]. Strains from Greece were obtained from the alreadyexisting "Leishmania cryobank collection" of the Medical School of the University of Crete. Informed written consent was obtained from all patients. Both the study and protocols used were approved by the Ethical Scientific Committee of the Medical School (University of Crete, Greece). Ethical considerations are described in Kuhls et al. [14,15,36,52,55,61,68]. Strains from Turkey were collected in the context of previous studies, as indicated above. Human strains were obtained from an already-existing "Leishmania Group cryobank collection". Informed written consent was obtained from all patients. The study was approved by Local Animal Care and Ethics Committee of the School of Medicine Izmir, Turkey.
All samples were coded and anonymized, and all Leishmania strains were isolated as part of normal diagnosis. Dog samples were collected by veterinarians after the owner's consent, following protocols that were approved by the Ethics Committees of the respective Institutions in Turkey, Cyprus, Greece and Croatia. Questionnaires with personal, epidemiological, and clinical data for each dog were also completed.
Approval from the respective Institutions to use the samples for research purposes was received in all cases. The present study was also reviewed and approved by the Ethics Committee of the Hellenic Pasteur Institute, Athens, Greece.

K26-PCR and MLEE analysis
Of the total 128 strain-set, 125 strains were characterized at species level and zymodeme type by K26-PCR, and 122 of them also by MLEE, in the context of either this or previous studies [15,36] (Table 1). A single 626 bp K26 amplicon characteristic for L. infantum MON-1 was presented by a few Turkish and Greek strains as well as by all strains from Cyprus and Croatia, including two clones of strain CD44 (CD44cl.4 and CD44cl.5) isolated from a canine with a mixed L. infantum MON-1/L. donovani MON-37 infection in Cyprus [34,36]. The 940 bp K26 amplicon, which could correspond to either L. infantum MON-1 or MON-98 zymodeme and was observed so far only in Greek strains, was presented by all strains from Turkey, Bulgaria and Albania, except for two Turkish strains from regions other than the Aegean. These were strains EP88 and EP134 isolated at the Black Sea and the Mediterranean region close to Cyprus, respectively. Besides the main 940 bp K26 amplicon observed for most strains from Greece, two Greek strains (GD192 and GD223) gave multiple K26-PCR products (626/870/ 940 bp and 626/940 bp, respectively). In addition to seven previously typed L. infantum MON-98 strains from Crete, MLEE analysis identified 14 new strains of MON-98 from eight human and canine isolates from mainland regions of Greece and the island of Crete as well as six canine isolates from the Aegean region of Turkey (Kuşadasi) (Figure 1). Overall, our sample-set comprises 101 MON-1 strains, 21 MON-98 strains from Turkey and Greece, and 6 Turkish L. infantum strains of unknown zymodeme type. Interestingly, 21 out of 70 Turkish and Greek strains (30%) were identified as MON-98, which is the highest percentage of this zymodeme found so far in Turkey and Greece.

Analysis of microsatellite loci and correlation with zymodeme type
The inter-and intra-zymodeme genetic diversity among 128 L. infantum MON-1 and MON-98 strains from south Europe was assessed by a panel of previously described microsatellite markers located on nine different chromosomes [15,36,52,53,61]. Table 2 summarizes the characteristics and descriptive statistics per analyzed locus for the overall strain-set.
Different degrees of polymorphism were observed among the analyzed loci. As shown in Table 2, the highest allelic richness, expected and observed heterozygosities were demonstrated for markers Lm4TA, Lm2TG and Li22-35 (E), with the expected heterozygosity (mean H e = 0.29) always being much higher than the observed (mean H o = 0.02). Inbreeding coefficients (F IS ) were substantially high (mean F IS = 0.92), as four markers presented absolute heterozygote deficiency (F IS = 1), and marker Li46-67 (C) was monomorphic as previously observed [52,54].
Separate analysis was performed to investigate intrazymodeme diversity and avoid possibly biased results due to the larger number of MON-1 strains in the total sample-set. The measures of genetic variability of Turkish/Greek strains demonstrated extensive inbreeding in both zymodeme populations and just trivial differences between MON-1 and MON-98 zymodemes (Table 3). However, some alleles at markers Lm2TG (22 repeats), CS20 (22 repeats), Li45-24 (G) (18 repeats) and LIST7039 (16 repeats) were only found in MON-98 strains.

Population structure inference
The Bayesian model-based clustering algorithm implemented in the STRUCTURE software was run by increasing K values from K = 2 to K = 10 and with seven geographically-defined strain-groups, as strains from France (FR), Spain (ES) and Portugal (PT) ( Table 1) representing SW Europe were considered to constitute one population, according to earlier reports [15,36,52,55,61].
Calculation of ΔK indicated the existence of four distinct genetic populations ( Figure 2Α Crete (Greece) were assigned to this population. Population assignment per analyzed strain at K = 4 is designated also in Figure 1 and Table 1.
Population structure at K = 4 was well-established in all ten iterations and different genetic distance analyses (Figures 3, 4 and Additional file 1: Figure S1) corroborated the existence of the four populations determined by the model-based analysis. However, the assignment of HR strains to POP3 and POP4 by STRUCTURE analysis at K = 4 ( Figure 2Α) was found inconsistent with their placement at the mid-point rooted NJ tree (Additional file 1: Figure S1) as well as with the distribution of the respective splits in the phylogenetic network (NeighborNet) (Figure 3). The FCA plot also confirmed the existence of at least four main populations (Figure 4). In overall, the three-dimensional representation mirrored the genetic isolation of CY strains (POP1), the differentiation of SW European strains (POP3) from all other populations as well as between the Balearic Islands and FR/mainland ES/ PT, and the diversification of TR/HR strains from SW European strains.
The population structure defined by STRUCTURE was however, unstable for other K values and subdivision was detected by re-analyzing these four populations independently. Population 1 split further in two subpopulations (not shown), but subdivision was not supported by other genetic distance analyses (Figures 3, 4 and Additional file 1: Figure S1).
Population 2 was divided in three subpopulations (not shown); however there was no correlation between subpopulations and host specificity, as both human and canine isolates were found in every subcluster. Additionally, the assignment of strains to these three subpopulations was not related to their geographical origin as strains from both mainland and coastal areas of TR and GR never divided, and the GR strains always grouped with those from BG and AL. Moreover, no differentiation based on zymodeme type was detected between GR MON-98 strains and BG/GR/AL MON-1 strains. In addition, the subdivision of POP2 was not supported by the different distance-based analyses used (Figures 3, 4 and Additional file 1: Figure S1).
Population 3 was divided in three stable subpopulations ( Figure 2B). All but one strain from the islands of ES were assigned to a distinct subpopulation (sub3A), as previously observed [15]. Sub3B 1 enclosed all strains from FR, PT and the mainland of ES as well as single strains from TR, HR and the ES island of Majorca. Of note is that the three strains from TR, HR and Majorca that grouped with strains from the mainland ES in sub3B 1 (Figure 2B), were placed in intermediate positions at the Splitstree and FCA plots (Figures 3 and 4) and that the single strains EP88 and AL from TR and HR respectively, formed a separate subgroup distant from other strains of POP3 at the NJ tree (Additional file 1: Figure S1). Moreover, strains EP88 and AL presented large splits in the phylogenetic network ( Figure 3) and were found more closely related to a group of TR/GR isolates (strains EP134, GD18 and GD19) that were assigned at POP4 by STRUCTURE. This finding is also illustrated in the FCA plot ( Figure 4). Accordingly in the NJ tree, strains EP88 and AL branched out and clustered with the same group of TR/GR strains after the formation of POP1 (Additional file 1: Figure S1). All other strains from HR (Lu, Ne, Me, SO and C) grouped in a separate homogeneous subpopulation (sub3B 2 ) ( Figure 2B), which was supported by all analysis methods (Figure 3, 4 and Additional file 1: Figure S1).
Re-analysis of population 4, which included the majority of the MON-98 strains, revealed further sub-structuring (4 stable subpopulations) based on both zymodeme type and geographical origin ( Figure 2C). Sub4A included only MON-98 strains from Turkey and Greece, while the Splitstree plot (Figure 3) illustrated gene flow between them. The group of HR strains (Pa, Mr, GS and U) and one strain from Albania (LEM3962), which was assigned to sub4B, presented a reticulate pattern in the Splitstree network ( Figure 3) rather closely related to two TR strains (EP58 and EP68) which were assigned to POP2. This was also apparent in the NJ tree (Additional file 1: Figure S1) and the FCA plot ( Figure 4). Turkish and Greek MON-1 and MON-98 strains were found in sub4C 1 , whereas sub4C 2 enclosed only GR MON-1 strains, most of which presented the 626 bp K26 amplicon (Table 1).

Genotypic characterization and F-statistics for populations and subpopulations
The measures of genetic diversity calculated for each main population defined by STRUCTURE are summarized in Table 4. Out of the 128 analyzed strains, 50 strains presented unique microsatellite profiles and 78 strains shared the same MLM Type with at least another strain. Identical profiles were found only among strains that were assigned in the same main population (Tables 1 and 4). The genetic polymorphism exhibited per analyzed locus for each population is given in Additional file 2: Table S1.
All strain populations were polymorphic at markers Lm2TG, Lm4TA, Li22-35 (E), Li23-41 (F) and CS20. The highest degree of polymorphism and allelic richness in all markers was observed within POP3 and POP4, where   Table 4). Regarding the subpopulations defined by STRUC-TURE and largely supported by the genetic distance analyses, the lowest levels of polymorphism, heterozygosity and inbreeding were found in sub3B 2 (HR strains) and sub4B (AL/HR strains) (Additional file 3: Table S2). The fact that lower F IS values are observed at subpopulation level indicate that Wahlund effect (population subdivision) accounts at least partially to the high inbreeding coefficients found for the main populations.
The four populations identified by model-and distancebased analyses were further supported by F-statistics (Table 5) indicating strong differentiation between them. CY strains that grouped in POP1, demonstrated the greatest genetic isolation from all other populations (F ST= 0.546-0.655). As shown in Additional file 4: Table S3, the observed subpopulations were also strongly supported by F-statistics (p ≤ 0.006).  (MON-1, MON-98), is demonstrated by large and small pie-charts, respectively. A. The likelihood of population number was calculated by STRUCTURE HARVESTER [64] implementing the Evanno method [63]. The major peak at K = 4 at the ΔK graph indicates the existence of four populations (populations 1-4) for the total strain set, which were constantly observed in all independent simulations (10 runs). Population subdivision was investigated by re-analyzing separately each population at K = 4; subdivision of populations 1 and 2 is not shown. B. Three geographically determined subpopulations (subs 3A, 3B 1 and 3B 2 ) were observed upon re-analysis of population 3. Bar plots for both K = 2 and K = 3 are shown. C. Four stable subpopulations (subs 4A, 4B, 4C 1 and 4C 2 ) based on both geographical origin and zymodeme type, were identified for population 4. Bar plots for both K = 3 and K = 4 are shown.

Correlation of MLMT profiles with origin, zymodeme type and host
The strain-set was also grouped based on geographical source and their MLM Types were analyzed to compare the degree of polymorphism between and within each endemic country. At least one matching pair was detected within each endemic country (Table 1).
Notably, the same MLMT profile (No 29) was identified for 17 out of 23 canine isolates from different prefectures of Cyprus and a second one (No 33) was shared by two canines from Paphos prefecture. Seven distinct profiles (No 45-51) were identified for the ten analyzed strains from Croatia. Of these, three multilocus genotypes (No 45, 47 and 48) were observed in more than one Croatian strain. Two canine isolates from Albania were indistinguishable by MLMT (No 54). The same MLMT profile (No 37) was identified among MON-98 canine isolates from Kuşadasi district (Aegean region), whereas a second one (No 41) was shared by strains of unknown zymodeme (EP58 and EP68) that were isolated from dogs in coastal areas of the Aegean (Muğla and Kuşadasi, respectively).
Only 20 out of 54 strains analyzed from Greece presented distinct MLM Types ( Table 1). The same MLMT profile (No 2) was shared by MON-1 strains (human and canine isolates) from the capital of Athens, the island of Mytilene and Lasithi prefecture in Crete. An identical profile (No 22) was shared by MON-1 strains (canine isolates) from Rhodes island and the capital of North Greece Thessaloniki, where the second largest seaport is located. Another two MLMT profiles (No 5 and 19) were identical between either MON-98 strains or MON-1 strains from Crete. MLM Type No 3 was common between MON-1 strains (human isolates) from Athens and MON-98 strains (canine isolates) from a harbor area in Attica, northern Peloponnese and Crete.
Furthermore, identical MLM Types were identified between strains originating from different countries, except  (Table 1) are underlined. Strains are represented by coloured boxes that correspond to their population or subpopulation assignment according to STRUCTURE analysis ( Figure 2). for strains from Croatia and Cyprus. Notably, the single canine isolate from Bulgaria shared the same MLMT profile (No 11) with one Albanian and three Greek strains. Also, two MLM Types (No 13 and No 25) were common between strains from Greece and Turkey. One of them was found in canine isolates from the Aegean region (Cyclades-Greece and Bodrum-Turkey, respectively) and a human isolate from Central Greece. The other MLM Type was shared solely by MON-1 canine isolates from north inland regions of Greece and Turkey and one strain of unknown zymodeme from Muğla (Aegean region). As shown in the Splitstree plot, these strains were separated in a clonal subgroup (Figure 3).
Characteristic alleles only for CY isolates were identified at markers Lm2TG and Lm4TA. Also, HR strains presented distinct alleles at marker Li71-7 (R) and GR strains at markers Li22-35 (E) and Li23-41 (F). No correlation was observed between specific microsatellite profiles and human or canine hosts.  (Table 1)

Discussion
The aim of the present study was to evaluate epidemiological aspects of L. infantum causing VL in the SE edge of Europe focusing on Turkey, Cyprus and Greece but also including strains from Balkan countries. Initially, the K26-PCR assay was applied to confirm the presence of L. infantum spp. and simultaneously to distinguish MON-1 from non MON-1 zymodemes within the strain-set [14]. Application of the K26 typing tool, which discriminates between different L. infantum and L. donovani zymodemes, was also essential due to the recent detection of L. donovani MON-37 in Cyprus [35] as well as the identification of the novel L. donovani MON-308 and MON-309 zymodemes in Turkey [36]. Notably, strains belonging to these three L. donovani zymodemes were found by MLMT to form a novel L. donovani sensu lato group [36].
In consistence with previous findings [14], K26 typing revealed that all strains from SE Mediterranean and Balkan countries were L. infantum presenting single K26 amplicons of 626 bp or 940 bp, except for two Greek MON-1 strains showing both K26-PCR products. In line with K26 typing, MLEE analysis confirmed the predominance of zymodeme MON-1. Of the 70 strains isolated from the Aegean region of Turkey, mainland regions of Greece and the island of Crete, 21 (30%) were typed, however, as zymodeme MON-98. This is the highest percentage of MON-98 strains found so far [5,14,69] and it should be noticed that Greek MON-1 and MON-98 strains were similarly distributed in the country. Since these are the only zymodemes found so far in the Balkans, the K26-PCR assay is adequate for Leishmania typing in this endemic region.
The MLMT profiles obtained for 109 L. infantum MON-1 and MON-98 strains from six SE European countries were compared to those of 19 previously analyzed L. infantum MON-1 strains from SW Europe (France, Spain and Portugal) [15,52,61]. A small number of Greek MON-1 and MON-98 strains and two Turkish L. infantum strains of unknown zymodeme type (strains EP3 and EP16) ( Table 1) had been also analyzed earlier [15,36,52,55,61]. In this study we have investigated additional strains from VL foci in Turkey and Greece and also included L. infantum strains from Cyprus, Bulgaria, Albania and Croatia, for which genetic diversity and population structure had not been studied until now. By this we could confirm the existence of substantial differentiation between L. infantum strains from the eastern and western parts of south Europe, and that population structuring is geographically determined, as previously demonstrated [70]. These findings are also in agreement with previous data according to which, Cypriot and Iberian strains could be distinguished upon processing their K26-PCR products by RFLP [14].
Different types of population genetic analyses, including Bayesian inference (as implemented in STRUCTURE), distance-based (NJ and Neighbor Net in SplitsTree) and FCA revealed the existence of four putative populations. Population 1 was well-defined by all methods and consisted exclusively of the Cypriot strains of L. infantum and was clearly separated from the other populations (Figures 3  and 4). As expected, all SW European strains formed a distinct population (POP3). Bayesian statistics has also assigned one strain from the Black Sea area in Turkey and six strains from Croatia to this population that was, however, not well-supported by the distance-based and FCA analyses, which indicates a rather intermediate position of these strains. Three subpopulations became apparent by separate re-analysis. The differentiation between strains from the Balearic Islands (sub3A) and strains from the Iberian Peninsula and France (sub3B 1 ) confirms previous analysis [15,52,61]. All but one strain from Croatia that grouped in POP3, were assigned to the third subpopulation (sub3B 2 ). The Turkish and Greek strains split among two populations, POP2 and POP4 with no correlation to their place and year of isolation, the zymodeme type, host specificity or pathology. POP2 was composed by most strains from Turkey, Greece, Albania and the single strain from Bulgaria. Nevertheless, distance and FCA analyses revealed that this population was rather homogenous and was not divided into subpopulations. In contrast, POP4 that enclosed some Turkish/Greek, one Albanian and four Croatian strains, appeared to be highly diverse and consisting of four subpopulations. Sub4A was formed only by Turkish/Greek MON-98 strains and sub4B enclosed only Albanian/Croatian strains. The other Turkish/Greek strains, whether of zymodeme MON-1 or MON-98, were assigned to sub4C 1 and sub4C 2 .
In general, model-and distance-based algorithms identified the same population and subpopulation divisions. The statistical support for the midpoint-rooted NJ tree was assessed by bootstrap analysis. However, sufficient bootstrap values were only observed for nodes between groups of strains found within the main clusters or subclusters (Additional file 1: Figure S1), while a spider-net pattern and large splits were often observed in Splitstree (Figure 3), probably due to the existence of homoplasy, intermediate or mixed genotypes and considerable genetic migration. This was expected, since the analyzed strains belong to genetically close L. infantum zymodemes and were isolated from geographically close countries with frequent movements of host populations and could explain the discrepancies found in the assignment of some groups of strains or individual strains by different model-and distance-based approaches. Analogous congruity has been also observed elsewhere [15,52]. Compared to previous studies [15,36,52,55,61], we have now analyzed a much larger number of MON-1 and MON-98 strains from Turkey and Greece, but as before, we were unable to differentiate these strains according to their geographical origin, Turkey vs. Greece or Greek mainland vs. the island of Crete, or between MON-1 and MON-98 strains. Most MON-98 strains from both Turkey and Greece were assigned to the highly diverse POP4 together with strains of zymodeme MON-1 from the same areas, but clustered in the two subpopulations sub4A and sub4C 1 ( Figure 2C). A distinct subpopulation of the MON-98 strains that grouped in POP2 was not supported by distance-based analyses. The reticulate network in Splitstree (Figure 3) clearly points to extensive gene flow and/or hybridization occurring among Turkish and Greek strains of both zymodemes. Our results were not affected when a group of representative non MON-1 strains, other than those of zymodeme MON-98, was included in our analysis (not shown). This finding is in agreement with previous analysis [15,36,52,55,61], according to which a small group of analyzed MON-98 strains from Crete did not cluster with other non MON-1 strains from SW Europe but always remained part of the MON-1 population. Our observations reflect the high polymorphism observed among different groups of L. infantum MON-1 strains rather than indicate differentiation between MON-1 and MON-98 zymodemes.
Identical MLM Types were found between MON-1 strains from Turkey and Greece, which along with the observed shared membership of POP2 and POP4 and the significantly high inbreeding level within TR/GR MON-1 (F IS ≥ 0.873) indicate considerable recombination and gene flow among them as previously suggested [15]. These findings may be associated with intensive host migration over long periods of time and common eco-epidemiological features, including a similar sand fly fauna [71], in Turkey and Greece. The detection of identical MLMT genotypes between Greek strains of both zymodemes in harbor areas of the mainland and in Greek islands suggests parasite transfer between the mainland and the islands.
The assignment of Greek human and canine isolates to the same genetic populations indicates the absence of host specificity. Since we have included a much larger number of strains, this observation enhances previous findings [15].
The genetic structure of strains from Crete appears to be a micrograph of the overall epidemiology of VL in Greece, which is in contrast to the clear division between strains from the mainland and islands of Spain, Majorca [8] and Ibiza [15]. At high hierarchy, the overall high genetic diversity, significant gene flow and population subdivision have most probably resulted from successive recombination or even hybridization events. In parallel, the observed low variability measures and considerable inbreeding observed within small subgroups of MON-1 or MON-98 strains from either inland regions of Greece or the island of Crete (Table 3) can be explained by successive subdivision followed by genetic drift within each of these groups and could indicate a founder effect. The existence of related clones of MON-1 as well as MON-98 strains could indicate selection of specific microsatellite genotypes, as previously observed for MON-1 strains from Madrid or Majorca and Ibiza [72].
L. infantum MON-1 canine strains from Cyprus were found genetically isolated from all other genetic groups and formed a monophyletic group (Additional file 1: Figure S1) presenting the lowest degree of polymorphism, allelic deprivation and traces of heterozygosity, which indicate genetic drift. Seventeen out of the 23 strains shared the same MLMT genotype suggesting an epidemic outbreak due to the same strain which has spread to different parts of the island. As expected the clone of strain CD44 (CD44cl.5) shared common alleles and grouped with the other MON-1 strains from Cyprus. The identification of different allelic variants than those typically found for L. donovani MON-37 strains [35], responsible for both HVL and CL, suggests the absence of recombination or hybridization between MON-1 and MON-37 strains. The epidemiology of CanL in Cyprus has been doubtless affected by the eradication campaigns held against malaria between 1940 and 1950 which have massively reduced the vector populations on the island. Another campaign was performed against echinococcosis between 1970 and 1975 which diminished the dog populations and altogether eliminated CanL for over 20 years [11]. Most probably this has caused severe bottleneck events in Leishmania populations leading to a limited parasite genetic pool characterized by significant inbreeding and genetic homogeneity. However, after the sand fly fauna has gradually recovered and CanL re-emerged, some degree of genetic variation was observed in the parasites as a result of adaptation to selective pressure.
The partition of canine MON-1 strains from Croatia between both SE (POP4) and SW (POP3) European populations (Figures 2, 4 and Additional file 1: Figure S1) mirrors the geographical position of the country. Despite the relatively small number of strains, we could document substantial genetic diversity for these strains. Croatia is situated between East and West Europe and parasites can be introduced from its neighboring countries. The microclimate, the presence of vectors and reservoirs favor the establishment and circulation of parasites in the Dalmatian coastal regions. The exportation from Croatia to neighboring countries is equally possible, as adult VL cases acquired at the Dalmatian coast were registered in Austria and Switzerland [40,42]. The genetic differentiation observed between Croatian and SW European strains on one side and between Turkish, Greek and Croatian strains, on the other side, can be also seen by the distinct and intermediate placement of the Croatian strains in the Splitstree (Figure 3). The rather low level of inbreeding and the substantially low variability values observed within these subpopulations (Tables 4 and Additional file 4: Table S3) indicate genetic drift, most probably caused by population substructure (Wahlund effect). Other explanations could be the restriction of the parasites' gene pool due to mass insecticide spraying for malaria control or the diversification of autochthonous cases circulating in Croatia for a long time [42]. Considering the unexpected genetic diversity of Croatian strains isolated from a rather small territory, future investigations should use a larger sample set including human and sand fly isolates, and should include strains from the northern Adriatic region and Italy.
Strains from Albania did not form a distinct genetic group in our analyses, but grouped either with strains from Turkey/Greece or Croatia. In contrast, a considerable amount of genetic flow was detected obviously reflecting the frequent human migration and the unrestrained movement of canine and sand fly populations between these countries. HVL in Albania still remains an infantile disease and presents by far the highest morbidity rates in Europe [45] although relapses and resistance to drug treatment are only rarely reported [44][45][46]. It raises concerns that 57.1% of seropositive HVL cases hospitalized in the Ioannina prefecture (North Greece) are immigrants from the neighboring villages in south Albania and that the highest prevalence of clinically affected HVL (42.8%) and CanL (82.8%) in Greece is reported in the Thesprotia prefecture (Epirus) close to the border between Greece and Albania [31]. These findings are alarming for both countries and suggest that the dynamics and complexity of VL in Albania could be similar to those observed in Greece. CanL is endemic across the country with an estimated 16-17% seroprevalence, while sand fly species Phlebotomus neglectus, Phlebotomus perfiliewi and Phlebotomus tobbi are widespread. It should not be underestimated that the maximum incidence of P. neglectus, which is highly present in all Greek islands except Gavdos [5,71] and is a proven L. infantum vector in the island of Corfu neighboring the Albanian coast [44], is reported in endemic foci presenting the highest VL morbidity rates. It should also be noticed that P. tobbi spp. despite its overall low abundance, is solely found in Albanian districts where VL morbidity exceeds 1/1000 [44].
The lack of isolated Leishmania strains from Bulgaria has so far restricted the investigation of the population structure and dynamics. We have included in our study the first strain that was isolated from a dog in the Blagoevgrad province (southwestern Bulgaria) that borders the Macedonia region of Greece (Figure 1). Our findings confirm the existence of at least L. infantum MON-1 zymodeme in Bulgaria [48] and demonstrate that the MLM Type of this particular canine isolate is identical to that of three strains from Greece and one strain from Albania. This is not surprising as Bulgaria shares common geo-ecological characteristics with its neighboring countries. The presence of L. infantum sand fly vectors and stable parasite circulation in Stara Zagora [49] highlights the necessity to conduct a comprehensive entomological and epizootic survey for the identification of Leishmania vectors and the isolation of Leishmania parasites from sand flies, human and canine hosts. This would give the opportunity to thoroughly investigate the population structure of L. infantum spp. and disease dynamics in Bulgaria.

Conclusions
This is the first study that gives an insight into the population structure of L. infantum in Cyprus located at the SE edge of Europe and the Balkans (Turkey, Bulgaria, Greece, Albania and Croatia). L. infantum MON-1 was found to be predominant in this region, while 16.8% of strains were typed as MON-98. Hence, the K26 typing tool proves adequate for typing strains from SE Europe and the Balkans. Herewith, we have identified distinct genetic populations that reveal the clear differentiation between SE and SW European strains, the intermediate placement and further diversification of Croatian strains from SW strains and the formation of a monophyletic group enclosing all canine MON-1 isolates from Cyprus. The latter could be explained by the natural isolation of Cyprus being an island and a bottleneck effect caused by the subsequent eradication campaigns against malaria and echinococcosis that were held in Cyprus. Also, we have analysed the highest number of MON-98 strains found so far in Turkey and Greece and did not identify a distinct population including all MON-98 strains. The latter finding is strengthened by the separate analysis of MON-98 and MON-1 strains, which shows that the genetic distances observed between MON-98 and MON-1 subgroups are similar to those found between different MON-1 subgroups. Altogether our results confirm and extend the study by Kuhls et al. [15] using a larger strainset from Turkey and Greece and including strains from other endemic Balkan countries. Notably, the clear diversification between Spanish strains from the mainland and the islands of Majorca and Ibiza is contrary to the clustering of Greek strains from the mainland and the islands of