Tick-borne pathogens in Finland: comparison of Ixodes ricinus and I. persulcatus in sympatric and parapatric areas

Background Almost 3500 tick samples, originally collected via a nationwide citizen science campaign in 2015, were screened to reveal the prevalence and distribution of a wide spectrum of established and putative tick-borne pathogens vectored by Ixodes ricinus and I. persulcatus in Finland. The unique geographical distribution of these two tick species in Finland allowed us to compare pathogen occurrence between an I. ricinus-dominated area (southern Finland), an I. persulcatus-dominated area (northern Finland), and a sympatric area (central Finland). Results Of the analysed ticks, almost 30% carried at least one pathogen and 2% carried more than one pathogen. A higher overall prevalence of tick-borne pathogens was observed in I. ricinus than in I. persulcatus: 30.0% (604/2014) versus 24.0% (348/1451), respectively. In addition, I. ricinus were more frequently co-infected than I. persulcatus: 2.4% (49/2014) versus 0.8% (12/1451), respectively. Causative agents of Lyme borreliosis, i.e. bacterial genospecies in Borrelia burgdorferi (sensu lato) group, were the most prevalent pathogens (overall 17%). “Candidatus Rickettsia tarasevichiae” was found for the first time in I. ricinus ticks and in Finnish ticks in general. Moreover, Babesia divergens, B. venatorum and “Candidatus Neoehrlichia mikurensis” were reported for the first time from the Finnish mainland. Conclusions The present study provides valuable information on the prevalence and geographical distribution of various tick-borne pathogens in I. ricinus and I. persulcatus ticks in Finland. Moreover, this comprehensive subset of ticks revealed the presence of rare and potentially dangerous pathogens. The highest prevalence of infected ticks was in the I. ricinus-dominated area in southern Finland, while the prevalence was essentially equal in sympatric and I. persulcatus-dominated areas. However, the highest infection rates for both species were in areas of their dominance, either in south or north Finland. Electronic supplementary material The online version of this article (10.1186/s13071-018-3131-y) contains supplementary material, which is available to authorized users.


Background
Ticks are recognized as the primary vectors for several pathogenic viruses, bacteria and protozoa worldwide [1,2]. In northern Europe, most notable tick-borne pathogens are Borrelia burgdorferi (sensu lato) spirochetes, of which at least seven genospecies are responsible for causing Lyme borreliosis (LB) [3,4]. In Finland, the primary vectors for tick-borne pathogens are Ixodes ricinus (Linnaeus, 1758) and Ixodes persulcatus (Schulze, 1930). The nationwide distribution of these two tick species was studied recently [5] showing I. ricinus dominance in southern Finland, a sympatric area in central Finland and I. persulcatus dominance in northern Finland. Tick distribution patterns may have an important role in the distribution and diversity of tick-borne pathogens as well.
Despite the presence of many potential pathogens in ticks in Finland, only a few tick-borne infections other than LB and TBE have been reported in humans. A single case of fatal babesiosis was described in a man with a rudimentary spleen, detected at their autopsy in 2004 [19]. Underlying reasons for the apparent discrepancy may be, for example, low pathogenicity of the putative pathogens, unclear clinical manifestations and unestablished diagnostic criteria of human infections, lack of awareness among health-care professionals of the emerging tick-borne diseases, and unavailable laboratory tests. On the other hand, the co-occurrence of several pathogens in ticks can lead to co-infections with different tick-borne pathogens in humans and animals [20][21][22][23]. Co-infections can alter the dynamics of pathogen transmission and pathogen interactions within a host animal, and increase the severity of manifestations in humans [23][24][25]. For example, A. phagocytophilum infects human neutrophils, modulates the immune response of the host and thereby increases susceptibility to other pathogens, including B. burgdorferi (s.l.) [26]. Thus, co-infections of tick-borne pathogens can have a significant impact on the disease manifestations making the diagnostics of these infections more challenging.
The aim of the present study was to map the major tick-borne pathogens circulating in ticks in Finland. The distribution of I. ricinus and I. persulcatus within the country creates an exceptional opportunity to study the potential differences in the prevalence of tick-borne pathogens in a sympatric area compared to areas dominated by a single tick species.

Origin of the samples
In 2015, citizens were asked to send ticks via postal mail to the University of Turku as a part of a tick collection campaign. This collection resulted in nearly 20,000 individual ticks received from all around Finland, up to the Arctic Circle. Detailed information about the collection, acquisition of the samples and tick identification has been described in a previous study [5]. Because that study also presented the occurrence of TBEV and B. miyamotoi in a subset of 2000 tick samples, these pathogens were not included in the present study. The presence of B. burgdorferi (s.l.) at the genospecies level was also described in the previous study [5]. However, B. burgdorferi (s.l.) is included in the present study since it is now identified to the genospecies level.
A subset of 3465 ticks (2014 I. ricinus and 1451 I. persulcatus) out of a total of 20,000 ticks were screened for B. burgdorferi (s.l.), Rickettsia spp., Babesia spp., Bartonella spp., Anaplasma spp., F. tularensis and "Ca. N. mikurensis". Of the 3465 samples, 175 were nymphs, four were larvae and the remaining were adults. The samples were selected to roughly represent both tick species (I. ricinus and I. persulcatus), and also the major collection areas, tick life stages and sex distribution of the whole collection [5]. Samples were further divided into three different distribution regions: 1, I. ricinus-dominated area (southern Finland); 2, sympatric area (central Finland) and 3, I. persulcatus-dominated area (northern Finland) (Fig. 1a).
DNA was extracted from the tick samples using NucleoSpin® RNA kits and RNA/DNA buffer sets (Thermo Fisher Scientific, Waltham, USA), following the kit protocols (RNA Kit: Rev. 16 May 2014 and RNA/ DNA buffer set: Rev. 08 May 2014). DNA extracts were stored at -20°C.

Tick species identification using genetic methods
Tick species, if unknown after morphological identification (n = 146), was determined (decisively 95 I. ricinus and 51 I. persulcatus) in a species-specific duplex real-time quantitative PCR (qPCR) assay as previously described [18]. Briefly, IXO-I2-F4 and IXO-I2-R4 primers targeting a fragment of Ixodes spp. internal transcribed spacer 2 (ITS2) gene were used to amplify genus-specific segments, and Ipe-I2-P4 and Iri-I2-P4 probes were used to match the ITS2 region for either tick species (I. persulcatus or I. ricinus, respectively; Table 1). DNA samples from I. ricinus and I. persulcatus confirmed by sequencing in an earlier study [27] were used as positive controls, and double-distilled water (ddH 2 O) was used instead of sample DNA as a negative control in each assay.

Detection of pathogens in tick samples
Bbsl-ospA-F and Bbsl-ospA-R primers, and a Bbsl-ospA-P probe (Table 1) amplifying a fragment of the outer surface protein A (ospA) gene, were used to detect B. burgdorferi (s.l.) DNA, as previously described [5].
Positive and negative controls [B. burgdorferi (sensu stricto) strain B31 ATCC 35210 and ddH 2 O, respectively] were included in all runs.
For screening Rickettsia, Anaplasma, "Ca. N. mikurensis", Babesia, F. tularensis and Bartonella, aliquots of original DNA samples were first pooled (10 samples per pool, 5 μl of each sample) due to low expected prevalence. For Anaplasma, Babesia and "Ca. N. mikurensis" screening, multiplex qPCR assays were first used. Briefly, qPCRs were first performed in 11 μl reaction volumes using the primers and probes displayed in Table 1. For Bartonella and Rickettsia screening, a duplex qPCR assay was first used in 8 μl volume using the primers and probes displayed in Table 1. Primers and probe (Table 1) targeting the 23 KDa gene were used to detect Francisella tularensis DNA. Individual samples from positive pools were analysed in 5 μl reaction volume. Samples were analysed in three replicate reactions carried out on 384-well plates. Three positive and negative control reactions were used in each assay. Samples were considered positive when successful amplification was detected in at least two replicate reactions. The thermal cycling profile used for all qPCR assays was: 95°C for 5 min, followed by 50 cycles of 95°C for 10 s and 60°C for 30 s. Thermal cycling was carried out at the Finnish Microarray and Sequencing Centre (FMSC, Turku, Finland) using QuantStudio™ 12 K Flex Real-Time PCR

Sequencing
The 5S-23S rDNA (rrfA-rrlB) intergenic spacer region (IGS) was sequenced from samples found positive for B. burgdorferi (s.l.) to identify the bacteria to genospecies level as previously described [28]. Samples found positive for Rickettsia spp., Anaplasma spp. and Babesia spp. by qPCR were sequenced using PCR primers displayed in Table 1. Detailed protocols of all PCRs are presented in Additional file 1: Methods.
Electrophoresis was carried out to confirm amplification success by running 1 μl of PCR product on 1.5% agarose gel. PCR products were purified by mixing 1 μl EXO I enzyme, 1 μl SAP enzyme and 8 μl of PCR product, after which the samples were first incubated for 5 min at 37°C and then heated for 10 min at 80°C. Purified samples were sent to Macrogen Inc. Europe (Amsterdam, Netherlands) for sequencing. The sequences were trimmed using Geneious 11.1.2 and run through BLAST (www.ncbi.nlm.nih.gov/BLAST/) and compared with reference sequences listed in the Gen-Bank nucleotide sequence database (www.ncbi.nlm.nih.gov/genbank/).

Data management and statistical analyses
Some of the tick samples received in the collection campaign were delivered in one shipment, indicating that the samples were from one sender. Such ticks were not necessarily independent observations, as generally assumed by basic statistical tests, because they likely shared either collection location or host animal/human individual, or both. It was often the case that ticks in one letter represented only one tick species and developmental stage. Therefore, we conducted formal statistical testing only for a couple of specific cases (see below), but otherwise tabulated frequencies and percentage values in a descriptive manner only.
With statistical tests we specifically analysed whether there was a difference in total infection rate (all pathogens combined) and B. burgdorferi (s.l.) or Rickettsia spp. prevalence between adult I. persulcatus and I. ricinus. Other pathogen groups were not statistically analysed due to low infection prevalence. In the subsequent phase, we tested whether there was a difference in pathogen prevalence between I. ricinus-and I. persulcatus-dominated areas and the sympatric area. Larvae and nymphs were ignored due to their relatively low sample sizes (n = 179) and low numbers of positive pathogen detections. We modelled the probability of a tick sample being positive for a pathogen by running a generalized estimating equation (GEE), a specific type of generalized linear mixed model for clustered observations, with binomial error distribution and logit link function. The shipment ID was set as a clustering factor, while the species and collection area of the tick were fixed explanatory factors in consecutive tests but never entered as fixed factors to the same model. Statistical testing was run with the IBM SPSS Statistics software v.23 (Armonk, NY, USA).

Overall pathogen prevalence
A total of 3465 tick samples, consisting of 2014 I. ricinus and 1451 I. persulcatus samples from the I. ricinus-and I. persulcatus-dominated areas and from the sympatric area, were analysed for the presence of pathogens, including the most common and some putative tick-borne pathogens circulating in ticks in Europe. Both the total infection rate and the diversity of different tick-borne pathogens was higher for I. ricinus (30.0%, five pathogen groups) than for I. persulcatus (24.0%, three pathogen groups) ( Table 2). The GEE model conducted for the adult samples indicated significantly higher probability of finding infected I. ricinus When investigating the number of infected ticks on the different latitudes in steps of one degree (approximately 110 km), the highest infection rates were found from the area below the latitude of 60°N (38.5%, n = 104), between the latitudes of 60°N and 61°N (33.7%, n = 808) and from the area between the latitudes of 65°N and 66°N (32.4%, n = 148) (Fig. 1a). In contrast, the lowest rate was from the area between the latitudes of 64°N and 65°N (21.3%, n = 127).

B. burgdorferi (s.l.) positivity
The most prevalent pathogen group was B. burgdorferi (s.l.), which was detected in 17.0% of the screened tick samples ( Table 2). The prevalence was 16.2% for I. ricinus and 18.1% for I. persulcatus. Out of the B. burgdorferi (s.l.) positive samples (n = 590), 394 were further identified to genospecies level ( The distribution maps drawn from the positive B. burgdorferi (s.l.) samples are shown in Fig. 1b-e. The highest prevalence of infected I. ricinus was observed in southern Finland (I. ricinus-dominated area), while the highest prevalence of infected I. persulcatus was observed in northern Finland (I. persulcatus-dominated area) ( Table 2).

Anaplasma spp. positivity
In total, Anaplasma spp. was detected in 0.6% of the screened DNA samples (  (Table 3). Amplicons were at least 99% identical to the 16S of A. phagocytophilum (GenBank: KY114936 from Croatia). Neither of the two positive I. persulcatus samples could be identified to species level due to a poor DNA sequence trace. The distribution map drawn from the positive A. phagocytophilum samples is shown in Fig. 1h.
Other pathogens "Ca. N. mikurensis" was detected in 0.5% of the screened DNA samples and the prevalence was 0.8% for I. ricinus (Table 2). "Ca. N. mikurensis" DNA was not detected in I. persulcatus samples. The distribution of the positive "Ca. N. mikurensis" was rather aggregated as most of the positive samples were collected from urbanized areas near the cities of Helsinki, Tampere and Turku in southern Finland (Fig. 1i). Babesia spp. was detected in 0.3% of all the screened DNA samples (Table 2), and the prevalence was 0.5% for I. ricinus. No infected I. persulcatus ticks were found. Nine positive samples were successfully sequenced, of which seven were identified as B. venatorum (77.8%) and two were identified as B. divergens (22.2%). Sequences were at least 99% identical to the reference sequences obtained from GenBank (B. divergens: U16370, origin not given, KY242392 from Poland; B. venatorum: KM289158 from Spain). The distribution of the positive Babesia samples was also rather aggregated (Fig. 1j). All positive samples from southern Finland were B. venatorum while B. divergens were found only from samples collected in central Finland.
Francisella tularensis and Bartonella spp. were not detected in either of the tick species.

Pathogens in larvae and nymphs
From 179 samples of juvenile life stages (nymph and larva), 33 (18.4%) were infected with B. burgdorferi (s.l.), Rickettsia or "Ca. N. mikurensis". Borrelia burgdorferi (s.l.) and Rickettsia prevalences in juvenile life stages (13.4 and 5.6%, respectively) were lower than in adults (17.3 and 11.1%, respectively). In contrast, the prevalence of "Ca. N. mikurensis" in juvenile life stages was higher than in adults (2.2 vs 0.4%, respectively). However, three of the four "Ca. N. mikurensis"-positive nymphs were collected from the same location, and positive samples were therefore strongly correlated. One I. ricinus larva carried B. garinii, but the rest of the positive B. burgdorferi (s.l.) samples from nymphs that could be identified by sequencing were B. afzelii. The co-infection prevalence for juvenile ticks was 2.6% (4/ 152) for I. ricinus and 3.7% (1/27) for I. persulcatus. Three of the four positive "Ca. N. mikurensis" samples were co-infected with B. afzelii.

Discussion
The distribution of I. ricinus and I. persulcatus in Finland is fairly unique [5,29,30]. The southern area of the country below the latitude 61°N is an area dominated by I. ricinus, while the area above 65°N is dominated by I. persulcatus. The belt between these latitudes is an area of sympatric occurrence of both species. In this study, we mapped the prevalence and distribution of an array of established and putative tick-borne pathogens in ticks in these three areas.
Overall, almost 30% of ticks were infected with at least one pathogen and 2% with more than one pathogen. Tick-borne pathogen diversity was higher in I. ricinus than in I. persulcatus. Of seven studied pathogen groups, five were detected in I. ricinus and three in I. persulcatus. A higher diversity of tick-borne pathogens in I. ricinus has also been observed in a previous study by Movila et al. [6], in which they investigated the differences of tick-borne microorganism communities in I. ricinus and I. persulcatus in distinct geographical regions of eastern Europe and European Russia. However, in our study, I. ricinus were also more frequently mono-and co-infected than I. persulcatus, which is contrary to the observations made by Movila et al. [6]. Infection rates can be influenced by many factors, e.g. life stage, sex, collection site of the ticks and by the season of tick collection. In Finland and neighboring countries, I. persulcatus adults have only one activity peak in April or May and are found to be questing only until July, while I. ricinus usually has two activity peaks during summer with the latter occurring during late August or September [5,31,32]. In our dataset, 85% of the I. persulcatus samples were collected by the end of May, while 85% of the I. ricinus samples were not collected until the end of July. Of all the analysed ticks, 95% were adults and less than one third were males, with these ratios equal in both species. Some of the studied pathogens in our study were analysed with single primers and some were multiplexed. Moreover, some were analysed individually while some were first pooled. These different methods could potentially cause a small bias on the prevalence of the studied pathogens. Nevertheless, our pathogen prevalence results mostly correspond with the prevalence results observed in the neighboring countries, suggesting that sample pooling and PCR multiplexing did not cause any major differences in the observed prevalences. However, slightly fewer positive findings could have been due to our "majority rules" approach when analyzing PCR results. Samples were considered positive when successful amplification was detected in at least two replicate reactions. A positive sample in one out of three replicates could suggest a low level of DNA close to the detection level rather than contamination of a sample. However, there were only a few samples that had only one positive replicate.
Not surprisingly, the most prevalent tick-borne pathogen was B. burgdorferi (s.l.) (17%), the causative agent of LB, which is the most prevalent tick-borne disease in Finland (around 6000-7000 cases yearly;~120 cases per 100,000 inhabitants) [33]. The observed prevalence in our study is in accordance with the average prevalence rates found in Europe (17.8% for adult I. ricinus ticks) [34]. Even in the sympatric area, prevalence was significantly higher in I. persulcatus than I. ricinus (17.2 and 12.4%, respectively). A higher prevalence of B. burgdorferi (s.l.) in I. persulcatus than in I. ricinus ticks has also been observed in previous studies conducted in sympatric regions [30,35,36].
The most common B. burgdorferi (s.l.) genospecies detected in our study were B. garinii and B. afzelii. In comparison to reported B. garinii prevalence in I. ricinus ticks in Europe, a higher prevalence was observed in the present study [34]. Interestingly, a particularly high B. garinii prevalence (62.6% from Borrelia-positive ticks) was observed in I. persulcatus ticks. Migratory songbirds such as Turdus species are known as common B. garinii reservoir hosts while rodents are known as the main B. afzelii reservoir hosts [37,38]. There can be differences in the occurrence of these reservoir species between countries, as well as fluctuations in their abundance (especially in voles) between different years [39], which might affect the proportions of Borrelia genospecies. Moreover, the activity peak of I. persulcatus in Finland co-occurs with the spring migration of Turdus spp. birds which might affect the higher B. garinii prevalence observed in I. persulcatus. Borrelia afzelii was observed relatively more often in I. ricinus than in I. persulcatus  Tick samples that were not categorized into collection areas due to inaccurate collection information provided by citizens samples. The same observation was made in the study by Movila et al. [6]. While the evidence for human pathogenicity of B. valaisiana is poor, B. afzelii, B. garinii and B. burgdorferi (s.s.) are the genospecies that commonly infect people [4,40]. These genospecies are also often associated with different clinical manifestations. In Europe, B. garinii is the main cause of Lyme neuroborreliosis, while B. afzelii is mostly associated with skin manifestations [41,42]. Among the identified genospecies, B. valaisiana was detected in 6.4% of ticks, which is in accordance with the findings in the neighboring countries (6% in Norway, Sweden and Estonia) [30,43,44]. Even though, B. valaisiana is not often detected in I. persulcatus, it has been shown that in the sympatric areas of I. ricinus and I. persulcatus, B. valaisiana may exchange vectors and can also be found in I. persulcatus [30,45]. In our study, B. valaisiana prevalence among B. burgdorferi (s.l.)-positive I. persulcatus samples (2.5%) was similar to the prevalence observed in Estonia, where these two tick species live in sympatric areas as well [30]. B. burgdorferi (s.s.) was not detected in I. persulcatus in this study, even though it has previously been found in I. persulcatus in a coastal Finnish region of the Gulf of Bothnia around the city of Kokkola [46]. Our study also revealed the occurrence of some less-known pathogens present in ticks in Finland. In contrast to B. burgdorferi (s.l.), all such pathogens had a higher prevalence in I. ricinus than in I. persulcatus samples. The most prevalent pathogen after B. burgdorferi (s.l.) spirochetes was Rickettsia spp. (10.8%). The observed prevalence is slightly higher than the reported prevalence in the neighboring country, Estonia (5.1%) [47]. Rickettsial DNA was more frequently detected in I. ricinus (13.9%) than in I. persulcatus (6.5%). The majority of the positive samples were identified as R. helvetica (90.9%). Ixodes ricinus is regarded as the main vector of R. helvetica, while the role of I. persulcatus is less studied. Rickettsia helvetica has been detected in I. persulcatus ticks before but with lower prevalence than in I. ricinus [6]. In our study, the prevalence of R. helvetica in I. ricinus (9.9%) was over four times higher than in I. persulcatus (2.2%), even though R. helvetica did not have a significant difference in infection prevalence between the collection areas. In contrast, "Ca. R. tarasevichae" was detected almost exclusively in I. persulcatus (18/19 positive samples). This is the first report of "Ca. R. tarasevichae" in Finland, and also the most western report. Moreover, to our knowledge this is the first report of the pathogen in I. ricinus. Interestingly, the positive I. ricinus sample is from southern Finland, the area of parapatric occurrence of I. ricinus. "Ca. R. tarasevichae" has previously been reported from neighboring countries Russia and Estonia [47,48], where I. persulcatus is a common tick species. A minority (1.2%) of the Rickettsia-positive samples were identified as R. monacensis. This pathogen was first detected in I. ricinus ticks in 2013-2014 in southwestern Finland with similar prevalence to this study [18]. Patient cases from Spain and China suggest that "Ca. R. tarasevichae" and R. monacensis are both capable of human infection [49,50].
Anaplasma spp. was the most abundant pathogen of the family Anaplasmataceae and almost all of the positive samples were detected in I. ricinus. The prevalence for I. ricinus was 1.0%, which is lower than the prevalence detected previously in questing I. ricinus adults in southwestern Finland (9.2%) [17]. However, reports from all over Europe have observed pronounced differences in prevalence among countries, study localities and tick life stages, ranging from 0 to 67% [9,51,52]. Anaplasma phagocytophilum is frequently detected in I. ricinus in Europe and it was the only Anaplasma species we detected. The prevalence for I. persulcatus was only 0.1% and neither of the two positive samples could be identified to species level by sequencing. Anaplasma phagocytophilum is the agent of human granulocytic anaplasmosis (HGA), and confirmed cases of HGA have been reported since 1997 in Europe [51,[53][54][55]]. Another human pathogenic member of the family Anaplasmataceae, "Ca. N. mikurensis" was detected for the first time in ticks from mainland Finland, and only in I. ricinus (overall prevalence of 0.5%). The observed prevalence was lower than in previous observations made in Europe [56,57]. The presence of "Ca. N. mikurensis" in I. ricinus ticks was previously analysed in ticks collected in 2013-2014 in southwestern Finland, but no positivity was detected [18]. However, our recent, still unpublished tick sampling data (Sormunen et al., in press), reveal the occurrence of this pathogen in southwestern Finland since 2015, agreeing with observations of the current study. "Ca. N. mikurensis" was also detected in I. ricinus ticks, collected in the years 2006-2013 in Estonia, with prevalence rates ranging from 1 to 9.1% [58].
Previous studies about Babesia spp. prevalence in Finnish ticks are scarce. The prevalence observed for Babesia spp. in I. ricinus samples in the current study (0.5%) corresponds to results reported from neighboring countries [30,[59][60][61]. Babesia venatorum was the most prevalent species in our study. In Europe and China this pathogen has been involved in the several documented cases of human babesiosis [62,63]. Babesia microti was not detected, even though it is known to commonly infect rodents in Finland and was previously found in I. persulcatus tick in the Kokkola coastal region [46,64]. The absence of B. microti in our study is likely due to low prevalence in Finnish ticks. B. divergens was detected in two samples. However, B. divergens is genetically very similar to B. capreoli. Since sequences should be identical to GenBank reference sequence U16370 to be regarded as B. divergens [65] and one of our positive samples were only 99% identical to this reference, we cannot be entirely sure whether this is truly B. divergens or B. capreoli. In Finland, a single case of fatal babesiosis in 2004 was caused by B. divergens and believed to be transmitted by a tick bite [19].
Neither Bartonella spp. nor Francisella tularensis were found in the ticks, although some Bartonella species have recently been reported from mammals and arthropods in Finland [66,67]. The relevance of ticks as vectors for human bartonellosis is not yet verified. Bartonella have been reported in one I. persulcatus tick in Estonia (with a total prevalence of 0.2%) [6], while the prevalence in questing I. ricinus ticks in Europe has varied up to 48.2% in nymphs and 12% in adult ticks [68]. Although hard ticks are important vectors of F. tularensis in North America and central Europe [69][70][71], their role in the transmission of F. tularensis in northern Europe is poorly understood. In Fennoscandia, the primary route of human F. tularensis infection is probably through mosquito bites [69,[72][73][74].
In our study, about 7% of infected ticks were found to be co-infected. Most of the co-infections were between B. burgdorferi (s.l.) and Rickettsia spp., as expected due to their high prevalence. Interestingly, a particularly high number of co-infections were observed in positive "Ca. N. mikurensis" samples, of which over half (6/11) were co-infected with B. burgdorferi (s.l.) ( Table 4). All of these co-infected B. burgdorferi (s.l.) samples were identified as B. afzelii and half of them were in nymph samples. "Ca. N. mikurensis" and B. afzelii share the same natural hosts (voles and mice) and co-infection with these pathogens has recently been shown to occur in nymphs more often than expected under random co-occurrence [57,75,76], as a result of larvae receiving both pathogens by feeding on co-infected hosts. The amount of different co-infection combinations was higher in I. ricinus (six) than in I. persulcatus (only one), partly due to the lower pathogen diversity and infection rates in I. persulcatus in our study. The co-infection rate in I. ricinus was 2.4%, which roughly corresponds with the results by Movila et al. (3.4%) [6]. However, in the same study, the co-infection rate in I. persulcatus samples was much higher than in ours (4.3 vs 0.8%). This is likely due to lower prevalence of Rickettsia infected I. persulcatus in our study. We examined co-infections only at genus level, and did not consider double infections among different pathogen species, e.g. B. afzelii and B. garinii. On average, 68% of the microbial sequences could be identified to species level. The rest of the sequence results were of poor quality (blurry or poorly resolved signal peaks in the chromatogram), which could have resulted from impurities in the DNA samples.
Species of B. burgdorferi (s.l.) and Rickettsia were detected beyond the latitude of 65°N, while the members of family Anaplasmataceae and Babesia spp. were detected only south of 65°N, perhaps due to lower abundance of I. ricinus in higher latitudes. The prevalence of infected ticks did not correlate directly with latitude, since the highest prevalence of infected ticks were found from the area below 61°N and from the area between 65°N and 66°N. Other ecological factors of the area, such as rainfall, distance to coast or inland waters, vegetation and variety of host animals, may have a bigger influence on pathogen diversity or proportion of the infected ticks [77]. According to previous observations, pathogen prevalence is also expected to correlate with the density of questing ticks at the collection sites [78][79][80]. Since our tick samples were gathered by citizens from all over Finland in a positive correlation to human population density, we cannot straightforwardly consider tick densities in our analyses. Interestingly, the highest infection rates were observed in areas of parapatric occurrence of one species, I. ricinus in southern and I. persulcatus in northern Finland. The assumption that in the zone of sympatry, the pathogen prevalence in one vector species would either increase or decrease by the influence of other closely related species, was not supported. Pathogen prevalence was lower for both species in the zone of their sympatry, indicating that other environmental factors might explain the lower prevalence in that area. When comparing infection rates between the collection areas by combining all the analysed ticks (excluding the information of tick species), the highest infection rate was observed in southern Finland, but infection rates did not differ between central and northern Finland (I. ricinus-dominated area: 32.6%; sympatric area: 25.4%; and I. persulcatus-dominated area: 25.9%). It remains unclear whether the higher prevalence in southern Finland is related to I. ricinus dominance or other environmental factors. Rickettsia spp. was the only pathogen group with a descending trend in infection prevalence towards the higher latitudes (I. ricinus-dominated area: 12.6%; sympatric area: 10.7%; and I. persulcatus-dominated area: 6.0%). However, the differences between the collection areas were not significant.

Conclusions
The unique distribution of I. ricinus and I. persulcatus in Finland allowed us to compare the pathogen distributions in parapatric and sympatric areas of tick occurrence. The highest infection rate was observed in southern Finland, but infection rates did not differ between sympatric and I. persulcatus-dominated areas. However, the highest infection rates for both species were observed in areas of their dominance, either in southern or northern Finland. Furthermore, our comprehensive subset of Finnish ticks revealed the presence of rare and potentially pathogenic bacteria, such as "Ca. N. mikurensis" and "Ca. R. tarasevichae", for the first time in mainland Finland. Although only a few human infections caused by organisms other than B. burgdorferi (s.l.) and TBEV have been reported in Finland so far, and it is known that infections due to these emerging pathogens can be either asymptomatic or mild [81,82], the risk of these tick-borne pathogens for public health should not be neglected.

Acknowledgements
We thank Tuula Rantasalo (Institute of Biomedicine, University of Turku) for helping with the pathogen screenings, and Manoj Fonville (RIVM) with excellent technical support in the determination of Borrelia genospecies. We thank the company Biotop for professional help with lab consumables. We thank the Finnish Functional Genomics Centre (FFGC, Turku, Finland) for the laboratory services. We also thank the two anonymous for their constructive comments.

Funding
This work was supported by the Jane and Aatos Erkko Foundation; Sakari Alhopuro; Jenny and Antti Wihuri Foundation; Pfizer Inc. (Finland); Maj and Tor Nessling Foundation; The Varsinais-Suomi Regional Fund of the Finnish Cultural Foundation; and the Academy of Finland (grants 292701).

Availability of data and materials
The data supporting the conclusions of this article are included within the article and its additional files. The datasets used and analysed during the present study are available from the corresponding author upon reasonable request.
Authors' contributions RP, ML, TK, JJS, EJV, IV, JHä, IES and KR planned and conducted the collection campaign of the ticks and ML, TK, JHy and EJV designed the current study on tick-borne pathogens. ML performed morphological identification of tick samples and conducted the main laboratory analyses, interpreted the data and drafted the manuscript. TK planned statistical tests and provided major comments on the manuscript. EF, JJS and SM performed DNA extractions and pathogen screenings. AP performed morphological identification of tick samples and pathogen screenings. RP lead the collection campaign and performed morphological identification of tick samples. KR determined coordinates of the samples. JHä, JHy, ML, IES and IV acquired funding. HS analysed Borrelia genospecies and provided IGS PCR method. EJV and JHy provided major comments on the manuscript. All authors read and approved the final manuscript.
Ethics approval and consent to participate Not applicable.

Consent for publication
Not applicable.