Tick range expansion to higher elevations: does Borrelia burgdorferi sensu lato facilitate the colonisation of marginal habitats?

Parasites can alter host and vector phenotype and thereby affect ecological processes in natural populations. Laboratory studies have suggested that Borrelia burgdorferi sensu lato, the causative agent of human Lyme borreliosis, may induce physiological and behavioural alterations in its main tick vector in Europe, Ixodes ricinus, which increase the tick’s mobility and survival under challenging conditions. These phenotypic alterations may allow I. ricinus to colonise marginal habitats (‘facilitation hypothesis’), thereby fuelling the ongoing range expansion of I. ricinus towards higher elevations and latitudes induced by climate change. To explore the potential for such an effect under natural conditions, we studied the prevalence of B. burgdorferi s.l. in questing I. ricinus and its variation with elevation in the Swiss Alps. We screened for B. burgdorferi s.l. infection in questing nymphs of I. ricinus (N = 411) from 15 sites between 528 and 1774 m.a.s.l to test if B. burgdorferi s.l. prevalence is higher at high elevations (i.e. in marginal habitats). Opposite of what is predicted under the facilitation hypothesis, we found that B. burgdorferi s.l. prevalence in I. ricinus nymphs decreased with increasing elevation and that Borrelia prevalence was 12.6% lower in I. ricinus nymphs collected at the range margin compared to nymphs in the core range. But there was no association between Borrelia prevalence and elevation within the core range of I. ricinus. Therefore the observed pattern was more consistent with a sudden decrease in Borrelia prevalence above a certain elevation, rather than a gradual decline with increasing elevation across the entire tick range. In conclusion, we found no evidence that B. burgdorferi s.l.-induced alterations of I. ricinus phenotype observed in laboratory studies facilitate the colonisation of marginal habitats in the wild. Rather, ticks in marginal habitats are substantially less likely to harbour the pathogen. These findings have implications for a better understanding of eco-evolutionary processes in natural host-parasite systems, as well as the assessment of Lyme borreliosis risk in regions where I. ricinus is newly emerging.

manipulations [7,8]. Phenotypic changes in infected hosts and vectors can also be a side effect of infection that is of no adaptive value for the parasite.
Surprisingly little is known about the ecological consequences of parasite-induced phenotypic alterations in hosts or vectors. Recent studies have suggested that parasite-induced changes may have profound effects on host population dynamics, trophic niche specialisation as well as interactions between competitors and predators at all trophic levels [9][10][11]. Consequently, parasites can have a strong ecological impact that becomes especially apparent when they are missing or introduced, e.g. during biological invasions [11,12]. Invasive host populations often harbour fewer parasites than populations in their native range, which may enhance their population growth and competitive ability (i.e. 'enemy release hypothesis' [13]). Parasites may be 'lost' because of low host densities and founder effects during the three phases of an invasion: introduction, establishment and range expansion [14]. However, empirical studies investigating the role of parasites during range expansion from core populations, i.e. without bottlenecks occurring during the introduction and establishment phases, are still rare [15], and both increasing and decreasing parasite loads at the host range margins have been reported (e.g. [16][17][18]).
The distribution of Ixodes spp. ticks is strongly influenced by abiotic factors, such as temperature and humidity, and there is accumulating evidence that in Europe and North America Ixodes spp. ticks have been expanding northwards as well as towards higher elevations during the last decades due to climate warming (reviewed in [19,20]). In Europe, I. ricinus is the main vector of the bacterium Borrelia burgdorferi sensu lato, the causative agent of human Lyme borreliosis. Borrelia burgdorferi s.l. forms a complex comprising of at least 18 genospecies [21]. A diverse host community, which includes rodents, insectivores, birds and reptiles, acts as reservoir hosts for B. burgdorferi s.l.. In Europe, B. afzelii and B. garinii are the most common genospecies [22]. Borrelia afzelii is a rodent specialist, infecting mice (Apodemus sylvaticus and A. flavicollis) and voles (Myodes glareolus and Microtus agrestis), but also shrews (Sorex araneus and S. minutus; [23]), whereas B. garinii is a bird specialist, infecting a range of bird species such as Turdus spp., Sturnus vulgaris, Sylvia spp. and Parus major [21,24].
Ixodes ricinus is a three host tick. Each life stage (i.e. larva, nymph and adult) seeks a host, feeds, and drops off to develop in the natural environment. Although I. ricinus is a host generalist, the host distribution differs across tick stages [25,26]. The host selection is, among other things, driven by habitat type, host behaviour, microclimate conditions, and questing height of the different stages [26]. Adults quest higher in the vegetation than nymphs, and nymphs higher than larvae. Larvae and nymphs feed on vertebrate hosts of all size classes but particularly on rodents and birds, while adult females usually feed on larger mammals [26]. Unlike hematophagous insects, ticks typically feed for several days on their host, and the transmission of tick-borne pathogens such as B. burgdorferi s.l. to the host may start several hours after tick attachment [26].
Various studies have reported that bacteria and viruses may induce phenotypic alterations in Ixodes spp. [26,27]. Specifically, laboratory studies have suggested that B. burgdorferi s.l.-infected I. ricinus are bigger, walk faster and longer, quest longer for a host, have higher energy reserves (for the same body mass) and survive better under challenging thermohygrometric conditions than uninfected I. ricinus. Yet to date, studies investigating such vector manipulation by pathogens remain scarce and results are heterogeneous ( [26,27], and references therein). These alterations of tick behaviour and physiology may benefit pathogen transmission by extending vector lifespan, or by increasing host-finding efficacy [26]. Importantly, both may facilitate the colonisation of marginal habitats by ticks, and lead to particularly high prevalence of B. burgdorferi s.l. in I. ricinus populations at range margins [28]. Particularly, higher fat contents in Borrelia-infected nymphs may provide an advantage under colder or more variable temperatures, such as at higher elevations, although laboratory data were not conclusive [29]. Furthermore, increased host-finding efficacy may facilitate the establishment of tick populations in less diverse or more fragmented host communities. Finally by modulating questing behaviour, such as adjusting questing height or moving to shade, Borrelia infection may shape the interaction between ticks and host species [26,27]. Specifically, questing higher up in the vegetation may allow ticks to attach to larger animals. Since larger animals have larger home ranges [30], questing at higher heights may allow ticks to disperse further. All these possible mechanisms may result in higher B. burgdorferi s.l. prevalences in I. ricinus populations at range margins. Assessing the role of B. burgdorferi s.l. in influencing vector range expansion to marginal habitats is crucial for our understanding of eco-evolutionary processes in natural host-parasite systems, but also for quantifying public health threats in regions where ticks are newly emerging. Yet, to our knowledge no study has tested this hypothesis in I. ricinus populations under natural conditions to date.
Although exclusion experiments (i.e. where a taxon is excluded) are the gold standard to assess the ecological role of parasites, such manipulations remain an enormous practical and logistic challenge in I. ricinus and its pathogens. Moreover, experiments in the wild testing whether B. burgdorferi s.l. enhances survival and reproduction of I. ricinus under harsh conditions (i.e. physiological alterations) would only provide a partial picture because B. burgdorferi s.l. might affect dispersal abilities of I. ricinus through its ability to find a host (i.e. behavioural alterations). Because of these practical challenges, we used a correlational approach in naturally infected tick populations in the Swiss Alps to evaluate the potential of B. burgdorferi s.l. to facilitate the colonisation of marginal habitats by I. ricinus. Previous studies monitored B. burgdorferi s.l. and I. ricinus within their core range [31][32][33][34][35][36], at the range margins without a direct comparison with core populations [37] or when B. burgdorferi s.l. prevalence was low [38][39][40] (but see [41]). In our study we monitored tick populations at different elevations, including core range and range margin populations, allowing for a direct comparison of Borrelia prevalence across the tick distribution range. We predict that if B. burgdorferi s.l. alters the phenotype of I. ricinus, thereby making them better colonizers of unfavourable habitats, B. burgdorferi s.l. prevalence will be disproportionally high at the tick range margin (i.e. at high elevations).

Discussion
Laboratory experiments have suggested that B. burgdorferi s.l. may increase tick tolerance to desiccation and boost host-finding efficacy ( [26,27], and references therein). Such vector manipulations by the pathogen may facilitate the range expansion of ticks to marginal habitats, such as higher elevations or higher latitudes [28]. Here we used a correlational approach to investigate whether patterns of B. burgdorferi s.l. prevalence in questing I. ricinus nymphs, and its variation with elevation, provide support for this hypothesis.
Overall, questing I. ricinus abundance was high and B. burgdorferi s.l. infection common in the study area: I. ricinus were found at 14 out of 15 sites, and B. burgdorferi s.l. infection in questing I. ricinus nymphs was detected at 12 sites. Above 1400 m.a.s.l., however, I. ricinus abundance strongly decreased and no I. ricinus were found above 1700 m.a.s.l.. B. burgdorferi s.l. prevalence in questing I. ricinus nymphs across sites was similar to previous reports from central Europe (e.g. [33], and references therein, [42]), but variation across sites was high (0-47% compared to 1-20% across Europe according to [22]), and it exceeded 40% at two sites of contrasting elevations (630 and 1123 m.a.s.l.).
As expected, the most common B. burgdorferi s.l. genospecies was B. afzelii, accounting for 83.2% of I. ricinus infections in our study. This is higher than what has been reported previously (7-68% of I. ricinus infections caused by B. afzelii; [22]) across Europe but similar to what has recently been reported from the Italian Alps ( [41] where 73.1% of the nymphs were infected by B. afzelii). Differences across studies and locations are likely due to differences in the composition of the local host community, i.e.  nii and B. valaisiana are bird specialists, and B. burgdorferi s.s. is a generalist that can infect birds and rodents (reviewed in [21]). Alternatively, strong bottlenecks within sites could lead to the local dominance of specific genospecies that may vary over time (e.g. [43,44]).
We found no evidence that B. burgdorferi s.l. prevalence is higher at the range margins of I. ricinus (i.e. at higher elevations), indicating that B. burgdorferi s.l.-mediated modifications of I. ricinus physiology and behaviour observed in the laboratory ( [26,27], and references therein) plays a minor role in the ongoing colonisation process of marginal habitats in the wild. Rather, I. ricinus nymphs at the range margin had substantially lower B. burgdorferi s.l. prevalence. This finding is in line with what was found in the Italian Alps, where higher Borrelia prevalence was observed at locations < 1400 m.a.s.l compared to > 1400 m.a.s.l (note that ticks were found up to > 1800 m.a.s.l at these study sites) [41]. Our analytical approach allows us to describe changes in Borrelia prevalence with elevation in more detail and highlights that the patterns we observe are more consistent with a sudden decrease above a certain elevation, with reduced Borrelia prevalence specifically at the range margin, rather than a gradual decline with increasing elevation across the entire tick range. Indeed, analyses within the core range of ticks revealed no association between Borrelia prevalence and elevation.
The observed patterns may be due to different nonexclusive processes. The colder microclimate at high elevations may decrease the survival of B. burgdorferi s.l. in I. ricinus and slow down its multiplication and transmission. Effects of meteorological temperature on pathogen development (i.e. the development from the ingested infectious stage to the stage in the salivary glands) and transmission efficiency have been observed in malariamosquito [45,46] and trypanosome-tsetse fly systems (references in [47]). Although protozoan development is not directly comparable to bacterial multiplication, many bacteria similarly utilize the environmental temperature as a signal to determine their location and to regulate expression of a large set of proteins necessary for survival, multiplication, migration and transmission [48,49]. For example, the expressions of Outer Surface Protein A and C of B. burgdorferi s.l., which are involved in its survival in the tick midgut and its dissemination into the tick salivary glands or into its vertebrate host, appear to be mediated by differences in ambient temperature [50]. To our knowledge, no study evaluated the effect of meteorological temperature on the expression of vector-specific or host-specific proteins, or the expression of virulence genes in B. burgdorferi s.l.. Although temperature and other abiotic factors are often involved to explain the lower pathogen prevalence at higher elevations [51], the lack of a gradual decline in Borrelia prevalence with increasing elevation within the core range of ticks makes it unlikely that such factors are the sole explanation for the low Borrelia prevalence at range margins.
Alternatively, the pattern may also be affected by other processes such as changes in host community, host species abundance, habitat quality or population structure of hosts/vectors in marginal habitats. Indeed, Patot et al. (2010) found no effect of temperature on filamentous virus transmission efficiency, but observed a clear relationship between virus prevalence and the density of its host, a parasitoid wasp. Similarly, B. burgdorferi s.l. transmission may be lower in I. ricinus populations at the range margin that are more fragmented [52], less dense [17] or/and less genetically diverse [53] than populations at the core range. Similar processes may act on the reservoir host populations of B. burgdorferi s.l. (e.g. rodents and birds). However detailed knowledge about the composition of the host community, host abundance, I. ricinus infestation of these host populations as well as B. burgdorferi s.l. prevalence would be necessary to link host abundance and B. burgdorferi s.l. prevalence in questing I. ricinus.
Finally, because edge populations are usually founded by only few individuals, stochastic processes may lead to the loss of pathogens in ticks at the range margins. Thus, the combination of stochastic events and host and vector population structure that hinder B. burgdorferi s.l. transmission may lead to local extinctions of B. burgdorferi s.l. at the range edges (e.g. [18]). It is often assume that vertebrate hosts dispersing ticks on wide ranges are mainly birds and wild cervids (e.g. red deers and chamois). One may argue that because the latter are incompetent hosts for Borrelia spp., it may dilute infection prevalence in ticks at the range margins. The lower Borrelia infection prevalence at the range margins supports this hypothesis. However a sudden drop in Borrelia prevalence above a certain elevation would also be consistent with a higher mortality of ticks driven by abiotic factors. Moreover the dominance of B. afzelii in this study, and in the Italian Alps [41], suggests that small rodents may play a particularly important role in the maintenance of Borrelia infection along elevational clines in the Alps. Though, further investigations are needed to evaluate these hypotheses in detail.
Our study is correlational and has limitations. First, fewer ticks were collected at the range margin, which will decrease the precision of the prevalence estimate at these sites. Second, a time-based approach (which is less reliable than an area-based approach) was used to estimate questing tick abundance because of the rugged mountain terrain of the study area. Third, vegetation structure varies substantially across the elevation range and could affect tick detection despite a consistent sampling effort. However, the decrease of I. ricinus abundance above 1150 m.a.s.l. was observed in different life-stages (Lemoine pers. obs.), which have different questing behaviours and biological needs. Therefore, we are confident that the decrease of I. ricinus abundance above 1150 m.a.s.l. is not a methodological artefact but describes an ecological pattern. Fourth, our approach does not allow to disentangle between range-expanding ticks bringing the pathogen with them (i.e., dispersal to a novel habitat) and the enzootic cycle being less likely to become established (i.e., establishment in a novel habitat, for example infecting less dense host populations). Fifth, whereas numerous studies in Europe and North America have documented a range expansion of ticks to higher elevations and latitudes due to climate change (reviewed in [19]), we do not have the long-term monitoring data to directly demonstrate this range expansion at our study sites. Finally, the presence of questing I. ricinus ticks in the vegetation does not prove that I. ricinus populations are established locally, or will survive through the winter. High elevation habitats might thus represent a 'ragged edge' rather than an expansion front. Ultimately, an experimental approach, in which the B. burgdorferi s.l. infection status of I. ricinus is manipulated and the colonisation of marginal habitats is monitored would be required to test conclusively whether B. burgdorferi s.l. facilitates tick range expansion to marginal habitats in the wild. But as stated previously, this would be an enormous practical challenge.

Conclusions
Overall, we found no evidence that B. burgdorferi s.l.induced changes in I. ricinus behaviour or physiology that facilitate I. ricinus range expansion to higher elevations in the Swiss Alps. Rather, questing I. ricinus in marginal habitats are less likely to carry B. burgdorferi s.l.. These findings show that when bitten by a tick, the risk of human B. burgdorferi s.l. infection is lower, rather than higher, in regions where I. ricinus is newly emerging. Low B. burgdorferi s.l. prevalence at I. ricinus range margins may enhance population growth and competitive ability of hosts and vectors. Less infected, hosts may, for example, invest differently in immunity and reproduction than hosts in core populations [15,54], which can affect hostparasite interactions when the parasite finally invades host populations at the range margins [55]. A better understanding of eco-evolutionary processes between pathogens, vectors and hosts at range margins, and their effect on pathogen life-history and virulence evolution, will therefore be a fruitful next step (e.g. [55]), and will contribute to a better prediction of zoonotic disease risks in regions where vectors and pathogens are newly emerging due to climate change.

Study area
We sampled questing I. ricinus at fifteen sites across an elevational range from 528 to 1774 m.a.s.l. in the Swiss Alps (Kanton Graubünden, Fig. 1, Table 1). Each site was visited three times, once in June, July and August 2014 to obtain a measure of questing tick abundance per site. During each dragging event, four people dragged a white blanket (1 m × 1 m each) slowly over the ground vegetation during 15 min each. During these 15 min, each collector regularly checked for ticks on their blanket and transferred them to Eppendorf tubes containing 95% ethanol. The count of each collector was pooled to provide an estimate of questing tick abundance per dragging event per hour ( Table 2). Although an area-based estimation is often more accurate than a time-based approach, it was not possible to apply such a method in our study because of the rugged mountain terrain of the Swiss Alps. Tick life stage and species were verified with a dissection microscope in the laboratory on the basis of morphologic features following [56].
We randomly selected 34-35 questing I. ricinus nymphs per site for B. burgdorferi s.l. detection. When fewer than 35 nymphs were collected at a site (Table 1), additional questing nymphs were collected for Borrelia detection one day after estimating questing tick abundance. All collected nymphs were analysed when fewer than 35 nymphs were available (Table 1). DNA was extracted using the 'HotShot' method [57,58] with slight modifications. Each nymph was incubated with 120 μL of alkaline lysis reagent (25 mM NaOH, 0.2 mM EDTA) and a metal bead (2 mm in diameter) at 95 °C for 15 min. After homogenisation with a Retsch TissueLyser (Haan, Germany) for 2 min, the samples were incubated at 95 °C for 15 min. After cooling, 60 μL 40 mM Tris-HCl was added to each tube.
To assess the B. burgdorferi s.l. infection status of I. ricinus nymphs, we used a combination of two complimentary approaches. First, we specifically targeted B. afzelii, the most common Borrelia genospecies at the study sites [54], using a highly sensitive quantitative real-time PCR (qPCR) assay using the B. afzelii-specific primers In a second step, all samples that were found to be B. afzelii-negative using the qPCR approach described above were screened for the presence of other B. burgdorferi s.l. genospecies using 16S-LD primers designed to target all Borrelia genospecies associated with Lyme borreliosis [59]. Amplifications were performed in a total volume of 10 μL containing 0.2 μL JumpStart Taq DNA Polymerase (Sigma-Aldrich), 0.5 μL of each primer (300 nM; LD-F: 5′-ATG CAC ACT TGG TGT TAA CTA-3′ and LD-R: 5′-GAC TTA TCA CCG GCA GTC TTA-3′) and 2 μL of DNA template. The PCR amplification protocol consisted of an initial denaturation step at 94 °C for 1 min, followed by 42 cycles of denaturation at 94 °C for 30 s, annealing at 50 °C for 30 s, and extension at 72 °C for 90 s, with a final elongation step at 72 °C for 10 min. PCR products were visualized under UV light on 1% agarose gels that were stained with SYBR ® Safe DNA gel stain (Thermo Scientific). The length of the amplicon was 351 bp. We used DNA from six reference genospecies B31 and Pbre (B. burgdorferi s.s), Phei (B. garinii), Pko and PVPM (B. afzelii) and Pbi (B. bavariensis) to verify that they are detected with the 16S approach. All references were successfully amplified.
Each run included B. afzelii-positive samples identified with the qPCR (see above) as controls. All of these samples were found to be B. burgdorferi s.l.-positive using the 16S method. The repeatability based on 39 samples (12 negative and 27 positive samples) was 85%. Out of the six samples that were not repeatable with the 16S method, three were previously quantified with the qPCR method and were found to have very low B. afzelii DNA concentrations (i.e. more than 30 cycles on the qPCR) suggesting that the sensitivity of the 16S approach decreases at low infection intensities. During the development of the methods, some fragments were sequenced and all were found to be B. burgdorferi s.l. suggesting that both methods were specific to the Lyme borreliosis causing genospecies.

Statistical analyses
Questing I. ricinus abundance ( Table 2) was analysed using a generalized mixed effect model with the negative binomial distribution (with the parameterization σ = μ(1 + μ/k)) and site ID was included as a random effect. Preliminary analyses showed that the negative binomial distribution (with the parameterization σ = μ(1 + μ/k)) described the abundance of nymphs (or nymphs and adults combined) best among the following distributions: Poisson, negative binomial (with the parameterization σ = μ(1 + μ/k) or σ = φμ) with or without zero-inflation (results not shown). We first tested for an association between tick abundance across life stages (nymphs and adults) and elevation or season. To this end, we included elevation and elevation 2 in interaction with tick stage, or in interaction with sampling session as fixed effects (i.e. as explanatory variables to quantify an average effect, β) in the mixed model. Finally, we defined the range margin based on questing tick abundance (see Results) and included range type (range core (C) vs range margin (M)) of each site as a fixed effect.
Borrelia infection (i.e. B. burgdorferi s.l. or B. afzelii infection) in I. ricinus nymphs ( Table 2) was analysed using a generalized mixed effect model with a binomial error structure. The number of infected and non-infected nymphs per dragging event was used as the dependant variable and site ID was included as a random effect. In a first step, we tested for an association between B. burgdorferi s.l. or B. afzelii prevalence and elevation by including elevation and elevation 2 as fixed effects (i.e. as explanatory variables to quantify an average effect, β) in the mixed model. In a second step, we tested if B. burgdorferi s.l. or B. afzelii prevalence was higher at the range margin of I. ricinus (i.e. at high elevations). To this end, we compared the B. burgdorferi s.l. or B. afzelii prevalence in I. ricinus at the range margin (sites > 1400 m.a.s.l.) to B. burgdorferi s.l. or B. afzelii prevalence in the core range of I. ricinus (sites < 1150 m.a.s.l.) by including range type (range core (C) vs range margin (M)) of each site as a fixed effect.
Finally, we focused on the patterns within the core range of I. ricinus by investigating the association between B. burgdorferi s.l. or B. afzelii prevalence and elevation, and the relationship between B. burgdorferi s.l. or B. afzelii prevalence and questing I. ricinus abundance within the core range of I. ricinus. The effect of questing I. ricinus abundance was investigated by including questing I. ricinus abundance and questing I. ricinus abundance 2 as fixed effects. Site ID was again included as a random effect.
To ensure independence between main and quadratic terms, all continuous variables were standardized (by subtracting the population mean from each sample and dividing it by the standard deviation). The significance of explanatory variables was assessed using likelihood ratio tests. The significance of the linear terms was assessed after removing the non-significant quadratic terms. Statistical analyses were performed using the package glmmADMB [60,61] in R 3.2.3 [62].